深度学习（5）——强度非均匀性下图像分割的水平集方法及其在MRI中的应用（下）

C.数值实现

我们的方法的实现很简单。公式（22）和（26）中的水平集演化可以通过使用与[11]中提供的DRLSE相同的有限差分方案来实现。虽然我们使用一个简单的完整域实现来实现本文中提出的水平集方法，但值得指出的是，[11]中提供的DRLSE的窄带实现也可以用于实现所提出的方法，将大大降低计算成本并使算法明显快于完整域实现。

在数值实现中，Heaviside函数$H$被近似为$H$的平滑函数代替，称为平滑的Heaviside函数$H_{ϵ}$，其定义为:
$$H_{ϵ}(x)=\frac{1}{2}[1+\frac{2}{\pi }\arctan (\frac{x}{ϵ})]\text{\hspace{1em}(27)}$$
如[4]，[10]中的$ϵ=1$。因此，作为Heaviside函数H的导数的狄拉克$\delta$函数的$\delta$被$H_{ϵ}$的导数代替，其由$H_{ϵ}$计算:
$${\delta }_{ϵ}=H_{ϵ}^{\prime }=\frac{1}{\pi }\frac{ϵ}{{ϵ}^2+x^2}\text{\hspace{1em}(28)}$$

在每个时间步骤，根据（23）和（24）更新常数$\vec{c}=(c_1,\cdot \cdot \cdot ,c_N)$和偏置$b$字段，其中$u_i=M_i(\varphi )$在第IV节中定义。注意，对于所有$i=1,\cdot \cdot \cdot ,N$，用于计算$e_i$的（17）中的两个卷积$b\ast K$和$b^2\ast K$也出现在（23）中的$\widehat{c_i}$的计算中。对于偏置场$b$，在（24）中计算另外两个卷积$（I{J}^{(1)}）\ast K$和$J^{(2)}\ast K$.因此，在$\varphi$的演化期间，在每个时间步骤总共计算四个卷积。卷积核$K$被构造为$w\ast w$掩模，其中$w$是最小的奇数，使得$w\ge 4\ast \sigma +1$，当$K$是在（11）中定义的高斯核。例如，给定比例参数$\sigma =4$，掩模尺寸为$17\times 17$。

在我们的模型中选择参数很容易。其中一些，例如参数$\mu$和时间步长$\Delta t$，可以固定为$\mu =1.0$和$\Delta t=0.1$。 我们的模型对参数的选择不敏感。 对于大多数强度范围为[0,255]的数字图像，参数$\nu$通常设置为$0.001\ast 255^2$作为默认值。 如我们在第III-C节中提到的，对于具有更多局部强度不均匀性的图像，参数$\sigma$和邻域${\Omega }_y$的大小（由其半径$\rho$指定）应该相对较小。

5.实验结果

图1.利蒙（上排）图像和血管（下排）CT图像的分割。
左，中，右列分别显示初始轮廓（利蒙图像的三角形和血管图像的四边形），中间轮廓和最终轮廓。

我们首先在两相情况下证明我们的方法（即$N=2$）。 除非另有说明，否则本节中的实验参数$\sigma$设置为4。所有其他参数都设置为第IV-C节中提到的默认值。图1示出了利蒙的相机图像和血管的计算机断层摄影血管造影（CTA）图像的结果。曲线演化过程通过显示初始轮廓（左列），中间轮廓（中间列）来描述，以及图像上的最终轮廓（右栏）。在这两个图像中可以清楚地看到强度不均匀性。我们的方法可以为这些图像提供理想的分割结果。

图2.我们的方法在乳房的MR图像，骨骼的X射线图像和前列腺的超声图像中的应用。
第1列：原始图像上的初始轮廓
第2栏：最终轮廓
第3栏：估计的偏差场
第4列：偏差校正图像

通过我们的方法估计的偏置场$\hat{b}$可以用于强度不均匀性校正（或偏差校正）。给定估计的偏置场$\hat{b}$，偏差校正图像被计算为商$I/\hat{b}$。为了证明在同时分割和偏置场估计中的方法的有效性，我们将其应用于具有强度不均匀性的三个医学图像：乳房的MR图像，X-ray图像的骨骼和前列腺的超声图像。这些图像表现出明显的强度不均匀性。超声图像也被严重的斑点噪声破坏。我们应用具有高斯核的卷积来平滑超声图像作为预处理步骤。选择高斯核的比例参数为$2.0$以平滑该超声图像。初始轮廓绘制在图2的第1列中的原始图像上。分割，偏置场估计和偏差校正的相应结果分别在第2,3和4列中示出。这些结果证明了我们的方法在分割和偏差校正方面的理想性能。

A.绩效评估与方法比较

作为水平集方法，我们的方法提供轮廓作为分割结果。 因此，我们使用以下基于轮廓的度量来精确评估分割结果。令$C$为轮廓作为分割结果，$S$为真实对象边界，也作为轮廓给出。对于轮廓$C$上的每个点$P_i,i=1,\cdot \cdot \cdot ,N$，我们可以计算从点$P_i$到地面实况轮廓$S$的距离，由$dist(P_i,S)$表示。然后，我们定义从轮廓$C$到地面实况$S$的偏差:
$$e_{mean}(C)=\frac{1}{N}\sum _{n=1}^{N}dist(P_i,S)$$
这被称为轮廓$C$的平均误差。该基于轮廓的度量可用于评估由轮廓给出的分割结果的子像素精度。

图3.我们的方法对轮廓初始化的稳健性通过其在(a)中具有不同初始轮廓的合成图像的结果来证明。初始轮廓(白色轮廓)和相应的分割结果(黑色轮廓)显示在(b-d)中。

图4.不同初始化和不同尺度参数$\sigma$的方法分割精度。
(a)20次不同初始化的结果的平均误差
(b)12个不同尺度参数$\sigma ,\sigma =4,5,\cdot \cdot \cdot ,15$的结果的平均误差。

1) 轮廓初始化的稳健性

通过上述指标，我们能够使用不同的初始化和不同的参数设置来定量评估方法的性能。我们将我们的方法应用于图3中的合成图像，其中20个不同的轮廓初始化和常数C.例如，我们在图3中显示了20个初始轮廓中的三个（白色轮廓）和相应的结果（黑色轮廓）。在这三个不同的初始化中，初始轮廓包围了感兴趣的对象[图3(b)]，穿过[图3©]中的物体，完全在一个物体内部[图3(d)]。尽管这些初始轮廓有很大差异，但相应的结果几乎相同，都能准确捕捉物体边界。通过根据平均误差评估这些结果来定量验证分割准确性。这些结果的平均误差均在0.21和0.24像素之间，如图4(a)所示。这些实验证明了我们的模型对轮廓初始化的稳健性和子像素级的理想精度。

2)不同尺度参数的稳定性能

我们还用不同尺度参数$\sigma$测试了我们方法的性能，这是我们模型中最重要的参数。对于这个图像，我们使用12个不同的$\sigma$从4到15应用我们的方法。这12个结果的相应平均误差绘制在图4（b）中。虽然平均误差随着$\sigma$增加而增加，但对于本实验中使用的所有12个不同$\sigma$值,它低于0.5像素。

B.与分段光滑模型的比较

图5.我们的方法和PS模型在不同图像条件下的性能（例如，不同的噪声，强度不均匀性和弱物体边界）。
顶行：在原始图像上绘制的初始轮廓
中间行：我们方法的结果
底行：PS模型的结果

图6.我们的模型和PS模型在精度和CPU时间方面的比较。
(a)平均错误
(b)CPU时间

我们还可以在合成图像上定量地比较我们的方法和PS模型。 我们使用相同的对象生成了15个不同的图像，这些图像的边界已知并用作基础事实。这15幅图像是通过平滑理想的二值图像，增加不同轮廓的强度不均匀性和不同的噪声水平而生成的。图5显示了这些图像中的三个作为示例，分别具有我们的模型和中间行和底行中的PS模型的相应结果。 我们对两个模型和所有15个图像使用相同的初始轮廓（顶行中的圆圈）。很明显，我们的模型比PS模型产生更准确的分割结果。为了定量评估准确度，我们计算了所有15个图像的两个模型的平均误差,如图6（a）所示，其中x轴代表15个不同的图像。如图6（a）所示，我们模型的误差明显低于PS模型的误差。

另一方面，我们的模型比PS模型更有效。这可以从两个模型消耗的15个图像的CPU时间中看出[见图6(b)]。 在这个实验中，我们的模型比PS模型快得多，在我们的实现中平均加速因子为36.43。本实验中的CPU时间是在带有Intel®Core™2 Duo CPU，2.40 GHz，2 GB RAM的Lenovo ThinkPad笔记本上运行我们的Matlab程序，在Windows Vista上使用Matlab 7.4。

C.应用于MR图像分割和偏差校正

图7.我们的方法在3T MR图像中的应用。
第1栏：原始图像
第2列：（红色）和（蓝色）的最终零级轮廓，即分割结果
第3栏：估计的偏差场
第4栏：偏差校正图像
第5列：原始图像的直方图（左图）和偏差校正图像

图8.应用于7T MR图像。
(a)原始图像
(b)偏差校正图像
©计算偏差场

在本小节中，我们将重点放在所提出的方法在脑MR图像的分割和偏差校正中的应用。我们首先在图7的第一列中显示3T MR图像的结果。这些图像表现出明显的强度不均匀性。分割结果，计算的偏置场，偏置校正图像分别显示在第二，第三和第四列中。可以看出，在偏差校正图像中，每个组织内的强度变得非常均匀。通过比较原始图像和偏差校正图像的直方图，也可以证明在强度均匀性方面的图像质量的改善。在第五列中绘制原始图像（左）和偏差校正图像（右）的直方图。在偏差校正图像的直方图中存在三个明确定义且分离良好的峰，每个峰对应于图像中的组织或背景。相反，由于偏差引起的强度分布的混合，原始图像的直方图不具有这种良好分离的峰值。

我们的方法也在7T MR图像上进行了测试，结果令人满意。在7T，由于信噪比的增加，可以获得图像分辨率的显着增益。 然而，磁化率引起的梯度与主场成比例，而成像梯度目前限于与在较低场强（即3T）下使用的基本相同的强度。 这种效果在空气/组织界面处最明显，如图8(a)中额叶底部所见。这似乎是一种高度本地化和强烈的偏见，这对传统的偏差校正方法具有挑战性。 该图像的结果显示了我们的方法纠正这种偏差的能力，如图8(b)和©所示。

6.与分段常量和分段平滑模型的关系

值得指出的是，我们在（14）中的两相水平集公式中的模型是众所周知的Chan-Vese模型[4]的推广，它是一个代表性的分段常数模型。当偏置场$b$是常数$b=1$时，我们在（14）中提出的能量$\epsilon$减少到Chan-Vese模型中的数据拟合项。 为了表明这一点，我们需要$\int K(y-x)\mathrm{d}x=1$的事实，并回想起$M_1(\varphi )=H(\varphi )$和$M_2(\varphi )=1-H(\varphi )$。因此，对于$b=1$的情况，通过改变（14）中的求和和积分的顺序，能量$\epsilon$可以被重写为:
$$\begin{array}{cc}\epsilon & ={\sum }_{i=1}^2\int \left(\int K(y-x)\mid I(x)-c_i{\mid }^2{M}_i(\varphi (x))\mathrm{d}y\right)\mathrm{d}x\\ & ={\sum }_{i=1}^2\int \left(\mid I(x)-c_i{\mid }^2{M}_i(\varphi (x))\int K(y-x)\mathrm{d}y\right)\mathrm{d}x\\ & ={\sum }_{i=1}^2\int \mid I(x)-c_i{\mid }^2{M}_i(\varphi (x))\mathrm{d}x\\ & =\int \mid I(x)-c_1{\mid }^2{H}_i(\varphi (x))\mathrm{d}x+\int \mid I(x)-c_2{\mid }^2(1-H_i(\varphi (x)))\mathrm{d}x\end{array}$$
这正是Chan-Vese模型中的数据拟合项(2)。Chan-Vese模型是分段常数模型，其目的在于找到分别与区域 Ω 1 = { ϕ &gt; 0 } \Omega_1=\{\phi &gt;0\} Ω1​={ϕ>0}和 Ω 2 = { ϕ &lt; 0 } \Omega_2=\{\phi &lt; 0 \} Ω2​={ϕ<0}中的图像I近似的常数$c_1$和$c_2$。

我们的模型也与分段光滑的Mumford-Shah模型密切相关。 Mumford-Shah模型通过在第II节中简要描述的计算上昂贵的过程,分别通过寻找在不相交区域${\Omega }_1,\cdot \cdot \cdot ,{\Omega }_N\subset \Omega$上定义的N个平滑函数$u_1,\cdot \cdot \cdot ,u_N$来执行图像分割。

与Mumford-Shah模型不同，我们的模型旨在找到图像$I$的乘法分量：平滑函数$b$和分段常数函数$J$. 得到的$b$和$J$产生分段平滑函数$bJ$作为图像I的近似。从我们的方法中的能量最小化过程和如前所述的Mumford-Shah模型，显然前者获得分段平滑近似，通过以比后者更有效的方式产生图像分割结果。

7.结论

我们已经提出了一种变分水平集框架，用于对具有强度不均匀性的图像进行分割和偏差校正。基于普遍接受的具有强度不均匀性的图像模型和导出的局部强度聚类特性，我们定义了表示图像域的分区的水平集函数的能量和表示强度不均匀性的偏置场。因此，通过最小化所提出的能量函数来联合执行分段和偏置场估计。从我们的变分框架中的数据项自然地确保了从所提出的能量导出的偏置场的缓慢变化特性，而不需要在偏置场上施加明确的平滑项。我们的方法比分段平滑模型更稳健。实验结果表明，该方法在精度，效率和鲁棒性方面具有优越的性能。作为一种应用，我们的方法已应用于MR图像分割和偏差校正，结果令人满意。