周勤,王远军
上海理工大学 健康科学与工程学院(上海,200093)
医学图像配准是在一个或多个浮动图像(moving image)上应用各种几何变换的优化过程,使浮动图像和参考图像(fixed image)在相同的结构上达到空间位置一致。为了方便描述,我们将两个图像序列的配准称为成对(pair-wise)配准,一组图像序列间的配准称为群组(group-wise)配准[1]。对于一些分析研究来说,配准是一个至关重要的过程,包括旨在了解表型群体趋势的研究、测量纵向变化(例如监测肿瘤组织的大小)、执行引导手术、将个体解剖与标准空间系统(即图谱)相关联以及其他应用[2]。对自不同类型设备,例如磁共振成像(magnetic resonance imaging, MRI)和计算机断层扫描(computer tomography, CT)正电子发射断层扫描(positron emission tomography,PET)以及等组合扫描的图像进行配准对医疗疾病诊断治疗十分重要[3]。
近年来,研究人员对同时配准两个或两个以上的图像越来越感兴趣[4-7]。在之前工作中[8],我们使用变分推断算法来实现图像群组配准,尽管这种方法可以实现图像集的良好配准,但是容易导致构建的参考图像不平滑并且收敛速度慢的问题。Che等[9-10]提出一个由PCA构建的模板图像引导的无偏差的深度群组配准框架,适用于多光谱图像。类似于文献[9-10]的做法,Roland等[11]使用了鲁棒的主成分分析(principal component analysis, PCA)的方法。可见,PCA算法应用在群组配准领域具有很大的潜力。
针对文献[8]存在容易导致构建的参考图像不平滑并且网络收敛速度慢的问题,提出使用PCA算法来计算构建参考图像时每幅图像的权值。然后,利用改进后的权值计算方法重新构建参考图像。最后,实验结果表明提出的算法在配准精度上虽然不如文献[8]方法的精度高,但在模板构建的收敛速度以及所构建的模板平滑性上表现更好。主要贡献如下:
(1)为了使构建的参考图像更加平滑和网络更快的收敛,提出一种基于PCA的无监督群组配准算法,可以使构建的参考图像快速收敛到图像集拓扑中心,减少网络训练时间。
(2)提出了一种基于U-net的编解码网络,通过跳跃连接有效结合深层和浅层特征,快速优化网络参数。
给定一组移动图像,目标是将N幅移动图像配准到模板图像上,使所有配准后的图像都与模板图像相似。方法流程如图1所示,其中输入图像、形变场、配准后的图像和模板图像均为三维体数据,图中为其中心切片展示,具体的有:图像的输入有N幅图像和PCA生成的模板图像。神经网络的输出是速度场,经过积分和采样生成N个形变场;随后,空间变换网络使用双线性插值来构建密集的空间变换,用于通过优化图像相似性度量来配准3D图像。这里形变场的估计、损失函数和神经网络模型构建的细节,请参考文献[8]最大的不同,仅在模板图像的构建上。
图1 基于PCA的无监督群组配准算法流程图Fig.1 Flow chart of unsupervised group-wise registration algorithm based on PCA
在我们的方法中,每个图像的像素坐标被采样作为单独的观察点,不同的图像被作为对应于N维空间的不同变量。使用PCA将维度降低到一维子空间。与最大特征值相关联的特征向量V可以用作构造模板图像t的权重w。
我们使用 LONI LPBA40[14](以下简称 LPBA40)和MICCAI2012大挑战中的数据[15](以下简称MICCAI)这两个不同的数据集去验证提出的无监督配准算法。数据集分别覆盖了从青少年脑MRI图像到老年人脑MRI图像。两个数据集的所有图像都是MRI的TI加权图像,对所有的图像均进行偏移场矫正,颅骨剥离,图像重采样、归一化和图像预配准等预处理步骤进行处理,具体描述见表1。两个数据均具有对应的标签,用于评估配准算法的精度,如图2所示。
表1 LPBA40和MICCAI数据集描述Tab.1 LPBA40 and MICCAI dataset descriptions
图2 两个数据集体数据和对应标签的水平位中心切片结果展示Fig.2 Two data collective data and corresponding label horizontal bit center slice results are displayed
为了更好地评估本文提出的算法,与文献[8]算法在 Dice系数 (Dice coefficient)[17]、变形场的雅克比行列式[18]、EM1[19]、EM2[19]和时间五个度量标准来评估配准准确率、变形场的折叠程度、模板与总体之间的偏差、模板的锐度差和配准的快慢。
训练时,模板图像通过对原始数据使用PCA算法进行初始化和更新,将移动图像和参考图像在通道维度上合并, 然后送到CNN网络学习预测形变场, STN利用预测的形变场得出扭曲后的图像,使用Adam算法反向更新网络参数直至使目标函数最优。测试时, 将参考图像和测试图像输入训练好的模型, 将预测的形变场作用于测试图像标签, 得到形变后的标签图, 通过这种方式计算不同标签的dice、EM1、EM2和变形场的重叠程度。具体的使用pytorch在6-core Intel i7-8700K CPU和6 GB NVIDIA GeFore RTX 2060 GPU机器上实现模型的训练与测试。实验环境为DUDA10.0并行计算架构, 操作系统为Win10, 软件为pyCharm。算法在pytorch框架中实现, 使用学习率为4e-4的Adam优化器优化训练,迭代次数为15 000次。超参数λa和λb分别设置为0.9和0.09。由于实验环境内存大小的限制, 仅将批处理大小设为1。通过设计以下两个不同的实验,对提出算法的精度和鲁棒性进行综合性评估。
(1) LPBA40数据集配准实验,并同文献[8]比较训练期间模板的变化。
(2)为进一步验证算法的鲁棒性, 将LPBA40数据集训练好的网络模型直接应用于MICCAI数据集。
为了分析提出算法的性能,使用LPBA40训练网络。首先,平均图像和PCA构建的模板被单独用作模板,用于比较训练期间清晰度的变化,如图3所示。其次,针对单个图像、平均图像和该组中的PCA用作模板图像,研究图像算法的性能,如表2和图4所示。最后,将基于LPBA40数据集训练好的模型直接应用到MICCAI数据集上进行配准,验证构建模型的迁移能力,配准结果如图5所示。
图3 训练期间模板变化Fig.3 Template changes during training
图3给出群组配准实验训练期间模板图像变化情况。图3从三个方位展示了训练期间平均法和PCA构建的模板图像的变化情况。从图3可以看出,相比之下,本文算法构建的模板更清晰,更靠近群组中心。在6 000 epoch时,PCA模板图像的清晰度已经非常高,超过平均模板图像的清晰度。显然,将PCA模板图像与深度网络相结合的技术比使用平均模板图像收敛得更快。收敛比较快的原因是基于PCA的模板图像考虑一组图像的主成分信息,使模板图像更接近群组中心。
表2给出不同算法在测试集上的配准结果。从表2中可以看出, 在大多指标上性能表现优异(粗体表示最优结果)。同本文算法相比,文献[8]方法具有更高的Dice得分、GPU运行时间少, 产生更好的微分同胚形变场(具有更低的非负雅可比位置数对)和生成的平均模板更清晰( 相对较小)。综合来看, 本算法表现较差。图4、图5为部分结果展示,其中图5为将基于LPAB40数据集训练的模型应用MICCAI数据集进行配准的结果展示。图4、图5第一列表示移动图像,第二列为模板图像,第三列为配准后的图像,第一行表示脑横截面图,第二行和第三行为脑3D重建和可视化结果。从图4和图5可以看出,配准后的图像和模板图像较为接近,从箭头来看,配准结果与模板图像很相似。在3D重建效果图中,虚线标注了大脑皮层结构中配准效果比较好的区域。从图4、图5中可以看出,在中央后回及中央脑沟皮层结构中得到较好的配准结果。
图5 MICCAI配准结果图Fig.5 MICCAI registration result
表2 不同算法在LPBA40测试集上多个度量标准的群组配准结果Tab.2 Group registration results of multiple metrics with different algorithms on LPBA40 test set
图4 LPBA40配准结果Fig.4 LPBA40 registration result
基于移动图像和模板图像之间的结构相似性,提出在构建模板阶段可以更多地关注每幅图像对模板构建的影响,而深度学习中PCA算法可以很好地适用于模板构建。同时,提出的算法结合了文献[8]的神经网络框架和损失函数,从而更好地构建模板。实验结果表明,虽然所提算法的配准精度不如文献[8],但模板的收敛速度较快且模型也具有很好的迁移能力。后续会考虑使用更为复杂的网络或者尝试更改网络的超参数来提高本算法的性能。