医学影像预处理之裁减

背景

在数据集训练之前，为了减小数据计算的工作量，提高训练的速度，通常会针对感兴趣的特征部分进行原图片的裁减，这样每张图片既保留了待提取的特征集，又缩小了整体的尺寸，可有效缩短模型训练耗费的时间。

理解三维图像

DICOM坐标系是相对于病人的方向来确定的，如下图所示：

三视图三个截面，分别称为轴状位（Transverse/Axisplane）、冠状位（Coronal/Frontal plane）和矢状位（Sagittal plane），如下图所示：

裁减原理

我们把一个影像文件读到一个三维数组，shape为（x，y，z），分别沿着x、y和z三个坐标轴方向滑动读取每一个slice中的图像数据，从而获取期望的特征边界。沿X轴方向滑动得到的每一个切面为图二所示的失状位，同理，沿y轴滑动得到的每一个切面为图二所示的冠状位，沿z轴滑动得到的每一个切面为图二所示的轴状位。以沿x轴滑动为例，对应在读图软件上，就是在矢状位视图滚动鼠标的效果，这个视图中滚动的每一张切片就是沿X轴滑动得到的，此外，可以看到在另外两个视图中有一个竖线在沿X轴左右滑动，如下图，其它视图同理。

为了方便对照和理解，我把三视图与轴状位、冠状位和矢状位一一对应起来，如下图：

理解了这些之后也就理解了裁减的原理和方法，还是以沿X轴滑动为例，先沿X轴正向滑动，遍历每个slice，获取特征的边界对应的Xmin值，然后再沿X轴反向滑动，遍历每个slice，获取特征的边界对应的Xmax值，这样我们就获取了Xmin:Xmax之间的部分，同理我们可以获取Ymin:Ymax和Zmin:Zmax。

代码实现

# encoding = utf-8
import nibabel as nib
import numpy as np
import os

DST_PATH = "filled with your data path"

def get_boundary(volume):

    for i in range(volume.shape[0]):
        if np.max(volume[i]) != 0:
            volume_Xmin = i
            break
    for i in range(volume.shape[0])[::-1]:
        if np.max(volume[i]) != 0:
            volume_Xmax = i
            break
    for i in range(volume.shape[1]):
        if np.max(volume[:, i, :]) != 0:
            volume_Ymin = i
            break
    for i in range(volume.shape[1])[::-1]:
        if np.max(volume[:, i, :]) != 0:
            volume_Ymax = i
            break
    for i in range(volume.shape[2]):
        if np.max(volume[:, :, i]) != 0:
            volume_Zmin = i
            break
    for i in range(volume.shape[2])[::-1]:
        if np.max(volume[:, :, i]) != 0:
            volume_Zmax = i
            break
    return (volume_Xmin, volume_Xmax,volume_Ymin,volume_Ymax,volume_Zmin,volume_Zmax)

if __name__ == "__main__":
    if not os.path.exists(DST_PATH):
      os.mkdir(DST_PATH)
    img_data = nib.load("./Zhangsan_Arterial_Data.nii.gz")
    img_label = nib.load("./Zhangsan_Arterial_Label.nii.gz")
    img_data_npy = img_data.get_fdata()
    img_label_npy = img_label.get_fdata()
    print("=======origin shape=======")
    print(img_data_npy.shape)
    volume_Xmin, volume_Xmax,volume_Ymin,volume_Ymax,volume_Zmin,volume_Zmax = get_boundary(img_label_npy)
    data_clip = img_data_npy[volume_Xmin:volume_Xmax, volume_Ymin:volume_Ymax, volume_Zmin:volume_Zmax]
    label_clip = img_label_npy[volume_Xmin:volume_Xmax, volume_Ymin:volume_Ymax, volume_Zmin:volume_Zmax]
    print("=======clip shape=========")
    print(data_clip.shape)
    nib.save(nib.Nifti1Image(data_clip.astype('int16'), affine=img_data.affine), os.path.join(DST_PATH, 'Clip_Data.nii.gz'))
    nib.save(nib.Nifti1Image(label_clip.astype('uint8'), affine=img_label.affine), os.path.join(DST_PATH, 'Clip_Label.nii.gz'))

裁减效果

裁减前：shape为(512, 512, 231)

裁减后：shape为(308, 98, 112)