【医学图像分割】CT医学图像的预处理（重采样）

孟想成真

CT图像中存在两个基本概念，窗口（window width）和窗位（window center），用于选取感兴趣的CT值范围，因为人体各组织结构不同，对X线吸收程度各异，形成不同的CT值，因此可以利用CT值来鉴别组织的性质。CT值的单位为Hounsfield
   Hounsfield单位(HU)是计算机断层扫描(CT)中普遍使用的无量纲单位，用于标准、便捷的表达CT数值。Hounsfield单位是通过对测量得到的衰减系数进行线性变换得到的。这种转换是基于空气和的密度，其中纯水被定义为0 Hounsfield单元，空气被定义为-1000 Hounsfield单元。组织密度越大，x射线吸收越强，其值为正，呈亮信号;密度较低的组织，x射线吸收较少，显示负值，呈暗信号。Hounsfield单位是以1979年诺贝尔生理学或医学奖得主Godfrey Hounsfield爵士的名字命名的，以表彰他在CT的发明中所做的贡献。
  窗宽是CT图像上显示的CT值范围，在此CT值范围内的组织和病变均以不同的模拟灰度显示。而CT值高于此范围的组织和病变，无论高出程度有多少，均以白影显示，不再有灰度差异; 反之，低于此范围的组织结构，不论低的程度有多少，均以黑影显示，也无灰度差别。增大窗宽，则图像所示CT值范围加大，显示具有不同密度的组织结构增多，但各结构之间的灰度差别减少。减小窗宽，则显示的组织结构减少，然而各结构之间的灰度差别增加。如观察脑质的窗宽常为-15~+85H，即密度在-15 ~+85H范围内的各种结构如脑质和脑脊液间隙均以不同灰度显示。而高于+85H的组织结构如骨质几颅内钙化，其间虽有密度差，但均以白影显示，无灰度差别;而低于-15H组织结构如皮下脂肪及乳突内气体均以黑影显示，其间也无灰度差别。
  窗位是窗的中心位置，同样的窗宽，由于窗位不同，其所包括CT值范围的CT值也有差异。例如窗宽同为100H，当窗位为0H时，其CT值范围为-50 ~ +50H ; 如窗位为+35H时，则CT值范围为-15~+85H。通常，欲观察某以组织结构及发生的病变，应以该组织的CT值为窗位。例如脑质CT值约为+35H，则观察脑组织及其病变时，选择窗位以+35H为妥。
  常见的人体组织CT值

(1) 液体CT值：-10~10Hu，水的CT值为0Hu;
(2) 空气CT值：-1000Hu;
(3) 脂肪CT值：-10~90Hu;
(4) 肝脾肾、脑实质等软组织CT值：20~50Hu;
(5) 骨性组织CT值一般超过300Hu；骨皮质一般CT值超过+1000Hu;
(6) 腹部常用CT值：-160~240Hu
下边是对腹部窗口的图像预处理：

1. # 腹部窗口的选取
2. def pretreatmentImage(image):
3.     image[image < -160] = 0
4.     image[image > 240] = 0
6.     return image

复制代码

医疗图像重采样：

对于一张大小为128*128的彩色图像，其在计算机中可以表示为128*128*3的矩阵，其中每一个像素点的取值范围为0-255，不同的数值代表不同的亮度。但是对于医疗图像其是由若干个slice组成的，假设每一个slice的大小为512*512的单通道的图像，其中每一个像素点表示的是一个体素的取值，其范围可以-1000~2000之间。接下来通过以胰腺分割数据集中PANCREAS_0015.nii.gz为例，对医疗图像中体素这个概念进行讲解。Spacing(0.78125, 0.78125, 1.0)表示的是原始图像体素的大小，也可以将Spacing想象成大小为(0.78125, 0.78125, 1.0)的长方体。而原始图像的Size为 (512, 512, 247)，表示的是原始在X轴，Y轴，Z轴中体素的个数。原始图像的size*对应的Spacing既可以得到真实3D图像大小（512*0.78125，512*0.78125，247*1 ），在图像重采样只是修改体素的大小，而真实3D图像大小是保持不变的，因此假设我们将Spacing修改成(1.0, 1.0, 2.0)的时候，则修改之后其对应的size应该为（（512*0.78125）/ 1.0，（512*0.78125）/ 1.0，（247*1 ）/ 2.0）即（400, 400, 124）。

1. import pydicom
2. path = ""
3. dcm_data = pydicom.read_file(path)
4. print(dcm_data.PixelSpacing, dcm_data.SliceThickness) # 查看医学图像的体素大小和切片间距

复制代码

重采样代码：

1. def img_resmaple(path, new_spacing=[1.0, 1.0, 1.0]):
2.     data = STK.ReadImage(path)
3.     image = STK.GetArrayFromImage(data) # 获取图像
4.     # print(image.shape)
5.     original_spacing = data.GetSpacing()
6.     # print(original_spacing)
7.     original_size = data.GetSize()
8.     # print(data.GetOrigin(), data.GetDirection())
9.     new_shape = [
10.         int(np.round(original_spacing[0] * original_size[0] / new_spacing[0])),
11.         int(np.round(original_spacing[1] * original_size[1] / new_spacing[1])),
12.         int(np.round(original_spacing[2] * original_size[2] / new_spacing[2])),
13.     ]
14.     resmaple = STK.ResampleImageFilter()
15.     resmaple.SetInterpolator(STK.sitkLinear)
16.     resmaple.SetDefaultPixelValue(0)
17.     resmaple.SetOutputSpacing(new_spacing)
18.     resmaple.SetOutputOrigin(data.GetOrigin())
19.     resmaple.SetOutputDirection(data.GetDirection())
20.     resmaple.SetSize(new_shape)
21.     data = resmaple.Execute(data)
22.     image = STK.GetArrayFromImage(data)
23.     # print(image.shape)
24.     return image

复制代码

按照图像的中心裁剪：

1. def centerCrop(image, label, output_size):
2.     if image.shape[0] <= output_size[0] or image.shape[1] <= output_size[1] or image.shape[2] <= output_size[2]:
3.         pw = max((output_size[0] - image.shape[0]) // 2 + 3, 0)
4.         ph = max((output_size[1] - image.shape[1]) // 2 + 3, 0)
5.         pd = max((output_size[2] - image.shape[2]) // 2 + 3, 0)
6.         image = np.pad(image, [(pw, pw), (ph, ph), (pd, pd)], mode='constant', constant_values=0)
7.         label = np.pad(label, [(pw, pw), (ph, ph), (pd, pd)], mode='constant', constant_values=0)
9.     (w, h, d) = image.shape
11.     w1 = int(round((w - output_size[0]) / 2.))
12.     h1 = int(round((h - output_size[1]) / 2.))
13.     d1 = int(round((d - output_size[2]) / 2.))
15.     # print(image.shape, output_size, get_center(label), w1, h1, d1)
16.     image = image[w1:w1 + output_size[0], h1:h1 + output_size[1], d1:d1 + output_size[2]]
17.     label = label[w1:w1 + output_size[0], h1:h1 + output_size[1], d1:d1 + output_size[2]]
19.     return image, label

复制代码

两种常用的归一化方法

1）min-max标准化
2）Z-score标准化方法
（1）min-max标准化（Min-Max Normalization）（线性函数归一化）
定义：也称为离差标准化，是对原始数据的线性变换，使得结果映射到0-1之间。
本质：把数变为【0,1】之间的小数。
转换函数：（X-Min）/(Max-Min)
如果想要将数据映射到-1,1，则将公式换成：（X-Mean）/(Max-Min)
其中：max为样本数据的最大值，min为样本数据的最小值，Mean表示数据的均值。
缺陷：当有新数据加入时，可导致max和min的变化，需要重新定义。
（2）0均值标准化（Z-score standardization）
定义：这种方法给与原始数据的均值（mean）和标准差（standard deviation）进行数据的标准化。经过处理的数据符合标准正态分布，即均值为0，标准差为1.
本质：把有量纲表达式变成无量纲表达式。
转换函数：（X-Mean）/(Standard deviation)
其中，Mean为所有样本数据的均值。Standard deviation为所有样本数据的标准差。
两种归一化方法的使用场景：
（1）在分类、聚类算法中，需要使用距离来度量相似性的时候、或者使用PCA技术进行降维的时候，第二种方法(Z-score standardization)表现更好。
因为：第一种方法(线性变换后)，其协方差产生了倍数值的缩放，因此这种方式无法消除量纲对方差、协方差的影响，对PCA分析影响巨大；同时，由于量纲的存在，使用不同的量纲、距离的计算结果会不同。
（2）在不涉及距离度量、协方差计算、数据不符合正太分布的时候，可以使用第一种方法或其他归一化方法。比如图像处理中，将RGB图像转换为灰度图像后将其值限定在（0 ，255）的范围。
因为：第二种归一化方式中，新的数据由于对方差进行了归一化，这时候每个维度的量纲其实已经等价了，每个维度都服从均值为0、方差1的正态分布，在计算距离的时候，每个维度都是去量纲化的，避免了不同量纲的选取对距离计算产生的巨大影响。
代码实现：

1. def MinMax(image):
2.     image = image - np.mean(image)  # 零均值
3.     if np.max(image) - np.min(image) != 0:
4.         image = (image - np.min(image)) / (np.max(image) - np.min(image))  # 归一化
6.     return image
8. def zsorce(image):
9.     image = (image - np.mean(image)) / np.std(image)
10.     return image

复制代码

参考：
CT图像之Hu值变换与窗宽窗位调整_HuiYu-Li的博客-CSDN博客_ct图像调窗
CT医学影像的窗高窗位、CT值（Hu值）_土豆娃wkl的博客-CSDN博客_hu影像学单位
医学图像预处理之重采样_YYLin-AI的博客-CSDN博客_医学图像重采样
机器学习-数据归一化方法（Normalization Method）_Microstrong0305的博客-CSDN博客_数据归一化


来源：https://blog.csdn.net/qq_37804582/article/details/124098153
免责声明：如果侵犯了您的权益，请联系站长，我们会及时删除侵权内容，谢谢合作！