影像组学nii瘤周扩展Python代码
影像组学分析过程中，有时候会进行瘤周的分析，这时候就要用到ROI的扩展。
本文记录我的瘤周扩展的思路，供大家参考。
因为代码是我自己摸索写出，如果其中有什么疑问或者觉得代码有什么问题的，请联系xh542428798@hotmail.com，以便修正。
思路
需要精确扩展指定毫米数的ROI，首先需要进行重采样。把dcm或nii文件的spacing修改到1* 1 *1。
再使用Scipy.ndimage的图像矩阵处理函数，进行图像扩展。
具体操作如下：
1、图像Resample
# 定义原始图像和ROI文件
img_path = r'original.nii'
roi_path = r'mask.nii'
tumorImage = sitk.ReadImage(img_path)
# 设置重采样，保持层间距不变，进行重采样以减小工作量
newSpacing=[1.0,1.0,tumorImage.GetSpacing()[2]]
# 进行resample
Niimage_resample = resample_image_itk(img_path,preference_path=img_path,newSpacing = newSpacing, resamplemethod=sitk.sitkLinear) # resample #这里要注意：mask用最近邻插值sitk.sitkNearestNeighbor，CT图像用线性插值
# Resample之后的ROI只能和resample之后的原图一起提取组学特征，spacing不匹配的话就会报错，所以原图也要Resample
image_array = sitk.GetArrayFromImage(Niimage_resample)
mask_img_resample = resample_image_itk(roi_path,preference_path=img_path,newSpacing = newSpacing, resamplemethod=sitk.sitkNearestNeighbor)
2、进行扩展
mask_img_arr = sitk.GetArrayFromImage(mask_img_resample)
iteration = 3 # newSpacing是1mm,扩展3mm
mask_img_arr_expand = expand_dilation(mask_img_arr,iterations=iteration)
3、对扩展前后的肿瘤ROI取差集
因为扩展只是外扩，也包含内部区域，而我们需要边缘
mask_img_arr_border = diffROI(mask_img_arr,mask_img_arr_expand)
4、保存处理后的ROI和原图的nii
outNII_path = '{}_resample.nii.gz'.format('3mm_origin')
array2nii(image_array,outNII_path,Niimage_resample)
outROI_path = 'Ex{}.nii.gz'.format('_3mm')
array2nii(mask_img_arr_border,outROI_path,Niimage_resample)
相关函数定义区
# imagePath是图像nrrd,preference_path是参考nrrd的路径
def resample_image_itk(imagePath,preferenceNii, newSpacing=[3.0,3.0,3.0], resamplemethod=sitk.sitkLinear): #这里要注意：mask用最近邻插值，CT图像用线性插值，sitk.sitkNearestNeighbor
image = sitk.ReadImage(preferenceNii)
resample = sitk.ResampleImageFilter()
resample.SetInterpolator(resamplemethod) #这里要注意：mask用最近邻插值，CT图像用线性插值 ##
resample.SetOutputDirection(image.GetDirection()) ##
resample.SetOutputOrigin(image.GetOrigin()) ##
newSpacing = np.array(newSpacing, float)
newSize = image.GetSize() / newSpacing * image.GetSpacing()
newSize = np.around(newSize,decimals=0)
newSize = newSize.astype(np.int)
resample.SetSize(newSize.tolist()) ##
resample.SetOutputSpacing(newSpacing) ##
image_resample = sitk.ReadImage(imagePath)
newimage = resample.Execute(image_resample)
return newimage #返回sitk类型的数据
def expand_dilation(mask_img_arr,iterations = 5):
shape_nrrd = mask_img_arr.shape
mask_img_arr_expand = np.zeros(shape_nrrd)
for index in range(shape_nrrd[0]):
mask_img_arr_expand[index,:,:] = scipy.ndimage.binary_dilation(mask_img_arr[index,:,:],iterations = iterations).astype('uint16')
return mask_img_arr_expand
def diffROI(mask_img_arr,mask_img_arr_expand):
# 通过传入原始ROI mask_img_arr和扩展了5mm的ROI mask_img_arr_expand
# 返回扩展的5mm的ROI区域
# 首先把二值ROI都转为0和1
mask_img_arr[mask_img_arr != 0] = 1
mask_img_arr_expand[mask_img_arr_expand != 0] = 1
shape_nrrd = mask_img_arr.shape
mask_img_arr_border = np.zeros(shape_nrrd)
for index in range(shape_nrrd[0]):
for x in range(shape_nrrd[1]):
for y in range(shape_nrrd[2]):
if mask_img_arr[index,x,y] != mask_img_arr_expand[index,x,y]:
mask_img_arr_border[index,x,y] = 1
return mask_img_arr_border
def array2nii(image_array, out_path, NIIimage_resample):
## image_array是矩阵，out_path是带文件名的路径，NIIimage_resample是sitk_obj
# 1.构建nrrd阅读器
image2 = NIIimage_resample
# 2.将整合后的数据转为array，并获取dicom文件基本信息
origin = image2.GetOrigin() # x, y, z
spacing = image2.GetSpacing() # x, y, z
direction = image2.GetDirection() # x, y, z
# 3.将array转为img，并保存为.nii.gz
image3 = sitk.GetImageFromArray(image_array)
image3.SetSpacing(spacing)
image3.SetDirection(direction)
image3.SetOrigin(origin)
sitk.WriteImage(image3, out_path)
PS：如果有疑问和错误，可以联系xh542428798@hotmail.com，但是苦逼上班族996，如果需要帮忙解决小白报错，求各位老板付费咨询，万分感谢！
