CN103544702A - 一种基于先验形状的核磁共振图像分割方法 - Google Patents

Abstract

本发明涉及一种基于先验形状的核磁共振图像分割方法，利用形状先验信息，融入自适应区域信息与边界信息，能够更好地分割复杂背景和图像中的目标物体，有效地解决图像中感兴趣的目标物体像素灰度不均、目标物体边缘梯度不明显、感兴趣的目标物体被遮挡等问题的同时，大大减低计算的复杂度。

Description

一种基于先验形状的核磁共振图像分割方法

技术领域

本发明涉及一种基于先验形状的核磁共振图像分割方法。

背景技术

乳腺癌是中老年女性最常见的恶性肿瘤之一，在中国已经高居女性恶性肿瘤的首位，并有逐年上升以及发病年轻化的趋势。目前对乳腺癌的预防尚无良策，对乳腺影像的检测分析是成功检测、诊断、治疗肿瘤的重要依据，尤其是对乳腺区域的分割对后续肿瘤的研究分析具有重要的意义。然而，由于乳腺组织和身体部分软组织密度相似，要单独分割出乳腺区域部分比较困难。目前，活动轮廓模型和水平集方法已经广泛使用在图像分割中，活动轮廓模型主要分为两大类：基于区域和基于边界的活动轮廓模型。基于边界的模型如GAC模型，根据边界梯度信息驱动活动边界向目标边界移动，虽然充分利用了边界附近的局部信息，但对于噪声和初始轮廓很敏感且过度依赖边界检测方程，容易陷入局部最小值。基于区域的模型，经典的如C-V模型，根据不同区域的灰度统计信息用一个区域描述子来引导边界移动从而分割不同的区域。然而，当目标的边缘和背景灰度差异较小时，演化曲线将会出现泄漏或过分割，导致分割失败。且当图像存在噪声、杂波或遮挡时，普通的分割模型不能提取感兴趣的目标。

发明内容

针对上述技术问题，本发明所要解决的技术问题是提供一种针对核磁共振图像，利用形状先验信息，融入自适应区域信息与边界信息实现图像分割的基于先验形状的核磁共振图像分割方法。

本发明为了解决上述技术问题采用以下技术方案：本发明设计了一种基于先验形状的核磁共振图像分割方法，包括如下步骤：

步骤01.针对待分割目标体核磁共振图像内不同灰度的区域，获得该模型的水平集能量泛函E_RSF；

步骤02.根据

获得先验形状能量项E_shape，其中，H为Heavside函数，表示待分割图像中目标体区域闭合边界的水平集函数，φ表示待分割图像在分割过程中，针对目标体区域的分割边界，Ψ为最优先验形状；

步骤03.构建新的混杂活动轮廓模型E_total＝E_RSF+Q_sE_shape，并搜索该模型曲线运动方程最优值实现针对待分割目标体核磁共振图像的分割，获得待分割的目标体区域轮廓，其中，Q_s为先验形状能量项E_shape对应的权重系数；

其中，最优先验形状Ψ通过如下步骤获得：

步骤001.获得一组大小相同的目标体核磁共振训练图像，分别针对各幅训练图像中的目标体轮廓进行标注，并分别针对各幅训练图像进行二值化处理；

步骤002.针对N幅训练图像，选择一幅训练图像作为基准，分别与其它N-1幅训练图像进行配准，加上作为基准的训练图像，共获得N幅配准训练图像，N为该组大小相同的目标体核磁共振训练图像的数量；

步骤003.根据配准训练图像中目标体的轮廓，将N幅配准训练图像分别表示成符号距离函数

获得平均水平集函数

并根据各个符号距离函数

分别获得各幅配准训练图像中目标体轮廓的偏移量函数

其中，i=1、···、N；

步骤004.根据各幅配准训练图像中目标体轮廓的偏移量函数

构成N维形状变化矩阵M，并对相应的协方差矩阵进行奇异值分解，获得N组针对目标体形状一一对应的特征值和特征向量，选取降维后矩阵的前h个维度所对应的特征向量构成矩阵W_h=（w₁、···、w_h），实现对形状变化矩阵M进行降维，并分别获得与矩阵W_h中各个特征向量一一对应的权系数，构成向量X_pca=（x₁、···、x_h），其中M^T为形状变化矩阵M的转置矩阵,1≤h<N；

步骤005.根据Φ=W_h·X_pca获得目标体的特征形状集Φ，根据

$P\left(\mathrm{\&Phi;}\left(X_{\mathrm{pca}}\right)\right)=\frac{1}{{\left(2\mathrm{\&pi;}\right)}^{\frac{h}{2}}{\left|{\mathrm{\&Lambda;}}_h\right|}^{\frac{1}{2}}}\exp (-\frac{1}{2}{X}_{\mathrm{pca}}^T{\mathrm{\&Lambda;}}_h^{-1}{X}_{\mathrm{pca}})$ 获得Φ(X_pca)的概率P(Φ(X_pca))，构成E＝-logP(Φ(X_pca))方程，求解此方程，获得最优X_pca，进而获得最优先验形状Ψ；

其中，Λ_h表示一个对角矩阵，这个对角矩阵的对角元素包含前h个特征向量。

作为本发明的一种优选技术方案：所述步骤01中，针对待分割目标体核磁共振图像，采用自适应区域像素灰度统计活动轮廓模型表示该待分割图像内不同灰度的区域，并获得该模型的水平集能量泛函E_RSF： $E_{\mathrm{RSF}}(C,f_1\left(x\right),f_2\left(x\right))=\&Integral;{\mathrm{\&epsiv;}}_x^{\mathrm{Fit}}(C,f_1\left(x\right),f_2\left(x\right))\mathrm{dx}+v\left|C\right|,$ 其中， ${\mathrm{\&epsiv;}}_x^{\mathrm{Fit}}(C,f_1\left(x\right),f_2\left(x\right))=\underset{g=1}{\overset{2}{\mathrm{\&Sigma;}}}{\mathrm{\&lambda;}}_g{\&Integral;}_{{\mathrm{\&Omega;}}_g}{K}_{\mathrm{\&sigma;}}(x-y){|I\left(y\right)-f_g\left(x\right)|}^2\mathrm{dy},$ 其中，K_σ为一个高斯核函数，C为待分割目标体核磁共振图像内不同灰度分布得到的前景区域Ω₁和背景区域Ω₂的闭合边界，|C|表示闭合边界C的长度，f₁(x)、f₂(x)分别为待分割目标体核磁共振图像中前景区域Ω₁和背景区域Ω₂的灰度值的拟合函数，λ₁、λ₂分别表示前景区域Ω₁、和背景区域Ω₂的权重，I(y)表示在高斯核函数K_σ中一点像素的灰度值，v|C|表示控制闭合边界C的长度，v表示闭合边界C的长度系数。

作为本发明的一种优选技术方案：所述步骤002中，采用欧式几何相似性变换针对N幅训练图像进行配准。

作为本发明的一种优选技术方案：所述步骤002中，针对N幅训练图像，选择一幅训练图像作为基准，分别与其它N-1幅训练图像采用欧式几何相似性变换按照如下配准函数进行配准：

$E_{\mathrm{align}}={\mathrm{\&Sigma;}}_{i=1}^N{\mathrm{\&Sigma;}}_{\underset{j\&NotEqual;i}{j=1}}^N\frac{{\&Integral;\&Integral;}_{\mathrm{\&Omega;}}{({\tilde{I}}^i-{\tilde{I}}^j)}^2\mathrm{dxdy}}{{\&Integral;\&Integral;}_{\mathrm{\&Omega;}}{({\tilde{I}}^i+{\tilde{I}}^j)}^2\mathrm{dxdy}}$

其中，Iⁱ表示第i幅训练图像，I^j表示第j幅训练图像，

和

分别表示第i幅训练图像Iⁱ、第j幅训练图像I^j经过配准后对应的图像。

作为本发明的一种优选技术方案：所述步骤03中，针对所述新的混杂活动轮廓模型E_total＝E_RSF+Q_sE_shape进行最小化，采用变分法及梯度下降流方程，实现针对待分割目标体核磁共振图像的分割。

本发明所述一种基于先验形状的核磁共振图像分割方法采用以上技术方案与现有技术相比，具有以下技术效果：本发明设计的基于先验形状的核磁共振图像分割方法，利用形状先验信息，融入自适应区域信息与边界信息，能够更好地分割复杂背景和图像中的目标物体，有效地解决图像中感兴趣的目标物体像素灰度不均、目标物体边缘梯度不明显、感兴趣的目标物体被遮挡等问题的同时，大大减低计算的复杂度。

附图说明

图1是本发明设计基于先验形状的核磁共振图像分割方法的流程示意图；

图2是本发明应用于乳腺区域核磁共振图像分割过程中选取先验样本对应的二值图；

图3是本发明应用于乳腺区域核磁共振图像分割过程中乳腺区域对应的初始轮廓线；

图4是乳腺区域核磁共振图像运用GAC模型进行分割的效果图；

图5是乳腺区域核磁共振图像运用C-V模型进行分割的效果图；

图6是乳腺区域核磁共振图像运用RSF模型进行分割的效果图；

图7是乳腺区域核磁共振图像运用本发明设计基于先验形状的核磁共振图像分割方法进行分割的效果图；

图8是不同模型针对核磁共振图像分割效果的戴斯系数值比较示意图；

图9是不同模型针对核磁共振图像分割效果的峰值信噪比值比较示意图。

具体实施方式

下面结合说明书附图对本发明的具体实施方式作进一步详细的说明。

如图1所示，本发明设计了一种基于先验形状的核磁共振图像分割方法，包括如下步骤：

步骤001.获得N幅大小相同的目标体核磁共振训练图像，分别针对各幅训练图像中的目标体轮廓进行标注，并分别针对各幅训练图像中目标体区域和背景区域分别采用两种不同的数值进行二值化处理；

步骤002.针对N幅训练图像，选择一幅训练图像作为基准，分别与其它N-1幅训练图像采用欧式几何相似性变换进行配准，加上作为基准的训练图像，共获得N幅配准训练图像，其中，按照如下配准函数进行配准：

$E_{\mathrm{align}}={\mathrm{\&Sigma;}}_{i=1}^N{\mathrm{\&Sigma;}}_{\underset{j\&NotEqual;i}{j=1}}^N\frac{{\&Integral;\&Integral;}_{\mathrm{\&Omega;}}{({\tilde{I}}^i-{\tilde{I}}^j)}^2\mathrm{dxdy}}{{\&Integral;\&Integral;}_{\mathrm{\&Omega;}}{({\tilde{I}}^i+{\tilde{I}}^j)}^2\mathrm{dxdy}}$

其中，Iⁱ表示第i幅训练图像，I^j表示第j幅训练图像，

和

分别表示第i幅训练图像Iⁱ、第j幅训练图像I^j经过配准（平移、旋转、尺度拉伸转换）后对应的图像；

步骤003.根据配准训练图像中目标体的轮廓，将N幅配准训练图像分别表示成符号距离函数

根据

获得平均水平集函数

并根据各个符号距离函数

分别获得各幅配准训练图像中目标体轮廓的偏移量函数

其中，i=1、···、N；

步骤004.根据各幅配准训练图像中目标体轮廓的偏移量函数构成N维形状变化矩阵M，并对相应的协方差矩阵C＝(1/N)M·M^T进行奇异值分解，获得N组针对目标体形状一一对应的特征值和特征向量，选取降维后矩阵的前h个维度所对应的特征向量构成矩阵W_h=（w₁、···、w_h），实现对形状变化矩阵M进行降维，并分别获得与矩阵W_h中各个特征向量一一对应的权系数，构成向量X_pca=（x₁、···、x_h），其中M^T为形状变化矩阵M的转置矩阵,1≤h<N；

步骤005.根据Φ=W_h·X_pca获得目标体的特征形状集Φ，根据

$P\left(\mathrm{\&Phi;}\left(X_{\mathrm{pca}}\right)\right)=\frac{1}{{\left(2\mathrm{\&pi;}\right)}^{\frac{h}{2}}{\left|{\mathrm{\&Lambda;}}_h\right|}^{\frac{1}{2}}}\exp (-\frac{1}{2}{X}_{\mathrm{pca}}^T{\mathrm{\&Lambda;}}_h^{-1}{X}_{\mathrm{pca}})$ 获得Φ(X_pca)的概率P(Φ(X_pca))，构成E＝-logP(Φ(X_pca))方程，求解此方程，获得最优X_pca，进而获得最优先验形状Ψ；

其中，Λ_h表示一个对角矩阵，这个对角矩阵的对角元素包含前h个特征向量；

步骤01.针对待分割目标体核磁共振图像，采用自适应区域像素灰度统计活动轮廓模型表示该待分割图像内不同灰度的区域，并获得该模型的水平集能量泛函：

$E_{\mathrm{RSF}}(C,f_1\left(x\right),f_2\left(x\right))=\&Integral;{\mathrm{\&epsiv;}}_x^{\mathrm{Fit}}(C,f_1\left(x\right),f_2\left(x\right))\mathrm{dx}+v\left|C\right|,$ 其中，

${\mathrm{\&epsiv;}}_x^{\mathrm{Fit}}(C,f_1\left(x\right),f_2\left(x\right))=\underset{g=1}{\overset{2}{\mathrm{\&Sigma;}}}{\mathrm{\&lambda;}}_g{\&Integral;}_{{\mathrm{\&Omega;}}_g}{K}_{\mathrm{\&sigma;}}(x-y){|I\left(y\right)-f_g\left(x\right)|}^2\mathrm{dy},$ 其中，K_σ为一个高斯核函数，C为待分割目标体核磁共振图像内不同灰度分布得到的前景区域Ω₁和背景区域Ω₂的闭合边界，|C|表示闭合边界C的长度，f₁(x)、f₂(x)分别为待分割目标体核磁共振图像中前景区域Ω₁和背景区域Ω₂的灰度值的拟合函数，λ₁、λ₂分别表示前景区域Ω₁、和背景区域Ω₂的权重，I(y)表示在高斯核函数K_σ中一点像素的灰度值，v|C|表示控制闭合边界C的长度，v表示闭合边界C的长度系数。

步骤02.根据

获得先验形状能量项E_shape，其中，H为Heavside函数，表示待分割图像中目标体区域闭合边界的水平集函数，φ表示待分割图像在分割过程中，针对目标体区域的分割边界；

步骤03.构建新的混杂活动轮廓模型E_total＝E_RSF+Q_sE_shape，并采用变分法及梯度下降流方程，搜索该模型曲线运动方程最优值实现针对待分割目标体核磁共振图像的分割，获得待分割的目标体区域轮廓，其中，Q_s为先验形状能量项E_shape对应的权重系数。

其中，新的混杂活动轮廓模型具体E_total＝E_RSF+Q_sE_shape具体表示如下：

$E_{\mathrm{total}}(\mathrm{\&phi;},\mathrm{\&psi;},f_1\left(x\right),f_2\left(x\right))=\underset{i=1}{\overset{2}{\mathrm{\&Sigma;}}}{\mathrm{\&lambda;}}_i\&Integral;(\&Integral;K_{\mathrm{\&sigma;}}(x-y){|I\left(y\right)-f_i\left(x\right)|}^2{M}_i\left(\mathrm{\&phi;}\left(y\right)\right)\mathrm{dy})\mathrm{dx}+v\&Integral;|\&dtri;H\left(\mathrm{\&phi;}\left(x\right)\right)|\mathrm{dx}+Q_s{\&Integral;\&Integral;}_{\mathrm{\&Omega;}}\{H\left(\mathrm{\&phi;}\right)-H\left(\mathrm{\&psi;}\right){\}}^2\mathrm{dxdy}$

其中，采用变分法及梯度下降流方程构建演化模型对应的数值计算格式如下：

$\frac{\&PartialD;\mathrm{\&phi;}}{\&PartialD;t}=-\mathrm{\&delta;}\left(\mathrm{\&phi;}\right)({\mathrm{\&lambda;}}_1{e}_1-{\mathrm{\&lambda;}}_2{e}_2)+\mathrm{v\&delta;}\left(\mathrm{\&phi;}\right)\mathrm{div}\frac{\&dtri;\mathrm{\&phi;}}{|\&dtri;\mathrm{\&phi;}|}+2\mathrm{\&delta;}\left(\mathrm{\&phi;}\right)Q_s(H\left(\mathrm{\&phi;}\right)-H\left(\mathrm{\&psi;}\right))$

当上式收敛且演化曲线和形状模板之间满足 ${\&Integral;\&Integral;}_{\mathrm{\&Omega;}}\{H\left(\mathrm{\&phi;}\right)-H\left(\mathrm{\&psi;}\right){\}}^2\mathrm{dxdy}\&le;{\&Integral;\&Integral;}_{\mathrm{\&Omega;}}\mathrm{\&alpha;\&delta;}\left(\mathrm{\&psi;}\right)\mathrm{dxdy}$ 时，曲线演化停止，分割结束。

本发明所述一种基于先验形状的核磁共振图像分割方法采用以上技术方案与现有技术相比，具有以下技术效果：本发明设计的基于先验形状的核磁共振图像分割方法，利用形状先验信息，融入自适应区域信息与边界信息，能够更好地分割复杂背景和图像中的目标物体，有效地解决图像中感兴趣的目标物体像素灰度不均、目标物体边缘梯度不明显、感兴趣的目标物体被遮挡等问题的同时，大大减低计算的复杂度。

如图1所示，本发明设计的基于先验形状的核磁共振图像分割方法在实际应用过程中，针对乳腺核磁共振图像的分割按照如下步骤执行：

步骤001.获得16个病人相同位置拍摄的大小为512×512的乳腺核磁共振图像，分别针对各幅乳腺核磁共振图像，人工标注乳腺区域轮廓，并分别针对各幅乳腺核磁共振图像中乳腺区域和背景区域分别采用两种不同的数值进行二值化处理，其中，乳腺区域灰度值定义为1，背景部分灰度值定义为0，二值化后的16幅图像作为训练图像，如图2所示；

步骤002.针对16幅训练图像，选择一幅训练图像作为基准，分别与其它15幅训练图像采用欧式几何相似性变换进行配准，加上作为基准的训练图像，共获得16幅配准训练图像，配准过程按如下配准函数进行配准：

$E_{\mathrm{align}}={\mathrm{\&Sigma;}}_{i=1}^{16}{\mathrm{\&Sigma;}}_{\underset{j\&NotEqual;i}{j=1}}^{16}\frac{{\&Integral;\&Integral;}_{\mathrm{\&Omega;}}{({\tilde{I}}^i-{\tilde{I}}^j)}^2\mathrm{dxdy}}{{\&Integral;\&Integral;}_{\mathrm{\&Omega;}}{({\tilde{I}}^i+{\tilde{I}}^j)}^2\mathrm{dxdy}}$

其中，Iⁱ表示第i幅训练图像，I^j表示第j幅训练图像，

和

分别表示第i幅训练图像Iⁱ、第j幅训练图像I^j经过配准（平移、旋转、尺度拉伸转换）后对应的图像；

步骤003.根据配准训练图像中乳腺的轮廓，将16幅配准训练图像分别表示成符号距离函数根据

获得平均水平集函数并根据各个符号距离函数

分别获得各幅配准训练图像中乳腺轮廓的偏移量函数其中，i=1、···、16；

步骤004.根据各幅配准训练图像中乳腺轮廓的偏移量函数

构成16维形状变化矩阵M，并对相应的协方差矩阵C＝(1/16)M·M^T进行奇异值分解，获得16组针对乳腺形状一一对应的特征值和特征向量，选取降维后矩阵的前h个维度所对应的特征向量构成矩阵W_h=（w₁、···、w_h），对形状变化矩阵M进行降维，并分别获得与矩阵W_h中各个特征向量一一对应的权系数，构成向量X_pca=（x₁、···、x_h），其中M^T为形状变化矩阵M的转置矩阵,1≤h<N；

步骤005.根据Φ=W_h·X_pca获得乳腺的特征形状集Φ，根据

$P\left(\mathrm{\&Phi;}\left(X_{\mathrm{pca}}\right)\right)=\frac{1}{{\left(2\mathrm{\&pi;}\right)}^{\frac{h}{2}}{\left|{\mathrm{\&Lambda;}}_h\right|}^{\frac{1}{2}}}\exp (-\frac{1}{2}{X}_{\mathrm{pca}}^T{\mathrm{\&Lambda;}}_h^{-1}{X}_{\mathrm{pca}})$ 获得Φ(X_pca)的概率P(Φ(X_pca))，构成E＝-logP(Φ(X_pca))方程，求解此方程，获得最优X_pca，进而获得最优先验形状Ψ；

其中，Λ_h表示一个对角矩阵，这个对角矩阵的对角元素包含前h个特征向量；

步骤01.如图3所示，针对待分割乳腺核磁共振图像，采用自适应区域像素灰度统计活动轮廓模型表示该待分割乳腺图像内不同灰度的区域，并获得该模型的水平集能量泛函： $E_{\mathrm{RSF}}(C,f_1\left(x\right),f_2\left(x\right))=\&Integral;{\mathrm{\&epsiv;}}_x^{\mathrm{Fit}}(C,f_1\left(x\right),f_2\left(x\right))\mathrm{dx}+v\left|C\right|,$ 其中：

${\mathrm{\&epsiv;}}_x^{\mathrm{Fit}}(C,f_1\left(x\right),f_2\left(x\right))=\underset{g=1}{\overset{2}{\mathrm{\&Sigma;}}}{\mathrm{\&lambda;}}_g{\&Integral;}_{{\mathrm{\&Omega;}}_g}{K}_{\mathrm{\&sigma;}}(x-y){|I\left(y\right)-f_g\left(x\right)|}^2\mathrm{dy},$ 其中，K_σ为一个高斯核函数，C为待分割目标体核磁共振图像内不同灰度分布得到的前景区域Ω₁和背景区域Ω₂的闭合边界，|C|表示闭合边界C的长度，f₁(x)、f₂(x)分别为待分割目标体核磁共振图像中前景区域Ω₁和背景区域Ω₂的灰度值的拟合函数，λ₁、λ₂分别表示前景区域Ω₁、和背景区域Ω₂的权重，I(y)表示在高斯核函数K_σ中一点像素的灰度值，v|C|表示控制闭合边界C的长度，v表示闭合边界C的长度系数。

步骤02.根据

获得先验形状能量项E_shape，其中，H为Heavside函数，表示待分割图像中乳腺区域闭合边界的水平集函数，φ表示待分割图像在分割过程中，针对乳腺区域的分割边界；

步骤03.构建新的混杂活动轮廓模型E_total＝E_RSF+Q_sE_shape，并搜索该模型曲线运动方程最优值实现针对待分割目标乳腺体核磁共振图像的分割，其中，Q_s为先验形状能量项E_shape对应的权重系数；

其中，新的混杂活动轮廓模型具体E_total＝E_RSF+Q_sE_shape具体表示如下：

$E_{\mathrm{total}}(\mathrm{\&phi;},\mathrm{\&psi;},f_1\left(x\right),f_2\left(x\right))=\underset{i=1}{\overset{2}{\mathrm{\&Sigma;}}}{\mathrm{\&lambda;}}_i\&Integral;(\&Integral;K_{\mathrm{\&sigma;}}(x-y){|I\left(y\right)-f_i\left(x\right)|}^2{M}_i\left(\mathrm{\&phi;}\left(y\right)\right)\mathrm{dy})\mathrm{dx}+v\&Integral;|\&dtri;H\left(\mathrm{\&phi;}\left(x\right)\right)|\mathrm{dx}+Q_s{\&Integral;\&Integral;}_{\mathrm{\&Omega;}}\{H\left(\mathrm{\&phi;}\right)-H\left(\mathrm{\&psi;}\right){\}}^2\mathrm{dxdy}$

采用变分法及梯度下降流方程构建演化模型对应的数值计算格式如下：

$\frac{\&PartialD;\mathrm{\&phi;}}{\&PartialD;t}=-\mathrm{\&delta;}\left(\mathrm{\&phi;}\right)({\mathrm{\&lambda;}}_1{e}_1-{\mathrm{\&lambda;}}_2{e}_2)+\mathrm{v\&delta;}\left(\mathrm{\&phi;}\right)\mathrm{div}\frac{\&dtri;\mathrm{\&phi;}}{|\&dtri;\mathrm{\&phi;}|}+2\mathrm{\&delta;}\left(\mathrm{\&phi;}\right)Q_s(H\left(\mathrm{\&phi;}\right)-H\left(\mathrm{\&psi;}\right))$

具体操作包括如下：

根据 $\frac{\&PartialD;\mathrm{\&phi;}}{\&PartialD;t}=-\mathrm{\&delta;}\left(\mathrm{\&phi;}\right)({\mathrm{\&lambda;}}_1{e}_1-{\mathrm{\&lambda;}}_2{e}_2)+\mathrm{v\&delta;}\left(\mathrm{\&phi;}\right)\mathrm{div}\frac{\&dtri;\mathrm{\&phi;}}{|\&dtri;\mathrm{\&phi;}|}+2\mathrm{\&delta;}\left(\mathrm{\&phi;}\right)Q_s(H\left(\mathrm{\&phi;}\right)-H\left(\mathrm{\&psi;}\right))$ 对待分割乳腺核磁共振图像进行水平集迭代计算，给定水平集的迭代次数为200，先验形状能量项E_shape对应的权重系数Q_s设为1/3，RSF模型（自适应区域像素灰度统计活动轮廓模型）中，针对待分割乳腺核磁共振图像中前景区域Ω₁和背景区域Ω₂的核函数的参数σ设为3.0，且前景区域Ω₁和背景区域Ω₂的权系数分别为λ₁=1、λ₂=2，Ω₁和Ω₂的闭合边界的长度系数v=0.001*255*255，曲线演化过程中迭代步长Δt设为0.4，当 $\frac{\&PartialD;\mathrm{\&phi;}}{\&PartialD;t}=-\mathrm{\&delta;}\left(\mathrm{\&phi;}\right)({\mathrm{\&lambda;}}_1{e}_1-{\mathrm{\&lambda;}}_2{e}_2)+\mathrm{v\&delta;}\left(\mathrm{\&phi;}\right)\mathrm{div}\frac{\&dtri;\mathrm{\&phi;}}{|\&dtri;\mathrm{\&phi;}|}+2\mathrm{\&delta;}\left(\mathrm{\&phi;}\right)Q_s(H\left(\mathrm{\&phi;}\right)-H\left(\mathrm{\&psi;}\right))$ 收敛且演化曲线和形状模板之间满足 ${\&Integral;\&Integral;}_{\mathrm{\&Omega;}}\{H\left(\mathrm{\&phi;}\right)-H\left(\mathrm{\&psi;}\right){\}}^2\mathrm{dxdy}\&le;{\&Integral;\&Integral;}_{\mathrm{\&Omega;}}\mathrm{\&alpha;\&delta;}\left(\mathrm{\&psi;}\right)\mathrm{dxdy}$ 时，曲线演化停止，分割结束，获得待分割的乳腺区域轮廓，其中，α是恒定值，本实施例中设为0.2，如果在演化过程中不满足条件，说明形状项的影响权重不够，需要增大Q_s或者减小λ₁和λ₂。

本发明设计的基于先验形状的核磁共振图像分割方法在针对乳腺核磁共振图像分割过程中，如图8和图9所示，其中，戴斯系数值DC＝2(C₁∩C₂)/(C₁+C₂)是一个从0到1的数，DC越大表明图像分割结果越接近参照图像；峰值信噪比值PSNR＝10log((2ⁿ-1)²/MSE)，峰值信噪比值PSNR越大说明图像分割结果越理想，由如图8和图9可以看出，本发明设计方法的分割效果更好；对比图4至图7所示，可以看出本发明设计方法相对于GAC模型、C-V模型和RSF模型，分割效果更好，由图中可以看出，由于乳腺和身体连接部分组织灰度非常相似，GAC模型、C-V模型和RSF模型的分割结果都不是很理想，而由图7所示本发明设计方法的分割效果图，可以看出，该模型在低对比度、边界模糊、灰度不均的情况下仍然可以较好地分割出乳腺区域轮廓。

上面结合附图对本发明的实施方式作了详细说明，但是本发明并不限于上述实施方式，在本领域普通技术人员所具备的知识范围内，还可以在不脱离本发明宗旨的前提下做出各种变化。

Claims (5)

1.一种基于先验形状的核磁共振图像分割方法，其特征在于，包括如下步骤：

步骤01.针对待分割目标体核磁共振图像内不同灰度的区域，获得该模型的水平集能量泛函E_RSF；

步骤02.根据

获得先验形状能量项E_shape，其中，H为Heavside函数，表示待分割图像中目标体区域闭合边界的水平集函数，φ表示待分割图像在分割过程中，针对目标体区域的分割边界，Ψ为最优先验形状；

步骤03.构建新的混杂活动轮廓模型E_total＝E_RSF+Q_sE_shape，并搜索该模型曲线运动方程最优值实现针对待分割目标体核磁共振图像的分割，获得待分割的目标体区域轮廓，其中，Q_s为先验形状能量项E_shape对应的权重系数；

其中，最优先验形状Ψ通过如下步骤获得：

步骤001.获得一组大小相同的目标体核磁共振训练图像，分别针对各幅训练图像中的目标体轮廓进行标注，并分别针对各幅训练图像进行二值化处理；

步骤002.针对N幅训练图像，选择一幅训练图像作为基准，分别与其它N-1幅训练图像进行配准，加上作为基准的训练图像，共获得N幅配准训练图像，N为该组大小相同的目标体核磁共振训练图像的数量；

步骤003.根据配准训练图像中目标体的轮廓，将N幅配准训练图像分别表示成符号距离函数

获得平均水平集函数

并根据各个符号距离函数分别获得各幅配准训练图像中目标体轮廓的偏移量函数其中，i=1、···、N；

步骤004.根据各幅配准训练图像中目标体轮廓的偏移量函数

构成N维形状变化矩阵M，并对相应的协方差矩阵进行奇异值分解，获得N组针对目标体形状一一对应的特征值和特征向量，选取降维后矩阵的前h个维度所对应的特征向量构成矩阵W_h=（w₁、···、w_h），实现对形状变化矩阵M进行降维，并分别获得与矩阵W_h中各个特征向量一一对应的权系数，构成向量X_pca=（x₁、···、x_h），其中M^T为形状变化矩阵M的转置矩阵,1≤h<N；

步骤005.根据Φ=W_h·X_pca获得目标体的特征形状集Φ，根据 $P\left(\mathrm{\&Phi;}\left(X_{\mathrm{pca}}\right)\right)=\frac{1}{{\left(2\mathrm{\&pi;}\right)}^{\frac{h}{2}}{\left|{\mathrm{\&Lambda;}}_h\right|}^{\frac{1}{2}}}\exp (-\frac{1}{2}{X}_{\mathrm{pca}}^T{\mathrm{\&Lambda;}}_h^{-1}{X}_{\mathrm{pca}})$ 获得Φ(X_pca)的概率P(Φ(X_pca))，构成E＝-logP(Φ(X_pca))方程，求解此方程，获得最优X_pca，进而获得最优先验形状Ψ；其中，Λ_h表示一个对角矩阵，这个对角矩阵的对角元素包含前h个特征向量。

2.根据权利要求1所述一种基于先验形状的核磁共振图像分割方法，其特征在于：所述步骤01中，针对待分割目标体核磁共振图像，采用自适应区域像素灰度统计活动轮廓模型表示该待分割图像内不同灰度的区域，并获得该模型的水平集能量泛函E_RSF：

$E_{\mathrm{RSF}}(C,f_1\left(x\right),f_2\left(x\right))=\&Integral;{\mathrm{\&epsiv;}}_x^{\mathrm{Fit}}(C,f_1\left(x\right),f_2\left(x\right))\mathrm{dx}+v\left|C\right|,$ 其中，

${\mathrm{\&epsiv;}}_x^{\mathrm{Fit}}(C,f_1\left(x\right),f_2\left(x\right))=\underset{g=1}{\overset{2}{\mathrm{\&Sigma;}}}{\mathrm{\&lambda;}}_g{\&Integral;}_{{\mathrm{\&Omega;}}_g}{K}_{\mathrm{\&sigma;}}(x-y){|I\left(y\right)-f_g\left(x\right)|}^2\mathrm{dy},$ 其中，K_σ为一个高斯核函数，C为待分割目标体核磁共振图像内不同灰度分布得到的前景区域Ω₁和背景区域Ω₂的闭合边界，|C|表示闭合边界C的长度，f₁(x)、f₂(x)分别为待分割目标体核磁共振图像中前景区域Ω₁和背景区域Ω₂的灰度值的拟合函数，λ₁、λ₂分别表示前景区域Ω₁、和背景区域Ω₂的权重，I(y)表示在高斯核函数K_σ中一点像素的灰度值，v|C|表示控制闭合边界C的长度，v表示闭合边界C的长度系数。

3.根据权利要求1所述一种基于先验形状的核磁共振图像分割方法，其特征在于：所述步骤002中，采用欧式几何相似性变换针对N幅训练图像进行配准。

4.根据权利要求3所述一种基于先验形状的核磁共振图像分割方法，其特征在于：所述步骤002中，针对N幅训练图像，选择一幅训练图像作为基准，分别与其它N-1幅训练图像采用欧式几何相似性变换按照如下配准函数进行配准：

$E_{\mathrm{align}}={\mathrm{\&Sigma;}}_{i=1}^N{\mathrm{\&Sigma;}}_{\underset{j\&NotEqual;i}{j=1}}^N\frac{{\&Integral;\&Integral;}_{\mathrm{\&Omega;}}{({\tilde{I}}^i-{\tilde{I}}^j)}^2\mathrm{dxdy}}{{\&Integral;\&Integral;}_{\mathrm{\&Omega;}}{({\tilde{I}}^i+{\tilde{I}}^j)}^2\mathrm{dxdy}}$

其中，Iⁱ表示第i幅训练图像，I^j表示第j幅训练图像，

和

分别表示第i幅训练图像Iⁱ、第j幅训练图像I^j经过配准后对应的图像。

5.根据权利要求1所述一种基于先验形状的核磁共振图像分割方法，其特征在于：所述步骤03中，针对所述新的混杂活动轮廓模型E_total＝E_RSF+Q_sE_shape进行最小化，采用变分法及梯度下降流方程，实现针对待分割目标体核磁共振图像的分割。

Images