随着智能医疗概念的兴起和推广，使用计算机视觉辅助进行医疗诊断及治疗逐渐成为现代医学领域的研究发展趋势。其中计算机辅助医学图像分割技术旨在从医学影像中分割出目标组织或结构，其分割结果的准确性对后续的病理分析、术中影像定位导航等相关临床操作均起到了至关重要的作用，是智能医疗领域的研究热点及难点(Masood等，2015)。然而在临床中，这一工作主要由影像科的医务工作者手工勾勒目标边界完成，其精度受到操作者业务水平及主观工作态度等多方面因素影响，在准确度方面难以得到保障，同时浪费了大量的人力资源。

临床中，脑部医学成像主要依赖核磁共振成像(MRI)技术，该技术受成像原理和人体组织结构复杂等因素限制，通常存在强度分布不均、组织相互重叠等问题，直接应用以阈值、几何约束或统计学为主的传统手段取得的分割效果较差，促使研究者探索新的图像分割方法(刘宇等，2017)。Kass等人(1988)提出Snake模型，将图像分割过程中因迭代而不断变化的边界视为活动轮廓并在该轮廓上定义能量函数，通过能量函数最优化确定轮廓演化方向及最优位置，活动轮廓模型的思想一经提出便受到研究人员的广泛关注。Cohen等人(1991)通过增加气球力约束对Snake模型进行优化，以此为基础大量带参数的活动轮廓模型相继诞生。

Osher等人提出水平集思想，在图像上定义高维函数，绘制等高线，并将水平高度为零的轮廓集合视为目标轮廓，极大程度地优化了能量函数求解过程并在图像处理领域产生了深远影响。Caselles等人(1993)为进一步解决曲线拓扑变化效果不佳问题构建了几何活动轮廓分割模型(GACM)，继而启发诞生了测地线活动轮廓(GAC)模型(Caselles等，1997)、Chan-Vese(C-V)模型(Chan等，2001)、局部强度聚类(LIC)模型(Li等，2011)、距离正则化水平集演化(DRLSE)模型(Li等，2010)等一系列水平集分割模型。

尽管水平集方法在图像分割领域中得到了广泛应用，但对于医学图像尤其是脑部MRI图像而言，常存在灰度分布不均、组织结构复杂且边界模糊等问题，直接应用水平集方法不能取得良好的实验效果。为了解决上述问题，研究者继而提出使用混合变分水平集的方法，并将概率分布模型、频率变换、机器学习等思想与水平集方法相结合以提升分割精度(陈红等，2018)。为了更充分地利用图像数据，研究者进一步结合了目标组织的结构特征、几何性状(Feng等，2016)、形状参数(Gao等，2018)、颜色知识等先验信息，使得分割结果的准确性和鲁棒性获得了一定程度的提升。

然而正如图 1所示，脑部图像中的肿瘤区域内部经常会存在大小及数量未知的明暗区域，单纯使用水平集分割算法所得的分割轮廓常包含孔洞或存在过度收缩等问题。若以常见的肿瘤区域轮廓约束分割过程中目标轮廓的形状，则可以减少孔洞或过度收缩问题的发生，最终更好地实现准确分割。常见的形状先验约束方法包括引入单一或多个预定义几何形状方程(Yang等，2017)、建立图像形状的统计学概率模型(Mesadi等，2018)、建立形状映射模型、建立稀疏形状模型(姚劲草，2017)等方式。但由于脑部肿瘤形状通常较为复杂，同时其轮廓形状因患者的个体差异而存在较大的差异，因此，准确地描述脑部肿瘤的形状结构不仅需要较多的像素数据，还需要大量的脑肿瘤MRI图像样本。稀疏表达在复杂形状表示方面具有一定的优势，故本文提出一种基于稀疏形状先验的脑肿瘤图像分割算法，对肿瘤形状进行描述，构建稀疏形状约束模型，并结合水平集约束方法构建能量泛函，利用高层稀疏约束和底层能量函数之间的关系进行目标轮廓的演化，最终实现对图像中肿瘤区域的分割提取，以结合水平集方法在拓扑结构描述和稀疏表达方法在复杂形状表达方面的优势同时减少肿瘤内部明暗区域对分割结果造成的影响。

本节简要介绍混合变分水平集分割理论的通用分割模型以及形状的稀疏表示模型。

基于混合变分水平集理论的目标分割方法是将目标的形状、颜色、比例等先验信息与变分水平集方法相融合，共同构建底层能量函数，其通用能量函数框架为

$ E\left( \phi \right) = {E_{\rm{p}}}\left( \phi \right) + {E_{\rm{i}}}\left( \phi \right) $

(1)

式中，$\phi $为水平集函数，$ {E_{\rm{p}}}(\phi)$代表先验参考信息项，$ {E_{\rm{i}}}(\phi)$代表变分水平集项，结合实际应用背景可使用不同的模型。

C-V模型以分割轮廓为界将图像分为目标区域和背景区域并进行灰度约束和轮廓长度约束，在目标与背景灰度差异大且分布均匀时效果良好。

LIC模型进一步应用于灰度分布不均的图像分割场景并加入局部灰度聚类能量项，可以更好地实现偏移场矫正。

DRLSE模型将能量函数分解为正则化项${E_{\rm{r}}}(\phi) $与扩展项$ {E_{{\rm{ext}}}}(\phi)$，前者定义为

$ {E_{\rm{r}}}\left( \phi \right) = \frac{1}{2}\int {{{\left( {\left| {\nabla \phi } \right| - 1} \right)}^2}} {\rm{d}}x $

后者因实际问题而异，有效避免了演化过程中$ \phi $重新初始化所带来的问题(Khadidos等，2017)。

本文研究基于形状先验的图像分割，研究难点除根据图像特征选择合适的变分水平集模型外，还包括选择合适的形状先验模型以适应脑肿瘤复杂多变的形状。

本文基于稀疏表示模型进行形状建模，该方法构建稀疏字典并使用尽可能少的原子进行组合以重构原始信号。常见的稀疏字典构成法包括预定义、K-奇异值分解(K-SVD)(Rubinstein等，2012)等方法，其中，当信号数据类型相对稳定时，使用预定义稀疏字典更有助于提升算法的运行效率。因此本文使用预定义的稀疏字典并基于经典的稀疏形状组合(SSC)表示方法(Cremers等，2003)进行形状表示。

SSC方法使用不同类型的形状信息作为原子构建过完备字典$\mathit{\boldsymbol{A}} \in {{\bf{R}}^{M \times K}} $，其中$M $表示原子维度，$ K$表示形状的类型数目。待观测信号$ \mathit{\boldsymbol{y}}$的重构公式表示为

$ \min {\left\| \mathit{\boldsymbol{s}} \right\|_0}{\rm{ s}}{\rm{. t}}{\rm{. }}\;\mathit{\boldsymbol{As}} = \mathit{\boldsymbol{y}} $

(2)

式中，$\mathit{\boldsymbol{s}} $表示稀疏系数，为了尽可能少地使用原子来重构信号以获得更为稀疏的表示形式，式中第1项约束了稀疏系数中非零项的数量。

然而L₀范数的非凸特性使得式(2)变为NP难问题，因此在实际操作中通常使用图松弛算法将L₀范数替换为高阶L₁范数，并使用拉格朗日乘子法以能量泛函的形式将稀疏系数求解过程表示为求解最小化能量泛函$ E$，即

$ \mathop {\min }\limits_{\mathit{\boldsymbol{s}} \in {{\bf{R}}^{K \times 1}}} \left\{ {E = \left\| {\mathit{\boldsymbol{As}} - \mathit{\boldsymbol{y}}} \right\|_2^2 + \lambda {{\left\| \mathit{\boldsymbol{s}} \right\|}_1}} \right\} $

(3)

式中，$ \lambda $为稀疏约束项的权重系数。

使用稀疏形状组合方法进行稀疏表示的关键在于过完备字典构成及形状特征的选择和提取，以实现在最大程度减少运算数据量的同时保证并提升形状拟合的准确率。

本文基于稀疏形状先验的脑肿瘤图像分割算法，结合肿瘤形状的稀疏表示模型与变分水平集方法共同构建能量表达式并给出能量函数优化过程。算法框架流程如图 2所示。

使用矩阵$ \boldsymbol{X}=\left[x_{1}, x_{2}, \cdots, x_{N}\right]$表示训练形状样本集，$ \mathit{\boldsymbol{X}} \in {{\bf{R}}^{K \times V}}$，其中$ K$表示单幅图像的像素总数量，$ N$表示训练样本中的图像数量，${{x_i}} $表示第$ i$个训练样本，$\mathit{\boldsymbol{X}} $中各样本均值使用$ \bar x = \frac{1}{N}\sum\limits_{i = 1}^N {{x_i}} $表示。

使用$ {m_i}$表示第$i $幅图像经中心化操作后的输出结果，具体为

$ {m_i} = {x_i} - \bar x $

(4)

由此可得到中心化矩阵$\mathit{\boldsymbol{M}} = {\left[ {{m_1}, {m_2}, \cdots, {m_N}} \right]^{\rm{T}}} $，将该矩阵的协方差矩阵表示为$ \boldsymbol{C}_{x}=\mathrm{E}\left\{\boldsymbol{M} \boldsymbol{M}^{\mathrm{\tau}}\right\}$，其中$E\left\{ \cdot \right\} $表示矩阵期望。由矩阵运算性质可知$ \boldsymbol{C}_{x}$为实对称矩阵。

将$ \boldsymbol{C}_{x}$的特征值及其相对应的特征向量分别表示为$\lambda_{1}, \lambda_{2}, \cdots, \lambda_{N} $和$\boldsymbol{e}_{1}, \boldsymbol{e}_{2}, \cdots, \boldsymbol{e}_{N} $，以此构建特征向量矩阵$\boldsymbol{U}=\left[\boldsymbol{e}_{1}, \boldsymbol{e}_{2}, \cdots, \boldsymbol{e}_{N}\right]^{\mathrm{T}} $。对于$ \mathit{\boldsymbol{X}}$中的任意样本${{x_i}} $，均可通过特征向量矩阵$\boldsymbol{U}=\left[\boldsymbol{e}_{1}, \boldsymbol{e}_{2}, \cdots, \boldsymbol{e}_{N}\right]^{\mathrm{T}} $映射为变量${y_i} $，具体为

$ {y_i} = \mathit{\boldsymbol{U}}\left( {{x_i} - \bar x} \right) $

(5)

对$ \mathit{\boldsymbol{X}}$中各样本依次进行映射可得到映射矩阵$\boldsymbol{Y}=\left[y_{1}, y_{2}, \cdots, y_{N}\right] $。由特征向量矩阵的性质可知，矩阵${\mathit{\boldsymbol{U}}^{\rm{T}}} = {\mathit{\boldsymbol{U}}^{ - 1}} $中的各行向量之间相互正交且满足$ {\mathit{\boldsymbol{U}}^{\rm{T}}} = {\mathit{\boldsymbol{U}}^{ - 1}}$，故样本${{x_i}} $可由

$ {x_i} = {\mathit{\boldsymbol{U}}^{\rm{T}}}{y_i} + \bar x $

(6)

恢复重建，将此重构形状样本集合记为$\mathit{\boldsymbol{\hat X}} = {\left[ {{{\hat x}_1}, {{\hat x}_2}, \cdots, {{\hat x}_N}} \right]^{\rm{T}}} $。

取$ n < N$个特征向量重新构建矩阵$ {\mathit{\boldsymbol{U}}_n}$的方式减少$ U$的维度是一种行而有效的降低运算复杂度操作。由于该操作将导致重建信号与原始信号间存在大小为${e_{ms}} = \sum\limits_{j = n + 1}^N {{\lambda _j}} $的损失。为使该损失项数值达到最小，取$ {\mathit{\boldsymbol{C}}_x}$的特征值中最大的前$ n$项对应的特征向量构成矩阵$ {\mathit{\boldsymbol{U}}_n}$，并将${{x_i}} $推广至任意目标形状$ x$。则式(6)推广为

$ \hat x = \mathit{\boldsymbol{U}}_n^{\rm{T}}\mathit{\boldsymbol{y}} + \hat x $

(7)

由于样本均值$ {\bar x}$仅与样本集$\mathit{\boldsymbol{X}} $相关，其值固定，故本文对中心化后的形状样本进行稀疏近似表示使其接近真实值，即对于任意形状样本$ x$，求其解稀疏组合$ \hat x - \bar x = \mathit{\boldsymbol{As}}$以使得$ \hat{x}-\bar{x} \approx x-\bar{x}$，并与式(7)对比，即

$ \left\{ {\begin{array}{*{20}{l}} {\hat x = \mathit{\boldsymbol{As}} + \bar x}\\ {\hat x = \mathit{\boldsymbol{U}}_n^{\rm{T}}\mathit{\boldsymbol{y}} + \bar x} \end{array}} \right. $

(8)

由于$ \mathit{\boldsymbol{U}}_n^{\rm{T}}$同样不随能量函数的迭代演化而发生变化，本文将问题进一步转化为求解稀疏组合$\mathit{\boldsymbol{As}} $以满足$ \mathit{\boldsymbol{y}} = \mathit{\boldsymbol{As}}$。构建稀疏约束模型

$ \min \left\| \mathit{\boldsymbol{s}} \right\|_p^p{\rm{ s}}{\rm{. t}}{\rm{. }}\;\mathit{\boldsymbol{As}} = \mathit{\boldsymbol{y}} $

(9)

结合式(5)对$y $进行变量替换

$ \min \left\| \mathit{\boldsymbol{s}} \right\|_p^p{\rm{ s}}{\rm{. t}}{\rm{. }}\;\mathit{\boldsymbol{As}} = \mathit{\boldsymbol{U}}\left( {x - \bar x} \right) $

(10)

本文使用K-means方法将映射矩阵$\mathit{\boldsymbol{Y}} $中的数据聚为$M $类并分别求取平均形状，得到集合$ \hat Y = \left[ {{{\bar y}_i}} \right.,\left. {{{\bar y}_2}, \cdots ,{{\bar y}_M}} \right]$，并使用该集合作为预定义的稀疏字典$ \mathit{\boldsymbol{A}}$，通过正交匹配追踪更新稀疏系数$ \mathit{\boldsymbol{s}}$。

实际情况中，由于待匹配的图像形状可能因存在位移、旋转和缩放而难以与模板实现完全匹配，因此在操作前需要对图像进行配准。使用人工标记点、图像形状特征点、角点等人工或人工特征点集合来进行图像间的模式匹配(Van等，1993)是一种基本的图像配准操作，然而该方法对特征点的选择提出了较高要求，同时其使用范围受到了限制。研究者进一步提出在图像配准过程中引入傅里叶变换(Luo等，2000)、互相关信息、小波变换等技术(Sharman等，2000)以避免特征点的筛选过程，在提高图像配准精度的同时提升了算法的灵活性。

本文基于Fourier-Melli(Chen等，1994)算法对图像进行配准，该方法利用笛卡尔坐标系中图像的缩放和旋转在相对应的对数极坐标系中分别表现为横向与纵向平移的原理，对形变图像与配准图像分别进行极坐标变换，并求其互能量谱以求解图像间的旋转和伸缩量，并在对形变图像进行变换后继续求解变换后形变图像与配准图像之间的互功率谱以求解平移变换量，从而对形变图像进行平移变换，以实现两幅图像之间的对齐。如图 3所示，原始数据集中脑肿瘤的位置和朝向不一，经过配准可获得整齐的图像数据集。

为方便描述，定义一种符号函数$ P_{\kappa, \theta, \tau}(q)$以表示待匹配的图像形状$ q$在进行配准操作化后所得的图像形状，$P_{\kappa, \theta, \tau}^{-1}(q) $表示其反变换，其中$\kappa $表示图像的放缩变换量，$ \theta $表示图像的旋转变换量，$ \tau $表示平移变换量，将训练集经对数极坐标变换后的均值表示为符号$ {\bar P}$。则构建稀疏能量函数项可以重新定义为

$ {E_{\rm{s}}} = \left\| {\mathit{\boldsymbol{As}} - \mathit{\boldsymbol{U}}{\mathit{\boldsymbol{P}}_{\kappa ,\theta ,\tau }}\left( q \right) + \mathit{\boldsymbol{U}}\bar P} \right\|_2^2 + \lambda {\left\| \mathit{\boldsymbol{s}} \right\|_1} $

(11)

本文使用数值为0和$ k$的两条活动水平集，定义水平函数在图像中肿瘤边沿的数值为0，在肿瘤外的脑部组织与纯黑背景交界处数值为$ k$。由此将MRI图像分为3个区域：1)肿瘤边沿内部区域，水平集函数数值为负；2)两水平集轮廓之间的区域，水平集函数数值保持在0~$ k$之间；3)大脑外部区域，水平集函数数值大于$ k$。

如图 4所示，本文使用绿色曲线表示肿瘤轮廓，红色曲线表示人脑部组织轮廓，将MRI图像分为3个不同区域。水平集函数在上述轮廓和各区域中的取值如图所示。对于图像中的任一像素$y $，取其高斯邻域窗口，如图中黑色曲线所示，与MRI图像的3个区域分别重叠于区域$ f_{1}(y)、f_{2}(y)$和$f_{3}(y) $。此外，为便于阅读，本文对原始的MRI图像进行了适当的灰度变换处理，此灰度变换处理仅用于示例展示用途，不应用于后续的实验部分。

使用区域能量函数项$ {E_{{\rm{re}}}}(\phi)$约束图像中的像素点在所属窗内进行灰度聚类，具体为

$ \begin{array}{*{20}{c}} {{E_{{\rm{re}}}}(\phi ) = \sum\limits_{i = 1}^3 {{\mu _i}} \int {\int {{M_i}(\phi (x)){K_\sigma }(x - y)} } }\\ {{{\left| {\mathit{\boldsymbol{I}}(x) - {\mathit{\boldsymbol{f}}_i}(y)} \right|}^2}{\rm{d}}x{\rm{d}}y} \end{array} $

(12)

式中，$K_{\sigma}(u)=\frac{1}{(2 \pi)^{n / 2} \sigma^{n}} \mathrm{e}^{-|u|^{2} / 2 \sigma^{2}} $定义了点$y $的邻域窗口，$\mu_{\mathrm{i}}, i=1, 2, 3 $分别代表各区域的权重系数项，此外，各区域使用函数

$ {M_i}(\phi (x)) = \left\{ {\begin{array}{*{20}{l}} {1 - H(\phi (x))}&{i = 1}\\ {H(\phi (x)) - H(\phi (x) - k)}&{i = 2}\\ {H(\phi (x) - k)}&{i = 3} \end{array}} \right. $

进行分区，$H(·) $代表单位阶跃函数。$ f_{i}(y), i=1, 2, 3$由梯度下降法求解，具体为

$ \left\{ {\begin{array}{*{20}{l}} {{f_1} = \frac{{{K_\sigma }*(I - H(\phi )I)}}{{{K_\sigma }*(I - H(\phi ))}}}\\ {{f_2} = \frac{{{K_\sigma }*(H(\phi )I - H(\phi - k)I)}}{{{K_\sigma }*(H(\phi ) - H(\phi - k))}}}\\ {{f_3} = \frac{{{K_\sigma }*(H(\phi - k)I)}}{{{K_\sigma }*H(\phi - k)}}} \end{array}} \right. $

(13)

进一步地，为了避免能量函数迭代过程中水平集函数$ \phi $重新初始化所带来的不良影响，本文加入正则化项

$ {E_{\rm{r}}}(\phi ) = \frac{1}{2}\int {{{(|\nabla \phi | - 1)}^2}} {\rm{d}}x $

(14)

以约束$ \phi $的距离函数特性。

由于$(1-H(\phi)) $可以将图像分割成为肿瘤区域为1、非肿瘤区域为0的二值图像，故本文将此图像作为稀疏能量函数中的目标形状数据，即

$ \begin{array}{*{20}{c}} {{E_{\rm{s}}}(\phi ) = \left\| {\mathit{\boldsymbol{As}} - \mathit{\boldsymbol{U}}{P_{\kappa ,\theta }}(1 - H(\phi )) + \mathit{\boldsymbol{U}}\bar P} \right\|_2^2 + }\\ {\lambda {{\left\| \mathit{\boldsymbol{s}} \right\|}_1}} \end{array} $

(15)

由此，本文将总体能量函数定义为稀疏能量函数项、区域能量函数项以及正则化项之和，表示为

$ E(\phi ) = {\omega _1}{E_{\rm{r}}}(\phi ) + {\omega _2}{E_{{\rm{re}}}}(\phi ) + {\omega _3}{E_{\rm{s}}}(\phi ) $

(16)

式中，各能量函数项权重系数$ \omega_{1}, \omega_{2}, \omega_{3}>0$。

能量函数优化算法包括参数初始化和迭代优化两部分，首先固定稀疏系数$ s$，更新水平集函数$ \phi $；然后固定$ \phi $，更新稀疏系数$ s$。当轮廓收敛或达到迭代终止条件时，停止函数迭代并输出分割结果。具体步骤为：

1) 目标属于字典中各种类型形状的初始概率相同，设置稀疏系数$ \mathit{\boldsymbol{s}} = \left[ {\frac{1}{N}, \frac{1}{N}, \cdots, \frac{1}{N}} \right]$；

2) 设置能量函数中各项系数；

3) 设置初始轮廓；

4) 使用Fourier-Melli算法更新$\tau, \kappa, \theta $；

5) 使用梯度下降法求解水平集函数$\phi $的梯度值$ \frac{\partial \phi}{\partial t}$；

6) 更新水平集函数$\phi=\phi_{\mathrm{pre}}+\Delta t \times \frac{\partial \phi}{\partial t} $，其中，$ \phi_{\mathrm{pre}}$代表更新前的水平集函数，$ \Delta t$代表时间步长；

7) 使用正交匹配追踪算法(Pati等，1993)更新$ \mathit{\boldsymbol{s}}$；

8) 返回步骤4)，直至轮廓收敛或已达到迭代次数限制；

9) 输出分割结果，结束。

在步骤5)中，水平集函数的梯度值$\frac{\partial \phi}{\partial t}$可由梯度下降法计算获得，即

$ \begin{array}{*{20}{c}} {\frac{{\partial \phi }}{{\partial t}} = {\omega _1}\left[ {{\nabla ^2}\phi - div\left( {\frac{{\nabla \phi }}{{\left| {\nabla \phi } \right|}}} \right)} \right] + }\\ {{\omega _2}\left[ {{e_1}\delta (\phi ) - {\mu _1}{e_2}(\delta (\phi ) - \delta (\phi - k)) - } \right.}\\ {\left. {{\mu _2}{e_3}\delta (\phi - k)} \right] + {\omega _3}\delta (\phi )P_{\kappa ,\theta ,\tau }^{ - 1}\left( {{\mathit{\boldsymbol{U}}^{\rm{T}}}\left( {\mathit{\boldsymbol{As}} - } \right.} \right.}\\ {\left. {\left. {\mathit{\boldsymbol{U}}{P_{\kappa ,\theta ,\tau }}(H(\phi )) - \mathit{\boldsymbol{U}}\bar x} \right)} \right)} \end{array} $

(17)

式中，$e_{i}=\left(K_{\sigma} * 1\right) I^{2}-2\left(K_{\sigma} * f_{i}\right) \cdot I+K_{\sigma} * f_{i}^{2} $。

本文使用与2017年医学影像计算与计算机辅助大会(MICCAI)同时举办的多模态脑胶质瘤分割挑战赛(BraTS)提供的多模态脑胶质瘤公开数据集(Menze等，2015)进行算法验证，该数据集包含脑部神经胶质瘤患者的脑部MRI图像，同时提供由专家团队标注的标准肿瘤区域数据。本文首先使用快速包围盒(FBB)算法(Saha等，2012)得到脑部肿瘤的初始矩形轮廓区域，并使用该轮廓重心作为区域增长法的种子点，继而进行区域增长得到水平集函数的初始轮廓位置。设置时间步长$ \Delta t=0.1$，各项系数为$ {w_1} = 0.2/\Delta t, {w_2} = 0.25, {w_3} = 2, {\mu _1} = 0.5, {\mu _1} = - 0.05$，设置算法的迭代次数上限为50。

为验证基于稀疏形状先验的脑肿瘤图像分割算法可行性，本文将该算法与部分现有图像分割算法进行实验对比：算法1(Joshi A等，2015)使用小波变换对脑部肿瘤进行预处理，继而定义了一种基于轮廓的水平集方法进行图像分割，最终使用软阈值的方法对分割结果的形状和大小进行筛选，而得到肿瘤目标；算法2(Zabir等，2015)使用K-means的方法寻找初始的肿瘤位置点并使用区域增长的方法获得水平集初始轮廓，并进行演化分割出目标的肿瘤区域；算法3(Kermi等，2018)前期处理阶段与本文类似，即利用脑部对称性寻找存在脑肿瘤的区域范围，并使用区域增长方式获得初始轮廓，其不同之处在于后续使用测地线活动轮廓模型方法进行脑部肿瘤分割；算法4(Mojtabavi等，2017)结合区域和边缘方法定义水平集方程并利用人工勾勒出脑部肿瘤的初始轮廓，使用fast marching方法(Sethian等，1996)进行迭代和优化。此外，为了进一步验证形状约束项对分割结果的影响，本文在测试本文算法过程中屏蔽了形状约束项，即将$ {w_3}$置0以进行对比实验。

为进一步衡量脑部肿瘤分割结果的准确性，选用相似度(DSC)、灵敏度(SEN)以及阳性预测率(PPV)作为技术指标。令$P $代表经图像分割算法处理后所得图像中的肿瘤区域，$T $表示经由专家标定后的标准肿瘤区域，| · |表示求取区域内像素总数目操作，$|P \wedge T| $表示$P、T $两区域彼此重叠部分的像素数量，则上述指标可具体定义为

$ {f_{{\rm{DSC}}}} = \frac{{\left| {P \wedge T} \right|}}{{\left( {\left| P \right| + \left| T \right|} \right)/2}} $

(18)

$ {f_{{\rm{SEN}}}} = \frac{{\left| {P \wedge T} \right|}}{{\left| T \right|}} $

(19)

$ {f_{{\rm{PPV}}}} = \frac{{\left| {P \wedge T} \right|}}{{\left| P \right|}} $

(20)

式中，$ {{f_{{\rm{DSC}}}}}$表示分割结果与专家数据之间重叠的比重, 其值越高, 代表分割结果越为准确；$ {{f_{{\rm{SEV}}}}}$表示分割结果中的正确部分占真值的比重, 其值越高, 证明分割结果的漏检率越低；$ {{f_{{\rm{PPV}}}}}$表示分割结果中正确部分占分割结果的比重, 其值越高, 证明分割结果的误检率越低。

图 5展示了数据库中部分随机样本的分割效果，其每一行代表使用一个随机样本进行图像分割后的分割结果。其中，图 5(a)中的红色线轮廓为专家手动标注的脑肿瘤轮廓实际值，图 5(b)-(g)中的红色线轮廓为使用各算法所得的脑部肿瘤轮廓线。使用各种类型图像分割算法的实验结果对比数据如表 1所示。

从实验结果对比图和对比数据可以看出，使用本文算法所得的脑部肿瘤分割结果优于其他算法。屏蔽形状约束项之后的实验结果准确性明显低于本文算法，且误判率较高，分割轮廓中包含了较多非肿瘤区域。算法1的准确率相对较低，同时SEN值相对较低，其分割轮廓容易存在过度收缩的问题，漏判率较高，同时该方法对于目标轮廓边沿的清晰度较为敏感，当目标轮廓边沿模糊时，分割结果容易产生较大的偏差；算法2的PPV值较低，其分割轮廓容易发生过度收缩的问题，造成实验结果的误判，同时对于目标的灰度均匀程度较为敏感，容易将目标内部较为明显的灰度团块误判成错误区域；算法3的分割结果整体表现较为稳定，但各方面表现相较于本文算法均有不足，对目标形状与预定义的几何形状是否接近较为敏感，当目标具有复杂轮廓边缘时所取得的分割效果不理想；算法4的SEN值也相对较低，且当肿瘤图像内部存在暗影区域时，其分割轮廓内部容易存在孔洞，增加了实验结果的漏判率。

实验结果表明，本文算法的分割结果与真实数据之间的平均相似度达到93.97%，灵敏度达到91.3%，阳性预测率达到95.9%。与其他方法相比，本文算法由于添加了稀疏形状约束而具有更好的轮廓形状适应性，能够避免因为轮廓内部存在不均匀团块而对分割造成的不良影响，具有更好的鲁棒性，分割结果也更加稳定。

本文提出一种基于稀疏形状先验的脑部肿瘤图像分割算法，通过对脑部肿瘤轮廓进行形状描述，构建其稀疏形状约束模型，并结合水平集约束方法构建能量泛函，而后利用高层稀疏约束和底层能量函数之间的关系进行目标轮廓的演化。通过在脑胶质瘤公开数据Brats2017上进行实验验证，可以发现当图像中的肿瘤区域内存在不均匀的明暗团块时，本文算法的实验结果可以避免分割轮廓过度收缩及内部存在孔洞等问题。同时，与其他对比算法相比，本文算法在结构复杂的图像上也取得了相对更好的实验结果。本文算法使用了固定字典进行系数建模，该方法在适用范围方面具有一定的局限性，本课题组在未来的研究中，将进一步对形状建模及训练部分展开研究，以提升脑部肿瘤形状模型的普适性和稳定性。此外，由于本文算法的研究重点在于2维图像的脑部肿瘤的分割，而在临床实际中，医生通常可以获得3维的脑部扫描图片。因此，为了更进一步地利用图像数据信息，本课题组未来也将研究多模态脑部肿瘤分割问题，从而实现更为精准的计算机辅助医疗诊断。