一种基于泊松模型的X光胸片图像肋骨抑制方法

摘要

本发明公开了一种基于泊松模型的X光胸片图像肋骨抑制方法，首先对图像进行Contourlet变换，利用变换系数在不同尺度间的相关性，对图像的纹理进行增强处理；然后根据图像的海森矩阵特征值，设计条形检测滤波器，检测出肋骨所在的区域；最后结合增强后的纹理和肋骨区域信息，建立并求解肋骨抑制的泊松模型，进而实现图像中的肋骨抑制。本发明利用Contourlet变换的各向异性和多方向特点，结合系数的尺度和方向信息，提高了纹理和噪声之间的区分度，增强纹理，同时抑制了噪声；通过求解泊松模型实现抑制肋骨，不需要对肋骨进行精确分割，避免了显式抑制肋骨时造成的边缘过渡不自然问题，有效地抑制肋骨，提高X光胸片图像的观测效果。

说明书

技术领域

本发明涉及图像处理技术领域，特别涉及一种基于泊松模型的X光胸片图像肋骨抑制方法。

背景技术

X光胸片图像被广泛应用于肺部疾病的诊断，如肺炎、肺结核、肺癌等。其中肺癌死亡人数为所有癌症死亡人数之首，若能早期检测出肺部疾病，能够提高疾病的治疗效果，并降低肺癌的死亡率。尽管在检测性能上X光胸片不如CT，但是其具有的高经济性和低辐射剂量，X光胸片仍然是目前最主要的检测手段。

X光胸片成像是基于肺部中不同的组织在密度和厚度上存在差异，使得X射线穿透过程中被吸收的程度不同，从而在射线接收面上形成不同强度的影像。由于射线是投影到一个二维平面上，因此在射线方向上的不同组织会重叠显现在图像上，造成局部病变区域难以观测与识别。叠影是X光胸片所固有的问题，也是很多肺部疾病诊断决策失误的主要原因。研究表明，医生对胸片诊断的误诊率达30％，而其中82％-95％的漏诊是由于肋骨对肺部软组织的叠加或遮挡。获取仅含肺部软组织的图像的常规方法是双能量差技术，这种技术利用专用设备对被检者进行两次不同能量的曝光透视，将肋骨和软组织进行分离，且需要更大的辐射剂量，仅局限于医院等机构。

近年来，研究人员利用数字图像处理技术解决胸片图像中的叠影问题，提出了一些解决方法。如Giger等人提出通过一种图像差分技术，提高肺部结节的检测率(Medical Physics.，vol.17，pp.861–865，1990.)。Keserci等设计一种滤波器对图像进行滤波，抑制图像中长条形目标来消除肋骨的影响，并增强圆形目标以突出图像中的肺结节(Medical Image Analysis.，vol.6，pp.431–447，2002)。Loog等提出一种基于回归的通用滤波器框架来抑制骨骼结构，这种方法依赖训练图像对滤波器进行学习，通过对测试图像的重建获取软组织图像(Computer Vision Approaches to Medical Image Analysis.，New York：Springer，2006，vol.4241，Lecture Notes in ComputerScience，pp.166–177.)。Suzuki等利用人工神经网络抑制肋骨并增强肺结节(IEEE Trans.Med.Imaging，2006，25(4)：406-416.)。Lee等人通过分割出肋骨，并采用遗传算法优化一个对照模型，进而消除肋骨(Computers&Mathematics with Applications,，2012，64(5)：1390-1399.)。概括起来，这些方法可以分为隐式和显式的肋骨抑制，显式方法定位出肋骨的位置后加以抑制，需要对肋骨进行精确分割，而隐式方法需要大量的样本进行训练，以建立回归模型，因此不易于临床应用。

发明内容

本发明所要解决的技术问题是提供一种基于泊松模型的X光胸片图像肋骨抑制方法，解决数字胸片中肋骨对软组织的叠加遮挡而造成难以识别的问题，满足阅片人或是计算机对数字胸片中病灶检测识别的需要。

为解决上述技术问题，本发明采用的技术方案是：

一种基于泊松模型的X光胸片图像肋骨抑制方法，包括以下步骤：

步骤一、读取原始的胸片图像I，对胸片图像I进行Contourlet变换；

步骤二、对Contourlet变换后的低频系数进行对比度拉升；

步骤三、对Contourlet变换后的高频系数进行增益变换，重构出纹理增强后的胸片图像E；

步骤四、根据原始的胸片图像I，在多个尺度下进行高斯滤波，计算滤波后图像的海森矩阵及其特征值，进而建立条状检测滤波器，利用条状滤波器检测定位肋骨，获得肋骨区域的二值掩膜图像；

步骤五、根据步骤三得到的胸片图像E和步骤四得到的肋骨区域的二值掩膜图像，得到输出图像I^*，即对胸片图像E采用拉普拉斯模板进行卷积运算，获得该图像梯度场的散度值建立泊松方程其中，G是输出图像I^*中肋骨区域内像素的灰度值，I是原始的胸片图像，是肋骨区域边界，即步骤四的二值掩膜图像中的肋骨区域边界。

根据上述方案，所述步骤一中，在进行Contourlet变换之前，在原始的胸片图像I的行和列方向上分别进行一定量的平移，以克服伪吉布斯问题。

根据上述方案，所述对比度拉升采用的公式为其中，为调整后的系数，h1，h2分别为高频和低频增益，M表示低频系数的均值，k∈[0,1]是对比度调节因子。

根据上述方案，所述增益变换公式为其中为增强后的变换域系数，p∈(1,∞)为增益因子，T_g是增益阈值，W_k为增益系数；所述增益系数W_k具体表示为其中，C′_j,k为尺度j，方向k上的变换域系数矩阵C_j,k通过插值后得到的矩阵，m,n为行列索引，J为尺度分解级数，0＜j≤J。

根据上述方案，步骤四中所述尺度为3个。

根据上述方案，所述条状检测滤波器根据滤波器F的输出结果，采用Otsu算法进行自适应阈值分割，获得肋骨区域的二值掩膜图像。

与现有技术相比，本发明的有益效果是：采用Contourlet变换，结合系数的尺度和方向信息，提高了纹理和噪声之间的区分度，使得增强纹理的同时抑制噪声。通过求解泊松模型抑制肋骨，在检测得到的肋骨区域中进行抑制。本发明不需要精确分割出肋骨，仅根据肋骨所在的大致区域即可实现，降低了处理问题的复杂性，同时解决了直接减除肋骨灰度而形成的过渡不自然问题。本发明方法仅需要很少的经验参数，整个处理过程无需人工干预，方便实用，易于集成至专用的数字胸片阅片设备或图像处理系统中。

附图说明

图1为本发明方法的实现流程示意图。

图2为本发明中原始图像与增强图像对比图。

图3为本发明中肋骨区域检测效果图。

图4为本发明中原始图像与肋骨抑制图像对比图。

具体实施方式

下面结合附图和具体实施方式对本发明作进一步详细的说明。如图1所示，一种基于泊松模型的X光胸片图像肋骨抑制方法，包括下述步骤：

步骤1，读取原始的胸片图像I，对原始的胸片图像I进行Contourlet变换。由于在Contourlet变换中，会对信号进行下采用和上采用处理，导致系数发生显著变化，这种情况称为伪吉布斯现象。也即是说，Contourlet变换不具备平移不变性，增强后重构图像中边缘点附近会产生较大的幅值振荡。为了克服伪吉布斯问题，在进行Contourlet变换前，在图像的行和列方向上分别进行一定量的平移。

步骤2，对Contourlet变换后的低频系数按照公式(1)进行对比度拉升，

${\hat{C}}_0(m,n)=(h1-h2)\left(k\right(C_0\left(m,n\right)-M)+M)---\left(1\right)$

其中，为调整后的系数，h1，h2分别为高频和低频增益，M表示低频系数的均值，k∈[0,1]是对比度调节因子，在保持图像原貌的基础上增强局部对比度。

步骤3，对Contourlet变换后的高频系数按照公式(2)进行增益变换后，重构出纹理增强后的胸片图像E。

${\hat{C}}_{j,k}(m,n)=\left(\begin{array}{cc}{\left(W_k\right(m,n)/T_g)}^p{C}_{j,k}(m,n)& \left|C_{j,k}(m,n)\right|\ge T_g\\ (1-{\left(\right(1-W_k(m,n))/(1-T_g\left)\right)}^p)C_{j,k}(m,n)& \left|C_{j,k}(m,n)\right|<T_g\end{array}\right)---\left(2\right)$

其中为增强后的变换域系数，p∈(1,∞)为增益因子，T_g是增益阈值，W_k为增益系数。综合考虑信号和噪声在尺度间的相关性计算增益系数，先对每个分解方向上的系数模值，通过二次插值变换到最细尺度上，然后进行求和并归一化。

$W_k(m,n)=\underset{j=1}{\overset{J}{\Sigma }}\left|C_{j,k}^{\prime }(m,n)\right|/max\left\{\underset{j=1}{\overset{J}{\Sigma }}\right|C_{j,k}^{\prime }(m,n)\left|\right\}---\left(3\right)$

其中，C′_j,k为尺度j，方向k上的变换域系数矩阵C_j,k通过插值后得到的矩阵，m,n为行列索引，J为尺度分解级数，这里只考虑高频系数，因此0＜j≤J。W_k是归一化后的系数矩阵，元素值范围为[0,1]。

步骤4，读取原始的胸片图像，在多个尺度下对图像进行高斯滤波得到I_s。尺度数量越多算法精度越高，但是计算量也会增大，综合考虑，尺度的数量一般选择3个为宜。尺度大小由高斯核的标准差确定，可根据肋骨在图像中的比例进行估计。在尺度为s的图像I_s中，计算所有像素点处的海森矩阵H_i。计算H_i的特征值λ₁,λ₂，构造条状检测滤波器F：

$F=e^{{({\lambda }_2/{\lambda }_1)}^2}(1-e^{{({\lambda }_1+{\lambda }_2)}^2})---\left(4\right)$

利用该滤波器增强多个尺度下的肋骨，取各个尺度下的最大值作为输出结果。对输出结果利用Otsu算法进行自适应阈值分割，获得肋骨区域的二值掩膜图像。

步骤5，根据纹理增强图像E的梯度信息和肋骨区域的边界信息，通过求解解泊松方程，获得最终的输出图像I^*。先对步骤3得到的纹理增强图像E，利用拉普拉斯模板进行卷积运算，获得该图像梯度场的散度值建立泊松方程：

$\begin{array}{ccc}\Delta G=DIV\left(E\right)& s.t.& G{|}_{\partial \Omega }=I{|}_{\partial \Omega }\end{array}---\left(5\right)$

其中，G是输出图像I^*中肋骨区域内像素的灰度值，I是原始图像，是由步骤4得到的肋骨区域边界。

图2给出了原始图像与纹理增强图像对比图，在图中，原始图像局部给出了一例肺部局部图像，增强后的图像突出了边缘细节特征，图像对比度也有了提升。图3给出了肋骨检测得到的二值掩膜图像，其中掩膜包含了肋骨所在的大致区域。图4给出了原始图像与肋骨抑制图像在肺野区域中的对比图，可以看出在肋骨抑制后，其重叠区域的纹理更加清晰可见。