定量超声系统中基于经验模态分解的组织微观结构检测方法

摘要

本发明公开了一种定量超声系统中基于经验模态分解的组织微观结构检测方法，其检测方法为：步骤一、获取生物组织超声背散射信号；步骤二、截取DOI区域的背散射信号作为输入数据；步骤三、对输入数据进行EMD分解，得到三阶本征模态函数和一个残差项；步骤四、提取三阶本征模态函数能量谱，并与原信号能量谱进行相似度分析，确定包含组织微结构信息的BIMF；步骤五、对确定出的BIMF进行估计，计算组织散射子间距信息。本发明具有如下优势：(1)提取模态函数项BIMF，去除非相干模态函数项中的无用信息；(2)消除弥散散射元信号及噪声干扰，抑制倒谱中的干扰峰值，提高MSS估计精度，估计精度可达90～95％。

说明书

技术领域

本发明涉及一种利用超声背散射技术对表征组织微观结构的平 均散射子间距进行检测的方法。

背景技术

医学超声广泛应用于人体组织器官成像中，虽然该技术在临床应 用中已经具有相当成熟的应用基础，并能提供实时的高精度人体器官 图像，但其不能提供由于组织钙化、变异等导致的微观结构变化。

目前一些等效散射子模型如二元混合模型、Faran圆柱模型以及 弱散射子模型等常用来表征组织微观结构变化。在这些模型中，组织 的检测转化为对超声衰减、背散射系数、散射子密度、散射子直径和 平均散射子间距(MSS)等的检测。其中MSS是一种最常用的方法， 该方法主要利用了生物组织在声学特性分布上具有一定的准周期性， 这种生物周期性来自组织结构中的规则散射元和弥散散射元，会随着 组织健康状况改变而发生变化，通过利用MSS测量该周期性，可以 对生物组织微观结构进行测量，并对如血管动脉粥样硬化、骨质疏松 等健康状况做出诊断。

在平均散射子间距估计中目前常用的方法主要有倒谱方法、功率 谱方法和瞬态相关图等方法。这些方法在一定程度上能够有效检测平 均散射子间距，但受背散射超声信号序列短、衰减大和噪声大的影响， 检测结果不稳定。一些自相关谱分析方法、周期性小波变换方法、奇 异谱分析方法等提高了检测鲁棒性，但计算量大，在去除噪声和弥散 散射信号干扰方面有待提高。

发明内容

为了解决定量超声系统中利用现有技术方法无法有效消除超声 背散射信号中的噪声和弥散散射信号干扰的问题，本发明提供了一种 定量超声系统中基于经验模态分解(EMD)的组织微观结构检测方 法。

本发明的目的是通过以下技术方案实现的：

一种定量超声系统中基于经验模态分解的组织微观结构检测方 法，包括如下步骤：

步骤一、获取生物组织超声背散射信号；

步骤二、截取DOI(Domain ofInterest)区域的背散射信号作为 输入数据；

步骤三、对输入数据进行EMD分解，得到三阶本征模态函数和 一个残差项；

步骤四、执行能量谱筛选：提取三阶本征模态函数能量谱，并与 原信号能量谱进行相似度分析，确定包含组织微结构信息的最优模态 函数项(BIMF)；

步骤五、对确定出的BIMF进行估计，计算组织散射子间距信息。

所述步骤一中，超声背散射信号获取有仿真和实验两种获取方 式。其中实验获取方式中，超声背散射信号收发由2MHz高频单阵元 聚焦探头实现，探头聚焦距离为25mm，探测方式为单一探头背散射 作用方式。

所述步骤二中，截取DOI区域信号时，采用矩形窗函数。仿真 数据截取信号序列时间长度为4μs，实验数据截取信号序列长度为 8.125μs。

所述步骤三中，基于经验模态分解的组织微观结构检测方法具体 通过以下迭代进行检测：

1)输入DOI区域背散射信号x(t)(t＝1，2，...，N，N为采样点数) 为原始函数；

2)确定输入的原始函数极值点，采用三次样条函数确定输入信 号上下包络x_L(t)x_U(t)；

3)计算上下包络均值定义h_i(t)为原始函数与 该包络均值的差值函数h_i(t)＝x_i(t)-m_i(t)；(i表示第i次分解，分解结束得 到第i个本征模态函数)

4)迭代计算判断差值函数h₁(t)是否满足IMF条件，判断决策过 程采用标准偏差(Standard Deviation，SD)准则：

$H_{\mathrm{SD}}={(h_{(k-1)}\left(t\right)-h_k\left(t\right))}^2/h_{(k-1)}^2\left(t\right)$(h_k(t)为第k次迭代差值函数)。若

H_SD≤0.15，则第i个本征模态函数成分IMF i定义为

c_i(t)＝h_i(t)(i＝1，2，...，N)，继续执行第5步；否则差值函数h_i(t)＝x_i(t)-m_i(t) 作为输入的原始函数x_i(t)返回执行第2步；

5)存储IMF分解次数为N。判断若N≥3，继续执行第6步， 否则剩余函数x_i(t)＝x_i-1(t)-c_i(t)返回执行第2步(c_i(t)为第i个本征模态 函数)；

6)经过上述迭代计算，获得原始函数3个本征模态函数成分 IMF1，IMF2，IMF3：{c₁(t)，c₂(t)，c₃(t)}。

所述步骤四中，执行能量谱筛选，分别求得三阶本征模态函数 IMF1-3能量谱，根据能量谱极值点值及其位置，确定包含组织微结 构信息的最优模态函数项(BIMF)。

所述步骤五中，利用倒谱估计方法对模态函数项BIMF进行估计， 由计算结果估计组织等效散射子间距，利用该散射子间距表征生物组 织特性。

本发明所述的一种定量超声系统中基于经验模态分解的组织微观 结构检测方法，可应用于生物医学检测、组织定征领域。与现有检测 技术及方法相比，具有如下优势：

(1)提取模态函数项BIMF，去除非相干模态函数项中的无用信 息；

(2)消除弥散散射元信号及噪声干扰，抑制倒谱中的干扰峰值， 提高MSS估计精度，估计精度可达90～95％。

附图说明

图1为本发明中技术方法流程图；

图2为本发明中超声组织背散射仿真数据获取示意图；

图3为本发明中生物组织目标区域(DOI区域)示意图；

图4为本发明中超声组织背散射实验数据获取示意图；

图5为本发明中仿真背散射信号；

图6为本发明中仿真背散射信号EMD分解所得三阶本征模态函 数(IMF1-3)；

图7为本发明中实验背散射信号；

图8为本发明中实验背散射信号EMD分解所得三阶本征模态函 数(IMF1-3)；

图9为本发明中仿真数据原始信号的平均散射子间距倒谱估计结 果；

图10为本发明中仿真数据最优模态函数项BIMF的平均散射子间 距倒谱估计结果；

图11为本发明中实验数据原始信号的平均散射子间距倒谱估计 结果；

图12为本发明中实验数据最优模态函数项BIMF的平均散射子 间距倒谱估计结果。

具体实施方式

下面结合附图对本发明的技术方案作进一步的说明，但并不局限 于此，凡是对本发明技术方案进行修改或者等同替换，而不脱离本发 明技术方案的精神和范围，均应涵盖在本发明的保护范围中。

本发明提供了一种基于有限元的超声理疗骨组织内部应力分布数 值模拟方法，其实现的过程如图1所示，具体包括以下步骤：

(1)获取生物组织超声背散射信号。

获取生物组织超声背散射信号的方式包括仿真和实验两种获取 方式。

仿真数据获取方式采用有限元数值求解方法。在有限元求解过程 中，对目标区域DOI进行离散化处理，离散单元采用2-D8节点单元 Plane-82，DOI区域大小为16×32mm，散射子大小为0.1±0.05mm， 散射子数量为32，散射子间距为0.5。采用自由边界条件，单元离散 化尺寸为超声波长的1/10。利用有限元方法求解弹性波动方程：

$\rho \frac{{\partial }^{2}u}{{\partial t}^2}=(\mu +\eta \frac{\partial }{\partial t}){▿}^{2}u+(\lambda +\mu +\xi \frac{\partial }{\partial t}+\frac{\eta }{3}\frac{\partial }{\partial t})▿▿·u;$

其中λ和μ分别为第一阶和第二阶拉梅常量，u＝u(x，y，z，t)为位移场矢 量，η和ξ为粘度，ρ为组织密度。

图2所示为当超声频率为1MHz，超声速度为2073m/s时弹性波 动方程的求解结果。

实验背散射信号的获取采用来自骨组织的超声背散射信号。图3 所示为选定的骨组织DOI区域，该区域位于股骨两端，区域上线边 界法线平行于受力方向。采用2MHz单阵元超声探头，焦距25mm， 由三维移动平台搭载进行DOI区域扫描，如图4所示。

(2)采用矩形窗函数截取DOI(Domain ofInterest)区域的背散 射信号作为输入数据。

获取超声背散射数据采用矩形窗函数截取来自DOI区域的组织 背散射信号。对于仿真超声背散射信号，矩形窗截取区间 7.5μs～11.5μs，采样数为805，如图5所示；对于实验超声背散射信 号，矩形窗截取区间为50μs～58.125μs，采样数为650，如图7。

(3)对输入数据进行EMD分解，得到三阶本征模态函数 (IMF1-3){c₁(t)，c₂(t)，c₃(t)}。

对输入信号进行EMD分解时，采用三次样条插值函数确定输入 信号上下包络。在迭代计算中标准偏差SD较大的区域，将该区域分 离，并对分离部分增加迭代次数。该操作通过定义阈值θ＝0.15和窗函 数w(t)实现：

$w\left(t\right)=\left(\begin{array}{cc}1& \left|\frac{m\left(t\right)}{\frac{1}{2}·(e_{\max }\left(t\right)-e_{\min }\left(t\right))}\right|>\theta \\ 0& \mathrm{others}\end{array}\right);$

则迭代计算时，步骤三第4步差值函数变为h_i(t)＝x_i(t)-w(t)·m_i(t)。

仿真和实验超声背散射数据经过EMD分解后得到的三阶模态函 数如图6和图8所示。

(4)分别求得3阶本征模态函数IMF1-3能量谱，根据能量谱极 值点值及其位置，筛选包含组织微结构信息的最优模态函数项 (BIMF)。

经窗函数截取的初始信号极值点值序列为V＝{v₁，v₂，v₃，...，v_n}，极值 点位置序列为P＝{p₁，p₂，p₃，...，p_n}。经EMD分解的三阶模态函数极值 点值序列分别为V₁，V₂，V₃，计算其与原信号极值点值序列V相关性为 C_V1，C_V2，C_V3；经EMD分解的三阶模态函数极值点位置序列分别为 P₁，P₂，P₃，计算其与原信号极值点值序列V相关性为C_P1，C_P2，C_P3。模态 函数项BIMF定义为C_Vi+C_Pi(i＝1，2，3)取最大值时的模态函数项。

(5)对确定出的BIMF做倒谱估计，计算组织散射子间距信息。

采用倒谱散射子估计方法，对选定的模态函数项BIMF进行倒谱 估计。仿真超声背散射数据原信号倒谱估计结果及BIMF散射子间距 估计结果如图9和图10所示。实验超声背散射数据原信号倒谱估计 结果及BIMF散射子间距估计结果如图11和图12所示。可以看出， 干扰峰值得到抑制，反应MSS的正确倒谱峰值得到增强。

倒谱峰值所对应时间为t_cmax，则MSS值可用下式进行计算：

$\mathrm{MSS}=\frac{t_{c\max }·V_s}{2};$

式中，V_s为超声在介质中的传播速度。