CN106361278A - 一种单次激励的感应式磁声快速成像方法 - Google Patents

Abstract

本发明涉及生物医学成像领域，公开了一种单次激励的感应式磁声快速成像方法，包括压缩感知算法设计模块和声源重建模块；压缩感知算法设计模块用于对曲波正交基的稀疏化、变换点扩散函数的压缩感知矩阵评估、Dantzig选择器的约束优化问题求解和曲波变换的图像重建；声源重建模块针对单次激励产生的磁声信号进行傅立叶变换、应用压缩感知算法的稀疏信号重建、通过曲波变换和傅立叶反变换的声源图像重建。本发明通过稀疏采样，选择最优传感器数目，减少数据采集量和处理时间，同时利用压缩感知算法的特性，实现单次激励模式的数据重建，提高数据的可靠性和稳定性，简化成像系统，具有很高的产业化前景。

Description

一种单次激励的感应式磁声快速成像方法

技术领域

本发明涉及生物医学成像领域，特别是涉及一种单次激励的感应式磁声快速成像方法。

背景技术

目前，公知的感应式磁声成像方法是一种无创超声及其电阻抗成像，其成像原理是将作为成像目标的生物体放置于永磁体产生的静态磁场当中，由外界施加的瞬变磁场在生物体内生成感应电流，感应电流在静态磁场作用下产生时变洛伦兹力，由洛仑兹力引发生物体振动并向外传播。超声传感器采集包含生物体内部信息的声波信号并应用相应算法重建生物体内的电导率分布图像。感应式磁声成像方法还处于实验室研究向实用转变的初级阶段。

公知的，现有的感应式磁声成像系统采用传感器环形扫描的多次激励(Multi-shot)方式，即步进电机每次驱动传感器移动一个角度，电压或者电流驱动线圈产生瞬变磁场并激发磁声信号，再由多次采集的磁声信号重建声源和电导率图像。环形扫描成像方式的理论基础是每次产生的瞬变磁场完全相同，从而对多次激励产生的磁声信号可以认为是源自同一声源。然而，在实验过程中，产生完全相同的瞬变磁场难度很大，并且声源对瞬变磁场的磁感应强度及其关于时间的变化率极度敏感。因此，在不能保证瞬变磁场完全相同的情况下，使用多次激励扫描成像方式重建的图像将与实际存在较大误差。

公知的，现有感应式磁声成像系统主要采用式(1)进行声源重建。

在实际计算时，需要对检测面Σ做空间离散、对声压信号p做时间离散，将式(5)转换成式(6)，

$\▿\·(J\×B)\≈\frac{-1}{2{\mathrm{\πc}}_s^3}\·\underset{i=1}{\overset{N_d}{\mathrm{\Σ}}}\left\{\frac{n\·(r_{di}-r)}{|r_{di}-r{|}^2}grad\left(grad\right(p_0\left(r_{di},\left|\right|r_{di}-r\left|\right|/c_s\right)\left)\right)\right\}---\left(2\right)$

其中，grad(grad(F))表示对F的二阶微分，Nd表示检测面上均匀分布的超声传感器数量。但是，这种重建方法存在一些技术缺陷：(a)随着点声源数量的增多，重建精度急剧下降；(b)针对Nd的最小取值缺乏相应的理论指导，即多少个采样点才能满足MAT-MI声源的精确重建。因而，在现有重建算法中，时间和空间域的离散将严重影响声源和电导率图像的重建精度。

公知的，压缩感知是Candès、Romberg、Tao、Donoho等学者提出的一套通过发掘信号的稀疏特性，利用低维空间、低分辨率、欠采样数据的非相关观测来实现高维、高分辨率信号的感知，实现信号无失真重建的算法。压缩感知因其在减少数据采集总量方面的潜在能力而成为信号稀疏采样和重建领域的研究热点，并在医学成像(MRI\CT等)、光学成像、图像处理等诸多领域获得应用。在医学成像领域，如果存在一个观测系统Ф生成关于成像目标f的感知信号Y＝Фf，以及一组用于稀疏表达f的基Ψ，其中，f＝xΨ，x是稀疏序列；同时，Ф和Ψ不相关。那么，应用压缩感知技术可以尽可能少地采集数据，并尽可能多地提取有用信息，从而缩短成像时间，获取较高质量的断层图像。因此，应用压缩感知方法设计感应式磁声成像的图像重建算法，不仅可以减少成像系统中的超声传感器数目，减少数据采集时间，而且可以采用单次激励方式，规避多次激励扫描方式的技术性缺陷。应用压缩感知算法设计的感应式磁声单次激励快速成像方法将推动感应式磁声成像系统的产业化。

发明内容

本发明针对现有技术中成像时间长、操作复杂的缺点，提供了一种单次激励的感应式磁声快速成像方法。

为了解决上述技术问题，本发明通过下述技术方案得以解决。

一种单次激励的感应式磁声快速成像方法，包括压缩感知算法设计模块和声源重建模块；压缩感知算法设计模块用于对曲波正交基的稀疏化、变换点扩散函数的压缩感知矩阵评估、Dantzig选择器的约束优化问题求解和曲波变换的图像重建；声源重建模块针对单次激励产生的磁声信号进行傅立叶变换、应用压缩感知算法的稀疏信号重建、通过曲波变换和傅立叶反变换的声源图像重建。

作为优选，压缩感知算法设计模块包括以下步骤：

(1)应用曲波正交基进行原始数据的稀疏化；

(2)根据变换点扩散函数判定压缩感知矩阵；

(3)通过Dantzig选择器求解约束优化问题，得到稀疏信号；

(4)采用曲波变换将稀疏信号转换成声源信号；

(5)应用自然梯度算法搜寻最优传感器组合，确定传感器数目m和相应位置，以及观测矩阵K。

作为优选，步骤(5)中，确定传感器数目m包括传感器数目初始值的选取、基于自然梯度法的组合优化。

作为优选，声源重建模块包括以下步骤：

(1)采集单次激励产生的时域感应式磁声信号，应用傅立叶变换将时域信号转换成频域信号；

(2)提取每一频率的超声信号组成矩阵，根据成像网格中心点、步骤(3)设计的传感器数目和位置、波数确定观测矩阵；

(3)应用压缩感知算法，求解约束优化问题，计算得到最稀疏解；

(4)将稀疏信号转换成频域的声源信号，应用傅立叶反变换将频率信号转换成时域信号序列；

(5)找到电磁激励时刻对应的时域信号，提取所有网格中心点的时域信号组建成一个矩阵，根据坐标提取相应数据生成声源的断层图像。

作为优选，步骤(1)中，采集单次激励产生的超声信号可重建声源的断层图像。

本发明由于采用了以上技术方案，具有显著的技术效果：通过采用曲波变换、变换点扩散函数、Dantzig选择器和感应式磁声成像波动方程设计压缩感知算法，通过单次激励采集的超声信号即可重建声源图像；通过自然梯度算法对传感器数目和位置进行寻优，确定最少传感器数目和相应位置，从而减少成像所需的传感器数目；通过频域信号的压缩感知重建，既能提高抗干扰能力，又能提高重建精度。

与现有感应式磁声成像技术相比，本发明方法具有显著的优势：

(1)通过采集单次激励产生的磁声信号进行图像重建，避免了现有的多次激励扫描成像方法对非同源的声源进行重建的技术性缺陷和风险，同时减少成像所需时间；

(2)通过压缩感知算法的稀疏采样，减少超声传感器数目，降低成像系统的成本，提高系统的性价比和紧凑度；

(3)通过频域信号的压缩感知重建，避免现有重建算法在空间和时间域的离散化导致的重建误差，提高重建图像的精度、抗干扰能力和稳定性；

(4)本发明方法将系统设计的复杂性限定在前期的压缩感知框架和算法构建，一旦算法确定，则只需要对时域信号进行傅立叶变换并应用压缩感知算法进行稀疏重建，最终重建声源图像。因此，本发明方法在设计阶段较为复杂，但是在应用阶段，操作相对简单，实现速度快。

综合上述四个方面的优势，本发明方法具有广泛的应用前景和潜在的商业价值。

附图说明

图1是本发明一种单次激励的感应式磁声快速成像方法的运行示意图；

图2是本发明一种单次激励的感应式磁声快速成像方法的原理示意图。

具体实施方式

下面结合附图与实施例对本发明作进一步详细描述。

如图1至图2所示，一种单次激励的感应式磁声快速成像方法，包括压缩感知算法设计模块和声源重建模块；压缩感知算法设计模块用于对曲波正交基的稀疏化、变换点扩散函数的压缩感知矩阵评估、Dantzig选择器的约束优化问题求解和曲波变换的图像重建；声源重建模块针对单次激励产生的磁声信号进行傅立叶变换、应用压缩感知算法的稀疏信号重建、通过曲波变换和傅立叶反变换的声源图像重建。

压缩感知算法设计模块包括以下步骤：

(1)应用曲波正交基进行原始数据的稀疏化；

(2)根据变换点扩散函数判定压缩感知矩阵；

(3)通过Dantzig选择器求解约束优化问题，得到稀疏信号；

(4)采用曲波变换将稀疏信号转换成声源信号；

(5)应用自然梯度算法搜寻最优传感器组合，确定传感器数目m和相应位置，以及观测矩阵K。

步骤(5)中，确定传感器数目m包括传感器数目初始值的选取、基于自然梯度法的组合优化。

声源重建模块包括以下步骤：

(1)采集单次激励产生的时域感应式磁声信号，应用傅立叶变换将时域信号转换成频域信号；

(2)提取每一频率的超声信号组成矩阵，根据成像网格中心点、步骤(3)设计的传感器数目和位置、波数确定观测矩阵；

(3)应用压缩感知算法，求解约束优化问题，计算得到最稀疏解；

(4)将稀疏信号转换成频域的声源信号，应用傅立叶反变换将频率信号转换成时域信号序列；

(5)找到电磁激励时刻对应的时域信号，提取所有网格中心点的时域信号组建成一个矩阵，根据坐标提取相应数据生成声源的断层图像。

步骤(1)中，采集单次激励产生的超声信号可重建声源的断层图像。

本发明的设计步骤如下：

步骤1、构建压缩感知框架。感应式磁声效应的频域声波方程如式(3)所示。

其中，c为声速，k为波速。

在医学成像中，式(3)可以用线性方程组表示，如式(4)所示。

Y＝KX+e (4)

其中，X为关于频率kf的MAT-MI声源强度(单位体积内边界声源密度乘以面积与内部声源乘以体积之后的总和)，Y为关于频率kf的声压测量信号，r_i,r_j分别为测量点和声源坐标，e为白噪声。

选择一组曲波正交基Φ，使得X＝ΦS且S是一个稀疏表达的向量，于是，Y＝KΦS+e。通过变换点扩散函数评估KΦ的不相关性，并判断KΦ是否为压缩感知矩阵。当i≠j时，TPSF(i,j)＜＜1且其值的分布类似于噪声，则可判定该矩阵为压缩感知矩阵。

$TPSF(i,j)=\frac{({\mathrm{K\Φ}}_i,{\mathrm{K\Φ}}_j)}{\left|\right|{\mathrm{K\Φ}}_i\left|\right|\left|\right|{\mathrm{K\Φ}}_j\left|\right|}---\left(5\right)$

其中，i,j分别表示Φ的第i,j列。

步骤2、确定压缩感知算法。如果KΦ是一个压缩感知矩阵，那么可以通过求解式(6)来重建稀疏向量S。

$\begin{array}{cc}min\left|\right|S|{|}_{L_1}& s.t.\left|\right|Y-K\mathrm{\&Phi;}S|{|}_{L_2}<\mathrm{\&epsiv;}\end{array}---\left(6\right)$

其中，表示误差。

应用Dantzig选择器算法(Candès,Tao,2007)，通过将形如式(6)的约束优化问题转换成线性规划的凸优化问题，可以用迭代收缩算法求解得到最稀疏解

步骤3、确定最少传感器数目及观测矩阵K。

经过步骤1的正交变换得到稀疏信号S，用M表示稀疏信号S的非零值个数，选取传感器数目m＝3M～5M，应用步骤2计算测量数目相关的稀疏信号重建精度，应用自然梯度算法搜寻确保重建精度的感应式磁声单次激励成像系统的最少测量数目和相应位置。根据传感器的数目和位置，确定观测矩阵K。

步骤4应用压缩感知算法处理超声信号。针对传感器采集的时域超声信号，应用傅立叶变换将时域信号转换成频域信号，提取每一频率的超声信号组成矩阵Y_f，根据成像网格中心点、步骤3设计的传感器数目和位置、波数k_f＝f/c确定观测矩阵K_f。应用步骤1和2设计的压缩感知算法，求解约束优化问题计算得到最稀疏解

步骤5重建声源的断层图像。通过将稀疏信号转换成频域的声源信号。收集成像网格中心点对应的各频率声源信号其中f_τ表示第τ个频率，r_i表示第i个网格的中心点位置。将这些频域信号按照频率值从小到大组成一个序列，应用傅立叶反变换将频率信号转换成时域信号序列找到电磁激励时刻对应的时域信号提取所有网格中心点的时域信号组建成一个矩阵此矩阵即是感应式磁声效应的三维声源矩阵，只需根据坐标提取相应数据即可生成断层图像。

本发明工作过程如下步骤：

步骤(1)构建压缩感知框架，如图2所示。将成像区域的空间离散成N个网格，成像系统变成形如Y＝KX+e的线性方程组。在此基础上，设计压缩感知框架的过程如下：

a.根据生理解剖和MRI、CT等成像数据，通过数值仿真计算感应式磁声的声源分布，应用曲波正交基将原始声源数据转换成尺度、方向和位置三个维度的数据，选取恰当的参数阈值，将转换后的数据转变成稀疏数据；

b.随机选取N个传感器的位置，构建矩阵KΦ，应用变换点扩散函数判定KΦ是否为压缩感知矩阵，通过改变传感器位置和曲波正交基的参数来改进KΦ矩阵，直到该矩阵满足压缩感知矩阵条件为止。

步骤(2)确定压缩感知算法，如图2所示。一旦确定KΦ是压缩感知矩阵，即可根据频域声源分布计算传感器采集的频域超声信号，得到矩阵Y。构建优化目标应用Dantzig选择器将约束优化问题转换成线性规划问题，如(7)式所示：

$\begin{array}{cc}\underset{U,S}{min}{1}^{T}U& s.t.\left\{\begin{array}{c}-U\&le;X\&le;U\\ -\mathrm{\&lambda;}\mathrm{\&sigma;}1\&le;{\left(K\mathrm{\&Phi;}\right)}^T(Y-K\mathrm{\&Phi;}S)\&le;\mathrm{\&lambda;}\mathrm{\&sigma;}1\end{array}\right\}\end{array}---\left(7\right)$

通过迭代收缩算法求解得到稀疏解并进一步

步骤(3)确定最少传感器数目及观测矩阵K，如图2所示。

经过步骤(1)的正交变换得到稀疏信号S，用M表示稀疏信号S的非零值个数，选取传感器数目m＝3M～5M，在矩阵K中随机选取m行，并在矩阵Y中选取对应的m行数据，约束优化问题成为应用步骤(2)的压缩感知算法重建图像并评估重建精度，再根据重建精度应用自然梯度法对m进行组合优化。根据m的选取，确定传感器的数目和位置，从而确定观测矩阵K。

上述三个步骤完成压缩感知算法的设计，包括成像网格中心点位置、传感器数目和相应位置、Dantzig选择器优化算法。

步骤(4)应用压缩感知算法处理超声信号，如图1的声源重建模块。根据步骤(3)确定的传感器数目和位置，采集单次激励产生的时域感应式磁声信号，应用傅立叶变换将时域信号转换成频域信号，提取每一频率的超声信号组成矩阵Y_f。根据成像网格中心点位置、传感器数目和位置、以及波数k_f＝f/c确定观测矩阵K_f。应用步骤(2)设计的压缩感知算法，求解约束优化问题计算得到最稀疏解

步骤(5)重建声源的断层图像，如图1的声源重建模块。通过将稀疏信号转换成频域的声源信号。收集成像网格中心点对应的各频率声源信号其中f_τ表示第τ个频率，r_i表示第i个网格的中心点位置。将这些频域信号按照频率值从小到大组成一个序列，应用傅立叶反变换将频率信号转换成时域信号序列找到电磁激励时刻对应的时域信号提取所有网格中心点的时域信号组建成一个矩阵此矩阵即是感应式磁声效应的三维声源矩阵，只需根据坐标提取相应数据即可生成断层图像。

综上所述，本发明提供一种通过采集单次激励产生的感应式磁声信号进行声源重建的方法。该方法通过稀疏采样，选择最优传感器数目，减少数据采集量和处理时间，同时利用压缩感知算法的特性，实现单次激励模式的数据重建，提高数据的可靠性和稳定性，简化成像系统，具有很高的产业化前景。

总之，以上所述仅为本发明的较佳实施例，凡依本发明申请专利范围所作的均等变化与修饰，皆应属本发明专利的涵盖范围。

Claims (5)

1.一种单次激励的感应式磁声快速成像方法，其特征在于：包括压缩感知算法设计模块和声源重建模块；压缩感知算法设计模块用于对曲波正交基的稀疏化、变换点扩散函数的压缩感知矩阵评估、Dantzig选择器的约束优化问题求解和曲波变换的图像重建；声源重建模块针对单次激励产生的磁声信号进行傅立叶变换、应用压缩感知算法的稀疏信号重建、通过曲波变换和傅立叶反变换的声源图像重建。

2.根据权利要求1所述的一种单次激励的感应式磁声快速成像方法，其特征在于：压缩感知算法设计模块包括以下步骤：

(1)应用曲波正交基进行原始数据的稀疏化；

(2)根据变换点扩散函数判定压缩感知矩阵；

(3)通过Dantzig选择器求解约束优化问题，得到稀疏信号；

(4)采用曲波变换将稀疏信号转换成声源信号；

(5)应用自然梯度算法搜寻最优传感器组合，确定传感器数目m和相应位置，以及观测矩阵K。

3.根据权利要求2所述的一种单次激励的感应式磁声快速成像方法，其特征在于：步骤(5)中，确定传感器数目m包括传感器数目初始值的选取、基于自然梯度法的组合优化。

4.根据权利要求1所述的一种单次激励的感应式磁声快速成像方法，其特征在于：声源重建模块包括以下步骤：

(1)采集单次激励产生的时域感应式磁声信号，应用傅立叶变换将时域信号转换成频域信号；

(2)提取每一频率的超声信号组成矩阵Y_f，根据成像网格中心点、步骤(3)设计的传感器数目和位置、波数k_f＝f/c确定观测矩阵K_f；

(3)应用压缩感知算法，求解约束优化问题，计算得到最稀疏解

(4)将稀疏信号转换成频域的声源信号，应用傅立叶反变换将频率信号转换成时域信号序列；

(5)找到电磁激励时刻对应的时域信号，提取所有网格中心点的时域信号组建成一个矩阵，根据坐标提取相应数据生成声源的断层图像。

5.根据权利要求4所述的一种单次激励的感应式磁声快速成像方法，其特征在于：步骤(1)中，采集单次激励产生的超声信号可重建声源的断层图像。