CN103295236B - 马尔可夫多特征随机场模型构建方法及其脑部mr图像分割技术 - Google Patents

Abstract

本发明公开了一种马尔可夫多特征随机场(MMFRF)模型构建方法及基于该模型的脑部MR图像分割技术。MMFRF模型构建方法中，本发明将传统MRF模型中单个特征随机场扩展为多个特征随机场情形，建立了对应的概率估计准则及优化求解方式。基于MMFRF模型的脑部MR图像分割技术包括：(a)将整个脑部图像域分成多个子域并将其建模成多个局部MMFRF模型，通过这些局部模型构建面向整个脑部域的组合MMFRF模型；(b)将组合MMFRF模型应用于脑部MR图像分割。本发明提供的MMFRF模型，适用于多种计算机视觉和图像分析应用；本发明提供的基于MMFRF模型的脑部MR图像分割技术实验结果准确、稳定。

Description

马尔可夫多特征随机场模型构建方法及其脑部MR图像分割技术

技术领域

本发明涉及一种使用马尔可夫随机场实现特征提取的模型构建方法和脑部MR图像分割技术，更具体地说，涉及用于纹理、灰度等多种特征提取的新颖的马尔可夫模型构建方法和基于该模型方法的脑部MR图像分割技术。背景技术

MRF理论提供了一种便利的、一致性的方法来建模诸如图像像素和相关性特征等物理现象间的上下文相关性，它从统计学的角度出发，对图像进行建模，把图像中的像素灰度值看成满足一定概率分布的随机变量。因为从条件MRF分布很难得到联合分布，Hammersley和Clifford提出了MRF和Gibbs分布的等价性理论，他们指出MRF的联合分布是一个具有简单形式的Gibbs分布，从而通过Gibbs势能的概念，利用规模较小的局部计算获得了全局性的结果，大大降低了MRF模型的计算复杂度。为了通过已有的优化理论公式化目标函数，MRF理论往往与统计决策和估计理论联合使用。其中，MAP概率是MRF建模中使用最为广泛的统计准则，而MAP-MRF框架的目标则可以表示成MRF标号的联合后验概率形式。基于MRF模型的贝叶斯框架包括似然项和先验项，前者称之为特征随机场，往往被用于描述诸如图像像素灰度的高斯分布等观测信息，后者称之为标号随机场，通常被用来捕捉诸如图像像素等站点的MRF先验知识。

一种特殊的MRF——隐马尔可夫随机场(HMRF)，是一个由MRF产生的随机过程，它的状态序列不能直接观测到，但能通过观测值随机场推断。基于MRF模型或者HMRF模型的优良品质，国内外许多专家和学者将其应用于求解模式识别和图像分析问题，如脑部MR图像分割问题。为方便记，将已有的MRF模型或者HMRF模型统称为传统的MRF模型。尽管基于传统MRF模型在计算机视觉和图像分析问题中取得了不错的研究成果，但是它们均存在一个明显的不足，模型只考虑了目标对象的灰度观测信息和上下文先验信息，忽略了目标中其它有用的观测特征信息。而实际上，一个目标对象，例如一幅图像特别是脑部MR图像，存在多种具有不同形式的特征模式，除像素灰度特征外，还存在纹理特征、形状特征等等。之后的研究工作并没有融合多种特征信息的新颖的MRF模型理论成果出现。

发明内容

本发明的一个目的是提供一种新颖的马尔可夫多特征随机场(MMFRF)模型构建方法，是为了更有效地捕获诸如图像等物理对象中存在的不同模式的特征信息，从而更加有效地应用于图像存储和分割、图像重建、边缘检测等计算机视觉和图像分析应用。

本发明的另一个目的是提供一种基于MMFRF模型的脑部MR图像分割技术，通过具体的构建纹理特征随机场、灰度特征随机场和标号随机场来获取脑部MR图像中存在的纹理、灰度和空间相关性等特征信息，有效、稳定的实现脑部MR图像分割。

为了达成上述目的，本发明提供一种MMFRF模型构建方法，包括如下步骤：

(1)定义，包括标号随机场、多特征随机场的定义；

(2)模型理论构建，给定多特征随机场集合Г情形下，通过Bayesian理论，得到标号随机场X的组态的后验概率计算方式；

(3)通过最大后验(MAP)概率估计准则对模型中标号随机场X的组态的后验概率进行估计；

(4)通过迭代条件模式(ICM)优化算法对的MAP估计进行优化求解，得到模型中标号随机场X的最终的标号值；

得到MMFRF模型。

其中，步骤(1)所述的定义，包括以下步骤：

1)标号随机场：

随机场X为集合S上关于邻域***N的马尔可夫随机场(MRF)，记为标号随机场；

2)多特征随机场：

记A＝{1，...，n}表示随机场个数，S＝{1，2，...，z}定义为z个站点索引位置的集合，随机变量组为定义在集合S上的一个特征随机场，其中每个随机变量取其状态空间中的一个值 表示随机场F_a的一个实现组态，F_a用于描述诸如图像像素和图像特征等目标对象的观测信息，记作特征随机场，Г＝{F_a，a∈A}表示集合S上的多个特征随机场的集合，称之为多特征随机场。

其中，步骤(2)所述的X的组态的后验概率计算方式：

给定特征随机场F_a组态f_a前提下组态x的后验概率P(x|f_a)为

$P\left(x\right|f_a)\∝P\left(x\right)\underset{i\∈S}{\mathrm{\Π}}{p}^a\left(y_i^a\right|x_i)$

其中，f_a中随机变量为条件独立的且每个变量都服从统一的依赖于x_i的条件密度函数P(x)为组态x的先验概率；在已知Г前提下，标号随机场X的一个组态的后验概率计算公式如下

$P\left(x\right|f_1,...,f_n)=\frac{P\left(f_1\right|x)...P\left(f_n\right|x)P\left(x\right)}{P(f_1,...,f_n)}\∝\left[\underset{a=1}{\overset{n}{\mathrm{\Π}}}\underset{i=1}{\overset{z}{\mathrm{\Π}}}{p}^a\left(y_i^a\right|x_i)\right]P\left(x\right).$

其中，步骤(3)所述的标号随机场X的MAP估计方式为：

在给定多个特征随机场情形下估计标号随机场X的真实标号x^*的计算公式为

$x^{*}=\underset{x\∈L^n}{\arg \max }\left\{\right[\underset{a=1}{\overset{n}{\mathrm{\Π}}}\underset{i=1}{\overset{z}{\mathrm{\Π}}}{p}^a\left(y_i^a\right|x_i)\left]P\left(x\right)\right\}$

其中，L＝{1，...，G}表示站点标号的集合。

为了达成上述目的，本发明提供一种基于MMFRF模型的脑部MR图像分割技术，包括如下步骤：

(1)使用子卷概率图谱(SVPAs)将整个脑部图像域分成多个子域的集合，将每个子域建模成一个局部MMFRF模型，最后通过这些局部模型构建统一的面向整个脑部域的组合MMFRF模型；

(2)计算组合MMFRF模型中的灰度特征随机场；

(3)计算组合MMFRF模型中的纹理特征随机场；

(4)计算组合MMFRF模型中的标号随机场；

(5)根据计算的多个特征随机场和标号随机场，使用ICM算法通过MAP估计实现图像分割。

其中，步骤(1)所述的构建面向整个脑部域的组合MMFRF，包括以下步骤：

1)MMFRF模型：

标号随机场X的真实标号x^*的MAP估计计算如下

$x^{*}=\underset{x\∈L^n}{\arg \max }\left\{\right[\underset{a=1}{\overset{n}{\mathrm{\Π}}}\underset{i=1}{\overset{z}{\mathrm{\Π}}}{p}^a\left(y_i^a\right|x_i)\left]P\left(x\right)\right\}$

其中，L＝{1，...，G}表示图像体素标号的集合，A＝{1，...n}为特征随机场个数，S＝{1，...z}为脑部MR图像体素集合，F_a表示一个特征随机场；表示特征随机场似然概率项，P(x)表示标号随机场先验概率项，标号随机场能量函数为

$U\left(x\right)=\underset{i=1}{\overset{z}{\mathrm{\Σ}}}(O\left(x_i\right|i)+\frac{1}{2}\underset{j\∈N_i}{\mathrm{\Σ}}Q(x_i,x_j|i,j))$

其中，N_i为体素i的邻域，O(x_i|i)和Q(x_i，x_j|i，j)分别表示两个分离的first-order和second-order先验能量函数；O(x_i|i)表示体素i为单一组织类别的概率，Q(x_i，x_j|i，j)表示给定某一体素其邻域体素组织类型的前提下该体素的组织类别概率；对于所有i和j，假设Q(x_i，x_j|i，j)＝Q(x_j，x_i|j，i)；

2)组合MMFRF模型：

将整个脑部图像域通过SVPAs划分成小脑、皮层下结构等12个脑部子结构；对于每个子域R_m，m＝1，...，M，定义特征随机场的区域项，P_m(x)为标号随机场P(x)的区域项；分别将所有区域特征随机场项和所有区域标号随机场项P_m(x)(m＝1，...M)合并成组合特征随机场和组合标号随机场；

组合特征随机场 $p^a\left(y_i^a\right|x_i)(a\∈A,A=\{1,...,n\left\}\right)$ 的计算如下：

$p^a\left(y_i^a\right|x_i)=\underset{m=1}{\overset{M}{\mathrm{\Σ}}}{r}_i\left(m\right)p_m^a\left(y_i^a\right|x_i)$

组合标号随机场的先验能量函数计算如下：

$U\left(x\right)=-\underset{i=1}{\overset{z}{\mathrm{\Σ}}}\{\ln \left(\underset{m=1}{\overset{M}{\mathrm{\Σ}}}{r}_i\left(m\right)\exp \right[-O_m\left(x_i\right|i)\left]\right)-\frac{1}{2}\underset{j\∈N_i}{\mathrm{\Σ}}Q(x_i,x_j)\}$

其中，O_m(x_i|i)和Q_m(x_i，x_j|i，j)分别表示先验能量函数O(x_i|i)和Q(x_i，x_j|i，j)的区域项。对于不同子区域，假设每个区域拥有相同的二阶MRF项，即对于所有区域R_m，m＝1，...，M，Q_m(x_i，x_j|i，j)＝Q(x_i，x_j)；

脑部MR图像分割的组合MMFRF模型组态x的后验概率计算如下：

$P\left(x\right|f_1,...,f_n)\∝\underset{i=1}{\overset{z}{\mathrm{\Π}}}\underset{a=1}{\overset{n}{\mathrm{\Π}}}\left[\underset{m=1}{\overset{M}{\mathrm{\Σ}}}{r}_i\left(m\right)p_m^a\left(y_i^a\right|x_i)\right]\×\left(\underset{m=1}{\overset{M}{\mathrm{\Σ}}}{r}_i\left(m\right)\exp \right[-O_m\left(x_i\right|i)-\frac{1}{2}\underset{j\∈N_j}{\mathrm{\Σ}}Q(x_i,x_j)\left]\right).$

其中，步骤(2)所述的计算组合MMFRF模型中的灰度特征随机场，包括以下步骤：

1)计算似然项包括(a)、(b)两种情形：

(a)纯组织类型体素，区域灰度似然函数计算如下

$p_m^1\left(y_i^1\right|{\mathrm{\θ}}_{\mathrm{lm}},x_i=l)=[\frac{1}{\sqrt{2\mathrm{\π}{\mathrm{\σ}}_{\mathrm{lm}}^2}}\exp \left(\frac{{(y_i^1-{\mathrm{\μ}}_{\mathrm{lm}})}^2}{{2\mathrm{\σ}}_{\mathrm{lm}}^2}\right)$

其中，表示在区域R_m中类别为l时的Gaussian参数集合，μ_m和分别表示均值和方差；

(b)混合组织类型体素，计算如下

$p_m^1\left(y_i^1\right|{\mathrm{\θ}}_{l_{1}m},{\mathrm{\θ}}_{l_{2}m},x_i=l)={\∫}_0^1\frac{\exp (-\frac{{[y_i^1-w{\mathrm{\μ}}_{l_{1}m}-(1-w){\mathrm{\μ}}_{l_{2}m}]}^2}{2[w^2{\mathrm{\σ}}_{l_{1}m}^2+{(1-w)}^2{\mathrm{\σ}}_{l_{2}m}^2]})}{\sqrt{2\mathrm{\π}[w^2{\mathrm{\σ}}_{l_{1}m}^2+{(1-w)}^2{\mathrm{\σ}}_{l_{2}m}^2]}}\mathrm{dw}$

其中，和分别表示在区域R_m中类别为l₁和l₂时的包含均值和方差的参数集合，l为图像体素标号类别；

2)参数估计：采用基于加权似然(WL)方法对图像参数进行估计，计算公式如下

${\mathrm{\Θ}}_m^{*}=\arg \underset{{\mathrm{\Θ}}_m}{\max }\underset{i=1}{\overset{z}{\mathrm{\Σ}}}{r}_i\left(m\right)\ln \left[\underset{l\∈L}{\mathrm{\Σ}}{\mathrm{\α}}_{\mathrm{lm}}{p}_m^1\left(y_i^1\right|{\mathrm{\θ}}_{\mathrm{lm}},l)\right]$

其中，表示图像灰度参数集合，α_lm表示组织类别l在区域R_m中的成分概率，满足μ_lm和分别表示组织类别为l时的均值和方差。

其中，步骤(3)所述的计算组合MMFRF模型中的纹理特征随机场，假设对于所有的子域R_m，m＝1，...，M，纹理似然能量函数是均一的，即 $p_m^2\left(y_i^a\right|x_i)=p^2\left(y_i^a\right|x_i),$ 其计算步骤如下：

1)特征提取，包括：

(a)对于图像中的每个体素i∈S，设置其感兴趣域(VOI)来提取基于GLCM的三维纹理特征；

(b)对于每个VOI，选择13个方向和两种距离(d＝1，2)，共26组不同偏移参数进行计算，共得到26个不同的共生矩阵；

(c)为每个矩阵计算角二阶矩、对比度和熵等12个基于GLCM的统计特征值，取26个统计值的平均值作为每个体素的最终特征值；

2)纹理特征随机场计算，包括以下步骤：

(a)为每个体素设定3×3×3的邻域；

(b)分别计算体素i与其邻接体素提取的纹理特征向量间的欧式距离，并将其作为衡量体素i与其邻接体素间相似性的评价尺度；

(c)对于每个体素i，根据其邻域窗口W_i内邻接体素的组织类型标号，利用计算的欧式距离，通过K最近邻(KNN)算法对体素i属于哪种组织类型进行投票，计算公式如下：

$p^2\left(y_i^2\right|x_i,y_{W_i}^2)=p^2(y_i^2={\mathrm{tf}}_i|x_i=l,y_{W_i}^2={\mathrm{tf}}_{W_i})$

$=\frac{\mathrm{ofvoxelslabeledaslofknearestneighbor}\sin W_i}{k}$

其中，表示体素i在给定其标号为l以及其邻域窗口内邻接体素纹理特征向量集为的前提下，它的纹理特征向量为tr_i的概率。k表示在邻域窗口W_i内的最邻接体素个数。

其中，步骤(4)所述的计算组合MMFRF模型中的标号随机场，包括以下步骤：

1)O_m函数计算如下

$O_m\left(x_i\right|i)=-\ln {\mathrm{\α}}_{\mathrm{lm}}^{*},x_i=l$

其中，表示α_lm的加权似然估计；

2)Q函数计算如下

$Q(x_i,x_j)=\frac{{\mathrm{\β\τ}}_{\mathrm{ij}}}{d(i,j)}$

其中，β是一个可调的实验参数，通过经验选取。d(i，j)表示体素i和体素j中心之间的距离。MRF的second-order相关性值τ_ij计算如下

${\mathrm{\τ}}_{\mathrm{ij}}=\left\{\begin{array}{c}-1,\mathrm{if}{x}_i=x_j,\\ 0,\mathrm{if}{x}_i\mathrm{and}{x}_j\mathrm{shareacomponent},.\\ 1,\mathrm{otherwise}.\end{array}\right.$

本发明的有益效果如下：(a)本发明构建的新颖的基于多个特征随机场的MMFRF模型，有效的融合了目标对象中的多种观测特征信息，并将其通过贝叶斯框架组合成统一的概率决策模型，建立了与MMFRF模型对应的MAP估计准则及ICM算法优化求解模式。传统的MRF模型成为本发明MMFRF模型取一个特征场(n＝1)时的特例，MMFRF模型是其扩展和提升。本发明提出的MMFRF模型适用于图像存储和分割、图像重建、边缘检测等多种计算机视觉和图像分析应用。(b)本发明构建的基于MMFRF模型的脑部MR图像分割技术，充分利用了MMFRF模型捕获图像中多种特征信息的优势。进一步地，利用SVPAs将整个脑部域分区成多个子域的集合，并通过MMFRF模型定义在这些局部的子域上组合成统一的基于MMFRF的全局模型，充分利用了图像的局部信息。结合以上两者，使得本发明的分割技术获得了更加全面、更加丰富的图像特征信息和图像局部细节信息，判别能力更强，分割效果更好、更稳定。

附图说明

图1是本发明实施例中所述MMFRF模型及其在脑部MR图像分割应用中的流程图。

图2是本发明实施例中所述一种新颖的MMFRF模型流程图。

图3是本发明实施例中所述基于MMFRF模型的脑部MR图像分割技术流程图。

具体实施方式

下面结合附图对本发明的具体实施方式进行详细的说明：

实施例

如图1所示，为MMFRF模型理论及其在脑部MR图像分割应用中的流程图。

如图2所示，MMFRF模型构建方法，包括如下步骤：

步骤201、定义，包括以下内容：

1)标号随机场：

随机场X为集合S上关于邻域***N的马尔可夫随机场(MRF)，记为标号随机场；

2)多特征随机场：

记A＝{1，...，n}表示随机场个数，S＝{1，2，...，z}定义为z个站点索引位置的集合，随机变量组为定义在集合S上的一个特征随机场，其中每个随机变量取其状态空间中的一个值 表示随机场F_a的一个实现组态，F_a用于描述诸如图像像素和图像特征等目标对象的观测信息，记作特征随机场，Г＝{F_a，a∈A}表示集合S上的多个特征随机场的集合，称之为多特征随机场。

步骤202、MMFRF模型的理论构建：

给定特征随机场F_a组态f_a前提下组态x的后验概率P(x|f_a)为

$P\left(x\right|f_a)\∝P\left(x\right)\underset{i\∈S}{\mathrm{\Π}}{p}^a\left(y_i^a\right|x_i)$

其中，f_a中随机变量为条件独立的且每个变量都服从统一的依赖于x_i的条件密度函数P(x)为组态x的先验概率。在已知Г前提下，标号随机场X的一个组态的后验概率计算公式如下

$P\left(x\right|f_1,...,f_n)=\frac{P\left(f_1\right|x)...P\left(f_n\right|x)P\left(x\right)}{P(f_1,...,f_n)}\∝\left[\underset{a=1}{\overset{n}{\mathrm{\Π}}}\underset{i=1}{\overset{z}{\mathrm{\Π}}}{p}^a\left(y_i^a\right|x_i)\right]P\left(x\right).$

步骤203、标号随机场X的MAP估计方式：

在给定多个特征随机场情形下估计标号随机场X的真实标号x^*的计算公式为

$x^{*}=\underset{x\∈L^n}{\arg \max }\left\{\right[\underset{a=1}{\overset{n}{\mathrm{\Π}}}\underset{i=1}{\overset{z}{\mathrm{\Π}}}{p}^a\left(y_i^a\right|x_i)\left]P\left(x\right)\right\}$

其中，L＝{1，...，G}表示站点标号的集合。

步骤204、ICM算法对X的MAP估计进行优化求解：

每次迭代过程中，均采用“贪婪”策略进行局部优化搜索，经过次数不多的迭代之后算法达到收敛。

如图3所示，基于MMFRF模型的脑部MR图像分割技术，包括如下步骤：

步骤301、图像预处理，将待分割的脑部MR图像去除脑壳，另外，将子卷概率图谱(SVPAs)配准到待分割图像的原始图像空间；

步骤302、使用SVPAs脑部子卷概率图谱将整个脑部图像域分成多个子域的集合，将每个子域建模成一个局部MMFRF模型，最后通过这些局部模型构建统一的面向整个脑部域的组合MMFRF模型：

1)MMFRF模型：

标号随机场X的真实标号x^*的MAP估计计算如下

$x^{*}=\underset{x\∈L^n}{\arg \max }\left\{\right[\underset{a=1}{\overset{n}{\mathrm{\Π}}}\underset{i=1}{\overset{z}{\mathrm{\Π}}}{p}^a\left(y_i^a\right|x_i)\left]P\left(x\right)\right\}$

其中，L＝{1，...，G}表示图像体素标号的集合，A＝{1，...n}为特征随机场个数，S＝{1，...z}为脑部MR图像体素集合，F_a表示一个特征随机场。表示特征随机场似然概率项，P(x)表示标号随机场先验概率项。标号随机场能量函数为

$U\left(x\right)=\underset{i=1}{\overset{z}{\mathrm{\Σ}}}(O\left(x_i\right|i)+\frac{1}{2}\underset{j\∈N_i}{\mathrm{\Σ}}Q(x_i,x_j|i,j))$

其中，N_i为体素i的邻域，O(x_i|i)和Q(x_i，x_j|i，j)分别表示两个分离的first-order和second-order先验能量函数。O(x_i|i)表示体素i为单一组织类别的概率，Q(x_i，x_j|i，j)表示给定某一体素其邻域体素组织类型的前提下该体素的组织类别概率。对于所有i和j，假设Q(x_i，x_j|i，j)＝Q(x_j，x_i|j，i)。

2)组合MMFRF模型：

将整个脑部图像域通过SVPAs划分成小脑、皮层下结构等12个脑部子结构，即将脑部图像域划分成12个子域R_m，m＝1，...，M的集合。

对于每个子域R_m，m＝1，...，M，定义特征随机场的区域项，P_m(x)为标号随机场P(x)的区域项。分别将所有区域特征随机场项和所有区域标号随机场项P_m(x)(m＝1，...M)合并成组合特征随机场和组合标号随机场。

组合特征随机场 $p^a\left(y_i^a\right|x_i)(a\∈A,A=\{1,...,n\left\}\right)$ 的计算公式如下：

$p^a\left(y_i^a\right|x_i)=\underset{m=1}{\overset{M}{\mathrm{\Σ}}}{r}_i\left(m\right)p_m^a\left(y_i^a\right|x_i)$

组合标号随机场的先验能量函数计算如下

$U\left(x\right)=-\underset{i=1}{\overset{z}{\mathrm{\Σ}}}\{\ln \left(\underset{m=1}{\overset{M}{\mathrm{\Σ}}}{r}_i\left(m\right)\exp \right[-O_m\left(x_i\right|i)\left]\right)-\frac{1}{2}\underset{j\∈N_i}{\mathrm{\Σ}}Q(x_i,x_j)\}$

其中，O_m(x_i|i)和Q_m(x_i，x_j|i，j)分别表示先验能量函数O(x_i|i)和Q(x_i，x_j|i，j)的区域项。对于不同子区域，假设每个区域拥有相同的二阶MRF项，即对于所有区域R_m，m＝1，...，M，Q_m(x_i，x_j|i，j)＝Q(x_i，x_j)。

脑部MR图像分割的组合MMFRF模型组态x的后验概率计算如下：

$P\left(x\right|f_1,...,f_n)\∝\underset{i=1}{\overset{z}{\mathrm{\Π}}}\underset{a=1}{\overset{n}{\mathrm{\Π}}}\left[\underset{m=1}{\overset{M}{\mathrm{\Σ}}}{r}_i\left(m\right)p_m^a\left(y_i^a\right|x_i)\right]\×\left(\underset{m=1}{\overset{M}{\mathrm{\Σ}}}{r}_i\left(m\right)\exp \right[-O_m\left(x_i\right|i)-\frac{1}{2}\underset{j\∈N_j}{\mathrm{\Σ}}Q(x_i,x_j)\left]\right)$

步骤303、计算组合MMFRF模型中的灰度特征随机场，包括以下步骤：

1)计算似然项包括(a)、(b)两种情形：

(a)纯组织类型体素，区域灰度似然函数计算如下：

$p_m^1\left(y_i^1\right|{\mathrm{\θ}}_{\mathrm{lm}},x_i=l)=[\frac{1}{\sqrt{2\mathrm{\π}{\mathrm{\σ}}_{\mathrm{lm}}^2}}\exp \left(\frac{{(y_i^1-{\mathrm{\μ}}_{\mathrm{lm}})}^2}{{2\mathrm{\σ}}_{\mathrm{lm}}^2}\right)$

其中，表示在区域R_m中类别为l时的Gaussian参数集合，μ_lm和分别表示均值和方差。

(b)混合组织类型体素，计算如下：

$p_m^1\left(y_i^1\right|{\mathrm{\θ}}_{l_{1}m},{\mathrm{\θ}}_{l_{2}m},x_i=l)={\∫}_0^1\frac{\exp (-\frac{{[y_i^1-w{\mathrm{\μ}}_{l_{1}m}-(1-w){\mathrm{\μ}}_{l_{2}m}]}^2}{2[w^2{\mathrm{\σ}}_{l_{1}m}^2+{(1-w)}^2{\mathrm{\σ}}_{l_{2}m}^2]})}{\sqrt{2\mathrm{\π}[w^2{\mathrm{\σ}}_{l_{1}m}^2+{(1-w)}^2{\mathrm{\σ}}_{l_{2}m}^2]}}\mathrm{dw}$

其中，和分别表示在区域R_m中类别为l₁和l₂时的包含均值和方差的参数集合，l为图像体素标号类别。

2)参数估计：将图像参数估计问题拟合成有限混合模型(FMM)参数优化问题。采用基于加权似然(WL)方法对图像参数进行估计，计算公式如下

${\mathrm{\Θ}}_m^{*}=\arg \underset{{\mathrm{\Θ}}_m}{\max }\underset{i=1}{\overset{z}{\mathrm{\Σ}}}{r}_i\left(m\right)\ln \left[\underset{l\∈L}{\mathrm{\Σ}}{\mathrm{\α}}_{\mathrm{lm}}{p}_m^1\left(y_i^1\right|{\mathrm{\θ}}_{\mathrm{lm}},l)\right]$

其中，表示图像灰度参数集合，α_lm表示组织类别l在区域R_m中的成分概率，满足μ_lm和分别表示组织类别为l时的均值和方差。

步骤304、计算组合MMFRF模型中的纹理特征随机场，假设对于所有的子域R_m，m＝1，...，M，纹理似然能量函数是均一的，即其计算包括以下步骤：

1)特征提取，包括：

(a)对于图像中的每个体素i∈S，设置其感兴趣域(VOI)来提取基于GLCM的三维纹理特征；

(b)对于每个VOI，选择13个方向和两种距离(d＝1，2)，共26组不同偏移参数进行计算，共得到26个不同的共生矩阵；

(c)为每个矩阵计算角二阶矩、对比度和熵等12个基于GLCM的统计特征值，取26个统计值的平均值作为每个体素的最终特征值。

2)纹理特征随机场计算，包括以下步骤：

(a)为每个体素设定3×3×3的邻域；

(b)分别计算体素i与其邻接体素提取的纹理特征向量间的欧式距离，并将其作为衡量体素i与其邻接体素间相似性的评价尺度；

(c)对于每个体素i，根据其邻域窗口W_i内邻接体素的组织类型标号，利用计算的欧式距离，通过KNN算法对体素i属于哪种组织类型进行投票，计算公式如下：

$p^2\left(y_i^2\right|x_i,y_{W_i}^2)=p^2(y_i^2={\mathrm{tf}}_i|x_i=l,y_{W_i}^2={\mathrm{tf}}_{W_i})$

$=\frac{\mathrm{ofvoxelslabeledaslofknearestneighbor}\sin W_i}{k}$

其中，表示体素i在给定其标号为l以及其邻域窗口内邻接体素纹理特征向量集为的前提下，它的纹理特征向量为tf_i的概率。k表示在邻域窗口W_i内的最邻接体素个数，实验中取值为15。

步骤305、计算组合MMFRF模型中的标号随机场，包括以下步骤：

1)O_m函数计算如下

$O_m\left(x_i\right|i)=-\ln {\mathrm{\α}}_{\mathrm{lm}}^{*},x_i=l$

其中，表示α_lm的加权似然估计；

2)Q函数计算如下

$Q(x_i,x_j)=\frac{{\mathrm{\β\τ}}_{\mathrm{ij}}}{d(i,j)}$

其中，β是一个可调的实验参数，通过经验选取。d(i，j)表示体素i和体素j中心之间的距离。MRF的second-order相关性值τ_ij计算如下

${\mathrm{\τ}}_{\mathrm{ij}}=\left\{\begin{array}{c}-1,\mathrm{if}{x}_i=x_j,\\ 0,\mathrm{if}{x}_i\mathrm{and}{x}_j\mathrm{shareacomponent},.\\ 1,\mathrm{otherwise}.\end{array}\right.$

步骤306、根据计算的多个特征随机场和标号随机场，利用ICM算法通过最大化如下的公式来实现图像分割：

$x^{*}=\underset{x\∈L^n}{\arg \max }\left\{\right[\underset{a=1}{\overset{n}{\mathrm{\Π}}}\underset{l=1}{\overset{z}{\mathrm{\Π}}}{p}^a\left(y_i^a\right|x_i)\left]P\left(x\right)\right\}.$

Claims (2)

1.一种马尔可夫多特征随机场模型构建方法，其特征在于，包括如下步骤：

(1)定义，包括标号随机场、多特征随机场的定义；

(2)模型理论构建，给定多特征随机场集合Γ情形下，通过Bayesian理论，得到模型中标号随机场X的组态的后验概率计算方式；

(3)通过最大后验概率MAP估计准则对模型中标号随机场X的组态的后验概率进行估计；

(4)通过迭代条件模式ICM优化算法对MAP估计进行优化求解，得到模型中标号随机场X的最终的标号值；

得到马尔可夫多特征随机场MMFRF模型；

其中步骤(1)所述的定义，包括以下步骤：

1)标号随机场：

随机场X为集合S上关于邻域***N的马尔可夫随机场MRF，记为标号随机场；

2)多特征随机场：

记A＝{1，...，n}表示随机场个数，S＝{1，2，...，z}定义为z个站点索引位置的集合，随机变量组为定义在集合S上的一个特征随机场，其中每个随机变量取其状态空间中的一个值 表示随机场F_a的一个实现组态，F_a用于描述目标对象的观测信息：图像像素和图像特征，记作特征随机场，Γ＝{F_a，a∈A}表示集合S上的多个特征随机场的集合，称之为多特征随机场；

其中步骤(2)所述的X的组态的后验概率计算方式：

给定特征随机场F_a组态f_a前提下组态x的后验概率P(x|f_a)为

$P\left(x\right|f_a)\∝P\left(x\right)\underset{i\∈S}{\Π}{p}^a\left(y_i^a\right|x_i)$

其中，f_a中随机变量为条件独立的且每个变量都服从统一的依赖于x_i的条件密度函数P(x)为组态x的先验概率；在已知Γ前提下，标号随机场X的一个组态的后验概率计算公式如下

$P\left(x\right|f_1,\mathrm{...},f_n)\frac{P\left(f_1\right|x)\mathrm{...}P\left(f_n\right|x)P\left(x\right)}{P(f_1,\mathrm{...},f_n)}\∝\[\underset{a=1}{\overset{n}{\Π}}\underset{i=1}{\overset{2}{\Π}}{p}^a\left(y_i^a\right|x_i)\]P\left(x\right);$

其中步骤(3)所述的标号随机场X的MAP估计方式为：

在给定多个特征随机场情形下估计标号随机场X的真实标号x^*的计算公式为

$x^{*}=\underset{x\∈L^n}{\arg \max }\{\[\underset{a=1}{\overset{n}{\Π}}\underset{i=1}{\overset{z}{\Π}}{p}^a\left(y_i^a\right|x_i)\]P\left(x\right)\}$

其中，L＝{1，...，G}表示站点标号的集合，G表示不同标号类别的个数。

2.一种基于权利要求1所述的MMFRF模型的脑部MR图像分割方法，其特征在于，包括如下步骤：

(1)使用子卷概率图谱SVPAs将整个脑部图像域分成多个子域的集合，将每个子域建模成一个局部MMFRF模型，最后通过这些局部模型构建统一的面向整个脑部域的组合MMFRF模型；

(2)计算组合MMFRF模型中的灰度特征随机场；

(3)计算组合MMFRF模型中的纹理特征随机场；

(4)计算组合MMFRF模型中的标号随机场；

(5)根据计算的多个特征随机场和标号随机场，使用ICM算法通过MAP估计实现图像分割；

其中步骤(1)所述的构建面向整个脑部域的组合MMFRF，包括以下步骤：

1)MMFRF模型：

标号随机场X的真实标号x^*的MAP估计计算如下

$x^{*}=\underset{x\∈L^n}{\arg \max }\{\[\underset{a=1}{\overset{n}{\Π}}\underset{i=1}{\overset{z}{\Π}}{p}^a\left(y_i^a\right|x_i)\]P\left(x\right)\}$

其中，L＝{1，...，G}表示图像体素标号的集合，A＝{1，...n}为特征随机场个数，S＝{1，...z}为脑部MR图像体素集合，F_a表示一个特征随机场；表示特征随机场似然概率项，P(x)表示标号随机场先验概率项，标号随机场能量函数为

$U\left(x\right)=\underset{i=1}{\overset{z}{\Σ}}\left(O\right(x_i|i)+\frac{1}{2}\underset{j\∈N_i}{\Σ}Q(x_i,x_j|i,j\left)\right)$

其中，N_i为体素i的邻域，O(x_i|i)和Q(x_i，x_j|i，j)分别表示两个分离的first-order和second-order先验能量函数；O(x_i|i)表示体素i为单一组织类别的概率，Q(x_i，x_j|i，j)表示给定某一体素其邻域体素组织类型的前提下该体素的组织类别概率；对于所有i和j，假设Q(x_i，x_j|i，j)＝Q(x_j，x_i|j，i)；

2)组合MMFRF模型：

将整个脑部图像域通过SVPAs划分成12个脑部子结构；对于每个子域R_m，m＝1，...，M，定义特征随机场的区域项，P_m(x)为标号随机场P(x)的区域项；分别将所有区域特征随机场项和所有区域标号随机场项P_m(x)(m＝1，...M)合并成组合特征随机场和组合标号随机场；

组合特征随机场 $p^a\left(y_i^a\right|x_i)(a\∈A,A=\{1,\mathrm{...},n\left\}\right)$ 的计算如下：

$p^a\left(y_i^a\right|x_i)=\underset{m=1}{\overset{M}{\Σ}}{r}_i\left(m\right)p_m^a\left(y_i^a\right|x_i)$

组合标号随机场的先验能量函数计算如下：

$U\left(x\right)=-\underset{i=1}{\overset{z}{\Σ}}\{ln\left(\underset{m=1}{\overset{M}{\Σ}}{r}_i\right(m)\exp \[-O_m(x_i|i)\])-\frac{1}{2}\underset{j\∈N_i}{\Σ}Q(x_i,x_j\left)\right\}$

其中，O_m(x_i|i)和Q_m(x_i，x_j|i，j)分别表示先验能量函数O(x_i|i)和Q(x_i，x_j|i，j)的区域项；对于不同子区域，假设每个区域拥有相同的二阶MRF项，即对于所有区域P_m，m＝1，...，M，r_i(m)表示站点i属于子域R_m的概率，Q_m(x_i，x_j|i，j)＝Q(x_i，x_j)；

脑部MR图像分割的组合MMFRF模型组态x的后验概率计算如下：

$P\left(x\right|f_1,\mathrm{...},f_n)\∝\underset{i=1}{\overset{z}{\Π}}\underset{a=1}{\overset{n}{\Π}}\[\underset{m=1}{\overset{M}{\Σ}}{r}_i\left(m\right)p_m^a\left(y_i^a\right|x_i)\]\×\left(\underset{m=1}{\overset{M}{\Σ}}{r}_i\right(m)\exp \[-O_m\left(x_i\right|i)-\frac{1}{2}\underset{j\∈N_j}{\Σ}Q(x_i,x_j)\]);$

其中步骤(2)所述的计算组合MMFRF模型中的灰度特征随机场，包括以下步骤：

1)计算似然项包括(a)、(b)两种情形：

(a)纯组织类型体素，区域灰度似然函数计算如下

$p_m^1\left(y_i^1\right|{\mathrm{\θ}}_{lm},x_i=l)=\[\frac{1}{\sqrt{2{\mathrm{\π\σ}}_{lm}^2}}\exp (-\frac{{(y_i^1-{\mathrm{\μ}}_{lm})}^2}{2{\mathrm{\σ}}_{lm}^2})\]$

其中，表示在区域R_m中类别为l时的Gaussian参数集合，μ_lm和分别表示均值和方差；

(b)混合组织类型体素，计算如下

$p_m^1\left(y_i^1\right|{\mathrm{\θ}}_{l_{1}m},{\mathrm{\θ}}_{l_{2}m},x_i=l)={\∫}_0^1\frac{\exp (-\frac{{\[y_i^1-{\mathrm{w\μ}}_{l_{1}m}-(1-w){\mathrm{\μ}}_{l_{2}m}\]}^2}{2\[w^2{\mathrm{\σ}}_{l_{1}m}^2+{(1-w)}^2{\mathrm{\σ}}_{l_{2}m}^2\]})}{\sqrt{2\mathrm{\π}\[w^2{\mathrm{\σ}}_{l_{1}m}^2+{(1-w)}^2{\mathrm{\σ}}_{l_{2}m}^2\]}}dw$

其中，和分别表示在区域R_m中类别为l₁和l₂时的包含均值和方差的参数集合，l为图像体素标号类别；

2)参数估计：采用基于加权似然WL方法对图像参数进行估计，计算公式如下

${\mathrm{\Θ}}_m^{*}=\arg \underset{{\mathrm{\Θ}}_m}{\max }\underset{i=1}{\overset{z}{\Σ}}{r}_i\left(m\right)ln\[\underset{l\∈L}{\Σ}{\mathrm{\α}}_{lm}{p}_m^1\left(y_i^1\right|{\mathrm{\θ}}_{lm},l)\]$

其中，表示图像灰度参数集合，α_lm表示组织类别l在区域R_m中的成分概率，满足μ_lm和分别表示组织类别为l时的均值和方差；

其中步骤(3)所述的计算组合MMFRF模型中的纹理特征随机场，假设对于所有的子域R_m，m＝1，...，M，纹理似然能量函数是均一的，即 $p_m^2\left(y_i^a\right|x_i)=p^2\left(y_i^a\right|x_i),$ 其计算步骤如下：

1)特征提取，包括：

(a)对于图像中的每个体素i∈S，设置其感兴趣域VOI来提取基于灰度共生矩阵GLCM的三维纹理特征；

(b)对于每个VOI，选择13个方向和两种距离d＝1，2，共26组不同偏移参数进行计算，共得到26个不同的共生矩阵；

(c)为每个矩阵计算基于GLCM的角二阶矩、对比度、熵、相关性、方差、和平均、和熵、和方差、差分熵、自相关性、相关性信息度量1和相关性信息度量2的12个统计特征值，取26个共生矩阵的统计特征值的平均值生成每个体素的最终的统计特征向量；

2)纹理特征随机场计算，包括以下步骤：

(a)为每个体素设定3×3×3的邻域；

(b)分别计算体素i与其邻接体素提取的纹理特征向量间的欧式距离，并将其作为衡量体素i与其邻接体素间相似性的评价尺度；

(c)对于每个体素i，根据其邻域窗口W_i内邻接体素的组织类型标号，利用计算的欧式距离，通过K最近邻KNN算法对体素i属于哪种组织类型进行投票，计算公式如下：

其中，表示体素i在给定其标号为l以及其邻域窗口内邻接体素纹理特征向量集为的前提下，它的纹理特征向量为tf_i的概率；k表示在邻域窗口W_i内的最邻接体素个数；

其中步骤(4)所述的计算组合MMFRF模型中的标号随机场，包括以下步骤：

1)O_m函数计算如下

$O_m\left(x_i\right|i)=-{\mathrm{ln\α}}_{lm}^{*},x_i=l$

其中，表示α_lm的加权似然估计；

2)Q函数计算如下

$Q(x_i,x_j)=\frac{{\mathrm{\β\τ}}_{ij}}{d(i,j)}$

其中，β是一个可调的实验参数，通过经验选取；d(i，j)表示体素i和体素j中心之间的距离；MRF的second-order相关性值τ_ij计算如下