脑胶质瘤是最常见的原发性脑肿瘤，其中一半以上为恶性度最高的胶质母细胞瘤（glioblastoma multiforme，GBM）^[1]. 无论是良性还是恶性GBM，均会使颅内压升高，压迫脑组织，损害中枢神经，危及患者生命. GBM病变组织的定位和量化计算（如肿瘤体积等）对于GBM的诊断、手术规划、术后分析等至关重要. 在临床上，放射科医生通常通过多模态磁共振（magnetic resonance，MR）图像手动分割肿瘤，工作繁琐耗时. 自动分割GBM肿瘤技术可减轻医生工作负担，有效避免医生因过度疲劳而产生错误. 因此，GBM肿瘤的计算机自动分割技术对GBM肿瘤的辅助治疗具有极其重要的意义^[2].

实现GBM肿瘤的精确自动分割存在较大困难，但因其具有巨大的辅助医疗意义，近年来众多科研人员致力于该项技术的研究. 已有的GBM多模态MR图像分割方法大致可以分为2类：一类是基于像素或体素的分割方法，另一类是基于图形的分割方法^[3]. 基于像素或体素的分割方法的基本思想是根据各个像素在多模态图像上亮度、纹理等信息将该像素分类至相应类别. Clark等^[4]通过模糊K均值分割，从正常脑组织和肿瘤中提取判别性知识，利用以自适应强度为基础的多光谱直方图阈值法检测肿瘤区域；类似的方法被Fletcher-Heath等^[5]运用于非增强肿瘤的分割. 基于图形的分割算法则用图的顶点来表征图像像素，用图的边表征2个像素的相似性，由此形成网络图，通过解决能量最小化问题将图分割成子网络图，使不同子网络图之间的差异和同一子网络图内部的相似性达到最大^[6-7]. 除了以上2类算法，基于水平集的分割算法也广泛应用于GBM多模态MR图像分割^[8-9]. 上述研究虽取得了一定的研究成果，但绝大多数集中在搜索和分割整个肿瘤区域，并未将整个肿瘤进一步细分为肿瘤不同子区域. 分离出肿瘤不同子区域具有十分重要的意义，因为通过肿瘤不同子区域可以获取对GBM预后具有重要预测作用的高层次影像特征（如坏死区比例等）^[10].

本研究提出基于混合特征和先验知识的GBM多模态MR图像分割算法，能够在适应临床实践需求的前提下将整个肿瘤细分为肿瘤不同子区域. 为了将头部区域方位调整到未旋转位置，对GBM多模态MR图像配准，并应用N4ITK法进行偏置场校正. 提取GBM多模态MR图像的局部位置特征、强度特征、纹理特征、对称特征和上下文特征，利用随机森林分类器将GBM多模态MR图像初步分割为脑部正常组织区、坏死区、肿瘤强化区、不包括坏死区和肿瘤强化区的T1异常区以及不包括T1异常区的液体衰减反转恢复（fluid attenuated inversion recovery，FLAIR）成像异常区. 考虑GBM肿瘤解剖结构先验知识，移除小区域和中值滤波后得到最终结果.

将肿瘤部分进一步细分为不同子区域，临床医生从而能够计算出肿瘤体积、增强区比重以及坏死区比重等对GBM预后具有重要预测作用的高层次MR影像特征. 因此，GBM多模态MR图像的分割是建立GBM预后预测模型中基本而又关键的问题，是由图像处理到图像分析的关键步骤. 提出的基于混合特征和先验知识的GBM多模态MR图像分割算法利用GBM 3种模态MR图像，分别为T1加权像（T1 weighted image，T1WI）在注射对比剂之前的图像（T1WI-Pre）、T1WI在注射对比剂之后的图像（T1WI-Post）以及FLAIR图像. 之所以选择3种模态MR图像，是因为GBM肿瘤组织区域在这些模态图像中能够呈现不同的特征：肿瘤强化区在T1WI-Post模态图像中呈现为高信号，坏死区在T1WI-Post模态图像中呈现为低信号，而肿瘤非强化区在T1WI-Pre模态图像中呈现为低信号，肿瘤区域在FLAIR模态图像中呈现为高信号.

算法框架如图1所示，包含4个步骤：1）图像预处理，包括GBM多模态MR图像的配准以及利用N4ITK法进行偏置场校正；2）底层特征提取，提取GBM多模态MR图像的局部位置特征、强度特征、纹理特征、对称特征和上下文特征；3）粗分割，结合提取的特征应用随机森林分类器将GBM多模态MR图像初步分割为脑部正常组织区、坏死区、肿瘤强化区、不包括坏死区和肿瘤强化区的T1异常区以及不包括T1异常区的FLAIR异常区；4）图像后处理，考虑GBM肿瘤解剖结构的先验知识，移除小区域并进行中值滤波，得到肿瘤不同子区域的最终结果. 图中绿色区域（区域1）表示肿瘤强化区，蓝色区域（区域2）表示坏死区，红色区域（区域3）表示不包括坏死区和肿瘤强化区的T1异常区，黄色区域（区域4）表示不包括T1异常区的FLAIR异常区.

图像预处理包括GBM 3种模态MR图像的共同配准和偏置场校正. 为了保证图像中同一位置的体素对应脑部同一组织，需要对GBM 3种模态MR图像进行共同配准. 采用FSL软件^[11]的FLIRT配准工具完成对GBM 3种模态MR图像T1WI-Pre、T1WI-Post和FLAIR的共同配准. 对于每位GBM患者的3种模态MR图像，将切片数最少的作为配准参考图像；根据大脑及其反射计算对称平面^[12]，将头部区域的方位调整到近似未旋转的位置.

此外，为了消除MR图像中的场偏移，尽量使得同一幅MR图像中相同的绝对亮度代表相同的物理意义，在共同配准后采用N4ITK偏置场校正法对GBM 3种模态MR图像的每个切片进行不均匀场偏移校正，可直接利用现有的ANTS软件实现^[13]：

(1) ${{{i}}_{{\rm{cor}}}} = {\rm{N4ITK}}\,\,({{i}}).$

式中：i为偏置场校正前图像，i_cor为偏置场校正后图像. 如图2所示为利用N4ITK偏置场校正法对FLAIR模态MR图像其中1个切片进行不均匀场偏移校正的结果. 由图可以看出，偏置场校正前后GBM肿瘤组织区域的强度有所减弱（变暗），正常组织区域的强度有所增强（变亮）.

对于每种模态图像中位置（x，y，z）处的体素点，通过提取其局部位置特征、强度特征、纹理特征、对称特征和上下文特征来捕捉肿瘤不同子区域的判别依据信息，在此基础上利用基于随机森林分类器的监督机器学习算法将体素点初步归类为脑部正常组织区、坏死区、肿瘤强化区、不包括坏死区和肿瘤强化区的T1异常区或者不包括T1异常区的FLAIR异常区. 需要提取的底层特征如下：

1）局部位置特征f_loc. 描述被分析像素点位置信息的3个坐标x、y、z. 第i个像素点的局部位置特征的表达式为

(2) ${{{f}}_{{\rm{loc}}}}(i) = {[x,\; y,\; z]_i}.$

2）强度特征f_int. 不同模态MR图像设置的窗位不同，导致对应同一组织的体素的强度有差异. 为了减小其影响，利用K均值聚类算法估算出每个切片图像中灰质区域体素的平均强度；将同一切片图像中每个体素强度减去该切片图像中灰质区域体素的平均强度，将不同模态MR图像的灰度级尽可能映射至相似范围内，表达式如下：

(3) $\bar I(x,y,z) = I(x,y,z) - u{\rm{.}}$

式中： $I(x,y,z)$ 为（x，y，z）处的体素强度，u为灰质区域体素的平均强度.

GBM肿瘤5类子区域体素的强度通常在GBM各模态MR图像中互不相同. 因此，利用体素及其邻域体素的强度作为强度特征描述符f_int. 以T1WI-Pre模态MR图像中（x，y，z）处的体素为例，（x，y，z）、（x，y，z−1）、（x，y，z+1）、（x−1，y，z）、（x+1，y，z）、（x，y−1，z）、（x，y+1，z）共7个位置体素的强度被提取：

(4) $\begin{split}{{f}}_{{\rm{int}}}^{{\rm{T1WI - Pre}}} =& [I(x, y, z),\; I(x, y, z - {1}),\; I(x, y, z + {1}),\\ &I(x - {1},y, z),\; I(x + {1},y, z),\\ &I(x, y - {1},z),\; I(x, y + {1},z)].\end{split}$

此外，以这7个位置的体素分别为中心的3×3像素的7个图像块的平均强度被作为捕捉中心体素的邻域强度信息：

(5) $\begin{split} {\bar{ f}}_{{\rm{int}}}^{{\rm{T1WI - Pre}}} =& [\bar I{(x, y, z)_{{3} \times {3}}},\; \bar I{(x, y, z - {1})_{{3} \times {3}}},\; \bar I{(x, y, z + {1})_{{3} \times {3}}},\\ &\bar I{(x - {1},y, z)_{{3} \times {3}}},\; \bar I{(x + {1},y, z)_{{3} \times {3}}},\\ &\bar I{(x, y - {1},z)_{{3} \times {3}}},\; \bar I{(x, y + {1},z)_{{3} \times {3}}}]. \end{split} $

式中： $\bar I{( \cdot )_{{3} \times {3}}}$ 表示以某一位置处的体素为中心的3×3像素图像块的平均强度. 因此，在T1WI-Pre模态MR图像中提取了 ${{f}}_{{\rm{int}}}^{{\rm{T1WI - Pre}}}$ 和 ${\bar{ f}}_{{\rm{int}}}^{{\rm{T1WI - Pre}}}$ 共14个强度特征. 同样，对于T1WI-Post和FLAIR模态MR图像可以提取 ${{f}}_{{\rm{int}}}^{{\rm{T1WI - Post}}}$、 ${\bar{ f}}_{{\rm{int}}}^{{\rm{T1WI - Post}}}$ 以及 ${{f}}_{{\rm{int}}}^{{\rm{FLAIR}}}$、 ${\bar{ f}}_{{\rm{int}}}^{{\rm{FLAIR}}}$ 共28个强度特征. 从而对于（x，y，z）处的体素，总共提取42个强度特征：

(6) $\begin{split}{{{f}}_{{\rm{int}}}} =& \left[{{f}}_{{\rm{int}}}^{{\rm{T1WI - Pre}}},\; {\bar{ f}}_{{\rm{int}}}^{{\rm{T1WI - Pre}}},\; {{f}}_{{\rm{int}}}^{{\rm{T1WI - Post}}},\right.\\ &\left.{\bar{ f}}_{{\rm{int}}}^{{\rm{T1WI - Post}}},\; {{f}}_{{\rm{int}}}^{{\rm{FLAIR}}},\; {\bar{ f}}_{{\rm{int}}}^{{\rm{FLAIR}}}\right].\end{split}$

3）纹理特征f_tex. 由于GBM 3种模态MR图像上的异常脑组织外观通常比正常脑组织外观更不规则和不均匀，有研究人员将基于纹理分析的方法用于MR图像分割和脑组织表征. 本研究采用Leung-Malik（LM）滤波器组提取GBM 3种模态MR图像T1WI-Pre、T1WI-Post和FLAIR的纹理特征. LM由48个多尺度、多方向的滤波器组成，分别为4个高斯滤波器、8个拉普拉斯高斯滤波器、18个1阶高斯导数滤波器（3个尺度、6个方向）以及18个2阶高斯导数滤波器（3个尺度、6个方向）^[14]. 高斯滤波器是传统的平滑方法，可提供区域信息；拉普拉斯高斯滤波器是传统的基于2阶导数过零点的边缘检测方法；1阶高斯导数滤波器可用于检测图像边缘；2阶高斯导数滤波器可用于检测图像脊. 若g为高斯核函数，那么拉普拉斯高斯滤波器的表达式为

(7) ${f_{{\rm{LoG}}}} = {\nabla ^2}g = {g_{xx}} + {g_{yy}}.$

1阶高斯导数滤波器的表达式为

(8) ${g'} = {g_x}{\rm{cos}}\;\theta + {g_y}{\rm{sin}}\;\theta .$

2阶高斯导数滤波器的表达式为

(9) ${g^{''}} = {g_{xx}}{\rm{co}}{{\rm{s}}^2}\;\theta + {g_{yy}}{\rm{si}}{{\rm{n}}^2}\;\theta - 2{g_{xy}}{\rm{cos}}\;\theta {\rm{sin}}\;\theta .$

式中： ${g_x} = {{\partial g}}/{{\partial x}},\;$ ${g_y} = {{\partial g}}/{{\partial y}},$ ${g_{xy}} ={{{\partial ^2}g}}/{{\partial x\partial y}},\; {g_{xx}} =$ ${{{\partial ^2}g}}/{{\partial {x^2}}}, $ ${g_{yy}} ={{{\partial ^2}g}}/{{\partial {y^2}}},$ $\theta $为高斯核函数g的方向.

对于（x，y，z）处的体素，利用LM滤波器提取该体素及其49×49像素邻域内（滤波器的大小为49×49像素）的特征值来代表该像素点的局部纹理特征，可以较好地反映图像的局部纹理特征. 因此，从GBM 3种模态MR图像中总共提取到144个纹理特征，即

(10) ${{{f}}_{{\rm{tex}}}} = [g{({{{i}}_{{\rm{cor}}}})_{{\rm{12}}}}, {f_{{\rm{LoG}}}}{({{{i}}_{{\rm{cor}}}})_{{\rm{24}}}}, {g'}{({{{i}}_{{\rm{cor}}}})_{{\rm{54}}}}, {g^{''}}{({{{i}}_{{\rm{cor}}}})_{{\rm{54}}}}].$

式中：下标12、24、54表示提取的特征个数. 如图3所示为利用LM滤波器组对FLAIR模态MR图像中1个切片进行滤波处理后得到的48个响应.

4）对称特征f_sym. 根据脑部解剖结构先验知识可知，正常大脑沿头部中心线大致对称. 因此，定位不对称区域可用于检测脑部异常组织. 利用文献[15]提出的方法计算对称平面，并利用K均值聚类算法生成聚类图像，提取原始图像和聚类图像的对称特征. 对于每种模态图像中（x，y，z）处的体素及其在大脑另一侧的对应体素，计算以下6个特征来描述对称性差异. 以FLAIR模态MR图像中（x，y，z）处的体素为例，其表达式为

(11) $\begin{split}{{f}}_{{\rm{sym}}}^{{\rm{FLAIR}}} =& \left[d_{{\rm{int}}}^{{\rm{FLAIR}}},\; d_{{\rm{seg}}}^{{\rm{FLAIR}}},\; d_{{\rm{cp}}}^{{\rm{FLAIR}}},\right.\\ &\left.d_{{\rm{sc}}}^{{\rm{FLAIR}}},\; d_{{\rm{sl}}}^{{\rm{FLAIR}}}{\rm{,\; }}d_{{\rm{clu}}}^{{\rm{FLAIR}}}\right].\end{split}$

式中： $d_{{\rm{int}}}^{{\rm{FLAIR}}}$ 为原图像中该体素与对称平面另一侧对应体素的强度差异， $d_{{\rm{seg}}}^{{\rm{FLAIR}}}$ 为聚类图像中该体素与对称平面另一侧对应体素的强度差异， $d_{{\rm{cp}}}^{{\rm{FLAIR}}}$ 为原图像中该体素沿着x轴的强度分布的对称性， $d_{{\rm{sc}}}^{{\rm{FLAIR}}}$ 为聚类图像中该体素沿着x轴的强度分布的对称性， $d_{{\rm{sl}}}^{{\rm{FLAIR}}}$ 为该像素位置所属的聚类对称性， $d_{{\rm{clu}}}^{{\rm{FLAIR}}}$ 为其聚类面积百分比. 同样，对于T1WI-Pre和T1WI-Post模态MR图像中 （x，y，z）处的体素可以提取到对称特征 ${{f}}_{{\rm{sym}}}^{{\rm{T1WI - Pre}}}$ 和 ${{f}}_{{\rm{sym}}}^{{\rm{T1WI - Post}}}$，表达式为

(12) $\begin{split}{{f}}_{{\rm{sym}}}^{{\rm{T1WI - Pre}}} =& \left[d_{{\rm{int}}}^{{\rm{T1WI - Pre}}},\; d_{{\rm{seg}}}^{{\rm{T1WI - Pre}}},\; d_{{\rm{cp}}}^{{\rm{T1WI - Pre}}},\right.\\ &\left.d_{{\rm{sc}}}^{{\rm{T1WI - Pre}}},\; d_{{\rm{sl}}}^{{\rm{T1WI - Pre}}}{\rm{,\; }}d_{{\rm{clu}}}^{{\rm{T1WI - Pre}}}\right],\end{split}$

(13) $\begin{split}{{f}}_{{\rm{sym}}}^{{\rm{T1WI - Post}}} =& \left[d_{{\rm{int}}}^{{\rm{T1WI - Post}}},\; d_{{\rm{seg}}}^{{\rm{T1WI - Post}}},\; d_{{\rm{cp}}}^{{\rm{T1WI - Post}}},\right.\\& \left.d_{{\rm{sc}}}^{{\rm{T1WI - Post}}},\; d_{{\rm{sl}}}^{{\rm{T1WI - Post}}}{\rm{,\; }}d_{{\rm{clu}}}^{{\rm{T1WI - Post}}}\right].\end{split}$

因此，对于GBM 3种模态的MR图像，在（x，y，z）处计算得到的对称特征f_sym共18个，即

(14) ${{{f}}_{{\rm{sym}}}} = \left[{{f}}_{{\rm{sym}}}^{{\rm{T1WI - Pre}}},\; {{f}}_{{\rm{sym}}}^{{\rm{T1WI - Post}}},\; {{f}}_{{\rm{sym}}}^{{\rm{FLAIR}}}\right].$

如图4所示，w为图像的长度，h为体素的强度. 如图4（a）所示为某对应的体素点对（*点）以及脑部的对称中心线（虚线），如图4（b）所示为聚类图像，如图4（c）所示为沿着图4（a）中水平实线的体素点强度分布情况.

5）上下文特征f_con. 上下文特征通常可以为医学图像的分割任务提供重要的邻域信息. 基于聚类图像，通过1个体素点垂直和水平方向邻近集群的4个体素来捕捉上下文信息. 如图5所示， *点表示坏死区域的1个体素，强度小于与其邻近集群的4个体素强度，表明坏死区域的体素强度比邻近区域低. 因此，可将每种模态MR图像中（x，y，z）处的体素与其邻近集群的4个体素之间的强度差异作为上下文特征. 以T1WI-Post模态图像中（x，y，z）处的体素为例，其上下文特征可通过式（15）~（18）计算得出：

(15) $\Delta I_{{\rm{abo}}}^{{\rm{T1WI - Post}}}{\rm{ = }}I_{{\rm{LBI}}}^{{\rm{T1WI - Post}}}(x, y) - I_{{\rm{LBI}}}^{{\rm{T1WI - Post}}}({x_{{\rm{las}}}}, y),$

(16) $\Delta I_{{\rm{bot}}}^{{\rm{T1WI - Post}}}{\rm{ = }}I_{{\rm{LBI}}}^{{\rm{T1WI - Post}}}(x, y) - I_{{\rm{LBI}}}^{{\rm{T1WI - Post}}}({x_{{\rm{nex}}}}, y),$

(17) $\Delta I_{{\rm{lef}}}^{{\rm{T1WI - Post}}}{\rm{ = }}I_{{\rm{LBI}}}^{{\rm{T1WI - Post}}}(x, y) - I_{{\rm{LBI}}}^{{\rm{T1WI - Post}}}(x, {y_{{\rm{las}}}}),$

(18) $\Delta I_{{\rm{rig}}}^{{\rm{T1WI - Post}}}{\rm{ = }}I_{{\rm{LBI}}}^{{\rm{T1WI - Post}}}(x, y) - I_{{\rm{LBI}}}^{{\rm{T1WI - Post}}}(x, {y_{{\rm{nex}}}}).$

即该模态MR图像的上下文特征 ${{f}}_{{\rm{con}}}^{{\rm{T1WI - Post}}}$ 的表达式为

(19) $\begin{split}{{f}}_{{\rm{con}}}^{{\rm{T1WI - Post}}} =& \left[\Delta I_{{\rm{abo}}}^{{\rm{T1WI - Post}}},\; \Delta I_{{\rm{bot}}}^{{\rm{T1WI - Post}}},\right.\\ & \left.\Delta I_{{\rm{lef}}}^{{\rm{T1WI - Post}}},\; \Delta I_{{\rm{rig}}}^{{\rm{T1WI - Post}}}\right].\end{split}$

式中： $I_{{\rm{LBI}}}^{{\rm{T1WI - Post}}}(x, y)$ 为T1WI-Post模态聚类图像中（x，y，z）处的体素强度； $I_{{\rm{LBI}}}^{{\rm{T1WI - Post}}}({x_{{\rm{las}}}}, y)$、 $I_{{\rm{LBI}}}^{{\rm{T1WI - Post}}}({x_{{\rm{nex}}}}, y)$、 $I_{{\rm{LBI}}}^{{\rm{T1WI - Post}}}(x, {y_{{\rm{las}}}})$、 $I_{{\rm{LBI}}}^{{\rm{T1WI - Post}}}(x, {y_{{\rm{nex}}}})$ 分别为该体素上、下、左、右邻近集群的4个体素强度. 同理，对于T1WI-Pre和FLAIR模态MR图像中（x，y，z）处的体素，亦可提取到上下文特征 ${{f}}_{{\rm{con}}}^{{\rm{T1WI - Pre}}}$ 和 ${{f}}_{{\rm{con}}}^{{\rm{FLAIR}}}$，表达式分别为

(20) $\begin{split}{{f}}_{{\rm{con}}}^{{\rm{T1WI - Pre}}} =& \left[\Delta I_{{\rm{abo}}}^{{\rm{T1WI - Pre}}},\; \Delta I_{{\rm{bot}}}^{{\rm{T1WI - Pre}}},\right.\\ & \left.\Delta I_{{\rm{lef}}}^{{\rm{T1WI - Pre}}},\;\Delta I_{{\rm{rig}}}^{{\rm{T1WI - Pre}}}\right],\end{split}$

(21) ${{f}}_{{\rm{con}}}^{{\rm{FLAIR}}} = \left[\Delta I_{{\rm{abo}}}^{{\rm{FLAIR}}},\; \Delta I_{{\rm{bot}}}^{{\rm{FLAIR}}},\; \Delta I_{{\rm{lef}}}^{{\rm{FLAIR}}},\; \Delta I_{{\rm{rig}}}^{{\rm{FLAIR}}}\right].$

根据GBM 3种模态MR图像可以计算得到 ${{f}}_{{\rm{con}}}^{{\rm{T1WI - Pre}}}$、 ${{f}}_{{\rm{con}}}^{{\rm{T1WI - Post}}}$、 ${{f}}_{{\rm{con}}}^{{\rm{FLAIR}}}$ 共12个上下文特征.

综上所述，对于位置（x，y，z）处的体素，总共提取了219个底层特征，分别为3个局部位置特征、42个强度特征、144个纹理特征、18个对称特征以及12个上下文特征. 随机森林分类器本身有特征优化的功能，能自动随机选取重要的特征，故不再单独实施特征优化步骤. 另外，在提取当前切片强度特征时使用了与其相邻的前后2个切片的强度信息，因此不再提取第1个切片和最后1个切片的底层特征. 底层特征的表达式为

(22) ${{{f}}_{{\rm{low}}}} = {\rm{[(}}{{{f}}_{{\rm{loc}}}}{)_{3}},({{{f}}_{{\rm{int}}}}{)_{{\rm{42}}}},({{{f}}_{{\rm{tex}}}}{)_{{\rm{144}}}},({{{f}}_{{\rm{sym}}}}{)_{{\rm{18}}}},({{{f}}_{{\rm{con}}}}{)_{{\rm{12}}}}].$

式中：3、42、144、18、12表示各类特征的个数.

在完成底层特征提取后，应用随机森林（random forest，RF）分类器对每个位置的体素进行分类，分为脑部正常组织区、坏死区、肿瘤强化区、不包括坏死区和肿瘤强化区的T1异常区或者不包括T1异常区的FLAIR异常区. RF分类器是综合考虑多个决策树而形成的集成分类器模型，在RF生成之后，取1个新样本进行测试，森林中每棵树都会独立判定分类结果，最终决策取相同判定最多的分类类别，在运算量没有显著增大的前提下能够提高预测精度. 投票决策过程如下：

(23) $H(x) = \mathop {{\text{argmax}}}\limits_Y \sum\limits_{i = 1}^k {I({h_i}(x) = Y)} .$

式中：H（x）为组合分类模型；h_i（x）为单个分类树模型；Y为输出变量； $I( \cdot )$ 为示性函数，即当集合内有此数时，函数值为1，否则为0.

为了进一步改善GBM 3种模态MR图像的分割性能，考虑GBM肿瘤解剖结构的先验知识，对初步分割结果进行精分割. 通常情况下肿瘤区域的体积大于200体素，因此所有小于200体素的肿瘤区域将被视为假阳性，予以移除. 运用5×5像素的中值滤波器平滑每类区域的边界. 只考虑没有切除后腔的GBM 3种模态MR图像，根据GBM肿瘤的典型构型，将由肿瘤强化区包围的所有体素重新分类为坏死区. 经过上述后处理所得结果为最终的肿瘤不同子区域分割结果.

利用国际上公开的GBM多模态MR图像数据集TCGA，建立国外人群GBM实验测试库（TCGA-GBM），使用了其中76位GBM患者的3种模态MR图像. TCGA由美国国立卫生院下属美国国家肿瘤研究院和美国国家人类基因组研究院共同发起，全部数据免费对公众开放（http://www.cancerimagingarchive.net/）. 另外，为了验证该算法分割中国人群GBM患者多模态MR图像的有效性，通过本课题组具有丰富GBM临床影像诊断经验的成员整理和收集的89位GBM患者的3种模态MR图像，建立中国人群GBM实验测试库（CH-GBM）. 实验数据的详细信息如表1所示，由于噪声和光照条件的影响，图像数据的品质有所差别. 为了获得2个实验测试库的真值图，本课题组在前期基于Matlab平台开发了GBM多模态MR图像标注软件，由课题组具有丰富GBM临床影像诊断经验的成员对每位GBM患者的3种模态MR图像数据进行手动标注.

在扫描整个GBM多模态MR图像序列时，脑部正常组织区通常比脑部肿瘤组织区大许多. 为了尽可能保持数据平衡，在训练RF分类器时，从脑部正常组织区随机选取50 000个体素，从坏死区、增强肿瘤区、不包括坏死区和增强肿瘤区的T1异常区以及不包括T1异常区的FLAIR异常区各随机选取5 000个体素. RF分类器由100棵决策树构成，每棵决策树随机选取5个变量（即从所有特征中随机选取5个特征）作为输入.

采用留一交叉验证法（leave-one-out cross-validation）训练和测试RF分类器，即每次从样本数据集中抽取1个不同的样本作为测试集，剩余的样本作为训练集，重复此过程直至每个样本均被测试1次. 对于TCGA-GBM，在进行留一交叉验证法时每次使用75位患者的GBM 3种模态MR图像作为训练集，1位患者的GBM 3种模态MR图像作为测试集；训练和测试过程重复76次，直到完成数据库中每位GBM患者3种模态MR图像所有特征的预测. 同样的，对CH-GBM也进行类似操作，在进行留一交叉验证法时每次均使用88位患者的GBM 3种模态MR图像作为训练集，1位患者的GBM 3种模态MR图像作为测试集.

采用Dice相似性系数D_dic、敏感度S_sen和特异性S_spe作为评价指标，表达式分别为

(24) ${D_{{\rm{dic}}}}({A} {\rm{,}}\;{M} ) = 2|{A} \cap {M} |/(|{A} | + |{M} |),$

(25) ${S_{{\rm{sen}}}} = {\rm{TP/(TP}} + {\rm{FN)}},\quad\quad\quad\quad\quad$

(26) ${S_{{\rm{spe}}}} = {\rm{TN/(TN}} + {\rm{FP)}}.\quad\quad\quad\quad\quad$

式中：A为自动分割轮廓区域所包围的点集；M为手动分割轮廓区域所包围的点集；TP表示真阳性，FP表示假阳性，TN表示真阴性，FN表示假阴性. Dice相似性系数越大，表明分割结果越好.

相关课题的前期研究表明当K=8时分类效果最佳，因此在利用K均值聚类算法生成聚类图像时，K=8. 运行平台的配置及Matlab版本如下：运行平台的处理器为Intel （R）Xeon（R）Bronze 3106 CPU @ 1.7 OGHz，操作系统为Windows Server 2008 R2企业版，开发测试软件平台为Matlab R 2014 b.

如图6、7所示为所提出算法分别对TCGA-GBM数据库和CH-GBM数据库中3个示例进行自动分割所得结果. 第1列为T1WI-Pre原图像，第2列为T1WI-Post原图像，第3列为FLAIR原图像，第4列为将所提算法分割结果叠加至FLAIR原切片上的图像，第5列为将专家手动分割结果叠加至FLAIR原切片上的图像. 其中，绿色区域（区域1）表示肿瘤强化区，蓝色区域（区域2）表示坏死区，红色区域（区域3）表示不包括坏死区和肿瘤强化区的T1异常区，黄色区域（区域4）表示不包括T1异常区的FLAIR异常区.

如表2所示为利用所提算法分别对TCGA-GBM和CH-GBM数据库进行分割，并使用Dice相似性系数D_dic、敏感度S_sen、特异性S_spe进行性能评估的结果. 表中，C区为整个肿瘤区域，C₁区为肿瘤强化区，C₂区为坏死区，C₃区为不包括坏死区和肿瘤强化区的T1异常区，C₄区为不包括T1异常区的FLAIR异常区. 可以看出，利用所提算法对CH-GBM数据库的测试结果略优于对TCGA-GBM数据库的测试结果，但分割结果差异不大，几乎对算法性能没有影响. 此外，表中还列出了利用所提算法对TCGA-GBM和CH-GBM数据库合并成的TCGA-CH-GBM数据库中的3种模态MR图像进行分割的结果.

近几年一些研究人员在试图将GBM肿瘤分割成不同子区域时使用了模态数量不同的MR图像和数据集，很难将所提算法与这些方法进行直接比较. 文献[16]、[17]不仅使用了本研究所提的3种模态的MR图像，还使用了T2模态MR图像，并且均只对10个示例进行测试. 本研究所提算法的平均Dice相似性系数稍大于文献[16]、[17]的结果. 文献[18]使用了包含40位患者MR图像的数据集，所得Dice相似性系数比所提算法得出结果稍高，但它除了要使用T1WI-Pre、T1WI-Post和FLAIR 3种模态MR图像外，还需利用T2-TSE、DTI-p和DTI-q 3种模态MR图像，即该方法需要借助6种模态MR图像，而这对分割精度具有重要影响，因此很难对本研究所提算法和文献[18]所提算法的性能进行直接比较.

本研究提出了基于混合特征和先验知识的GBM多模态MR图像分割算法，能够在适应临床实践需求的前提下将整个肿瘤分为脑部正常组织区、坏死区、肿瘤强化区、不包括坏死区和肿瘤强化区的T1异常区以及不包括T1异常区的FLAIR异常区，并对TCGA-GBM和CH-GBM数据库中3种模态MR图像进行分割. 以Dice相似性系数作为评价指标，利用所提算法分割TCGA-GBM数据库的整个肿瘤及其不同子区域的平均Dice相似性系数分别为0.871、0.863、0.759、0.775、0.686；分割CH-GBM数据库的整个肿瘤及其不同子区域的平均Dice相似性系数分别为0.882、0.872、0.776、0.783、0.692，在一定程度上优于现有的一些典型算法^[16-17]. 所提算法能以较高的准确率对GBM多模态MR图像进行分割，适用于临床GBM多模态MR图像分割任务. 如何进一步提高对GBM多模态MR图像坏死区的分割精度是今后的主要研究方向.