[发明专利]一种基于超声图像纹理特征定征生物组织的方法

权利要求书

1.一种基于超声图像纹理特征定征生物组织的方法，其特征在于，包括以下步骤：

步骤1：由超声探头扫描目标区域后接收到的回波信号进行波束成形处理，得到原始的M行N列的超声背向散射信号矩阵；N列：生成的矩阵中有N条扫描线，M行：每条扫描线中包含M个采样点；两个相邻采样点之间的距离为Intaxi米，两个相邻扫描线之间的距离为Intlat米；

步骤2：对步骤1中得到的背向散射信号矩阵进行包络检测处理得到M行N列的超声的包络信号矩阵；

步骤3：对步骤2中得到的包络信号矩阵进行对数压缩/扫描变换，同时进行尺度缩放后转化为360x360大小的矩阵，再通过256级灰度映射为B-mode图像；

步骤4：对步骤2中得到的M行N列的包络信号矩阵，采用滑动窗口法来估算其大小同样为M行N列的Nakagami-m参数矩阵；

采用一个大小为M_wхN_w的矩形滑动框从包络矩阵的左上角第一个元素开始从左向右、从上到下以一个元素为步长逐行滑动遍历整个包络矩阵；每个矩形窗包含M_wхN_w个数据点，M_w＝3хpulselength/Intlat，N_w＝3хpulselength/Intaxi，其中pulselength为超声脉冲长度，表示向上取整，每滑动一次窗口位置，则对窗口内的包络信号的概率密度函数建模为Nakagami分布：

其中Γ(·)表示伽马函数，U(·)表示单位阶跃函数，r表示每一个数据变量；该分布的尺度参数Ω和形状参数m分别通过下式计算：

Ω＝E(R²) (2)

其中，E(·)表示统计均值，R表示大小为1xM_w*N_w大小的包络数据向量；由此得到一个大小为(M-M_w)х(N-N_w)的m参数矩阵，再对此矩阵进行三样条插值，得到一个大小为MхN的m参数矩阵；

步骤5：对步骤4中得到大小为MхN的m参数矩阵进行伪彩色的参数映射，最终得到一幅Nakagami图像；

步骤6：将步骤5中得到的Nakagami伪彩色成像进行256级数的灰度映射得到360x360大小的Nakagami灰度成像；

步骤7：从步骤3获得的B-mode图像和步骤6获得的Nakagami灰度图像中提取纹理特征；

其具体操作过程如下：

7.1手动的在B-mode图像和Nakagami灰度图像中的高亮区域分别圈定出变性区域，该区域以外的区域则视为正常组织区域；

7.2如图3所示，用像素大小为33x33的矩形框分别在每幅B-mode图像和Nakagami灰度图像的变性区域和正常组织区域提取出4组样本数据矩阵；

7.3针对每组数据矩阵，分别计算其灰度共生矩阵GLCM、灰度游程矩阵GLRLM中的相应特征；

下面先分别介绍一下GLCM与GLRLM：

GLCM：是一种描述图像中灰度对的联合概率密度函数PDF估计值的矩阵；设该图像具有H个灰度级，则GLCM用公式(4)表述：

P_θ(i,j)(i,j＝0,1,2,...,H-1) (4)

其中i和j是特定位置、距离下两个像素的灰度值；θ为两个像素间的位置关系；不同θ的决定了两像素之间的距离、角度关系；选择0°、45°、90°和135°四种角度；

能量表示图像中像素对的重复性，用于度量图像灰度均匀性的特征，由公式(5)定义：

对比度定义见式(6)：

同质性也称“逆差矩”，同质性定义见公式(7)：

相关性用于度量相邻像素、特定位置像素之间灰度的线性相关性，相关性定义见公式(8):

其中，为均值，

为标准差；

GLRLM的元素ρ(I,J|θ)可被定义为在特定角度θ下某一灰度级数I和相应的最大游程数目J出现的次数；其中θ分别可选0°、45°，90°和135°；在GLRLM中的0°矩阵中对应位置元素为1；同理，方向为45°、游程数目为2的有两组；对应在GLRLM的45°矩阵中对应位置元素为2；

通过对GLRLM进行不同的运算进一步提取各种纹理特征；假定图像的灰度级为N_g，出现的不同游程的数量定义为N_r，定义图像中的像素总数为N_p；

所需的GLRLM图像纹理特征，计算方法见公式9至12：

1)长游程增强LRE特征

2)短游程低灰度增强SRLGE特征

3)长游程低灰度增强LRLGLE特征

4)长游程高灰度增强LRHGLE特征

7.4对从每幅B-mode图像中提取出的33x33样本矩阵，由7.3计算出的纹理特征中选择Contrast、Homogeneity、LRE、LRLGLE和LRHGLE这5个特征为作为最终特征参量；对每幅Nakagami灰度图像中提取出的33x33样本矩阵，由7.3计算出的纹理特征中选择Contrast、Correlation、Energy、Homogeneity和SRLGE这5个特征作为最终特征参量；由此，从每幅B-mode图像中的变性组织区域和正常组织区域中提取的特征向量为[Contrast,Homogeneity,LRE,LRLGLE,LRHGLE]；从每幅Nakagami灰度图像中的变性组织区域和正常组织区域中提取的特征向量为[Contrast,Correlation,Energy,Homogeneity,SRLGE]；

步骤8：将步骤7的7.1-7.4中得到的B-mode图像对应的组织变性区域特征向量与Nakagami灰度图像对应的变性组织区域特征向量组合在一起作为最终的变性组织区特征向量矩阵A；同样，将步骤7的7.1-7.4中得到的B-mode图像对应的正常组织区域特征向量与Nakagami灰度图像对应的正常组织区域特征向量组合在一起作为最终的正常组织特征向量矩阵Nor；

步骤9：为步骤8中生成的变性组织区域特征向量矩阵A创建对应标签向量ALabel，标签为-1；为步骤8中生成的正常组织区域特征向量矩阵D创建对应标签向量DLabel，标签为1；

步骤10：将步骤9中得到的特征向量矩阵A和D，以及对应的标签向量ALabel和DLabel作为支持向量机SVM的输入参量，利用SVM的有监督学习训练机制构建二元分类器模型，以便区分正常组织和变性组织；

通过SVM来训练出二分类模型，而训练样本由长度为N的提取了特征向量的数据集和相应标签构成，其形式为(x_k,y_k)，k＝1,2,…,N为下标，x属于实数空间，y_k∈(1,-1)，每一个样本由特征向量和对应标签组成；这里x_k为特征向量，y_k的值定义了x_k的类别，这里用y_k的值来进行区分，即如果y_k＝-1，x_k属于变性组织，如果y_k＝1，x_k属于正常组织；若样本数据线性可分，其超平面即二维空间里的一条直线，它的代数表达式为：

g(x)＝w^T·x_k+b (13)

其中w代表权系数向量，b是常数项；训练的目的就是构造一个决策函数，将测试数据尽可能的正确分类；将分类超平面进行归一化处理来建立一个比较规范的SVM数学模型；使其满足：

使两类中距离最近的两个点的间隙越大，说明把两个类的点分的越清晰，这样就更加准确的进行分类；所以需要的就是使距离M达到最大：

寻找最优超平面的问题就转化为一个经典的最优化约束的求解问题；而对于最基本的线性可分的情况，需要解决的是如下这个二次规化的问题：

它的约束条件为：

y_k(w^Tx_k+b)≥1,k＝1,2,3,...,N (17)

对于解决这个数学问题的通常方法是引入拉格朗日乘子α^[1]，这样公式(16)变为：

其中，k为下标1,2，…,N；约束优化的求解方法取决于公式(18)的鞍点，因此：

把公式(19)和(20)带入公式(18)，得出该优化问题的对偶形式：

其中，k和d表示下标1,…,N；其约束公式为：

于是问题转化为找到最大的Q，只要确定了α_k^*表示Q取得最大值时的α_k，k＝1,…,N；的一组唯一解；随即确定即最优的超平面的权系数向量w^*是样本向量的线性组合；根据KKT条件，该优化问题的解必须满足：

a_k{y_k[(w·x_k)+b]-1}＝0,k＝1,...,N (23)

因此，对多数样本对应的α_k将为0，把α_k≠0对应于使公式(14)中等号成立的样本称为支持向量，求解完上述问题后，就求得分类函数的判别式：

由于非支持向量对应的α_k均为0，因此上式的求和实际上只对支持向量进行；b^*是Q取得最大值时的最佳分类阈值，由任意一个支持向量用式(14)求得；

引入非负的松弛变量ξ_k≥0(k＝1,…,N为下标)和误差惩罚因子C^[3]，用以缓解点被分错的情况，以保证正常的判决结果，由此公式(16)和(17)变成：

y_k(w^Tx_k+b)≥1-ξ_k,k＝1,2,...,N (26)

对于这个问题的唯一改进是对拉格朗日乘子的限制，它表示为：

0≤α_k≤C,k＝1,2,...,N (27)

假设有非线性映射Φ：R_n→Z将输入空间的样本映射到高维特征空间Z中，当在特征空间中构造最优超平面时，训练算法及使用特征空间中的点积，即Φ(x_k)·Φ(x_d)；所以，若能找到一个函数K使得K(x_k,x_dd)＝Φ(x_k)·Φ(x_d)，k,d＝1,..,N为下标；

根据泛函数的有关理论，只要一种函数K(x_k,x_d)满足Mercer条件，它就对应某一变换空间中的内积；因此，在最优分类面中采用满足Mercer条件的内积函数K(x_k,x_d)就实现某一非线性变换后的线性分类；

设训练样本输入空间的维数为e，表示从输入空间到特征空间的一个非线性变换的集合，q是特征空间的维数；定义一个特征空间中充当决策面的超平面：

若令w＝[b,w₁,w₂,...,w_q]^T，Φ(x)＝[1,φ₁(x),φ₂(x),...,φ_q(x)]^T，T表示转置，则决策超平面进一步化简为如下形式：

w^TΦ(x)＝0 (29)

此时，式(21)所示的对偶形式的目标函数变为：

其中，为Q取得最大值时的最优解，下标k和d表示所有数据点中可用作支持向量的点的下标，总数量s＜N,最优的超平面的权系数向量w^*表示为

相应的决策函数就成为：

式中，sgn表示符号函数，b^*代表Q取得最大值时的最佳分类阈值；算法的其他条件不变，由于最终的判别函数中实际只包含与支持向量的内积以及求和，因此判别时计算复杂度取决与支持向量的个数；在公式(30)和(32)中，满足Mercer条件的内积函数K(x_k,x_d)称为核函数，选择高斯径向基核函数，其公式为：

K(x_k,x_d)＝exp(-||x_k-x_d||²/2g) (33)

其中，g表示作用宽度因子，S确定C和g分别取值为2.6和0.5；

步骤11：当利用步骤10中的方法构建好SVM分类识别模型后，对于新生成的任意一组B-mode图像和Nakagami灰度图像以33x33大小的矩阵进行逐点滑动识别；具体而言，首先对B-mode图像和Nakagami灰度图像以0进行边界填充，将填充后的图像以33x33像素大小的像素块进行逐点块分割处理；原始B-mode图像和Nakagami灰度图像大小为每幅360x360，对其进行边界填充并块处理后共分割为129600个子图像；对每个子图像采用步骤7.3-7.4的方法提取特征向量，随后采用步骤8的方法构建SVM输入特征向量矩阵，送入步骤10中构建好的SVM分类器中进行分类判别，得到凝固区检测的0-1二值图，此时0表示分类结果为变性组织，1表示分类结果为正常组织；

步骤12：对步骤11中得到的初步识别的结果图进行形态学处理得到形态学处理后的变性区识别结果；

步骤13：对步骤12中得到形态学处理后的结果图分别沿每行和每列进行4阶最小二乘多项式拟合处理，并且以曲线勾勒出拟合后的-6dB等高线，此等高线以内部分被视为组织变性区；

步骤14：将步骤13中得到的多项式拟合后的-6dB等高线叠加到对应的B-mode图像上，则得到了生物组织变性区域的最终识别结果图。