如何具体实现registration

0. 心得

一开始其实根本没理解配准是具体需要做什么，只是一味的vibe coding，什么重采样，刚性配准什么的，最后生成出来的图像全是噪声，第二周感觉把这些放一起就没问题了，结果惊奇的发现人体直接错位了，锁骨来到了人体的中间地方，经过重新认真严谨思考，才意识到配准的真正含义，接下来是对具体配准过程的学习笔记。

完整代码见 GitHub：registration-scripts

1. 任务目标

当前任务是多期相 CT 图像配准。

对于同一个病例，可能存在多个期相或序列，例如：

  • 平扫期

  • 动脉期

  • 门静脉期

  • 延迟期

  • std5mm

  • lung5mm

这些序列的原始 DICOM 可能存在：

  • 切片数量不同

  • 扫描起始位置不同

  • spacing 不同

  • origin 不同

  • 扫描范围不同

  • 呼吸状态导致器官位置轻微变化

因此，不能直接认为不同期相中相同编号的 DICOM 切片对应同一个身体位置。

当前目标是：

以一个期相作为 fixed image，将其他期相作为 moving image，通过配准和重采样，使不同期相在统一空间中具有一一对应的切片关系。

最终输出为 nnU-Net 多通道输入格式。

2. 大体流程

当前采用的整体流程为：

  1. 读取同一病例下的 DICOM 文件；

  2. 根据 SeriesInstanceUID 区分不同序列；

  3. 排除 localizer、scout、切片数过少、不同部位的序列；

  4. 选择一个序列作为 fixed image；

  5. 其他序列作为 moving image；

  6. 对 moving image 进行配准；

  7. 将 moving image 重采样到 fixed image 的空间；

  8. 根据共同有效扫描范围进行统一裁剪；

  9. 检查同编号切片是否对应同一部位；

  10. 保存为 nnU-Net 多通道 NIfTI 格式。

简化表示为：

fixed image

↓

提供统一空间标准

moving image

↓

配准到 fixed image

↓

重采样到 fixed image 空间

↓

裁剪共同覆盖区域

↓

输出 nnU-Net 多通道图像

3. 多期相 CT 配准脚本实现与流程说明

3.1 参数设置部分

脚本开头需要设置几个关键参数。

fixed_dir 表示基准期相的 DICOM 文件夹路径，例如门静脉期、std5mm，或者老师指定的标准期相。

moving_dir 表示待配准期相的 DICOM 文件夹路径，例如动脉期、延迟期、lung5mm 等。

output_dir 表示输出文件夹。

case_id 表示 nnU-Net 中的病例编号，例如 case_001。

最终会输出：

  • case_001_0000.nii.gz：fixed 通道

  • case_001_0001.nii.gz：moving 配准后的通道

此外，还有三个重要参数：

  • default_pixel_value = -1024.0

  • body_threshold = -500

  • crop_margin = 10

其中，default_pixel_value 是 moving 重采样到 fixed 空间后，moving 原本没有覆盖到的区域所填充的默认值。对于 CT 图像，空气区域通常接近 -1000 HU，因此使用 -1024 比使用 0 更合理。

body_threshold 用于粗略区分身体区域和空气背景。CT 中空气大约是 -1000 HU，而人体软组织、骨骼等通常大于 -500 HU，因此可以用 image > -500 粗略提取人体区域。

crop_margin 表示裁剪时额外保留的边缘范围，避免裁剪过紧导致身体边缘被切掉。

3.2 读取 DICOM 序列

脚本通过 read_dicom_series() 函数读取一个 DICOM 序列文件夹。

需要注意的是，一个 DICOM 序列通常不是一个文件，而是一组 .dcm 切片。SimpleITK 会把这一组 DICOM 文件读取成一个三维图像。

读取过程大致为：

  1. 使用 ImageSeriesReader 创建读取器；

  2. 在文件夹中查找 DICOM series ID；

  3. 根据 series ID 获取该序列下的所有 DICOM 文件；

  4. 将多张 DICOM 切片读取为一个 3D 图像。

如果一个文件夹中存在多个 series，脚本会给出提示，并默认读取第一个 series。实际处理时，最好保证一个文件夹中只放一个 series，或者在后续批量化脚本中手动选择正确的 series。

3.3 打印图像空间信息

脚本中使用 print_image_info() 函数打印图像的空间信息，包括：

  • Size

  • Spacing

  • Origin

  • Direction

  • PixelType

这些信息分别表示：

  • Size：图像尺寸，例如 512 × 512 × 140

  • Spacing：体素间距，例如 0.97 × 0.97 × 5.0

  • Origin：图像原点坐标

  • Direction：图像方向矩阵

  • PixelType：像素类型，例如 32-bit float

打印这些信息的作用是：

  1. 检查 fixed 和 moving 原始空间信息是否一致；

  2. 检查 moving 配准重采样后是否已经进入 fixed 空间；

  3. 检查裁剪后的 fixed 和 moving 是否仍然保持空间一致；

  4. 判断结果是否符合 nnU-Net 多通道输入要求。

3.4 刚性配准：求 moving 到 fixed 的空间变换

脚本中的 rigid_registration() 是核心函数之一。

刚性配准的目标是求出 moving image 到 fixed image 的空间变换，即 moving 应该如何平移和旋转，才能和 fixed 对齐。

刚性配准允许：

  • x、y、z 方向平移

  • 绕 x、y、z 轴旋转

不允许：

  • 缩放

  • 拉伸

  • 局部形变

这对于医学分割任务比较安全，因为它不会改变器官和病灶的真实形态。

在配准前，脚本会先将 fixed 和 moving 转换为 Float32 类型。因为 DICOM CT 原图常常是 int16，而 SimpleITK 的配准算法通常需要浮点图像，否则可能报像素类型不支持的错误。

接着，脚本使用 Euler3DTransform 初始化一个 3D 刚性变换，并用 CenteredTransformInitializer 将 fixed 和 moving 的几何中心先大致对齐。这样可以为后续优化提供一个较好的初始位置。

3.5 相似性指标：互信息

配准时，算法需要判断 fixed 和 moving 当前是否对齐，这就需要相似性指标。

脚本中使用的是 Mattes Mutual Information，也就是互信息。

使用互信息的原因是：当前任务是多期相 CT 配准。不同期相虽然都是 CT，但由于增强时间不同，灰度分布可能不同。例如：

  • 动脉期血管更亮

  • 门静脉期肝实质强化不同

  • 延迟期病灶表现不同

如果直接使用 MSE 这类灰度差异指标，可能会因为灰度变化而误判配准效果。

互信息不要求两个图像灰度值完全一致，而是衡量两幅图像之间的统计相关性，因此更适合多期相 CT 配准。

3.6 优化器与多分辨率策略

脚本使用梯度下降优化器来寻找最优的刚性变换参数。

优化器的作用是不断调整 transform 的参数，使相似性指标逐渐变好。简单来说，算法会不断尝试不同的平移和旋转参数，并判断哪一种变换能让 moving 更接近 fixed。

脚本中的关键参数包括：

  • learningRate：每次更新参数的步长

  • numberOfIterations：最大迭代次数

  • convergenceMinimumValue：收敛判断阈值

  • convergenceWindowSize：用于判断是否收敛的窗口大小

此外，脚本使用了多分辨率策略：

  • 先在低分辨率图像上粗配准

  • 再在中等分辨率图像上细化

  • 最后在原始分辨率图像上微调

对应参数为：

  • shrinkFactors = [4, 2, 1]

  • smoothingSigmas = [2, 1, 0]

这样做的好处是可以减少噪声干扰，提高配准稳定性，也能让优化过程更容易收敛。

3.7 重采样到 fixed 空间

刚性配准得到的结果是一个 transform。这个 transform 只是描述 moving 应该如何变换，还没有真正生成最终图像。

因此，脚本下一步使用 resample_moving_to_fixed() 将 moving image 真正重采样到 fixed image 的空间中。

重采样后，moving_registered 会继承 fixed image 的空间网格，包括：

  • size

  • spacing

  • origin

  • direction

也就是说，从空间坐标上看，moving_registered 已经进入 fixed 的空间。

这一步完成后，fixed 和 moving_registered 的同一个 slice index 才可能对应同一个物理位置。

需要注意的是，重采样过程中，如果 moving 原本没有覆盖 fixed 的某些区域，这些区域会使用默认值填充。对于 CT 图像，脚本使用 -1024 作为默认填充值，以接近空气背景。

3.8 计算共同有效区域

配准和重采样后，moving_registered 会被放到 fixed 的完整空间中。但如果两个期相的原始扫描范围不同，就可能出现某些区域 fixed 有内容，而 moving 没有内容。

例如：

  • fixed 从锁骨扫到腹部

  • moving 从肺部扫到腹部

此时，moving_registered 在 fixed 前面锁骨区域可能是空白或默认填充值。

因此，脚本使用 get_common_body_bbox() 计算 fixed 和 moving_registered 的共同有效区域。

具体做法是：

  1. 对 fixed 使用阈值 fixed > -500 得到 fixed 的身体区域；

  2. 对 moving_registered 使用阈值 moving_registered > -500 得到 moving 的身体区域；

  3. 对两个 mask 取交集，得到两个期相都具有身体内容的共同区域；

  4. 使用连通域分析保留最大连通区域；

  5. 计算该区域的 bounding box；

  6. 在 bounding box 外额外加 margin，避免裁剪过紧。

这样得到的 bbox 就是两个期相共同覆盖的有效扫描区域。

3.9 统一裁剪与 nnU-Net 输出

脚本使用 crop_image() 对 fixed 和 moving_registered 进行裁剪。

关键点是：

fixed 和 moving_registered 必须使用同一个 start 和 crop_size。

也就是说：

  • fixed_cropped 使用共同 bbox 裁剪；

  • moving_cropped 使用同一个共同 bbox 裁剪。

不能让 fixed 和 moving 分别自己计算 bbox 后单独裁剪。否则两个图像的起点和范围可能不同，导致多通道之间空间错位。

裁剪完成后，脚本将结果保存为 nnU-Net 双通道格式：

  • case_001_0000.nii.gz

  • case_001_0001.nii.gz

其中：

  • 0000 是 fixed 图像；

  • 0001 是 moving 配准并裁剪后的图像。

最后脚本会检查两个输出文件的空间信息是否一致，包括：

  • size

  • spacing

  • origin

  • direction

如果这些全部一致，说明从空间网格上看，它们已经可以作为 nnU-Net 的多通道输入。

但最终仍然需要用 ITK-SNAP 或 3D Slicer 进行 overlay 检查，确认图像内容上也真正对齐。

3.10 脚本解决的问题总结

这个脚本主要解决了多期相 CT 配准中的几个实际问题。

第一，解决了两个期相切片数不同的问题。

moving 可以被重采样到 fixed 的空间中，从而获得与 fixed 一致的空间网格。

第二，解决了同编号切片不对应的问题。

通过配准和重采样，使 fixed 与 moving_registered 的同一个 slice index 对应相同物理位置。

第三，解决了扫描范围不同的问题。

通过共同有效区域裁剪，去掉某一期相原本没有扫描到的空白区域。

第四，解决了 nnU-Net 多通道输入的空间一致性问题。

最终输出的 case_001_0000.nii.gz 和 case_001_0001.nii.gz 具有相同的 size、spacing、origin 和 direction。

4. 当前方法的不足与后续改进方向

虽然当前脚本已经实现了“刚性配准 + 重采样 + 共同有效区域裁剪”的基本流程，并且能够输出 nnU-Net 所需的多通道 NIfTI 文件，但该方法仍然存在一些不足。

4.1 目前只使用刚性配准，无法处理明显局部形变

当前脚本采用的是刚性配准，只允许 moving image 发生整体平移和旋转。

这种方法的优点是稳定、安全，不会改变器官和病灶形态。但它的缺点是无法处理局部形变。

在多期相 CT 中，不同期相之间可能存在呼吸状态差异。例如，肝脏、膈肌、胃肠道等腹部结构可能发生局部位置变化。此时，仅靠刚性配准可能无法让所有局部结构都完全对齐。

因此，当前方法更适合处理整体位置偏移，而不适合处理明显的局部器官变形。

后续可以考虑：

  • 在刚性配准基础上尝试仿射配准；

  • 对明显局部形变的样本谨慎尝试非刚性配准；

  • 但非刚性配准可能改变病灶形态，因此需要严格验证。

4.2 互信息指标虽然稳妥，但不一定适合所有情况

当前脚本使用 Mutual Information 作为相似性指标。

互信息适合多期相 CT，因为不同期相的灰度分布可能不同。但是互信息也不是万能的。

如果两个序列重叠区域较少、扫描范围差异较大，或者图像中存在大量空白区域，互信息可能受到干扰，导致配准结果不稳定。

后续可以尝试比较不同 metric 的效果，例如：

  • Mutual Information

  • Normalized Mutual Information

  • NCC

  • MSE

对于多期相 CT，互信息通常更稳妥，但仍然需要通过可视化检查验证。

4.3 共同有效区域裁剪依赖阈值，可能不够精确

当前脚本使用 image > -500 来粗略判断身体区域，并据此计算 fixed 和 moving_registered 的共同有效区域。

这种方法简单有效，但它本质上是一个经验阈值方法。

可能的问题包括：

  • 阈值设置不合适时，可能漏掉部分低密度组织；

  • 床板、伪影或增强区域可能被错误纳入身体区域；

  • 对不同扫描协议、不同窗宽窗位、不同部位的数据适应性有限；

  • 如果图像噪声较多，mask 可能不稳定。

后续可以改进为：

  • 根据 CT HU 值范围做更稳定的身体 mask；

  • 结合形态学操作去除小噪声区域；

  • 只保留最大连通区域；

  • 对不同部位设置不同阈值；

  • 或结合已有器官分割结果进行更精确裁剪。

4.4 当前方法默认两个期相具有足够重叠区域

当前流程的前提是 fixed 和 moving 至少有较大的共同扫描区域。

如果两个序列本身不是同一部位，或者共同区域太少，即使运行成功，配准结果也没有实际意义。

例如：

  • fixed 是胸腹部；

  • moving 是盆腔；

  • 或者 moving 只是 localizer/scout；

  • 或者两个期相扫描范围重叠很少。

这种情况下，脚本可能仍然输出文件，但结果不适合作为 nnU-Net 输入。

因此，后续批量化时必须加入序列筛选机制：

  • 排除 localizer / scout；

  • 排除切片数过少的序列；

  • 检查 BodyPartExamined、SeriesDescription、ProtocolName；

  • 检查 ImageOrientationPatient 是否一致；

  • 检查 z 方向扫描范围是否有足够重叠。

4.5 当前脚本只处理两个期相，还没有完全扩展到多期相批量处理

当前脚本主要验证了一个 fixed image 和一个 moving image 的配准流程。

但实际数据中，一个病例可能包含多个期相或多个序列，例如：

  • 平扫期；

  • 动脉期；

  • 门静脉期；

  • 延迟期；

  • std5mm；

  • lung5mm。

后续还需要扩展为：

  • 一个 fixed image；

  • 多个 moving image；

  • 每个 moving image 分别配准到 fixed；

  • 最后统一裁剪共同有效区域；

  • 输出 case_001_0000.nii.gz、case_001_0001.nii.gz、case_001_0002.nii.gz 等多通道文件。

也就是说，目前只是单病例、双期相流程，后续还需要做多期相和多病例批量化。

4.6 目前主要依赖可视化检查，缺少定量评价指标

当前配准效果主要通过 overlay 或人工观察判断，例如：

  • 肺是否对肺；

  • 肝是否对肝；

  • 骨骼轮廓是否重合；

  • 同编号切片是否对应同一部位。

这种方式直观，但主观性较强。

后续可以增加定量评价指标，例如：

  • Mutual Information 配准前后变化；

  • NCC；

  • Dice；

  • Hausdorff Distance；

  • 器官边界重合程度。

如果有器官或病灶标签，可以用标签重叠程度评价配准效果。但需要注意，标签重采样必须使用最近邻插值。

4.7 重采样仍然会引入插值误差

无论是刚性配准还是仿射配准，最终都需要进行重采样。

重采样过程中使用线性插值，会重新计算体素灰度值，可能导致：

  • 边缘轻微模糊；

  • 灰度分布轻微变化；

  • 小结构细节损失；

  • 病灶边界不如原图锐利。

因此，重采样次数应尽量减少。

当前方法中，应避免对同一图像反复重采样。比较合理的做法是：

  • 先求最终 transform；

  • 再一次性 resample 到 fixed 空间；

  • 不要多次中间保存再重复 resample。

4.8 当前没有专门处理标签文件

当前脚本主要处理 CT 原图。

如果后续需要把分割标签也变换到同一空间，需要额外处理标签文件。

标签文件不能使用线性插值，而必须使用最近邻插值。

否则标签值可能从整数类别变成小数，导致类别错误。

例如：

  • 原标签 0 = background；

  • 1 = liver_cyst；

  • 2 = liver_hemangioma；

  • 3 = tumor。

如果线性插值，可能出现 1.3、2.6 这种无意义标签值。

因此，后续处理标签时应使用：

  • 图像：线性插值；

  • 标签：最近邻插值。

4.9 配准失败时缺少自动报警机制

目前脚本可以输出结果，但如果某个样本配准失败，脚本不一定能自动判断。

例如：

  • 图像被转歪；

  • 共同区域过小；

  • 输出大面积空白；

  • fixed 和 moving 实际不是同一部位；

  • optimizer 收敛到错误位置。

这些情况可能需要人工打开图像才能发现。

后续批量化时可以加入自动检查机制，例如：

  • 检查共同有效区域体积是否过小；

  • 检查裁剪后 z 层数是否太少；

  • 检查 moving_registered 中有效体素比例；

  • 检查配准前后 metric 是否明显改善；

  • 记录每个 case 的日志。

4.10 当前流程还没有完全封装成 nnU-Net 数据集格式

当前脚本已经可以输出：

  • case_001_0000.nii.gz

  • case_001_0001.nii.gz

但完整 nnU-Net 数据集还需要：

  • imagesTr/

  • labelsTr/

  • dataset.json

  • 正确的通道名称；

  • 正确的标签定义；

  • 训练集/测试集划分。

因此，配准只是 nnU-Net 输入准备的一部分。

后续还需要把配准结果进一步整理成标准 nnU-Net 数据集结构。

4.11 当前不足总结

当前方法已经完成了多期相 CT 配准的基本流程，但仍然存在以下不足：

  1. 只使用刚性配准，难以处理局部形变；

  2. 互信息指标不一定适合所有异常情况；

  3. 共同区域裁剪依赖阈值，鲁棒性有限；

  4. 默认两个期相有足够重叠区域，不适合不同部位数据；

  5. 当前只验证了双期相，还未完全扩展到多期相批量处理；

  6. 配准效果主要依赖人工可视化检查，缺少定量评价；

  7. 重采样会引入插值误差；

  8. 标签文件需要单独使用最近邻插值处理；

  9. 缺少配准失败的自动报警机制；

  10. 还需要进一步整理为完整 nnU-Net 数据集格式。

后续改进方向是：

在当前刚性配准流程稳定的基础上，进一步实现多期相批量处理、自动质量检查、共同区域裁剪优化，以及 nnU-Net 数据集格式自动生成。