Python医学影像分析完整入门教程：DICOM读取、预处理、病灶分割实战

文章目录摘要一、前期准备环境搭建与核心库说明1.1 依赖库安装1.2 公开测试数据集获取二、基础实操1使用pydicom读取单张DICOM影像三、基础实操2SimpleITK读取DICOM序列整套CT扫描四、医学影像核心预处理HU值转换与窗宽窗位调节五、图像基础处理滤波降噪、阈值分割六、进阶实战MONAI搭建基础医学图像分割流程摘要医学影像分析是AI医疗落地最成熟的方向之一CT、MRI、X光、病理切片等影像数据大多以DICOM标准格式存储。Python凭借pydicom、SimpleITK、MONAI等开源工具链成为科研人员、算法工程师开展医学影像处理的首选编程语言。本篇教程从环境搭建、DICOM基础读取、CT窗宽窗位调整、影像预处理延伸至深度学习分割基础案例所有代码均经过本地运行验证适合零基础入门也可作为项目开发参考。一、前期准备环境搭建与核心库说明1.1 依赖库安装建议新建独立Python虚拟环境避免包版本冲突推荐Python3.9~3.11版本。pipinstallpydicom simpleitk numpy matplotlib nibabel monai torch torchvision各库核心用途pydicom轻量DICOM读写工具适合读取单张dcm文件、调取患者元信息SimpleITK工业级医学影像处理库支持DICOM序列、NIfTI格式擅长3D体积处理、图像配准、滤波numpy/matplotlib数组运算与影像可视化nibabel专门读取nii.gz(NIfTI)标准医学影像MONAI基于PyTorch的医疗影像深度学习框架替代原生PyTorch完成数据增强、3D模型训练。新手常见踩坑Windows环境安装SimpleITK如果报错可以升级pip后重试读取医院原始DICOM文件遇到元数据缺失读取时增加force参数。1.2 公开测试数据集获取本地练习不建议直接使用医院真实患者影像涉及隐私合规推荐免费公开数据集MedMNIST轻量化标准数据集包含胸片、眼底、病理图像一键加载无需处理DICOMTCIA包含肺部CT、脑部MRI等完整DICOM序列BraTS脑肿瘤MRI多模态NIfTI数据集常用于3D分割练习。二、基础实操1使用pydicom读取单张DICOM影像DICOM文件不止保存像素灰度还存储患者姓名、扫描层厚、像素间距、设备参数等几十项临床元数据。很多初学者直接读取像素数组忽略元信息后续三维重建、病灶测距都会出现误差。importpydicomimportmatplotlib.pyplotasplt# 读取单张dcm文件dcm_path./data/test.dcmdspydicom.dcmread(dcm_path)# 打印基础元信息print(f患者ID{ds.PatientID})print(f扫描层厚{ds.SliceThickness})print(f像素间距{ds.PixelSpacing})print(f图像尺寸{ds.Rows}×{ds.Columns})# 获取像素数组img_arrayds.pixel_array# 显示影像plt.figure(figsize(6,6))plt.imshow(img_array,cmapgray)plt.axis(off)plt.title(单张DICOM影像)plt.show()常见报错解决部分导出的DICOM缺少文件头直接读取会抛出InvalidDicomError修改读取代码强制读取dspydicom.dcmread(dcm_path,forceTrue)三、基础实操2SimpleITK读取DICOM序列整套CT扫描一次CT扫描通常数十至数百张dcm切片存放在同一个文件夹称为DICOM序列。pydicom不擅长批量自动排序切片优先选择SimpleITK处理三维CT/MRI体积数据。importSimpleITKassitkimportmatplotlib.pyplotasplt dicom_dir./data/ct_series/# 自动读取文件夹内完整DICOM序列readersitk.ImageSeriesReader()series_idsreader.GetGDCMSeriesIDs(dicom_dir)file_namesreader.GetGDCMSeriesFileNames(dicom_dir,series_ids[0])reader.SetFileNames(file_names)# 生成3D影像volumevolume_3dreader.Execute()print(f三维影像尺寸(W,H,Z):{volume_3d.GetSize()})print(f空间间距{volume_3d.GetSpacing()})# 提取中间切片可视化slice_idxvolume_3d.GetSize()[2]//2slice_imgvolume_3d[:,:,slice_idx]slice_arrsitk.GetArrayFromImage(slice_img)plt.figure()plt.imshow(slice_arr,cmapgray)plt.axis(off)plt.title(fCT中间切片 第{slice_idx}层)plt.show()四、医学影像核心预处理HU值转换与窗宽窗位调节CT影像原始像素值并非直接可用存储的是原始探测器数值需要转换成HU亨氏单位。临床医生观察肺部、骨骼、软组织会使用不同窗宽窗位这是医学影像区别普通图像最关键的知识点。肺窗窗宽1500窗位-600观察肺结节纵隔窗窗宽350窗位50观察淋巴结、软组织骨窗窗宽2000窗位500用于骨骼观察。importSimpleITKassitkimportnumpyasnpimportmatplotlib.pyplotaspltdefapply_window(image_array,window_width,window_level):通用窗宽窗位映射函数min_valwindow_level-window_width/2max_valwindow_levelwindow_width/2image_arraynp.clip(image_array,min_val,max_val)image_array(image_array-min_val)/(max_val-min_val)*255returnimage_array.astype(np.uint8)# 读取volume转换numpy数组dicom_dir./data/ct_series/readersitk.ImageSeriesReader()series_idsreader.GetGDCMSeriesIDs(dicom_dir)file_namesreader.GetGDCMSeriesFileNames(dicom_dir,series_ids[0])reader.SetFileNames(file_names)volumereader.Execute()vol_arrsitk.GetArrayFromImage(volume)slice_datavol_arr[vol_arr.shape[0]//2]# 设置肺窗与纵隔窗lung_windowapply_window(slice_data,window_width1500,window_level-600)mediastinum_windowapply_window(slice_data,window_width350,window_level50)fig,(ax1,ax2)plt.subplots(1,2,figsize(12,6))ax1.imshow(lung_window,cmapgray)ax1.set_title(肺窗)ax1.axis(off)ax2.imshow(mediastinum_window,cmapgray)ax2.set_title(纵隔窗)ax2.axis(off)plt.show()开发提醒很多开源数据集预处理缺失窗宽标准化直接送入模型训练模型很难学到病灶特征这是新手训练分割模型效果差的高频原因。五、图像基础处理滤波降噪、阈值分割在不使用深度学习的场景可采用传统图像处理方法完成简单组织分割适合快速提取胸腔、骨骼区域。importSimpleITKassitk# 高斯滤波降噪smooth_imgsitk.SmoothingRecursiveGaussian(volume,sigma0.8)# Otsu自动阈值分割seg_filtersitk.OtsuThresholdImageFilter()seg_filter.SetInsideValue(0)seg_filter.SetOutsideValue(1)seg_resultseg_filter.Execute(smooth_img)# 保存分割结果为nii.gz文件可使用3D Slicer打开查看sitk.WriteImage(seg_result,./output/seg_result.nii.gz)print(分割文件保存完成)生成的nii.gz文件推荐使用3D Slicer免费软件打开可视化三维分割效果是医学影像工程师必备工具。六、进阶实战MONAI搭建基础医学图像分割流程传统算法难以识别微小病灶行业主流方案使用U-Net系列网络完成病灶分割。MONAI封装了医学影像专用数据加载、空间增强、3D网络省去大量重复开发工作。下面展示标准数据预处理管道模板。frommonai.transformsimport(Compose,LoadImaged,AddChanneld,ScaleIntensityRanged,RandFlipd,Resized,ToTensord)frommonai.dataimportDataset,DataLoaderimportglobimportos# 数据集路径配置data_root./dataset/train_imagessorted(glob.glob(os.path.join(data_root,imagesTr,*.nii.gz)))train_labelssorted(glob.glob(os.path.join(data_root,labelsTr,*.nii.gz)))train_files[{image:img,label:lab}forimg,labinzip(train_images,train_labels)]# 构建训练数据增强流水线train_transformsCompose([LoadImaged(keys[image,label]),AddChanneld(keys[image,label]),ScaleIntensityRanged(keys[image],a_min-1000,a_max400,b_min0.0,b_max1.0,clipTrue),Resized(keys[image,label],spatial_size[128,128,64]),RandFlipd(keys[image,label],spatial_axis0,prob0.5),ToTensord(keys[image,label])])train_dsDataset(datatrain_files,transformtrain_transforms)train_loaderDataLoader(train_ds,batch_size1,shuffleTrue)print(数据集加载完成)基于该数据加载器可直接接入MONAI内置UNet、SegResNet等网络开展训练。