4. PET/CT空间配准基础

PET与CT在同一台设备上采集，但两套图像来自不同物理通道：CT记录的是X射线衰减系数，PET记录的是正电子湮灭产生的符合事件计数。即便患者躺在同一张检查床上，呼吸、肠道蠕动、扫描时间差以及重建参数差异，都会让两幅图像在空间上出现偏移。空间配准的任务，就是寻找一个空间变换，使PET体素与CT体素在解剖位置上对齐。本章聚焦刚性配准，即只允许平移和旋转，不允许缩放和剪切，这是PET/CT融合中最常用、最稳健的配准模型。

刚性配准的数学形式与自由度

三维刚性变换可写为 y = R·x + t，其中 x 是浮动图像（PET）中的点，y 是参考图像（CT）中的对应点，R 是3×3旋转矩阵，t 是3×1平移向量。旋转矩阵由三个欧拉角决定，加上三个平移分量，共6个自由度。若采用四元数表示旋转，则参数为4个四元数分量加3个平移分量，共7个参数，但四元数需满足单位范数约束，实际独立参数仍为6个。在ITK/SimpleITK中，通常使用Euler3DTransform或VersorRigid3DTransform，前者以欧拉角为参数，后者以四元数为参数，后者在优化过程中数值更稳定。

需要特别说明坐标系约定：ITK中图像坐标以体素索引为基准，物理坐标通过方向矩阵、间距和原点映射到世界坐标。旋转中心默认取图像几何中心，也可显式指定。若旋转中心设置不当，优化器会在平移和旋转之间产生耦合，导致收敛缓慢甚至失败。因此，在配准前应确认两幅图像的方向矩阵一致，必要时用SetDirection统一到同一解剖方向。

互信息与归一化互信息

PET与CT的灰度关系不是线性的，甚至不是单调的：高代谢区域在PET上为高计数，在CT上可能是软组织密度，也可能是骨密度。因此，基于灰度差平方和或相关性的相似性测度不适用。互信息（Mutual Information, MI）从信息论角度衡量两幅图像灰度分布的统计相关性，不要求灰度之间存在函数关系，因此成为多模态配准的首选测度。

设参考图像灰度 A 与浮动图像灰度 B，其联合概率 p(a,b) 和边缘概率 p(a)、p(b) 可由联合直方图估计。互信息定义为：

MI(A,B) = Σ p(a,b) · log( p(a,b) / (p(a)·p(b)) )

当两幅图像完全对齐时，联合直方图趋于集中，互信息达到最大。但MI对重叠区域大小敏感：当变换把浮动图像移出参考图像视野时，重叠区域减小，边缘熵下降，MI可能反而增大，造成“越移越优”的假象。归一化互信息（Normalized Mutual Information, NMI）通过归一化联合熵与边缘熵，削弱这种重叠依赖：

NMI(A,B) = ( H(A) + H(B) ) / H(A,B)

其中 H(A)、H(B) 为边缘熵，H(A,B) 为联合熵。NMI在重叠区域变化时更稳定，适合PET/CT这类视野可能不完全重合的场景。SimpleITK中对应NormalizedMutualInformation，ITK中对应NormalizedMutualInformationImageToImageMetric。

配准参数搜索与优化策略

配准本质是一个多维优化问题：在6维参数空间中寻找使相似性测度最大（或代价最小）的参数组合。由于MI/NMI曲面通常非凸，存在多个局部极值，直接使用梯度下降容易陷入局部最优。工程上常采用多分辨率策略：先在下采样图像上粗配准，捕获大尺度平移和旋转，再逐步提高分辨率精细调整。SimpleITK的ImageRegistrationMethod支持设置SetShrinkFactorsPerLevel和SetSmoothingSigmasPerLevel，典型配置为3层，收缩因子[4,2,1]，平滑sigma[2,1,0]。

优化器选择上，Powell法不需要梯度，对噪声和插值不连续更鲁棒，适合MI这类非光滑测度；L-BFGS-B收敛快，但需要测度可微，且对初始值敏感。实际PET/CT配准中，若初始偏移较小（同机架、同床、同次扫描），可用L-BFGS-B；若初始偏移较大或跨设备，建议先用Powell粗搜，再用L-BFGS-B精调。学习率、迭代次数和收敛容差需根据测度曲线调整：若测度在若干次迭代内变化小于1e-6，可认为收敛。

SimpleITK/ITK实战：PET-CT刚性配准

下面给出一段可直接运行的SimpleITK Python示例，假设已读取CT和PET的NIfTI文件，且已统一方向。代码中显式设置旋转中心为CT几何中心，采用NMI测度、3层多分辨率、Powell优化器，并输出最终变换参数和重采样后的PET。

import SimpleITK as sitk

fixed = sitk.ReadImage("ct.nii.gz", sitk.sitkFloat32)
moving = sitk.ReadImage("pet.nii.gz", sitk.sitkFloat32)

# 统一方向，避免方向矩阵不一致导致旋转中心偏移
fixed = sitk.DICOMOrient(fixed, "LPS")
moving = sitk.DICOMOrient(moving, "LPS")

initial_transform = sitk.CenteredTransformInitializer(
    fixed, moving,
    sitk.Euler3DTransform(),
    sitk.CenteredTransformInitializerFilter.GEOMETRY)

registration = sitk.ImageRegistrationMethod()
registration.SetMetricAsMattesMutualInformation(numberOfHistogramBins=50)
registration.SetMetricSamplingStrategy(registration.RANDOM)
registration.SetMetricSamplingPercentage(0.1)
registration.SetInterpolator(sitk.sitkLinear)
registration.SetOptimizerAsPowell(
    numberOfIterations=200,
    maximumLineIterations=20,
    stepLength=1.0,
    stepTolerance=1e-6,
    valueTolerance=1e-6)
registration.SetOptimizerScalesFromPhysicalShift()
registration.SetShrinkFactorsPerLevel([4, 2, 1])
registration.SetSmoothingSigmasPerLevel([2, 1, 0])
registration.SmoothingSigmasAreSpecifiedInPhysicalUnitsOn()
registration.SetInitialTransform(initial_transform, inPlace=False)

final_transform = registration.Execute(fixed, moving)
print("Final metric value:", registration.GetMetricValue())
print("Transform:", final_transform)

resampled = sitk.Resample(
    moving, fixed, final_transform,
    sitk.sitkLinear, 0.0, moving.GetPixelID())
sitk.WriteImage(resampled, "pet_registered.nii.gz")

运行后应关注三个输出：最终测度值、变换参数和重采样图像。测度值本身没有绝对意义，但可与初始测度比较，若提升不明显，说明初始对齐已较好或优化未充分收敛。变换参数中平移分量若超过约20 mm，需检查是否因旋转中心或方向矩阵设置错误导致。重采样图像可用ITK-SNAP或3D Slicer叠加到CT上，观察心脏、肝脏、脊柱等标志性结构是否对齐。若PET高代谢灶与CT对应解剖位置存在系统性偏移，可尝试增加直方图bin数或改用SetMetricAsNormalizedMutualInformation()。

结果判断与常见问题

判断配准是否成功，不能只看测度值。建议从三个层面检查：第一，视觉检查，在冠状面、矢状面和横断面分别叠加，观察骨性标志和器官轮廓；第二，定量检查，计算配准后PET与CT在感兴趣区内的互信息或归一化互信息，与初始值比较；第三，临床检查，确认高代谢灶是否落在CT对应的解剖位置。若发现平移方向正确但旋转角度异常，可尝试将初始变换改为VersorRigid3DTransform，并设置SetOptimizerScalesFromPhysicalShift()，使优化器在物理空间中等尺度搜索。

常见失败原因包括：直方图bin数过少导致测度曲面粗糙，采样比例过低导致测度估计噪声大，插值器选择不当导致测度不连续，以及优化器步长与物理尺度不匹配。针对这些，可分别调整numberOfHistogramBins至50~64，采样比例至0.1~0.2，插值器用线性或B样条，并启用SetOptimizerScalesFromPhysicalShift()。若仍不收敛，可先用sitk.CenteredTransformInitializer的MOMENTS模式替代GEOMETRY，利用图像矩估计初始平移。