吴香奕,曹锋,曹瑞芬,吴茜,董江宁,徐榭,裴曦,5
磁共振成像(Magnetic Resonance Imaging,MRI)具有优良的软组织对比度,在现代放射治疗中起着至关重要的作用。目前,临床上常采用MRI 和计算机断层扫描(Computer Tomography,CT)图像进行多模态图像配准,用于勾画靶区和危及器官。然而,MRI-CT 图像配准带来很大的不确定性,会影响整个治疗过程。例如,对于头部和前列腺,该配准不确定性能达到0.5~3.5 mm[1-2]。而在MRI-only 放疗计划中,由于不需要CT 图像以及MRI-CT 图像配准,因此就不存在配准导致的系统误差,同时还能简化放射治疗工作流程,减少患者CT 成像导致的剂量[3]。此外,MRI-only 放疗过程中,可以根据实时MRI图像调整放疗计划,为在线自适应放射治疗的临床实施提供一种新的解决方案,避免MRI 与计划CT 配准带来的不确定性。然而,在MRI-only 放疗临床实施中,MRI 图像的体素强度与放疗计划中计算剂量所需的电子密度信息之间缺乏直接关联,因此需要将MRI图像转换成含有电子密度信息的伪CT(Synthetic CT,sCT)图像。目前,已有方法包括体积密度分配方法[4-6]、基于体素方法[6-12]和基于图谱方法[13-22],以及最近发展的深度学习方法。研究显示使用卷积神经网络(Convolutional Neural Network,CNN)[23-24]或生成对抗网络(Generative Adversarial Network,GAN)[25-26]的模型优于使用体素和图谱的方法,但这些深度学习方法都设置了基于sCT与其对应真实CT之间体素级比较的损失项,因此必须采用同一患者的配对MRI 和CT 图像进行训练。然而利用配对图像训练主要存在以下两个问题:一方面,同一患者的配对MRI 和CT 图像在临床实践中难以获取;另一方面,基于体素级比较的损失项,容易受到MRI-CT 图像配准误差的影响[27]。另外,训练集中配对MRI和CT 图像也可能对GAN 模型本身的训练产生负面影响[26,28]。因此,不同于之前的研究,本文提出一种新的MRI和CT 图像转换方法[29],引入基于归一化互信息的相似性约束损失项,使用非配对患者的MRI 和CT 图像来训练循环一致性生成对抗网络(CycleGAN)模型,将盆腔MRI 图像转换为sCT 图像。使用平均绝对误差(Mean Absolute Error,MAE)和平均误差(Mean Error,ME)评估测试集中所有患者的sCT 精度,并计算sCT 与其对应CT 图像间的放疗剂量差异和伽马通过率,以验证CycleGAN 模型利用非配对训练数据生成sCT图像的可行性,以及基于生成的sCT图像做放疗计划的准确性。
从中国科学技术大学第一附属医院西区影像科直肠癌患者图像数据集中,随机抽取到35例患者的MRI图像,36例患者的CT图像,以及10例患者的MRI和CT图像。患者年龄32~83岁。由于均是诊断图像,扫描时不采用热塑模,因此不像放疗流程中所采集图像那样具有很好的固定性,这些诊断图像中显示出的身体轮廓并不规则,差异较大(图1)。T1-MRI(重复时间500~520 ms,回波时间6.9~8.0 ms,翻转角度90°)图像由3.0T 磁共振设备(Signa HDxt,GE Healthcare Technologies,Milwaukee,Wisconsin,USA)扫描获得。CT 图像是在CT 机上(Discovery 750,GE Healthcare Technologies,Milwaukee,Wisconsin,USA)以120 kV和176~420 mA扫描获得。MRI图像的分辨率从(6.000×0.625×0.625)mm3至(6.000×0.820×0.820)mm3不等,CT图像的分辨率均为(5.000×0.854×0.854)mm3。所有图像都重新采样到统一分辨率(5.000×1.000×1.000)mm3,横断位片层裁剪为324×424大小。采用最大类间方差法(Otsu)和形态学操作提取出CT和MRI图像中的盆腔区域,排除金属支架和空气。需要说明的是,训练集中MRI和CT图像无需任何形式的配准。
图1 从非配对训练集中随机选出的CT(a)和MRI(b)图像Fig.1 CT(left)and MRI(right)images randomly selected from the unpaired training set
MRI 和CT 图像是人体组织两种不同的医学成像方式,因此MRI 和CT 图像之间存在着潜在的相关性。这些相关性可以通过Zhu等[30]提出的CycleGAN模型来求得,并实现人体组织的MRI 和CT 图像之间的相互转换。该模型最大的优点是其训练不需要配对的MRI 和CT 图像,正好符合本研究基于非配对图像生成盆腔sCT的目的。图2展示了由两个循环和4个独立全卷积网络(Fully Convolutional Networks,FCNs)组成的CycleGAN 模型结构。两个循环分别为前向循环和后向循环。4 个FCNs包括两个生成器GMRI-CT和GCT-MRI,两个判别器DMRI和DCT。在前向循环中,训练生成器GMRI-CT将输入的MRI 图像转换为sCT 输出图像,GMRI-CT尽量让判别器DCT难以辨别sCT 和CT 图像。训练生成器GCT-MRI将上一步生成的sCT 图像再转换回MRI 图像,称为重建MRI 图像(rMRI),并且尽量让rMRI 与输入的MRI 图像一致。总的来说,前向循环实现MRI 到sCT 再到rMRI 的循环。后向循环与前向循环相反,实现CT 到伪MRI 图像(sMRI)再到重建的CT图像(rCT)的循环。两个生成器采用的FCNs 结构相同,包括两个跨步卷积层、9个残差网络和两个反卷积层[31]。该模型能够灵活地适应不同尺寸的输入图像,并得到与输入图像尺寸相同的输出图像。且模型中采用的残差网络能够简化MRI 和CT 图像之间复杂关系的学习。本文两个判别器的FCNs 采用PatchGAN 结构[32],能将输入图像以局部图像块为单位进行真假判别。这样只需训练较少的参数从而提高训练效率,同时基于局部图像块进行判别也可以提高生成图像的精度。
图2 CycleGAN模型的整体结构和两个循环Fig.2 The overall structure and two cycles of CycleGAN model
本文模型采用的总目标函数由对抗损失项(Adversarial Loss,Ladv)、循环一致性损失项(Cycleconsistency Loss,Lcycle)和相似性约束损失项(Similarity-constraint Loss,Lsc)3 部分组成,以保持伪图像与其对应输入图像之间的结构一致性。Ladv引导生成器GMRI-CT学习将输入的MRI图像转换为具有CT 图像分布特征的sCT 图像,以使判别器DCT也无法判别出sCT 是不是真实的CT 图像,即无法准确地将其判别为0。GMRI-CT同理。Ladv函数如下:
其中,Pdata(CT)是CT 图像的分布,Pdata(MRI)是MRI 图像的分布,GCT-MRI( CT )是给生成器输入训练集中的CT图像之后生成的sMRI 图像,DMRI(GCT-MRI( CT ))是将sMRI 图像输入判别器DMRI之后得到的实际判别值(取值范围为[0,1]),1 是理想状态下的判别值,是基于训练集中CT 图像生成的sMRI 输入判别器DMRI之后得到的实际判别值([0,1])与理想值1之间均方根误差的期望值。同理。
Lcycle是输入真实图像和输出伪图像之间的间接结构约束项。生成器GMRI-CT和GCT-MRI应互逆,两个映射都应该满足双射,以使在两个循环中都可以将输入图像转换成伪图像之后还能再转换回输入图像,比如前向循环中的MRI 到GMRI-CT生成的sCT 再到GCT-MRI生成的rMRI,此rMRI图像理想情况下应该与真实的输入MRI图像一样。该损失项的函数如下:
由于本文是研究基于非配对训练数据生成精确sCT 图像的可行性,训练过程中并没有对应的真实CT图像与生成的sCT图像做对比。因此为了进一步确保sCT图像的精度,需要对生成的sCT图像与输入的MRI 图像进行直接的相似性约束。本文引入了多模态图像配准中常用的相似性测度—归一化互信息(NMI)来作为相似性约束损失项Lsc。该损失项的函数由分别对应两个生成器的两部分组成:
其中,生成器GMRI-CT部分的损失项函数如下:
H(MRI) 是 MRI图像的信息熵,H(GMRI-CT( MRI ))是GMRI-CT生成的sCT 图像的信息熵,H(GMRI-CT( MRI ),MRI)是两图像的联合熵。当两图像完全独立时,联合熵最大,Lsc(GMRI-CT)取1;当两图像完全对齐时,联合熵最小,Lsc(GMRI-CT)取0。p(MRIi)和p((GMRI-CT(MRI))j)是MRI 和sCT 图像的像素值分布,p(MRIi,(GMRI-CT(MRI))j)是两图像像素值的联合分布。同理可得Lsc(GCT-MRI)。
判别器的目标函数引导其提高判别真实图像和伪图像的能力。DCT的目标函数如下:
DMRI的目标函数同理可得。
训练集包括35例患者的MRI图像和另外36例患者的CT图像,总计932张MRI和914张CT横断位图像。测试集中10例患者共有253张MRI图像和220张CT横断位图像。CT像素值分布范围约束为-600~1 500 HU,MRI约束为0~2 500。循环一致性损失项和相似性约束损失项的权重均设为10。为了增加训练样本的数量,输入图像在进入网络训练之前被随机裁剪为256×256大小。模型采用Adam优化器[33]进行200代优化训练,在前100代训练中,学习率固定为0.000 2,在后100代训练中,学习率从0.000 2线性降低至0。Batch的大小设为6。在单个NVIDIA Quadro M6000 24 GB上,模型训练耗时约25 h。测试集中10例患者的sCT图像生成耗时约9 s。
由于测试集中同一患者的MRI 和CT 图像扫描时间和扫描摆位不同,两者并非完全一一对应。因此本文先使用Elastix[34]将生成的sCT 图像配准到其对应的真实CT 图像,以便对两者进行对比,从而评估sCT精度。根据盆腔部位CT值特征,设置-600 HU为阈值对sCT 和真实CT 图像进行阈值分割,提取感兴趣的身体区域,排除背景和空气,并将空气和背景的值设置为常见的空气HU 值:-1 000 HU。利用2×2的棋盘格图像展示真实CT和sCT图像之间的结构一致性。为直接比较两图像的体素值差异,在两图像感兴趣区域内计算MAE 和ME 以定量评估sCT 精度。MAE和ME的计算如下:
其中,N表示真实CT 图像和sCT 图像感兴趣区域内的体素个数,CTi和sCTi分别表示真实CT 图像和sCT图像感兴趣区域内每个体素的像素值。
验证剂量分布的准确性对于预测临床应用的效果非常重要,常用的临床验证指标包括剂量差异和伽马通过率。国际上sCT的放射治疗验证研究中,也均采用这些指标。因此,针对测试集中所有直肠癌患者的sCT,本研究在10%、50%和90%的剂量阈值区域内分别计算剂量差异,以及分别计算使用3%/3 mm、2%/2 mm和1%/1 mm标准时的伽马通过率。对于测试集中的所有患者,均在Pinnacle放疗计划系统中,以(3×3×3)mm3计算网格,使用卷积叠加算法,基于真实CT进行调强放射治疗计划,再将得到的原始放疗计划以相同的算法和设定参数,加载到sCT上重新计算出基于sCT的剂量分布。处方剂量为58~71 Gy。
伽马通过率(Gamma pass rate)和每个体素的剂量差异(Dose difference)的计算如下:
其中,ΔDose 表示sCT 中当前体素与CT 中离当前体素相距Δdistance 的体素的剂量差异,Δdistancemax和ΔDosemax为计算伽马通过率设置的限制标准。例如3%/3 mm 标准下的伽马通过率计算时的Δdistancemax和ΔDosemax分别为3 mm和3%。
其中,DoseCT表示基于CT计算出的每个体素的剂量,DosesCT表示基于sCT 计算出的每个体素的剂量,Doseprescribe表示处方剂量。
图3采用棋盘格图像展示了随机挑选的4 例患者不同部位的CT 和sCT 图像的结构一致性,从中可以看到在这些病例的真实CT和sCT图像之间,轮廓、骨骼、肌肉和皮下脂肪等连接都很自然。但在第4例患者的sCT和真实CT图像中出现了内部空腔不一致的现象。图4是第4例患者的MRI图像、sCT图像、真实CT 图像和两者棋盘格图像的内部细节展示。从图中可以看到,病例4的MRI和sCT图像内部结构是保持一致的,但该病人的MRI 和CT 图像之间直肠中内容物存在差异。这很可能是由于这两种图像是在不同时期扫描得到的(该病例MRI 扫描时间为2018年7月27日,CT 扫描时间为2018年7月26日)。所以虽然sCT 与MRI图像的结构一致性保持很好,sCT空腔部位也无法与CT 对应。也就是说,同一病人的MRI和对应CT之间存在差异,导致sCT和CT图像之间的空腔部位差异变大。
图5显示了测试集中所有患者的真实CT 和sCT图像间的MAE 和ME 值箱型统计分布。平均MAE(±std)和ME(±std)值分别为35.537(±4.537)HU 和11.278(±7.425)HU。其中结果最好的一个病例的MAE 和ME 值分别为30.784 和3.713 HU。结果最差的一个病例的MAE 和ME 值分别为41.559 和25.645 HU。图6则显示了某测试病例的sCT 及其对应真实CT图像间的HU 值差异图。可以看出,HU 值差异比较大的主要区域是在皮质骨和皮肤边界,且生成的sCT图像相比真实CT图像要模糊一些。
图3 由4个患者的sCT和真实CT图像组成的棋盘格图像Fig.3 Checkboard images of sCT and real CT images of 4 patients from the test set
图4 第4例患者的MRI、sCT和真实CT图像的内部空腔细节图Fig.4 Detail image of air package in MRI,sCT and CT images of patient#4
图7显示了对测试集中某病例基于相同放疗计划和参数,在真实CT 及对应sCT 图像上计算得到的剂量分布,及两者间剂量差异。从图中可以发现,基于sCT及其CT图像计算的剂量分布几乎无肉眼可见差异。如剂量分布差异图7所示,绝大部分体素的剂量差异都小于1%,只有皮肤边界上少量体素的剂量差异大于3%。
图5 测试集中所有患者sCT与CT图像的MAE和ME 值的箱型统计分布(橙色实线表示中值,绿色虚线表示平均值)Fig.5 Box-wise statistical distribution of MAE and ME values of sCT and CT images of all patients in the test set(The solid orange line indicates the median value,and the green dotted line indicates the mean value)
图8显示了不同剂量阈值区域内(10%、50%和90%)相对于处方剂量的平均剂量差异,以及不同剂量阈值区域内不同标准下(3%/3 mm、2%/2 mm 和1%/1 mm)的伽马通过率。在不同的剂量阈值区域内,平均剂量差异均小于0.5%。3D 伽马分析也表明,对于测试集中所有患者,基于sCT 图像进行放疗计划的平均伽马通过率分别大于99%、98%和95%。
表1统计了测试集中所有直肠癌患者,在10%、50%和90%处方剂量为阈值的感兴趣区域内,相对于处方剂量的平均剂量差异和伽马通过率。结果表明,在测试集的所有患者中,真实CT和sCT图像间的剂量分布是一致的;在10%处方剂量阈值区域,病例的最佳平均剂量差异值为0.40%,最差值为0.65%。3%/3 mm 标准下伽马通过率最佳值为100%,最差值为99.3%。
图6 MRI,CT和sCT图像,以及CT与sCT的HU值差异图Fig.6 MRI,CT and sCT images and difference map of HU value difference between CT and sCT
图7 病例的剂量分布差异图Fig.7 Dose distribution difference map of an example
图8 不同剂量阈值区域内的平均剂量差异和平均伽马通过率Fig.8 Mean dose differences and mean gamma pass rates within different dose threshold regions
表1 测试集中所有直肠癌患者的剂量差异和伽马通过率(± s[最小值,最大值])Tab.1 Dose differences and gamma pass rates for all 10 test patients with rectal cancer(Mean±SD[min,max])
表1 测试集中所有直肠癌患者的剂量差异和伽马通过率(± s[最小值,最大值])Tab.1 Dose differences and gamma pass rates for all 10 test patients with rectal cancer(Mean±SD[min,max])
伽马通过率/%感兴趣区域剂量差异/%γ3%,3mm γ2%,2mm γ1%,1mm D>10%D>50%D>90%0.49±0.08[0.40,0.65]0.35±0.12[0.17,0.55]0.34±0.16[0.15,0.60]99.83±0.21[99.30,100.00]99.99±0.03[99.90,100.00]100.00±0.00[100.00,100.00]98.84±0.54[98.00,99.80]99.93±0.15[99.60,100.00]99.98±0.06[99.80,100.00]95.79±1.23[94.30,98.00]98.82±0.64[98.00,99.90]98.89±1.42[95.50,100.00]
与传统基于相同患者的配准CT 和MRI 图像进行训练的研究不同,本文采用不同患者的非配对数据来训练CycleGAN 网络,研究了从盆腔的T1-MRI生成sCT 图像的可行性和精确性。使用CycleGAN消除了MRI-CT 多模态图像配准的不确定性,也避免了Wolterink等[26]提到的相同患者的配对图像可能对网络训练带来的负面影响。此外,与真实CT 图像对比,sCT 精度和剂量分布分析结果都表明生成的sCT图像符合临床应用要求—在3%/3 mm 标准下,伽马通过率须达到90%,严格一些会要求达到95%;根据美国物理师协会(American Association of Physicists in Medicine,AAPM)TG-51 号报告[35],总剂量差应控制在±5%以内。
由于MRI 和CT 图像之间扫描日期和摆位不同,空腔位置和尺寸不可避免地出现了不一致。为了排除由于原始图像空腔不一致对于模型以及sCT 评估的影响,将像素值在-600 HU 以下的体素默认为背景和空腔,排除在感兴趣区域以外,不参与到MAE 和ME值的计算。这样能排除空腔影响,计算出的MAE和ME 值才更合理有效。最后得到的盆腔部位的平均MAE 小于36 HU,在此参考标准上优于之前其他学者报告的结果。虽然采用的测试数据不同,但是也具有一定的参考价值[36-38]。然而,体内空腔确实有可能对剂量分布的合理分析产生关键影响。
与基于CT 图像计算得到的剂量分布相比,基于sCT 图像计算得到的剂量分布的平均剂量差异在处方剂量的±0.5%以内,根据Persson 等[39]研究,这对整个放射治疗工作流程中的总不确定性的贡献可以忽略。从剂量差异图上也可以看出,空腔中的不匹配对剂量计算的影响很小,部分原因是沿射束路径的空腔数量少且尺寸小。此外,在皮肤边界上呈现了较大的剂量差异(3%~5%),这主要归因于MRI和CT图像是在不同摆位下进行扫描的,导致MRI 图像生成的sCT 图像显示出与对应CT 图像不同的皮肤边界。即使对sCT 和CT 图像进行配准,两者的皮肤边界也不可能完全匹配。并且在皮肤边界上呈现的较大剂量差异对整体平均剂量差异影响并不大。这一方面是由于呈现较大剂量差异的网格点数量很小;另一方面是由于光子放射治疗本身的特征—高剂量区域发生在皮肤以下1.5 cm处,而皮肤边界处所受剂量较小。所以尽管皮肤边界存在轻微的剂量误差,但肿瘤剂量几乎不受影响,整体的平均剂量差异依然小于0.5%。在1%/1 mm、2%/2 mm 和3%/3 mm 标准的平均伽马通过率分别大于95%、98%和99%,进一步证明了使用该模型生成sCT 图像进行放疗计划剂量计算时,其剂量精度达到了临床要求。
本研究只使用了不同患者非配对的单一传统序列MRI图像和CT图像来训练模型并生成精确的sCT图像。根据之前的研究[38],使用多序列MRI 图像作为深度学习网络的输入,可以提高sCT的质量。但多序列MRI 图像需要更长的采集时间,这会导致更严重的运动伪影。且临床中使用的其他MRI 序列的视野很小,这可能导致生成的sCT 不够完整。因此,未来验证多序列MRI 图像的使用对提高sCT 精度的效果时,还必须慎重考虑提高sCT准确性和上述问题之间的平衡。
该研究是基于诊断用的影像科MRI 和CT 图像进行的,扫描摆位差异、MRI 图像的视野不足以及较大的图像厚度都在一定程度上限制了sCT 生成模型的性能。在未来的工作中,将研究用放疗流程中扫描的图像来进一步研究该模型的性能,以及该模型在其他身体部位的通用性。
本文首次使用不同病人的非配对MRI 和CT 图像来研究CycleGAN 模型,验证了CycleGAN 模型从盆腔T1-MRI生成sCT 的可行性和准确性。训练集中不需要任何形式的配准,MRI 和CT 图像也都来自不同的患者,从而消除了多模态图像配准误差对模型性能和sCT 图像精度的影响以及相同患者的配对图像可能对模型训练造成的影响。研究结果也表明,该模型基于非配对数据训练可以从盆腔单一传统MRI图像生成精确的sCT图像,基于该sCT图像进行外部光子放疗的剂量计算也完全满足临床应用的要求。未来将展开该sCT 图像生成模型在其他身体部位的通用性研究,为临床MRI-only 放射治疗工作提供必需的基础。