-
X射线计算机断层扫描(Computed Tomography, CT)在医学检查和诊断中被广泛应用的同时,CT扫描时产生的电离辐射也给患者带来了潜在的健康风险。研究表明,过量的X射线照射会给患者带来诱发遗传病和癌症的风险[1-2]。因此,减少CT的照射剂量已成为医学影像方向的一个主要的研究热点。降低X射线照射强度的低剂量CT(Low-Dose CT, LDCT)成像方法逐步成为减少CT照射剂量的主要方法之一。该方法减少了用于X射线成像的光子数量,但增加了重建图像的噪声和伪影,降低了CT重建图像的质量。近年来,随着人工智能技术的发展,越来越多的学者致力于采用深度学习技术来减少LDCT成像方法在CT重建图像中引入的噪声和伪影[3-4]。
将深度学习技术应用于LDCT重建,显著地提高了LDCT图像的成像质量。但受深度神经网络黑箱特性的影响,深度学习技术与LDCT重建方法相结合仍然存在一些挑战,即重建网络的健壮性和泛化性[5]。究其原因,大多数的深度学习技术都是纯数据驱动的,很难解释其底层网络结构的工作机理。为了克服传统深度学习方法缺乏可解释性的缺点,文献[6-8]通过展开基于模型的迭代算法构造深度神经网络。通过对潜在的物理成像过程进行建模,传统的迭代算法通常具有较高的可解释性。因而,通过迭代展开方法构造神经网络结构,增加了所构造深度神经网络的可解释性,进而改善了深度学习技术的健壮性和泛化性。
将迭代重建算法中各次迭代的计算步骤映射为单个的浅层神经网络层,并将有限数量的神经网络层堆叠在一起,从而实现了迭代重建算法的深度神经网络迭代展开。由于迭代重建算法具有很强的可解释性,因而相应的迭代展开重建网络也具有一定的可解释性,提高了所构造深度神经网络的健壮性和泛化性。文献[9]采用近端原始–对偶优化方法求解优化问题得到CT迭代重建算法,并使用卷积神经网络替换迭代算法中的邻近算子,从而实现了迭代展开重建网络。为了提升网络的特征提取能力,该方法使用多通道的特征图实现相邻两次迭代之间的数据传递。文献[10]利用快速迭代收缩阈值算法求解CT重建的优化问题,并将相应的迭代算法展开为迭代重建网络。文献[11]将基于梯度下降法的迭代重建算法展开为LDCT重建网络。为了提高重建图像的质量,该迭代重建算法使用了专家领域正则项。由于可以通过卷积神经网络替换专家领域正则项[12],该算法利用多个并行的卷积神经网络来代替专家领域正则项。在该算法的基础上,文献[13]采用图卷神经网络增加正则项部分神经网络的特征提取能力。文献[3-4]使用近端向前–向后分裂算法求解LDCT优化问题得到迭代重建算法,并将迭代求解方法展开成LDCT迭代重建网络,以及使用卷积神经网络直接替换优化问题的正则项。文献[14]首先使用半二次分裂算法将CT重建优化问题分解为两个等价的子问题;然后使用共轭梯度算法求解子问题,从而得到迭代重建算法;最后对迭代重建算法进行神经网络展开得到迭代展开重建网络。
相比传统深度学习方法,利用迭代展开方法构造LDCT重建网络,可以有效地改善LDCT重建图像的视觉质量。基于模型的迭代重建算法的目标函数由数据保真项和正则项构成,数据保真项用来描述潜在的物理成像过程,正则项用来挖掘投影域和图像域的先验信息,常见正则项如低秩、稀疏、平滑和非局部自相似等。只有数据保真项与正则项相互作用,才能有效地提升CT迭代重建算法的计算能效。但是,在这些迭代展开方法[3-4, 9-11, 13-14]中,主要针对CT重建优化问题中的数据保真项进行了神经网络展开,使用神经网络直接替换正则项。但使用神经网络替换正则项,同样会降低CT重建网络的可解释性。针对该问题,文献[15-19]尝试采用神经网络对CT迭代重建算法中的数据保真项和正则项进行展开。但是,由于包含正则项的CT迭代重建算法过于复杂,研究人员在迭代求解算法中采用了近似计算,以便于进行神经网络展开。
在临床实践中,医学图像大都是非稀疏的,可以通过稀疏变换来获得其稀疏表示,总变分(Total Variation, TV)正则化约束是一种常用的稀疏表示方法。由于TV正则项对图像边缘敏感且算法复杂度较低,常将TV正则化应用到LDCT图像重建中以便于恢复图像中的解剖结构。本文首先采用原始–对偶算法求解基于TV的LDCT重建问题,得到易于神经网络展开的LDCT迭代重建算法。然后,对该迭代算法进行神经网络展开,得到了基于全变分展开的LDCT重建网络(A Total Variation Prior Unfoldment Based Low Dose CT Reconstruction Network, TV-CTNet)。该方法对LDCT迭代重建算法中的数据保真项和正则项进行神经网络展开,提高了所构造LDCT重建网络的可解释性,进一步改进了LDCT重建图像的视觉质量。
-
TV模型在去除图像噪声的同时,能较好地保持图像的边缘和纹理信息[20]。将TV先验模型应用于迭代重建算法中,能够取得理想的重建效果。基于TV模型的LDCT重建优化问题可以表示成:
$$ \mathop {\min }\limits_x \left\| {{\boldsymbol{Ax}} - {\boldsymbol{g}}} \right\|_2^2 + \lambda {\left\| {\nabla {\boldsymbol{x}}} \right\|_1} $$ (1) 式中,
$ {\boldsymbol{x}} $ 为待重建图像;$ {\boldsymbol{g}} $ 为投影数据;$ {\boldsymbol{A}} $ 为系统矩阵;$ \nabla $ 表示图像的梯度;$ \lambda $ 为权重参数;$ {\left\| {\nabla {\boldsymbol{x}}} \right\|_1} $ 为TV正则项;$ \left\| {{\boldsymbol{Ax}} - {\boldsymbol{g}}} \right\|_2^2 $ 为数据保真项。采用原始–对偶算法[21]求解优化问题(1),从而得到该优化问题的LDCT迭代重建算法:$$ {{\boldsymbol{y}}_{n + 1}} = \frac{{{{\boldsymbol{y}}_n} + \sigma \left( {{\boldsymbol{A}}{{\bar {\boldsymbol{u}}}_n} - {\boldsymbol{g}}} \right)}}{{1 + \sigma }} $$ (2) $$ {{\boldsymbol{z}}_{n + 1}} = \frac{{\lambda \left( {{{\boldsymbol{z}}_n} + \sigma \nabla {{\bar {\boldsymbol{u}}}_n}} \right)}}{{\max \left( {\lambda {{\mathbf{1}}_I},\left| {{{\boldsymbol{z}}_n} + \sigma \nabla \bar {\boldsymbol{u}}} \right|} \right)}} $$ (3) $$ {{\boldsymbol{u}}_{n + 1}} = {{\boldsymbol{u}}_n} - \tau {\boldsymbol{A}}^{\rm{T}}{{\boldsymbol{y}}_{n + 1}} - \tau {\nabla ^{\rm{T}}}{{\boldsymbol{z}}_{n + 1}} $$ (4) $$ {\bar {\boldsymbol{u}}_{n + 1}} = {{\boldsymbol{u}}_{n + 1}} + \theta \left( {{{\boldsymbol{u}}_{n + 1}} - {{\boldsymbol{u}}_n}} \right) $$ (5) 式中,
$ \sigma $ 、$ \tau $ 和$ \theta $ 为采用原始–对偶算法求解优化问题所引入的权重参数;$ n $ 为迭代循环变量;$ {{\boldsymbol{y}}_n} $ 和$ {{\boldsymbol{z}}_n} $ 为中间变量;$ {{\mathbf{1}}_I} $ 为元素值都为1的矩阵。进一步地推导上述迭代算法,以便得到易于神经网络展开的LDCT迭代重建算法。式(2)两端同时乘以
${\boldsymbol{A}}^{\rm{T}}$ ,令${{\boldsymbol{p}}_n} = {\boldsymbol{A}}^{\rm{T}}{{\boldsymbol{y}}_n}$ ,则式(2)可以改写为:$$ {{\boldsymbol{p}}_{n + 1}} = \frac{{{{\boldsymbol{p}}_n} + \sigma ( {{{\boldsymbol{A}}^{\rm{T}}}\boldsymbol A{{\bar {\boldsymbol{u}}}_n} - {{\boldsymbol{A}}^{\rm{T}}}{\boldsymbol{g}}} )}}{{1 + \sigma }} $$ (6) 从式(6)中可以看出,变量
$ {{\boldsymbol{p}}_{n + 1}} $ 为变量$ {{\boldsymbol{p}}_n} $ 、${{\boldsymbol{A}}^{\rm{T}}}{\boldsymbol{A}}{\bar {\boldsymbol{u}}_n}$ 和$ {{\boldsymbol{A}}^{\rm{T}}}\boldsymbol g $ 的加权和。式(3)两端同时乘以$ {\nabla ^{\rm{T}}} $ ,令$ {{\boldsymbol{q}}_n} = {\nabla ^{\rm{T}}}{{\boldsymbol{z}}_n} $ ,则式(3)表示为:$$ {{\boldsymbol{q}}_{n + 1}} = \frac{{\lambda ( {{{\boldsymbol{q}}_n} + {\nabla ^{\rm{T}}}\nabla {{\bar {\boldsymbol{u}}}_n}} )}}{{\max \Big( {\lambda {{\mathbf{1}}_I},\left| {{{\left( {{\nabla ^{\rm{T}}}} \right)}^{ - 1}}{{\boldsymbol{q}}_n} + \sigma \nabla {{\bar {\boldsymbol{u}}}_n}} \right|} \Big)}} $$ (7) 令
$\eta = \max \Big( {\lambda {{\mathbf{1}}_I},\Big| {{{( {{\nabla ^{\rm{T}}}} )}^{ - 1}}{{\boldsymbol{q}}_n} + \sigma \nabla {{\bar {\boldsymbol{u}}}_n}} \Big|} \Big)$ ,则式(7)可以改写为:$$ {{\boldsymbol{q}}_{n + 1}} = \lambda ( {{{\boldsymbol{q}}_n} + {\nabla ^{\rm{T}}}\nabla {{\bar {\boldsymbol{u}}}_n}} )/\eta $$ (8) 从式(8)中可以看出,变量
$ {{\boldsymbol{q}}_{n + 1}} $ 为变量$ {\nabla ^{\rm{T}}}\nabla {\bar {\boldsymbol{u}}_n} $ 和$ {{\boldsymbol{q}}_n} $ 的加权和。进行上述变量代换后,式(4)相应地变换为:$$ {{\boldsymbol{u}}_{n + 1}} = {{\boldsymbol{u}}_n} - \tau {{\boldsymbol{p}}_{n + 1}} - \tau {{\boldsymbol{q}}_{n + 1}} $$ (9) 上述公式中给出了式(1)中数据保真项和TV正则项的迭代计算步骤。从这些公式中可以看出,迭代算法的各计算步骤仅涉及变量的加权求和,便于采用神经网络进行展开。
-
在迭代重建算法中,优化问题的目标函数通常由数据保真项和正则项构成。数据保真项被用于对潜在的物理过程进行建模,而正则项被用来捕获图像域和投影域中数据的先验知识。为了获取先验信息,采用各种精心设计的正则项来估计投影数据或重建图像的统计分布特征,包括压缩感知、低秩、非局部块相似、TV、字典学习、稀疏变换学习等。在目标函数中使用了这些手动设计的正则项,显著地改善了迭代重建算法的重建效果。然而,现有的迭代展开重建算法大都只针对数据保真项进行了展开,直接使用神经网络替换了正则项部分,降低了这一部分神经网络的可解释性,不可避免地限制了迭代展开重建网络的性能。对TV正则项部分的迭代算法进行神经网络展开,可以进一步地改善整体神经网络的特征提取能力[22-23]。
对基于TV的LDCT迭代重建算法进行神经网络展开,所得LDCT重建网络的网络结构如图1所示。分别采用4个浅层的神经网络
$ {\rm{CN}}{{\rm{N}}_{\boldsymbol{p}}} $ 、$ {\rm{CN}}{{\rm{N}}_{\boldsymbol{q}}} $ 、$ {\rm{CN}}{{\rm{N}}_{\boldsymbol{u}}} $ 和$ {\mathrm{CN}}{{\mathrm{N}}_{\bar {\boldsymbol{u}}}} $ 替换迭代重建算法中的4个变量的计算步骤,即式(6)、式(8)、式(9)和式(5)。利用子网络$ {\rm{CN}}{{\rm{N}}_{\boldsymbol{p}}} $ 替换式(6)的计算步骤,子网络$ {\rm{CN}}{{\rm{N}}_{\boldsymbol{p}}} $ 的输入为$ {{\boldsymbol{p}}_n} $ 、$ {{\boldsymbol{A}}^{\rm{T}}}{\boldsymbol{A}}{\bar {\boldsymbol{u}}_n} $ 和$ {{\boldsymbol{A}}^{\rm{T}}}{\boldsymbol{g}} $ ,输出为$ {{\boldsymbol{p}}_{n + 1}} $ 。采用固定参数的网络层$ {\rm{CN}}{{\rm{N}}_{\boldsymbol{A}}} $ 和$ {\rm{CN}}{{\rm{N}}_{{{\boldsymbol{A}}^{\rm{T}}}}} $ 来实现正投影操作$ {\boldsymbol{A}} $ 和反投影操作$ {{\boldsymbol{A}}^{\rm{T}}} $ [24]。分别采用射线驱动和像素驱动的方法来计算这两个网络层中的正投影操作$ \boldsymbol A $ 和反投影操作$ {{\boldsymbol{A}}^{\rm{T}}} $ 。如图1所示,子网络$ {\rm{CN}}{{\rm{N}}_{\boldsymbol{p}}} $ 包含了3个卷积层,卷积模板大小为$ 5 \times 5 $ 。每个卷积操作后都有一个批量归一化层(Batch Normalization, BN)和一个Tanh激活函数以提高网络特征提取能力。利用子网络$ {\mathrm{CN}}{{\mathrm{N}}_{\boldsymbol{q}}} $ 替换式(8)的计算步骤。子网络$ {\rm{CN}}{{\rm{N}}_{\boldsymbol{q}}} $ 的输入为变量$ {{\boldsymbol{q}}_n} $ 和$ {\nabla ^{\rm{T}}}\nabla {\bar {\boldsymbol{u}}_n} $ ,输出为$ {{\boldsymbol{q}}_{n + 1}} $ 。分别采用固定参数的卷积层$ {\mathrm{CN}}{{\mathrm{N}}_\nabla } $ 和$ {\mathrm{CN}}{{\mathrm{N}}_{{\nabla ^{\rm{T}}}}} $ 来实现操作$ \nabla $ 和操作$ {\nabla ^{\rm{T}}} $ 。通过卷积层$ {\rm{CN}}{{\rm{N}}_\nabla } $ ,计算得到输入图像$ {\bar {\boldsymbol{u}}_n} $ 在$ x $ 方向和$ y $ 方向的梯度分布。将图像$ {\bar {\boldsymbol{u}}_n} $ 的梯度分布输入到卷积层$ {\rm{CN}}{{\rm{N}}_{{\nabla ^{\rm{T}}}}} $ ,通过该网络层的计算,得到图像$ {\bar {\boldsymbol{u}}_n} $ 的散度分布,即$ {\nabla ^{\rm{T}}}\nabla {\bar {\boldsymbol{u}}_n} $ 。如图1所示,子网络$ {\rm{CN}}{{\rm{N}}_{\boldsymbol{q}}} $ 包含了3个卷积层,卷积模板大小为$ 5 \times 5 $ 。每个卷积操作后都有一个BN层和一个Tanh激活函数。进一步地,通过深度学习从大数据中学习得到式(8)中参数$ \eta $ 的数值,改善LDCT重建网络的计算精度。利用子网络$ {\rm{CN}}{{\rm{N}}_{\boldsymbol{u}}} $ 替换式(9)的计算步骤,子网络$ {\rm{CN}}{{\rm{N}}_{\boldsymbol{u}}} $ 的输入为$ {{\boldsymbol{u}}_n} $ 、$ {{\boldsymbol{p}}_{n + 1}} $ 和$ {{\boldsymbol{q}}_{n + 1}} $ ,输出为$ {{\boldsymbol{u}}_{n + 1}} $ 。如图1所示,子网络$ {\rm{CN}}{{\rm{N}}_{\boldsymbol{u}}} $ 包含了3个卷积层,卷积模板大小为$ 5 \times 5 $ 。每个卷积操作后都有一个BN层和一个线性修正单元(Rectified Linear Unit, ReLU)作为激活函数。利用子网络$ {\rm{CN}}{{\rm{N}}_{\bar {\boldsymbol{u}}}} $ 替换式(5)的计算步骤,子网络$ {\rm{CN}}{{\rm{N}}_{\bar {\boldsymbol{u}}}} $ 的输入为$ {{\boldsymbol{u}}_n} $ 和$ {{\boldsymbol{u}}_{n + 1}} $ ,输出为$ {\bar {\boldsymbol{u}}_{n + 1}} $ 。如图1所示,子网络$ {\rm{CN}}{{\rm{N}}_{\bar {\boldsymbol{u}}}} $ 包含了3个卷积层,其中前两个卷积层的卷积模板大小为$ 5 \times 5 $ ,最后一个卷积层的卷积模板大小为$ 1 \times 1 $ 。每个卷积操作后都有一个BN层和ReLU激活函数。本文提出的LDCT重建网络为迭代展开方法。将迭代次数为
$ N = 10 $ 的LDCT迭代重建算法进行神经网络展开,得到LDCT重建网络。通过LDCT重建网络,将输入的投影数据$ {\boldsymbol{g}} $ 重建为断层图像$ {\bar {\boldsymbol{u}}_n} $ 。操作$ {\nabla ^{\rm{T}}}\nabla {\bar {\boldsymbol{u}}_n} $ 表示计算图像$ {\bar {\boldsymbol{u}}_n} $ 的散度,由式(8)可知向量$ {{\boldsymbol{q}}_n} $ 表示图像$ {\bar {\boldsymbol{u}}_n} $ 的散度信息。在每次迭代中使用散度信息$ {{\boldsymbol{q}}_n} $ 对变量$ {{\boldsymbol{u}}_n} $ 进行更新。在$ N $ 次迭代后,利用子网络$ {\mathrm{CN}}{{\mathrm{N}}_{\bar {\boldsymbol{q}}}} $ 输出图像$ {\bar {\boldsymbol{u}}_N} $ 的散度信息$ {\bar {\boldsymbol{q}}_N} $ 。如图1所示,子网络$ {\mathrm{CN}}{{\mathrm{N}}_{\bar {\boldsymbol{q}}}} $ 包含了3个卷积层,卷积模板大小为$ 5 \times 5 $ 。每个卷积操作后都有一个BN层和一个Tanh激活函数。为了提高LDCT重建网络的特征提取能力,相邻两次迭代之间变量$ {{\boldsymbol{p}}_n} $ 、$ {{\boldsymbol{q}}_n} $ 、$ {{\boldsymbol{u}}_n} $ 和$ {\bar {\boldsymbol{u}}_n} $ 的通道数分别被设置为10、10、10和1。在LDCT重建网络中,各子网络中特征图的数量如图1所示。在网络的训练过程中,采用
$ {\ell _2} $ 范数作为损失函数来减少LDCT重建网络输出的重建图像$ {\bar {\boldsymbol{u}}_N} $ 和标签图像$ \hat {\boldsymbol{u}} $ 之间的差异,以及减少重建图像$ {\bar {\boldsymbol{u}}_N} $ 的散度$ {\bar {\boldsymbol{q}}_N} $ 与标签图像$ \hat {\boldsymbol{u}} $ 的散度$ \hat {\boldsymbol{q}} $ 之间的差异。因此,本文所设计LDCT重建网络的损失函数可以写为:$$ {\mathrm{loss}} = \left\| {{{\bar {\boldsymbol{u}}}_N} - \hat {\boldsymbol{u}}} \right\|_2^2 + \xi \left\| {{{\bar {\boldsymbol{q}}}_N} - \hat {\boldsymbol{q}}} \right\|_2^2 $$ (10) 式中,
$ \xi $ 为权重参数。本文中,$ \xi $ 值被设置为0.5。 -
本文使用公开数据集“Low-Dose CT Image and Projection Data”[25]来验证提出的LDCT重建网络的有效性。该数据集通过螺旋扫描获得,不适合直接用于扇束CT重建。通过对该数据集提供的CT图像进行正投影,模拟得到扇束投影数据。从数据集中随机选取20位病人的CT图像进行网络的训练和测试,其中18位病人的数据作为训练集,剩余2位病人的数据作为测试集。在扇束投影数据集的生成过程中,射线源到线阵列探测器中心的距离为1068 mm;射线源到旋转中心的距离为595 mm;围绕病人进行等间隔扫描得到360个投影,扫描间隔为
$ {1^ \circ } $ ;线阵列探测器拥有768个等间隔分布的传感器,相邻两传感器间的距离为1 mm;LDCT重建网络生成的CT断层图像大小为$ 512 \times 512 $ ,像素尺寸为0.5859 mm$ \times $ 0.5859 mm。模拟所得投影数据为无噪声的投影数据$ \tilde {\boldsymbol{g}} $ 。为了模拟真实的临床环境,将泊松噪声和电子噪声加入到无噪声的投影数据集,生成不同X射线照射强度下的投影数据[13],其公式为:$$ {\boldsymbol{g}} = \ln \frac{{{I_0}}}{{{\mathrm{Poisson}}\left( {{I_0}\exp \left( { - \tilde {\boldsymbol{g}}} \right)} \right) + {\mathrm{Normal}}\left( {0,{\varepsilon ^2}} \right)}} $$ (11) 式中,
$ {I_0} $ 为入射X射线的光子数;$ {\mathrm{Poisson}} $ 为泊松函数;$ {\mathrm{Normal}} $ 为标准正态分布函数,其方差为$ \varepsilon $ 。依据文献[13]进行参数设置,生成不同剂量率下的投影数据。在正常剂量情况下,入射X射线的光子数被设置为$ 1.0 \times {10^6} $ 。在剂量率为15%、10%和5%这3种情况下,入射X射线的光子数等比例减少,分别被设置为$ 1.5 \times {10^5} $ 、$ 1.0 \times {10^5} $ 和$ 0.5 \times {10^5} $ 。在这3种低剂量率投影数据的模拟情况下,参数$ \varepsilon $ 的值均被设置为$ \sqrt {10} $ 。在深度学习 Pytorch 框架内进行本文算法的训练和测试,采用Adam优化算法训练本文所设计的LDCT重建网络,该网络的初始学习率设置为
$ {10^{ - 4}} $ ,批量大小设置为1。本文算法在工作站上进行训练和验证,其操作系统为Ubuntu20.04 LTS,其硬件参数配置如下:CPU型号为Intel Core i7-9700k 3.6 GHz,内存大小为64 GB,GPU为NVIDIA RTX 309024 GB。为了验证本文算法的有效性和可行性,本文算法与一种传统的CT解析重建算法(FBP[26])和5种最新流行的CT重建网络算法(REDCNN[27]、LPD[9]、Air-Net[28]、FISTA-Net[29]和PD-Net[24])进行了比较,并采用峰值信噪比[30](Peak Signal-to-Noise Ratio, PSNR)、结构相似性[30](Structural Similarity Index Measure, SSIM)和归一化均方误差[30](Normalized Mean Square Error, NMSE)定量地分析不同重建算法的计算结果。 -
在15%剂量的模拟数据集下,不同LDCT重建算法的计算结果如图2所示。为了更好地对比不同重建算法的噪声去除和细节恢复能力,选取参考图像中红色方框位置作为感兴趣区域(Region of Interest, ROI),并将ROI内的图像进行展示,如图3所示。重建图像和ROI内图像的图像质量评价指标均按PSNR、SSIM和NMSE的顺序列于相应图像下方。图4是不同重建算法计算结果与参考图像的残差图像。从图2b和图3b中可以清晰地观察到,FBP算法不能有效地滤除低剂量投影数据中的噪声,导致重建所得LDCT图像中出现了严重的伪影残留和高强度的散斑噪声分布。REDCNN算法是一种经典的残差卷积神经网络降噪算法,采用该方法对FBP算法的重建结果进行去噪,其降噪结果过于平滑(图2c和图3c)。与FBP和REDCNN算法相比,LPD算法提高了重建图像的视觉质量,如图2d和图3d所示。但是,其重建结果仍然存在一些纹理细节的丢失。Air-Net算法采用神经网络直接替换正则项的神经网络展开,导致该算法不能很好地保留图像的边缘特征(图2e和图3e)。FISTA-Net算法可以有效地抑制计算结果中的噪声和伪影。但是,该算法在抑制噪声的同时,会模糊一些局部结构特征(图2f和图3f)。TV-CTNet算法针对TV正则项进行了神经网络展开,计算结果中保留了较为完整的结构和纹理信息(图2和图3)。相比PD-Net算法,TV-CTNet算法保留了更多与参考图像相同的图像纹理特征。
图4为不同重建算法计算结果的残差图像,残差图像的数值大小可以反映出不同算法的重建效果。其中,残差图像中的信息量越大,则表明在重建过程中伪影抑制和细节保留的数量越少。从图4中可以看出,本文算法的伪影抑制和细节恢复数量最多,具有较强的伪影抑制能力和结构恢复能力。表1给出了15%剂量下不同重建方法计算结果的评价指标(PSNR、SSIM和NMSE)的平均值和方差。与重建图像的视觉效果一致,本文算法在所有对比算法中获得了最佳的评价指标。
表 1 15%剂量下不同重建算法计算结果的定量评价结果
评价指标 算法 FBP REDCNN LPD Air-Net FISTA-Net PD-Net TV-CTNet PSNR 14.6623±0.7277 33.6915±0.7931 33.8402±1.0074 33.2282±1.1199 33.8719±0.9114 35.7497±1.0559 35.7880±0.9594 SSIM 0.2384±0.0105 0.8010±0.0280 0.7951±0.0285 0.7903±0.0269 0.8120±0.0263 0.8575±0.0260 0.8632±0.0240 NMSE 0.8023±0.2166 0.0099±0.0022 0.0097±0.0031 0.0111±0.0034 0.0096±0.0026 0.0061±0.0011 0.0061±0.0014 -
在10%剂量的模拟数据集下,不同LDCT重建算法的重建结果如图5所示,其中ROI内的图像在图6中显示。如图5b和图6b所示,在FBP算法的重建结果中可以容易地观察到明显的条纹伪影。REDCNN和Air-Net算法的重建结果存在引入一些假的结构和丢失一些边缘细节的问题(图5c、图5e和图6c、图6e)。LPD算法重建结果的解剖结构保存完整,但细节部分较为模糊(图5d和图6d)。从图5f和图6f中可以观察到,FISTA-Net算法的重建结果中存在一些低频失真。PD-Net算法在纹理保存和噪声抑制方面具有较高的鲁棒性(图5g),但是重建结果中依然存在细微结构的微小变形,如图6中红色圆圈位置所示。与5种先进的深度学习算法相比,本文算法所生成重建图像的视觉质量更高,具有更强的纹理结构保存和噪声抑制能力(图5h和图6h)。
图7给出了不同重建算法在10%剂量下的残差图像。从图中可以观察到,本文提出算法的残差信息量最少,表示结构保留最为完整。表2给出了10%剂量下所有算法重建结果的评价指标(PSNR、SSIM和NMSE)的平均值和方差。与其他对比算法相比,本文TV-CTNet算法所得计算结果的PSNR值最高,NMSE值最低,SSIM值略低于PD-Net算法的评价指标。
表 2 10%剂量下不同重建算法计算结果的定量评价结果
评价指标 算法 FBP REDCNN LPD Air-Net FISTA-Net PD-Net TV-CTNet PSNR 14.4895±0.6726 33.2347±0.8024 33.5149±1.1182 32.9440±1.0957 33.4898±0.8269 35.4082±1.1821 35.4752±1.0151 SSIM 0.2249±0.0106 0.7977±0.0277 0.7931±0.0265 0.7836±0.0244 0.8019±0.0243 0.8590±0.0254 0.8522±0.0258 NMSE 0.8349±0.2304 0.0110±0.0028 0.0105±0.0031 0.0119±0.0037 0.0104±0.0026 0.0067±0.0016 0.0065±0.0012 -
在5%剂量的模拟数据集下,不同LDCT重建算法的计算结果如图8所示,其中ROI内的图像在图9中显示。从图8和图9中可以观察到,由于投影数据中存在过多的噪声信息,FBP算法重建得到的图像中出现了严重的图像退化,存在着许多斑点噪声(图8b和9b)。如图8c、图8e和图9c、图9e所示,REDCNN和Air-Net算法在一定程度上模糊了图像的边缘结构,导致LDCT重建结果过于平滑。LPD和FISTA-Net算法的重建图像(图8d、图8f和图9d、图9f)引入了在一些伪影,降低了临床诊断的准确性。从图8g、图8h和图9g、图9h中可以看出,PD-Net和TV-CTNet算法的计算结果在视觉效果上优于其他重建算法,尤其是本文提出的TV-CTNet算法在红色圆圈处的重建效果与参考图像更为相似。
图10给出了不同重建算法计算结果与参考图像的残差图像。从图10中可以很明显地观察到,本文TV-CTNet算法的重建结果与参考图像的差异最小,在抑制LDCT重建图像噪声的同时保留了更多的纹理细节。表3给出了5%剂量下,不同重建算法计算结果的评价指标(PSNR、SSIM和NMSE)的平均值和方差。从表3中可以观察到,本文TV-CTNet算法取得了最优的PSNR和NMSE评价指标。
表 3 5%剂量下不同重建算法计算结果的定量评价结果
评价指标 算法 FBP REDCNN LPD Air-Net FISTA-Net PD-Net TV-CTNet PSNR 14.0046±0.6856 32.8915±0.8240 33.2776±1.2278 32.6899±0.8237 32.9173±0.8055 34.5460±1.2544 34.9261±0.9706 SSIM 0.1952±0.0106 0.7888±0.0269 0.7770±0.0339 0.7665±0.0270 0.7861±0.0271 0.8450±0.0276 0.8413±0.0271 NMSE 0.9285±0.2417 0.0118±0.0025 0.0112±0.0042 0.0125±0.0033 0.0118±0.0026 0.0081±0.0018 0.0074±0.0016 -
在临床实践中,扫描剂量可能随着诊断要求的变化而发生改变。如果LDCT重建网络需要针对不同的扫描剂量水平重新进行网络训练,则LDCT重建网络在临床上的应用就会失去灵活性。为了验证本文TV-CTNet算法的健壮性,采用10%剂量数据训练的模型来测试不同剂量水平(5%和15%)下的LDCT重建结果,如图11所示,其中ROI内的重建结果如图12所示。同时,重建图像和ROI内图像的质量评价指标(PSNR、SSIM和NMSE)放在相应图像的下方。在图11和图12中,本文TV-CTNet算法分别与FBP和LPD算法进行了比较。从图中可以看出,FBP算法重建结果的视觉质量随着扫描剂量的减小明显变差(图11b1~图11b3和图12b1~图12b3)。LPD算法与本文TV-CTNet算法采用相同的优化问题求解方法,但是其仅利用卷积神经网络求解邻近算子问题,并未对正则项进行有效地展开,图像伪影去除能力有限且细节恢复较差(图11c1~图11c3和图12c1~图12c3)。与LPD算法相比,本文提出TV-CTNet算法有效地降低了重建图像的噪声,所生成的重建图像具有更准确的解剖细节(图11d1~图11d3和图12d1~图12d3)。重建图像的定量分析如表4所示。从表4可以看出,随着扫描剂量的增加,TV-CTNet算法所得重建图像中的质量评价指标逐渐提高。实验结果表明,TV-CTNet算法可以从与训练集投影数据分布略微不同的投影数据中重建得到质量较高的重建结果,进而表明TV-CTNet算法具有一定的健壮性。
表 4 采用10%剂量下训练所得模型分别对不同剂量数据重建结果的定量分析
算法 5%剂量 10%剂量 15%剂量 PSNR SSIM NMSE PSNR SSIM NMSE PSNR SSIM NMSE FBP 14.0046 0.1952 0.9285 14.4895 0.2249 0.8349 14.6623 0.2384 0.8023 LPD 32.2606 0.7455 0.0139 33.5149 0.7931 0.0105 33.6982 0.7838 0.0099 TV-CTNet 34.7350 0.8355 0.0077 35.4752 0.8522 0.0065 35.6091 0.8530 0.0063 -
为了验证迭代展开方法中不同神经网络模块对TV-CTNet算法性能的影响,构造了4种重建网络进行消融实验,即Net 1、Net 2、Net 3和Net 4。表5给出了不同重建网络中所包含的神经网络模块,其中“×”表示未采用神经网络模块替换迭代重建算法中变量的计算步骤,“√”表示采用神经网络模块替换相应变量的计算步骤。在15%剂量的模拟数据集下,对不同重建网络的性能进行了实验验证。图13分别给出了不同重建网络的计算结果。从图13中可以观察到,重建结果的视觉质量随着神经网络模块的增加而持续改善。不同重建网络计算结果的评价指标在表5中显示。与重建结果的视觉质量一致,采用神经网络模块替换迭代重建算法中全部变量的计算步骤所得的TV-CTNet算法取得了最佳的评价指标。消融实验中的实验结果验证了采用神经网络模块替换迭代重建算法中变量计算步骤的有效性。
表 5 消融实验中不同重建网络所得计算结果的定量分析
算法 重建网络 结果 ${\mathrm{CN}}{{\mathrm{N}}_{\boldsymbol{p}}}$ ${\mathrm{CN}}{{\mathrm{N}}_{\boldsymbol{q}}}$ ${\mathrm{CN}}{{\mathrm{N}}_{\boldsymbol{u}}}$ ${\mathrm{CN}}{{\mathrm{N}}_{\bar {\boldsymbol{u}}}}$ PSNR SSIM NMSE Net 1 × × × × 19.0029±0.9224 0.4468±0.0190 0.2952±0.0844 Net 2 √ × × × 27.9718±1.1683 0.6435±0.0251 0.0386±0.0142 Net 3 √ √ × × 30.2834±1.4191 0.7146±0.0266 0.0238±0.0108 Net 4 √ √ √ × 34.2744±1.1466 0.8094±0.0236 0.0089±0.0032 TV-CTNet √ √ √ √ 35.7880±0.9594 0.8632±0.0240 0.0061±0.0014
A Total Variation Prior Unfoldment Based Low Dose CT Reconstruction Network
-
摘要: 针对CT迭代展开重建网络仅对数据保真项进行神经网络展开降低了重建网络计算性能的问题,通过对基于全变分的CT迭代重建算法进行神经网络展开,提出一种对数据保真项和全变分正则项全部进行神经网络展开的重建网络,从而改善了CT重建图像的视觉质量。首先,采用原始–对偶算法求解基于全变分的CT重建问题,得到易于神经网络展开的迭代重建算法。然后,对该迭代重建算法进行神经网络展开,尤其是对正则项部分的算法进行神经网络展开,得到迭代展开CT重建网络。在模拟的低剂量CT数据集上验证了该算法的有效性。实验结果表明,与6种低剂量CT重建算法相比,该算法在抑制低剂量CT图像噪声的同时,很好地保留了图像中的结构和细节纹理。重建图像的定量评价分析显示,该算法取得了良好的峰值信噪比和归一化均方误差指标值,验证了提出的低剂量CT重建算法具有较好的噪声抑制能力和较强的鲁棒性。Abstract: Without unrolling the prior terms, most unrolling approaches for Computed Tomography (CT) reconstruction primarily unroll the fidelity term of iterative reconstruction methods to neural networks, which may reduce the computational efficiency of the reconstruction network. To overcome this drawback, a new CT reconstruction network is formed by unrolling a CT iterative reconstruction algorithm based on total variation, especially for the unfoldment of the total variation prior. The unfoldment of the prior improves the visual quality of CT reconstructed images. Firstly, the primal-dual algorithm is utilized to solve the CT reconstruction problem based on the total variation prior, to obtain an iterative reconstruction algorithm which can be easily unrolled to the neural network. Then, the unrolling approach for CT reconstruction is obtained by unrolling this iterative reconstruction algorithm. The effectiveness of the proposed algorithm is tested on a simulated low dose CT dataset. The experimental results show that, compared with six kinds of low-dose CT reconstruction algorithms, the new algorithm effectively preserves the structure and texture details of the image while removing noise in low-dose CT images. The quantitative analysis shows that the proposed algorithm scored the highest PSNR and the lowest NMSE, which indicates that the proposed algorithm is good at noise suppression in the low-dose CT reconstruction.
-
表 1 15%剂量下不同重建算法计算结果的定量评价结果
评价指标 算法 FBP REDCNN LPD Air-Net FISTA-Net PD-Net TV-CTNet PSNR 14.6623±0.7277 33.6915±0.7931 33.8402±1.0074 33.2282±1.1199 33.8719±0.9114 35.7497±1.0559 35.7880±0.9594 SSIM 0.2384±0.0105 0.8010±0.0280 0.7951±0.0285 0.7903±0.0269 0.8120±0.0263 0.8575±0.0260 0.8632±0.0240 NMSE 0.8023±0.2166 0.0099±0.0022 0.0097±0.0031 0.0111±0.0034 0.0096±0.0026 0.0061±0.0011 0.0061±0.0014 表 2 10%剂量下不同重建算法计算结果的定量评价结果
评价指标 算法 FBP REDCNN LPD Air-Net FISTA-Net PD-Net TV-CTNet PSNR 14.4895±0.6726 33.2347±0.8024 33.5149±1.1182 32.9440±1.0957 33.4898±0.8269 35.4082±1.1821 35.4752±1.0151 SSIM 0.2249±0.0106 0.7977±0.0277 0.7931±0.0265 0.7836±0.0244 0.8019±0.0243 0.8590±0.0254 0.8522±0.0258 NMSE 0.8349±0.2304 0.0110±0.0028 0.0105±0.0031 0.0119±0.0037 0.0104±0.0026 0.0067±0.0016 0.0065±0.0012 表 3 5%剂量下不同重建算法计算结果的定量评价结果
评价指标 算法 FBP REDCNN LPD Air-Net FISTA-Net PD-Net TV-CTNet PSNR 14.0046±0.6856 32.8915±0.8240 33.2776±1.2278 32.6899±0.8237 32.9173±0.8055 34.5460±1.2544 34.9261±0.9706 SSIM 0.1952±0.0106 0.7888±0.0269 0.7770±0.0339 0.7665±0.0270 0.7861±0.0271 0.8450±0.0276 0.8413±0.0271 NMSE 0.9285±0.2417 0.0118±0.0025 0.0112±0.0042 0.0125±0.0033 0.0118±0.0026 0.0081±0.0018 0.0074±0.0016 表 4 采用10%剂量下训练所得模型分别对不同剂量数据重建结果的定量分析
算法 5%剂量 10%剂量 15%剂量 PSNR SSIM NMSE PSNR SSIM NMSE PSNR SSIM NMSE FBP 14.0046 0.1952 0.9285 14.4895 0.2249 0.8349 14.6623 0.2384 0.8023 LPD 32.2606 0.7455 0.0139 33.5149 0.7931 0.0105 33.6982 0.7838 0.0099 TV-CTNet 34.7350 0.8355 0.0077 35.4752 0.8522 0.0065 35.6091 0.8530 0.0063 表 5 消融实验中不同重建网络所得计算结果的定量分析
算法 重建网络 结果 ${\mathrm{CN}}{{\mathrm{N}}_{\boldsymbol{p}}}$ ${\mathrm{CN}}{{\mathrm{N}}_{\boldsymbol{q}}}$ ${\mathrm{CN}}{{\mathrm{N}}_{\boldsymbol{u}}}$ ${\mathrm{CN}}{{\mathrm{N}}_{\bar {\boldsymbol{u}}}}$ PSNR SSIM NMSE Net 1 × × × × 19.0029±0.9224 0.4468±0.0190 0.2952±0.0844 Net 2 √ × × × 27.9718±1.1683 0.6435±0.0251 0.0386±0.0142 Net 3 √ √ × × 30.2834±1.4191 0.7146±0.0266 0.0238±0.0108 Net 4 √ √ √ × 34.2744±1.1466 0.8094±0.0236 0.0089±0.0032 TV-CTNet √ √ √ √ 35.7880±0.9594 0.8632±0.0240 0.0061±0.0014 -
[1] LI M, WANG J, CHEN Y, et al. Low-dose CT image synthesis for domain adaptation imaging using a generative adversarial network with noise encoding transfer learning[J]. IEEE Transactions on Medical Imaging, 2023, 42(9): 2616-2630. doi: 10.1109/TMI.2023.3261822 [2] ZHANG X, HAN Z, SHANGGUAN H, et al. Artifact and detail attention generative adversarial networks for low-dose CT denoising[J]. IEEE Transactions on Medical Imaging, 2023, 40(12): 3901-3918. [3] CHEN M, PU Y, BAI Y. Low-dose CT image denoising using residual convolutional network with fractional TV loss[J]. Neurocomputing, 2021, 452: 510-520. doi: 10.1016/j.neucom.2020.10.004 [4] DING Q, NAN Y, GAO H, et al. Deep learning with adaptive hyper-parameters for low-dose CT image reconstruction[J]. IEEE Transaction on Computational Imaging, 2021, 7: 648-660. doi: 10.1109/TCI.2021.3093003 [5] 王革. X射线成像和深度学习的交叉融合[J]. CT 理论与应用研究, 2022, 31(1): 1-12. WANG G. X-ray imaging meets deep learning[J]. CT Theory and Applications, 2022, 31(1): 1-12. [6] MONGA V, LI Y, ELDAR Y C. Algorithm unrolling: Interpretable, efficient deep learning for signal and image processing[J]. IEEE Signal Processing Magazine, 2021, 38(2): 18-44. doi: 10.1109/MSP.2020.3016905 [7] XIA W, SHAN H, WANG G, et al. Physics-/model-based and data-driven methods for low-dose computed tomography: A survey[J]. IEEE Signal Processing Magazine, 2023, 40(2): 89-100. doi: 10.1109/MSP.2022.3204407 [8] ZHANG J, CHEN B, XIONG R, et al. Physics-inspired compressive sensing: Beyond deep unrolling[J]. IEEE Signal Processing Magazine, 2023, 40(1): 58-72. doi: 10.1109/MSP.2022.3208394 [9] ALDER J, ÖKETM O. Learned primal-dual reconstruction[J]. IEEE Transactions on Medical Imaging, 2018, 37(6): 1322-1332. doi: 10.1109/TMI.2018.2799231 [10] XIANG J, DONG Y, YANG Y. FISTA-Net: Learning a fast iterative shrinkage thresholding network for inverse problems in imaging[J]. IEEE Transactions on Medical Imaging, 2021, 40(5): 1329-1339. doi: 10.1109/TMI.2021.3054167 [11] CHEN H, ZHANG Y, CHEN Y, et al. LEARN: Learned experts’ assessment-based reconstruction network for sparse-data CT[J]. IEEE Transactions on Medical Imaging, 2018, 37(6): 1333-1347. doi: 10.1109/TMI.2018.2805692 [12] YANG H, LU J, ZHANG H, et al. Field of experts regularized nonlocal low rank matrix approximation for image denoising[J]. Jouranl of Computational and Applied Mathematics, 2022, 412: 114244. doi: 10.1016/j.cam.2022.114244 [13] XIA W, LU Z, HUANG Y, et al. MAGIC: Manifold and graph integrative convolutional network for low-dose CT reconstruction[J]. IEEE Transactions on Medical Imaging, 2021, 40(12): 3459-3472. doi: 10.1109/TMI.2021.3088344 [14] ZHANG H, LIU B, YU H, et al. MetaInv-Net: Meta inversion network for sparse view CT image reconstruction[J]. IEEE Transactions on Medical Imaging, 2021, 40(2): 621-634. doi: 10.1109/TMI.2020.3033541 [15] WU W, HU D, NIU C, et al. DRONE: Dual-domain residual-based optimization network for sparse-view CT reconstruction[J]. IEEE Transactions on Medical Imaging, 2021, 40(11): 3002-3014. doi: 10.1109/TMI.2021.3078067 [16] WU W, GUO X, CHEN Y, et al. Deep embedding- attention-refinement for sparse-view CT reconstruction[J]. IEEE Transactions on Instrumentation and Measurement, 2023, 72: 1-11. [17] HU D, ZHANG Y, LIU J, et al. DIOR: Deep iterative optimization-based residual-learning for limited-angle CT reconstruction[J]. IEEE Transactions on Medical Imaging, 2022, 41(7): 1778-1790. doi: 10.1109/TMI.2022.3148110 [18] ZHANG Y, HU D, HAO S, et al. DREAM-Net: Deep residual error iterative minimization network for sparse-view CT reconstruction[J]. IEEE Journal of Biomedical and Health Informatics, 2023, 27(1): 480-491. doi: 10.1109/JBHI.2022.3225697 [19] SU T, CUI Z, YANG J, et al. Generalized deep iterative reconstruction for sparse-view CT imaging[J]. Physics in Medicine and Biology, 2022, 67(2): 025005. doi: 10.1088/1361-6560/ac3eae [20] ZHANG Z, CHEN B, XIA D, et al. Directional-TV algorithm for image reconstruction from limited-angular -range data[J]. Medical Image Analysis, 2021, 70: 102030. doi: 10.1016/j.media.2021.102030 [21] SIDKY E Y, JORGENSEN H J, PAN X C. Convex optimization problem prototyping for image reconstruction in computed tomography with the Chambolle–Pock algorithm[J]. Physics in Medicine and Biology, 2012, 57: 3065-3091. doi: 10.1088/0031-9155/57/10/3065 [22] JIU M, PUSTELNIK N. A deep primal-dual proximal network for image restoration[J]. IEEE Journal of Selected Topics in Signal Processing, 2021, 15(2): 190-203. doi: 10.1109/JSTSP.2021.3054506 [23] JIU M, PUSTELNIK N. Alternative design of DeepPDNet in the context of image restoration[J]. IEEE Signal Processing Letters, 2022, 29: 932-936. doi: 10.1109/LSP.2022.3153229 [24] ZHANG P, REN S, LIU Y, et al. A total variation prior unrolling approach for computed tomography reconstruction[J]. Med Phys, 2023, 50(5): 2816-2834. doi: 10.1002/mp.16307 [25] MOEN T R, CHEN B, HOLMES D R, et al. Low-dose CT image and projection dataset[J]. Med Phys, 2021, 48(2): 902-911. doi: 10.1002/mp.14594 [26] SHEPP L A, LOGAN B F. The Fourier reconstruction of a head section[J]. IEEE Trans Nucl Sci, 1974, 21(3): 218-227. [27] SIDKY E Y, PAN X. Image reconstruction in circular cone-beam computed tomography by constrained, total-variation minimization[J]. Physics in Medicine and Biology, 2008, 53(17): 4777-4807. doi: 10.1088/0031-9155/53/17/021 [28] CHEN H, ZHANG Y, KALRA M K, et al. Low-dose CT with a residual encoder-decoder convolutional neural network[J]. IEEE Transactions on Medical Imaging, 2017, 36(12): 2524-2535. doi: 10.1109/TMI.2017.2715284 [29] CHEN G, HONG X, DING Q, et al. AirNet: Fused analytical and iterative reconstruction with deep neural network regularization for sparse-data CT[J]. Med Phys, 2020, 47(7): 2916-2930. doi: 10.1002/mp.14170 [30] JIAO F, GUI Z, LI K, et al. A dual-domain CNN-based network for CT reconstruction[J]. IEEE Access, 2021, 9: 71091-71103.