面向多模态MRI脑胶质瘤区域三维分割与生存期预测的级联U-Net网络

余力，刘宵雪，闫朝阳，李建瑞，张志强，黄韫栀，徐军

1. 南京信息工程大学人工智能学院医学人工智能联合研究院, 南京 210044； 2.东部战区总医院影像科， 南京 210002

0 引 言

脑胶质瘤是最常见、死亡率最高的原发性脑肿瘤。根据世界卫生组织的分类，脑胶质瘤按照其恶性程度可以分为Ⅱ—Ⅳ级(Louis等，2016)。低级别脑胶质瘤，如星型细胞瘤和少突胶质瘤，10年生存率为57%(Glas等，2009)。然而高级别脑胶质瘤，即Ⅳ级的胶质母细胞瘤，5年生存率仅为5%(Saunders等，2021)。因此，早期发现以及良好的预后对于患者的健康是至关重要的。

磁共振图像具有软组织对比度高、空间分辨率高等优点，广泛应用于评价肿瘤的异质性(Leach等，2005)。脑胶质瘤具有不同的组织学子区域，包括肿瘤周围水肿/侵入组织、坏死核心以及增强和非增强肿瘤核心。脑胶质瘤的这种固有异质性也体现在其成像表型上，核磁共振扫描成像所显示的肿瘤亚结构是由不同强度分布组成的，反映了不同的肿瘤生物学特征。在对异常肿瘤组织进行准确检测和分割的基础上，进行诊断分级以及生存预测，进而制订不同的治疗方案，将会极大地改善患者的生活质量。然而，手动分割肿瘤亚区并对其进行分析诊断十分耗时，并且需要经验丰富的神经放射科医生，此外，不同医生间的主观差异同样不可避免。

弥散性脑胶质瘤的显著异质性意味着需要有效的风险分级以便对不同患者制订最为适合的治疗方案。许多研究开发了基于分子亚型的分类器，为患者的生存分析提供了有效办法(Patel等，2014)。放射组学是一个快速发展的领域，旨在通过提取广泛的定量特征来挖掘医学影像数据(Aerts等，2014)。基于放射组学的分析结果对于患者预后越来越重要。传统的基于机器学习的方法，如支持向量机(support vector machine，SVM)、K近邻算法(k-nearest neighbor，KNN)和随机森林(random forest, RF)通常用于脑肿瘤分析(Kaur等，2019；Reza等，2019)，然而这些方法在建模阶段都存在特征不充分、不完全的共同局限性。

基于深度学习的方法(张巧丽 等，2017；田萱 等，2019；李锵 等，2020；江宗康 等，2020)克服了手工特征设计的缺点，通过建立大规模可训练的模型以学习给定任务所需的最佳特征。Havaei 等人(2017)提出使用卷积神经网络(convolutional neural network, CNN)利用局部和全局特征进行鲁棒性脑肿瘤分割，然而，该算法以2D切片作为输入，因此忽略了核磁共振(magnetic resonance, MR)图像3D空间的上下文信息。Kamnitsas 等人(2017)提出的DeepMedic网络采用并行3D卷积进行多尺度处理，有效结合了局部和上下文信息，大大地提高了分割性能，获得两项竞赛的冠军。然而该网络仅能处理局部图像块，并且分割效率较低。

在BraTS2018(Brain Tumor Segmentation 2018)竞赛中，Myronenko(2019)提出了一种非对称U型网络结构，其中较大的编码器分支用于特征提取，较小的解码器用于标签重建，此外，该方法利用变分自编码器分支来正则化解码器分支，从而提高泛化性能。2019年的BraTS 竞赛中，Jiang等人(2019)提出了一个两阶段的网络，第1阶段使用不对称U型网络来进行粗分割，然后在第2阶段使用类似但更宽的网络来进行二次预测分割。同样地，该方法在第2阶段也增加了一个解码器分支，用于正则化编码器。Clough等人(2019)提出基于变分自编码器(variational auto-encoder, VAE)的分类算法，该算法利用编码器获取的低维隐空间特征进行图像重建和图像分类。可见，在解码器中加入变分自编码器分支不仅可以减少过拟合，提高模型鲁棒性，而且可以发掘数据的潜在有效表示，提高模型性能。

对于特定的肿瘤分割任务，一个具有挑战性的问题是不同肿瘤的发展阶段以及不同病例中肿瘤形态和位置具有可变性。为了提高预测精度，许多分割方法(Wang等，2017；Zhou等，2018)将任务分解为检测和分割两大任务，并使用附加的前向模型进行对象定位。例如，Wang等人(2017)根据肿瘤亚区层次顺序训练了3个网络。Oktay等人(2018)证明，在胰腺肿瘤分割任务中，通过将注意力门(attention gates, AGs)引入标准卷积神经网络框架，可以达到与先检测后分割方法同样的效果。

受到上述工作的启发，本文针对脑胶质瘤提出一套完整的计算机辅助诊断算法。首先，本文算法利用两阶段级联网络对脑胶质瘤图像进行分割以及特征提取，接下来根据上述特征进行总体生存期的预后工作。本文算法采用不对称的U型网络结构作为第1阶段网络，从而获得相对粗糙的分割结果。第2阶段中，将第1阶段网络的初步分割图与MR(magnetic resonance)图像的拼接作为输入，对肿瘤子区域进行更加精确的预测，此外，本文算法中第2阶段包含有变分自编码器分支，不仅可以防止分割模型过拟合，提高泛化性，还可以基于该分支的鲁棒性特征实现患者总体生存期的预后。

1 本文方法

图1展示了本文方法的整体流程。图1(a)包括4种原始MRI(magnetic resonance imaging)模态：T1模态(native), T1Gd模态(post-contrast T1 weighted), T2模态(T2-weighted)以及T2-FLAIR模态(T2 fluid attenuated inversion recovery)；图1(b)展示了多模态图像经过裁剪、z-score标准化等预处理操作后的结果；图1(c)展示了本文所提出的多任务网络，该网络可同时实现脑胶质瘤多个子区域的自动分割任务和分割结果的重建任务；图1(d)展示了该网络的分割结果；图1(e)展示了基于重构任务提取的隐藏层特征应用随机森林算法进行生存期预测的结果。

图1 本文方法整体工作流程Fig.1 Overall workflow of the proposed ((a)original MR images of four modes; (b) preprocessed images; (c) the overall network architecture; (d) segmentation results of glioma;(e) prediction results of survival)

1.1 两级级联网络结构

本文提出了一种带有VAE分支的两级级联的网络结构。第1阶段网络以多模态MR图像为输入，预测粗略分割图。初步分割结果与多模态MR图像进行拼接作为第2级网络的输入，从而改进分割效果。另外，第2阶段添加了VAE分支，起到防止过拟合作用，并且能够提取图像有效的潜在特征进行生存期预测。整体网络结构如图1(c)所示。

首先，第1阶段的网络(图2)包括用于语义特征提取的编码器以及用于肿瘤分割结果预测的解码器。

为了解决随着网络加深而产生梯度消失或梯度爆炸的问题，He等人(2016)提出的ResNet引入了跳跃连接的思想，在残差块的输入和输出之间添加一个恒等映射，有效地缓解深层网络产生的模型退化现象。本文方法同样采用了恒等映射思想。编码器部分由残差块组成，共计4个空间级别，分别包括1，2，2和4个残差块。每个残差块包括两个带有组归一化和ReLU激活函数的3×3×3大小、步长为1的3D卷积，并行添加一个加性跳跃连接。编码器部分的输入是尺寸为4×128×128×128的图像块，由4种模态的裁剪为128×128×128的MR图像块组成。经过32×3×3×3大小的卷积层以及一个随机丢失率为0.2的Dropout层后，再经过上述一系列的残差模块。这里，不同空间层次的残差块通道数分别为32，64，128，256。最终，编码器部分得到256×16×16×16大小的特征图。

解码器部分与编码器部分不同的是在残差块的基础上引入了SE-Net(squeeze-and-excitation networks)(Hu等，2018)的思想，利用3组SE-残差模块以及反卷积将编码器所提取的高维特征恢复至原始图像块分辨率。首先利用1×1×1大小的卷积块将特征图通道数减少一半，同时利用步长为2的3维反卷积将空间大小进行翻倍操作，然后将编码器部分中相同的空间级别的输出特征图经过AGs模块(注意力门)(Oktay等，2018)，进而作为下一SE残差模块的输入。AGs模块生成一个门控信号，用来控制不同空间位置的特征的重要性。在解码器的最后部分使用1×1×1大小的卷积将特征通道数从32映射为4，最后使用Softmax激活函数将特征映射转换为概率层面的分割结果。

第2阶段网络(图3)的输入则是由4种模态的MR图像与第1阶段网络生成的分割图进行拼接构建的，因此，第2阶段网络输入通道数为8(4组通道为4种模态的MR图像，4组通道为第1阶段网络生成的概率分割图)，编码器和解码器其他部分与第1阶段保持一致。

图3 第2阶段网络结构Fig.3 Network structure of the second stage

为了防止模型过拟合，提高模型的泛化性，本文方法在该阶段添加变分自编码器(VAE)分支进行图像重构。在该分支中，首先利用自适应平均池化层以及全连接层将编码器输出降低至256维，用其表示具有正态分布的128个平均值以及128个标准差。然后，从具有该均值和标准差的正态分布中提取样本，按照与解码器相同的策略逐步重构为第1阶段的分割结果的维数。需要注意的是，编码器与VAE分支之间不添加跳跃连接，并且上采样使用的是三线性插值方法。

1.1.1 SE-残差模块

由于不同模态对于不同的肿瘤子区域区分度是不同的，例如T1Gd模态能够更加清晰地反映出肿瘤核心区域以及坏死区域，而T2模态对于水肿区域有着更高的对比度，因此期望网络能够反映出不同的通道之间重要性的差异。传统卷积对于通道维度上的特征融合默认是对于输入特征图的所有通道赋予相同的权重进行融合，而SENet(Hu等，2018)提出了SE(squeeze-and-excitation)模块，能够学习到不同通道间的重要性。首先该模块对卷积得到的特征图进行压缩操作，得到通道级别的全局特征，然后对这个全局特征进行激励操作，学习不同通道间的关系，从而得到不同通道间的权重。这种通道上的注意力机制能够关注信息量更大的通道特征，抑制不重要的通道特征，与不同肿瘤子区域的预测需要关注不同的模态这个想法相吻合。因此本文将连续的两个卷积层替换为带有跳跃连接的SE模块，如图2中SE-残差模块所示。该模块由两个带有组归一化和ReLU激活层的3×3×3大小的卷积块、一个SE模块以及一个跳跃连接组成。由于本文数据以及网络较大，批大小仅能设为1，在批大小较小时，组归一化相较批归一化有更好的性能，因此选择组归一化替代了常用的批归一化。

1.1.2 attention gate(AG)模块

当预测目标在形状和大小方面表现出较大的患者间差异时，先检测再分割的多级级联卷积神经网络结构(Zhou等，2018)能够取得较好的效果。然而，这种方法会导致计算资源浪费和模型参数冗余。为了解决这个问题，attention U-Net(Oktay等，2018))提出了注意力门的方法,在不引入大量参数和计算量的前提下，可以抑制无关区域的特征来提高模型灵敏度和准确性。因此，本文在传统U-Net的加性跳跃连接处增加了注意力门机制(Oktay等，2018)。

1.2 损失函数

损失函数包含3个部分

L=L1+L2+0.1×L3

(1)

式中，L1、L2分别为第1阶段输出和第2阶段预测的分割与专家标记之间的损失函数，即

L1=Ld+Lf

(2)

L2=Ld+Lf

(3)

Ld表示解码器的输出与真实分割标记之间的广义Dice损失。

(4)

式中，rln代表类别l在第n个位置的真实像素类别，而pln表示相应的预测概率值，ωl表示每个类别的权重，计算公式为

(5)

Lf是在交叉熵损失函数的基础上进行了修改，解决了样本类别失衡的问题，公式为

Lf=(1-pt)γlog(pt)

(6)

式中，pt表示真实标记的预测概率，γ则是调节简单样本权重降低速率的参数，本文默认设为2。

L3=LL2+LKL

(7)

LL2是VAE分支的输出Ip与第1阶段预测的分割图像Ii之间的均方误差损失

(8)

LKL是KL(Kullback-Leibler)散度，作为VAE分支的惩罚项，使得隐藏层特征接近正态分布,即

(9)

式中，N表示图像块中体素的个数，μ为隐藏层中代表均值的128维向量，σ为隐藏层中代表标准差的128维向量。按照Myronenko(2019)的方法，将超参数权重设置为0.1，从而在各项损失之间达到良好的平衡。

1.3 总体生存期预测

患者的信息包括4种模态的MR影像以及患者的年龄。以往研究表明，仅凭患者的年龄就可以相对较好地预测患者的预后，并且优于综合更复杂的影像组学特征的方法(Weninger等，2018；Kofler等，2019)。鉴于Clough等人(2019)利用变分自编码器(VAE)提取的低维潜在向量进行分类任务，本文摒弃传统的影像组学特征，使用VAE分支所提取的高级语义特征进行总体生存期的预测。如图4所示，首先，利用VAE提取分割后的肿瘤特征以及患者的年龄这一临床特征，二者结合构建总特征集。接下来，利用皮尔森相关系数检验出高度相关的特征，并去除冗余特征。通过基于随机森林回归的递归特征消除方法进一步减少特征的数量。通过对特征的重要性评估，迭代消除不重要的特征。按照10个月与15个月的总生存期将患者生存期划分为短期、中期以及长期3大类。参考Agravat和Raval (2019)的方法，分别使用基于随机森林的回归器对总体生存期进行预测，使用基于随机森林的分类器对于总体生存期的类别进行预测。当回归模型的预测结果与分类模型的类别不一致时，如果分类模型预测该类别的概率超过80%时，将最终的预测结果按照分类模型的预测类别划分为一个固定的总体生存时间，即短期对应于10个月，中期对应于12.5个月，长期预测生存期对应于15个月。

图4 预测总体生存期流程Fig.4 Pipeline of predicting overall survival period

2 实验结果与分析

2.1 实验数据

2.1.1 数据介绍

本文使用来自MICCAI(Medical Image Computing and Computer Assisted Intervention Society)多模态脑肿瘤分割挑战赛(BraTS 2020) 的实验数据进行脑胶质瘤的分割以及预后工作。训练数据共计369例，包括76例低级别脑胶质瘤患者样本以及293例高级别脑胶质瘤患者，其中肿瘤完全切除状态的患者包括118例。验证集共计125例无标签患者样本，其中肿瘤完全切除状态的患者包括29例。所有模态图像大小均为240×240×155，分辨率重采样为1 mm×1 mm×1 mm。所有用于预测总体生存期的患者均提供了年龄这一临床指标。

2.1.2 数据预处理

由于原始MR图像尺寸较大，本文将其根据脑体积裁剪至128×160×192的大小，既可以使数据去除冗余信息(强度值为0的位置，即非脑部信息)，又可以减小图像块体积，满足计算机显存大小。其次，由于扫描仪中患者位置的不同，扫描仪本身以及许多未知因素将导致MR图像上的强度值差异，这称为偏置场的低频平滑的不良信号，因此在进行分割前，使用SimpleITK(https://simpleitk.org/)进行偏置场校正。最后将执行过偏置场校正的图像块进行z-score标准化处理。最终得到的图像如图1(b)所示。z-score标准化公式为

(10)

式中，x为原始数据，而u和σ则分别为所有原始数据的均值和标准差，计算方法分别为

(11)

(12)

2.1.3 数据增强方法

为了防止过拟合，本文采取了一系列数据增强方法，包括对每个轴应用随机翻转、随机旋转、伽马变换以及弹性形变。每种策略皆以50%的概率实施。

2.2 实验环境及设置

实验环境为：Ubuntu16.04 LST64位操作系统，128 GB内存，显卡为NVIDIA GeForce RTX 2080Ti，CPU E5-2630。

本文网络基于PyTorch进行实施。最大迭代轮数设为300。本文使用Adam优化器进行权重更新，批大小设置为1，初始学习率设为1E-4，并且进行更新

(13)

式中,lr表示当前学习率，lr0表示初始学习率，e是当前迭代轮数，Ne为总迭代轮数。本文还使用了1E-5的L2权重衰减进行正则化。由于显存受限，利用PyTorch中的checkpoint包有效减少了显存占用。

2.3 后处理

结果表明，当预测的增强肿瘤区域较小时，本文算法容易将非增强型肿瘤区域错误预测为增强型肿瘤区域。在后处理中，根据经验，当预测的增强型肿瘤区域小于500个体素时，使用坏死区域替代增强肿瘤区域。

2.4 评估指标

2.4.1 分割任务评估指标

脑胶质瘤分割结果的评估针对3个亚肿瘤区域，包括全肿瘤区域、肿瘤核心区域以及增强型肿瘤区域。肿瘤核心区域是由增强型肿瘤区域以及坏死区域共同组成，而全肿瘤区域则是由肿瘤核心区域以及水肿区域组成。

本文分别以Dice相似系数(Dice similarity coefficient, DSC)以及豪斯多夫距离(Hausdorff distance)来评估预测分割结果的好坏。

其中，DSC计算方式为

(14)

式中，P表示预测结果，T表示真实标记。

豪斯多夫距离则是两个点集之间最大不匹配程度的一种度量，定义为

H(P,T)=max(h(P,T),h(T,P))

(15)

(16)

(17)

式中，H(P,T)称为双向豪斯多夫距离，h(P,T)称为从点集P到点集T的单向豪斯多夫距离，h(T,P)则称为从点集T到点集P的豪斯多夫距离。

2.4.2 生存期预测评估指标

总体生存期以10个月和15个月为界限分为短期、中期以及长期3大类生存期。本文主要以这3大类的准确率以及生存期的均方误差来评估生存期预测的效果。

准确率(accuracy, ACC)是最常见的评价指标，即

(18)

式中，S、M、L代表短期、中期和长期类患者的个数，TS、TM、TL则代表被正确分类的短期、中期和长期类患者的个数。

均方误差(mean sqkare error, MSE)则是衡量各数据偏离真实值的距离平方和的平均数

(19)

2.5 结果与分析

2.5.1 脑肿瘤分割结果

为了验证本文算法的有效性，分别选择以下5个网络结构进行对比：

1)3D U-Net(Çiçek等，2016)。对称的编码解码结构，在相同的空间维度增加了跳跃连接。

2)UNet-VAE(Myronenko，2019)。增加了变分自编码分支，防止过拟合。

3)attention U-Net(Oktay等，2018)。在U-Net的跳跃连接处加入了注意力门机制。

4)attention U-Net-VAE。结合了2)和3)两种方法。

5)two stage attention U-Net: 两阶段级联注意力U-Net。将第1阶段的粗分割结果与输入的多模态MR图像拼接在一起作为第2阶段的输入，进而进行更加精细的分割。

图5展示了不同算法在BraTS 2020代表性病例的轴向面、冠状面和矢状面的肿瘤组织分割结果。其中黄色区域为增强型肿瘤区，红色区域为坏死和非增强型肿瘤区，绿色区域为水肿区域。可以看到各个算法都或多或少地存在将增强型肿瘤区域预测为坏死及非增强型肿瘤区域的现象。相较之下，本文算法很好地解决了这个问题。两阶段级联注意力U-Net出现将水肿区域错误地预测为增强型肿瘤区域，本文算法并未出现这种情况。

图5 不同算法在水平面、矢状面和冠状面的分割效果对比图Fig.5 The comparison of the segmentation effects of different algorithms on the horizontal plane, sagittal plane and coronal plane((a)FLAIR; (b)label; (c) 3D U-Net; (d) U-Net-VAE; (e) attention U-Net; (f) attention U-Net-VAE; (g) two stage attention U-Net; (h)ours)

表1展示了不同算法在BraTS2020验证集上的结果。相较于其他算法，本文算法在3个肿瘤子区域皆取得了最好的结果。可以看到在残差U-Net的基础上添加变分自编码器分支，相较于原始残差U-Net在增强型肿瘤区域的分割结果上并没有明显优势，反而Dice相似系数低了1%，但是在全肿瘤区域以及肿瘤核心区域上取得了更好的结果。由此可见，加入了该分支后，模型存在将增强型肿瘤区域预测为坏死区域的现象。而attention U-Net相较于残差U-Net则是在每个肿瘤子区域都取得了0.4%～0.7%的优势，相较于U-Net-VAE的效果更好。由此可见，在跳跃连接处加入注意力门机制的方法是有效的。通过对比attention U-Net-VAE与U-Net-VAE同样可以发现注意力机制对于分割结果的有效改善。本文算法采用两阶段级联结构，通过对粗略分割的细化，在每个肿瘤子区域都取得了最优的分割结果，对于增强型肿瘤区域及肿瘤核心区域的分割效果得到了较为显著的提升，相较于3D U-Net在这两个区域分别提升了4.3%和5.36%的Dice相似系数。各个分割算法在全肿瘤区域都取得了89%～90%的Dice相似系数，这是因为全肿瘤区域是由增强型肿瘤区域、坏死及非增强型肿瘤区域以及水肿区域共同组成，水肿区域具有明显的边界，因此各个算法在该区域都取得了较好的结果。

表1 各类算法在BraTS2020验证集上的分割效果对比表Table 1 Comparison table of the segmentation effect of various algorithms on the BraTS2020 verification set

2.5.2 总体生存期预测结果

表2展示了本文方法在BraTS2020竞赛验证集上的总体生存期预测结果。本文尝试了在VAE分支提取的特征后添加两层全连接层以进行端到端的预后任务，同样将年龄这一临床特征作为输入，具体结构如图6所示。相较于利用传统机器学习方法，直接使用端到端的方法进行预后，可以发现准确率有所下降。由于显存有限，批大小只能设为1，每次只能训练一个样本，损失振荡剧烈，因此网络较难拟合。并且本文肿瘤完全切除状态的患者(有预后信息)仅118例，在小样本数据集的情况下，传统机器学习有所优势。可以看到使用传统机器学习的方法进行预后时，当仅使用年龄这一临床特征时，准确率可以达到0.414，而结合了本文使用VAE分支所提取的潜在特征，准确率可以达到0.552。另外与竞赛排名第1的Agravat和Raval (2019)方法进行了对比，可以看到本文方法的准确率与其相比低了3.4%，然而从均方误差的角度来衡量，本文的方法略有优势。

表2 各类算法在BraTS2020验证集上的总体生存期预测效果对比表Table 2 Comparison table of the overall survival prediction effect of various algorithms on the BraTS2020 validation set

图6 基于神经网络的端到端预后方法Fig.6 End to end prognosis method based on neural network

以上结果说明本文算法提取的特征具有鲁棒性，并且不需要另外单独计算大量的影像组学特征，大大提高了生存期预测的效率。

3 结 论

本文将基于深度学习和传统机器学习的方法应用于多模态MR脑部影像，实现了脑肿瘤分析中的分割多肿瘤组织区域以及总体生存期预测两大任务。由于变分自编码器既可以提供额外的监督信息，增强模型的鲁棒性，又可以发掘数据的潜在有效表示，本文提出了一种两阶段级联的带有VAE分支的肿瘤分割方法。进而利用VAE分支所提取的潜在有效特征进行总体生存期预测的任务。具体来说，本文基于上述提取的特征，利用基于随机森林的递归特征消除方法进行特征选择，最后利用随机森林基于所选特征进行总体生存期的预测任务。本文的工作集多肿瘤子区域分割与总体生存期预测于一体，可以在临床医生对患者进行诊断时提供有力的参考依据。但是本文模型整体计算量较大，后续工作将在保持模型分割精确度以及所提取特征有效性的基础上，减小模型的计算量。