CN102999917A - 基于t2-mri和dw-mri的宫颈癌图像自动分割方法 - Google Patents

Abstract

一种基于T2加权的核磁共振成像T2-MRI和弥散加权的核磁共振成像DW-MRI的宫颈癌图像自动分割方法，包括利用非线性配准方法将DW-MR图像配准到T2-MR图像，并对配准后的DW-MR图像进行分类；采用非线性各向异性扩散滤波技术对T2-MR图像进行滤波，分割出膀胱和直肠，并利用膀胱和直肠的分割结果将感兴趣区分割出来；对T2-MR图像的感兴趣区和DW-MR图像采用联合最大后验概率CMAP的方法进行肿瘤的精确分割。本发明充分利用了T2-MR图像和DW-MR图像的有效信息，可以有效地克服T2-MR图像中的噪声、局部容积效应和强度重叠的影响，是一种精确有效的宫颈癌分割方法，对预防、诊断和治疗宫颈癌有着重要的临床意义和应用价值。

Description

基于T2-MRI和DW-MRI的宫颈癌图像自动分割方法

技术领域

本发明属于图像处理领域，具体涉及一种基于T2加权的核磁共振成像T2-MRI和弥散加权的核磁共振成像DW-MRI的宫颈癌图像自动分割方法。

背景技术

宫颈癌是严重威胁女性生命健康的常见恶性肿瘤之一。宫颈癌的准确分割，对预防、诊断和治疗宫颈癌有着重要的临床意义和应用价值。

随着影像学技术的发展，医学图像分割已经成为医学图像分析领域里关键和具有挑战性的问题。而宫颈癌分割又由于复杂的人体腹部组织结构变得尤为复杂，单一成像模式T2-MRI不能完全显示宫颈癌的有效信息。如图1(a)和(b)分别为宫颈癌病人腹部的T2-MR图像和DW-MR图像，宫颈位于膀胱和直肠的中间。从图1(a)可以看出，T2-MR图像具有较高空间分辨率，肿瘤边界比较清晰，但是宫颈的正常组织、肿瘤、膀胱壁和直肠相互之间都有较严重的强度重叠；从图1(b)可以看出在DW-MR图像中肿瘤具有明显较高的灰度值，但是其分辨率低，肿瘤边界较模糊。因此单一成像模式下的宫颈癌分割的自动实现是相当困难的，一些常规的方法例如区域增长和阈值等都无法很好的分割肿瘤。

发明内容

本发明的目的在于提供一种基于T2-MRI和DW-MRI的宫颈癌图像自动分割框架以及利用联合最大后验概率(CMAP)精确分割宫颈癌肿瘤区域的方法，从而进行准确的宫颈癌分割。

为达到上述目的，一种基于T2加权的核磁共振成像T2-MRI和弥散加权的核磁共振成像DW-MRI的宫颈癌图像自动分割方法，包括：

步骤1：利用非线性配准方法将DW-MR图像配准到T2-MR图像，并对配准后的DW-MR图像进行分类；

步骤2：采用非线性各向异性扩散滤波技术对T2-MR图像进行滤波，分割出膀胱和直肠，并利用膀胱和直肠的分割结果将感兴趣区分割出来；

步骤3：对T2-MR图像的感兴趣区和DW-MR图像采用联合最大后验概率CMAP的方法进行肿瘤的精确分割。

本发明充分利用了T2-MR图像和DW-MR图像的有效信息，可以有效地克服T2-MR图像中的噪声、局部容积效应和强度重叠的影响，是一种精确有效的宫颈癌分割方法，对预防、诊断和治疗宫颈癌有着重要的临床意义和应用价值。

附图说明

图1是宫颈癌病人的解剖结构图，(a)是T2-MR图像；(b)是DW-MR图像；

图2是基于T2-MRI和DW-MRI的自动分割框架图；

图3是最大后验概率方法(CMAP)的流程图；

图4是基于T2加权的核磁共振成像(T2-MRI)和弥散加权的核磁共振成像(DW-MRI)的宫颈癌图像自动分割框架实验效果图，(a)是原图；(b)红色轮廓线内为包含肿瘤和正常组织的感兴趣区；(c)配准到T2-MR图像的DW-MR图像；(d)只在T2-MR图像上采用MAP方法的宫颈癌分割结果(即β＝0时)；(e)在T2-MR图像和DW-MR图像上采用CMAP方法的宫颈癌分割结果(即β＝1时)；(f)专家手工分割结果。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚明白，以下结合具体实施例，并参照附图，对本发明进一步详细说明。

本发明的核心思想是一种基于T2加权的核磁共振成像(T2-MRI)和弥散加权的核磁共振成像(DW-MRI)的宫颈癌图像自动分割框架以及利用联合最大后验概率(CMAP)精确分割宫颈癌肿瘤区域的方法，具体步骤包括：首先，利用非线性配准方法将DW-MR图像配准到T2-MR图像(这里采用互信息配准方法作为示例)，并对配准后的DW-MR图像进行分类；然后采用非线性各向异性扩散滤波技术对T2-MR图像进行滤波(这里采用P-M非线性各向异性扩散滤波作为示例)，接着分割出膀胱和直肠，并利用膀胱和直肠的分割结果将感兴趣区(包含正常组织和肿瘤)分割出来；最后对T2-MR图像的感兴趣区和DW-MR图像采用联合最大后验概率(CMAP)的方法进行肿瘤的精确分割。

以下结合具体的实施例对根据本发明提供的这种基于T2加权的核磁共振成像(T2-MRI)和弥散加权的核磁共振成像(DW-MRI)的宫颈癌图像自动分割框架以及利用联合最大后验概率(CMAP)精确分割宫颈癌肿瘤区域的方法进行详细描述。如图2所示为本发明提供的基于T2加权的核磁共振成像(T2-MRI)和弥散加权的核磁共振成像(DW-MRI)的宫颈癌图像自动分割框架的流程图，该方法包括以下步骤：

步骤1：利用非线性配准方法(如互信息配准方法、Demons算法等)将DW-MR图像配准到T2-MR图像，并对配准后的DW-MR图像进行分类；

步骤2：采用非线性各向异性扩散滤波技术对T2-MR图像进行滤波，接着分割出膀胱和直肠，并利用膀胱和直肠的分割结果将感兴趣区(包含正常组织和肿瘤)分割出来；

步骤3：对T2-MR图像的感兴趣区和DW-MR图像采用联合最大后验概率(CMAP)的方法进行肿瘤的精确分割。

上述步骤1包括以下两小步：1)利用非线性配准方法将DW-MR图像配准到T2-MR图像，这里采用互信息配准方法作为示例；2)对配准后的DW-MR图像进行分类是采用的自动阈值分类方法，实现肿瘤的初步分割及定位。

上述步骤2包括以下四步：1)采用非线性各向异性扩散滤波技术对T2-MR图像进行滤波，在去除噪声的同时保持边缘信息，这里采用P-M非线性各向异性扩散滤波作为示例；2)膀胱分割；3)直肠分割；4)感兴趣区分割。

上述第一小步中P-M非线性扩散滤波器公式如下：

$\left\{\begin{array}{c}\frac{\∂I(i;t)}{\∂t}=\mathrm{div}(c(\▿I(i;t))\▿I(i;t))\\ I(i;0)=I_0\left(i\right)\end{array}\right.$

其中I₀(i)为图像在i处的像素强度值，i表示图像的位置信息为(m，n)，是梯度算子，t代表时间，c(·)是扩散系数，两种形式为：

$c(\▿I)=\exp (-\frac{{\left|\right|\▿I\left|\right|}^2}{k^2})$

$c(\▿I)=\frac{1}{1+{\left(\frac{\left|\right|\▿I\left|\right|}{k}\right)}^2}$

其中k是梯度门限。

上述第二小步中膀胱分割采用的方法是主动轮廓模型。该模型的建立是基于我们已有的先验知识——膀胱是腹腔T2-MR图像中灰度值最高的较均匀的整块区域。

上述第三小步中直肠分割利用了直肠位于宫颈的下方这一先验知识和步骤1中的肿瘤初步分割结果，在去除膀胱的T2-MR图像上采用了模糊C均值的算法，得到直肠分割结果。

上述第四小步中感兴趣区分割，在去除膀胱和直肠的T2-MR图像上，利用步骤1中肿瘤的初步分割结果，采用模糊C均值的算法，分割出包含肿瘤和正常组织的感兴趣区，结果如图4(b)。

上述步骤3为对T2-MR图像的感兴趣区和DW-MR图像采用联合最大后验概率(CMAP)的方法进行肿瘤的精确分割，结果如图4(e)。

下面具体介绍上述联合最大后验概率(CMAP)方法具体实施过程。

设y表示图像，y_i，i＝1，2，…N，表示图像在i处的灰度值，N表示图像的像素的总个数。假定图像要被分为K类，以x_i＝k，k＝1，2，…，K代表像素i被归为第k类。联合最大后验概率(CMAP)方法的流程图如图3所示，其具体步骤如下：1)计算T2-MR图像的能量函数U_T2(x)；2)计算DW-MR图像的能量函数U_DW(x)；3)计算T2-MR图像和DW-MR图像的联合能量函数U_T2(x)+βU_DW(x)；4)判断是否满足终止条件，若满足，根据能量最小原则判定肿瘤和正常组织的类别，从而输出精确分割的肿瘤区域；若不满足终止条件返回步骤1)。每一小步中对应的计算公式在如下的讲述中会展示出。

传统的MAP分割算法是求取分割结果

使得后验概率P(x|y)最大。即：

${\stackrel{\mathrm{\Λ}}{x}}_{\mathrm{MAP}}=\arg \underset{x\∈X}{\max }P\left(x\right|y)$

根据贝叶斯公式，P(x|y)可以表示为：

$P\left(x\right|y)=\frac{P\left(y\right|x)P\left(x\right)}{P\left(y\right)}\∝P\left(y\right|x)P\left(x\right)$

从上式可以看出，图像分割的目标转化为并求取分割结果x，使得P(y|x)P(x)最大，即：

${\stackrel{\mathrm{\Λ}}{x}}_{\mathrm{MAP}}=\arg \underset{x\∈X}{\max }\left\{P\left(y\right|x)P\left(x\right)\right\}$

其中P(y|x)和P(x)分别为在给定分割x条件下的概率密度和先验概率。下面计算P(y|x)和P(x)。

假定图像中噪声是高斯白噪声，则图像模型为：

y_i＝μ_ik+n_ik

其中，μ_ik是第k类组织在像素i处的灰度均值，n_ik第k类组织在像素i处对应的高斯白噪声，它的分布服从正态分布

故P(y|x)可以表示为：

$P\left(y\right|x)=\underset{k=1}{\overset{K}{\mathrm{\Π}}}\underset{i\∈R_k}{\mathrm{\Π}}{P}_k\left(y_i\right|x)=\underset{k=1}{\overset{K}{\mathrm{\Π}}}\underset{i\∈R_k}{\mathrm{\Π}}\frac{1}{2\mathrm{\π}{\mathrm{\σ}}_{\mathrm{ik}}}\exp [-\frac{1}{2}{\left(\frac{y_i-{\mathrm{\μ}}_{\mathrm{ik}}}{{\mathrm{\σ}}_{\mathrm{ik}}}\right)}^2]$

$=\frac{1}{{\left(2\mathrm{\π}\right)}^{1/2}}\exp \{-\underset{k=1}{\overset{K}{\mathrm{\Σ}}}\underset{i\∈R_K}{\mathrm{\Σ}}[\ln \left({\mathrm{\σ}}_{\mathrm{ik}}\right)+\frac{1}{2}{\left(\frac{y_i-{\mathrm{\μ}}_{\mathrm{ik}}}{{\mathrm{\σ}}_{\mathrm{ik}}}\right)}^2\left]\right\}$

P(x)由吉布斯分布表示为：

$P\left(x\right)=\frac{1}{a}\exp \{-U_G\left(x\right)\}$

其中，a是归一化常数，U_G(x)为能量函数：

$U_G\left(x\right)=\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{U}_{\mathrm{Gk}}\left(x\right)=\underset{k=1}{\overset{K}{\mathrm{\Σ}}}\underset{i\∈R_k}{\mathrm{\Σ}}\underset{j\∈N_i}{\mathrm{\Σ}}b\·V_k(x_i,x_j)$

其中，N_i为第i位置的邻域，V(x_i，x_j)为势函数，b为常数。

最终将后验概率公式转化为：

P(x|y)∝exp{-U(x)}

分割目标进一步转化为求取分割结果x，使得能量函数最小，即

${\stackrel{\mathrm{\Λ}}{x}}_{\mathrm{MAP}}=\arg \underset{x\∈X}{\max }U\left(x\right)$

$U\left(x\right)=\underset{k=1}{\overset{K}{\mathrm{\Σ}}}\underset{i\∈R_k}{\mathrm{\Σ}}[\frac{1}{2}{\left(\frac{y_i-{\mathrm{\μ}}_{\mathrm{ik}}}{{\mathrm{\σ}}_{\mathrm{ik}}}\right)}^2+\ln \left({\mathrm{\σ}}_{\mathrm{ik}}\right)]+\underset{k=1}{\overset{K}{\mathrm{\Σ}}}\underset{i\∈R_k}{\mathrm{\Σ}}\underset{j\∈N_i}{\mathrm{\Σ}}b\·V_k(x_i,x_j)$

上述第一小步中的T2-MR图像的能量函数U_T2(x)可以根据U(x)的计算公式计算出。

上述第二小步中的DW-MR图像的能量函数U_DW(x)可以根据U(x)的计算公式计算出。

由于T2-MR图像的后验概率P_T2(x|y_T2)和DW-MR图像的后验概率P_DW(x|y_DW)是独立的，因此CMAP的分割方法就是求取分割结果

使得联合后验概率最大，即：

${\stackrel{\mathrm{\Λ}}{x}}_{\mathrm{CMAP}}=\arg \underset{x\∈X}{\max }\left(P_{T2}\left(x\right|y_{T2})P_{\mathrm{DW}}\left(x\right|y_{\mathrm{DW}})\right)$

根据MAP算法中转化为能量函数的公式，将上式转化为

P_T2(x|y_T2)P_DW(x|y_DW)∝exp{-(U_T2(x)+βU_DW(x))}

其中β是权重系数，用来平衡T2-MR图像和DW-MR图像对分割结果的影响程度。最后CMAP分割算法的问题就转化为了求取能量函数最小问题：

${\stackrel{\mathrm{\Λ}}{x}}_{\mathrm{CMAP}}=\arg \underset{x\∈X}{\max }(U_{T2}\left(x\right)+\mathrm{\β}{U}_{\mathrm{DW}}\left(x\right))$

上述第三小步的T2-MR图像和DW-MR图像的联合能量函数U_T2(x)+βU_DW(x)就可以根据上述公式计算出来。

当β＝0时，CMAP方法就成为了MAP方法，此时的实验分割结果如图4(d)；当β取不同的值时，T2-MR图像和DW-MR图像对分割结果的影响程度也就不同。

运行结果

为了验证本发明方法，我们采用手动分割宫颈癌作为标准参考。

图4给出了本发明给出的这种基于T2加权的核磁共振成像(T2-MRI)和弥散加权的核磁共振成像(DW-MRI)的宫颈癌图像自动分割框架实验效果图。(a)原图；(b)红色轮廓线内为包含肿瘤和正常组织的感兴趣区；(c)配准到T2-MR图像的DW-MR图像；(d)只在T2-MR图像上采用MAP方法的宫颈癌分割结果(即β＝0时)；(e)在T2-MR图像和DW-MR图像上采用CMAP方法的宫颈癌分割结果(即β＝1时)；(f)专家手工分割结果。

实验表明，本发明方法-一种基于T2加权的核磁共振成像(T2-MRI)和弥散加权的核磁共振成像(DW-MRI)的宫颈癌图像自动分割框架以及利用联合最大后验概率(CMAP)精确分割宫颈癌肿瘤区域的方法-准确的分割出了宫颈癌图像的肿瘤区域，达到了分割肿瘤的目的。

以上所述，仅为本发明中的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉该技术的人在本发明所揭露的技术范围内，可理解想到的变换或替换，都应涵盖在本发明的包含范围之内，因此，本发明的保护范围应该以权利要求书的保护范围为准。

Claims (12)

1.一种基于T2加权的核磁共振成像T2-MRI和弥散加权的核磁共振成像DW-MRI的宫颈癌图像自动分割方法，包括：

步骤1：利用非线性配准方法将DW-MR图像配准到T2-MR图像，并对配准后的DW-MR图像进行分类；

步骤2：采用非线性各向异性扩散滤波技术对T2-MR图像进行滤波，分割出膀胱和直肠，并利用膀胱和直肠的分割结果将感兴趣区分割出来；

步骤3：对T2-MR图像的感兴趣区和DW-MR图像采用联合最大后验概率CMAP的方法进行肿瘤的精确分割。

2.根据权利要求1所述的方法，其特征在于所述非线性配准方法包括互信息配准方法或Demons算法。

3.根据权利要求1所述的方法，其特征在于所述对配准后的DW-MR图像进行分类采用自动阈值分类方法，实现肿瘤的初步分割及定位。

4.根据权利要求1所述的方法，其特征在于所述非线性各向异性扩散滤波技术包括P-M非线性各向异性扩散滤波。

5.根据权利要求4所述的方法，其特征在于由下式实现P-M非线性各向异性扩散滤波：

$\left\{\begin{array}{c}\frac{\∂I(i;t)}{\∂t}=\mathrm{div}(c(\▿I(i;t))\▿I(i;t))\\ I(i;0)=I_0\left(i\right)\end{array}\right.$

其中I₀(i)为图像在i处的像素强度值，

是梯度算子，div是散度算子，t代表时间，c(·)是扩散系数，两种形式为：

$c(\▿I)=\exp (-\frac{{\left|\right|\▿I\left|\right|}^2}{k^2})$

$c(\▿I)=\frac{1}{1+{\left(\frac{\left|\right|\▿I\left|\right|}{k}\right)}^2}$

其中k是梯度门限。

6.根据权利要求1所述的方法，其特征在于所述膀胱分割采用主动轮廓模型，其中，所述膀胱是腹腔T2-MR图像中灰度值最高的较均匀的整块区域。

7.根据权利要求1所述的方法，其特征在于，所述直肠分割利用直肠位于宫颈的下方这一先验知识和步骤1中的肿瘤初步分割结果，在去除膀胱的T2-MR图像上采用了模糊C均值的算法，得到直肠分割结果。

8.根据权利要求1所述的方法，其特征在于所述感兴趣区分割包括在去除膀胱和直肠的T2-MR图像上，利用步骤1中肿瘤的初步分割结果，采用模糊C均值的算法，分割出包含肿瘤和正常组织的感兴趣区。

9.根据权利要求1所述的方法，其特征在于所述步骤3包括：

1)计算T2-MR图像的能量函数U_T2(x)；

2)计算DW-MR图像的能量函数UD_W(x)；

3)计算T2-MR图像和DW-MR图像的联合能量函数U_T2(x)+βU_DW(x)；

4)判断是否满足终止条件，若满足，根据能量最小原则判定肿瘤和正常组织的类别，从而输出精确分割的肿瘤区域，若不满足终止条件返回步骤1)。

10.根据权利要求9所述的方法，其特征在于所述能量函数U_T2(x)由下式计算：

$U\left(x\right)=\underset{k=1}{\overset{K}{\mathrm{\Σ}}}\underset{i\∈R_k}{\mathrm{\Σ}}[\frac{1}{2}{\left(\frac{y_i-{\mathrm{\μ}}_{\mathrm{ik}}}{{\mathrm{\σ}}_{\mathrm{ik}}}\right)}^2+\ln \left({\mathrm{\σ}}_{\mathrm{ik}}\right)]+\underset{k=1}{\overset{K}{\mathrm{\Σ}}}\underset{i\∈R_k}{\mathrm{\Σ}}\underset{j\∈N_i}{\mathrm{\Σ}}b\·V_k(x_i,x_j)$

其中，y_i，i＝1，2，…N，表示图像在i处的灰度值，N表示图像的像素的总个数。假定图像要被分为K类，以x_i＝k，k＝1，2，…，K代表像素i被归为第k类。μ_ik是第k类组织在像素i处的灰度均值，

第k类组织在像素i处对应的高斯白噪声的方差，N_i为第i位置的邻域，V(x_i，x_j)为势函数，b为常数，R_k是所有属于第k类组织的像素的集合。

11.根据权利要求9所述的方法，其特征在于所述联合能量函数U_T2(x)+βU_DW(x)由下式计算：

其中，β是权重系数。

12.根据权利要求11所述的方法，其特征在于通过比较不同的β值下的分割结果，得到精确的肿瘤分割结果。

Images