一种基于稀疏表示的深度学习磁共振波谱重建方法与流程

本发明涉及磁共振波谱重建方法，尤其是涉及一种基于稀疏表示的深度学习磁共振波谱重建方法。

背景技术：

磁共振波谱能够提供分子结构的原子水平的信息，是医学、化学和生命科学领域的重要分析工具。磁共振实验中，信号的采样时间随分辨率和采样维度的增大而增大。非均匀采样技术通过获取部分数据的方法，被广泛地用来加速实验数据的采集，但它需要先进的波谱重建方法来得到完整的波谱。

在波谱重建中，一些研究人员利用磁共振时域和频域信号的特性来重建波谱。其中效果较好的一种利用了磁共振频域信号的稀疏特性。比如，qu等(xiaoboqu,xuecao,diguo,zhongchen,"compressedsensingforsparsemagneticresonancespectroscopy,"ininternationalsocietyformagneticresonanceinmedicine19thscientificmeeting,pp.3371,2010.；xiaoboqu,diguo,xuecao,shuhuicai,zhongchen,"reconstructionofself-sparse2dnmrspectrafromundersampleddatainindirectdimension,"sensors,vol.11,pp.8888-8909,2011.)指出磁共振波谱具有自稀疏特性，并从稀疏性和相干性的角度分析，建议通过最小化波谱的p范数(0<p≤1)来重建磁共振波谱。之后，利用同样的自稀疏特性，研究人员将压缩感知成功应用于蛋白质的二维甚至更高维磁共振波谱的重建中(krzysztofkazimierczukandvladislavorekhov,"acceleratednmrspectroscopybyusingcompressedsensing,"angewandtechemieinternationaledition,vol.50,pp.5556-5559,2011.；danielholland,markbostock,lynngladden,danielnietlispach,"fastmultidimensionalnmrspectroscopyusingcompressedsensing,"angewandtechemieinternationaledition,vol.50,pp.6548-6551,2011.；svenhyberts,alexandermilbradt,andreaswagner,haribabuarthanari,gerhardwagner,"applicationofiterativesoftthresholdingforfastreconstructionofnmrdatanon-uniformlysampledwithmultidimensionalpoissongapscheduling,"journalofbiomolecularnmr,vol.52,pp.315-327,2012.)。此外，qu等(xiaoboqu,maximmayzel,jian-fengcai,zhongchen,vladislavorekhov,"acceleratednmrspectroscopywithlow-rankreconstruction,"angewandtechemieinternationaledition,vol.54,pp.852-854,2015.)还针对宽谱峰重建效果不理想等问题，提出基于低秩汉克尔矩阵对欠采样的波谱进行高质量的重建方法。但是，以上方法受限于最优化算法的迭代计算过程，波谱重建时间较长。

随着人工智能的发展，深度学习为解决传统最优化重建方法存在的问题提供了一种新的思路。qu等(xiaoboqu,yihuihuang,hengfalu,tianyuqiu,diguo,tatianaagback,vladislavorekhov,zhongchen,"acceleratednuclearmagneticresonancespectroscopywithdeeplearning,"angewandtechemieinternationaledition,doi:10.1002/anie.201908162,2019.)提出了在磁共振频域信号上利用卷积神经网络和合成数据训练的深度学习网络来进行欠采样磁共振波谱的重建的新方法，明显降低了波谱重建时间。chen和wang等(dichengchen,ziwang,diguo,vladislavorekhov,xiaoboqu,"reviewandprospect:deeplearninginnuclearmagneticresonancespectroscopy,"chemistry-aeuropeanjournal,doi:10.1002/chem.202000246,2020.)综述了深度学习在磁共振波谱领域的历史和前沿进展。

但是，磁共振波谱重建中，现有方法的重建质量仍有待提高，尚无利用频域信号的稀疏特性，建立深度学习神经网络来实现快速且高质量的磁共振波谱重建的方法。

技术实现要素：

本发明目的在于提供重建速度快、重建质量高和可解释性强的一种基于稀疏表示的深度学习磁共振波谱重建方法。

本发明包括以下步骤：

1)利用磁共振波谱的时域信号的指数函数特性，仿真产生全采样的时域信号；

2)对时域信号进行欠采样，建立包含全采样时域信号对应的波谱、欠采样时域信号和对应欠采样模板的训练集；

3)设计基于稀疏表示的深度学习网络模型、网络的反馈功能及损失函数；

4)利用步骤2)获得的训练集，求解基于稀疏表示的深度学习网络的最优参数；

5)将待重建的欠采样的磁共振时域信号输入网络重建磁共振波谱。

在步骤1)中，所述利用磁共振波谱的时域信号的指数函数特性，仿真产生全采样的时域信号的具体方法为：根据指数函数生成全采样的磁共振波谱的时域信号，该时域信号可以被表示为有限个指数函数的叠加。

在步骤2)中，所述对时域信号进行欠采样，建立包含全采样时域信号对应的波谱、欠采样时域信号和对应欠采样模板的训练集的具体方法为：定义为时域中的欠采样且在未采样点填零的算子，其对应的欠采样模板为u，对全采样时域信号s进行欠采样操作，得到欠采样时域信号y：

定义yn与un中的下标n表示欠采样时域信号y与欠采样模板u的第n个元素，算子u的具体操作为：若un＝1，全采样时域信号s的第n个元素保留到欠采样时域信号y，表示全采样时域信号s的第n个元素有采样到；若un＝0，则对应的yn＝0，表示全采样时域信号s的第n个元素没有采样到；这一过程表示为先对全采样时域信号s进行傅里叶变换得到对应的全采样波谱x，再将欠采样时域信号y，全采样波谱x和欠采样模板u共同组成训练集其中欠采样时域信号y和全采样波谱x分别作为神经网络的输入数据和标签。

在步骤3)中，所述基于稀疏表示的深度学习网络模型以迭代块为核心，通过叠加若干个迭代块作为整体的网络结构，每个迭代块包含三个子块；单个迭代块的网络结构如下：

a)子块d用于更新网络中间量d，它由ld层卷积神经网络密集连接构成，卷积核大小均为id×id；第一层输入为和欠采样时域信号y的并集，其中f表示傅里叶变换，上标h表示复共轭转置，表示上一个迭代块的最终输出；第2～ld-1层之间以线性整流函数(rectifiedlinearunit,relu)连接，且每层输入均为之前网络层输出的并集；最后一层输出子块d的计算结果d；完整的子块d用如下非线性映射函数表示：

其中，θ^d表示该迭代块中子块d的内部参数；fd(·)表示子块d所训练的非线性映射；

b)子块r用于更新网络中间量r；它由lr层卷积神经网络密集连接构成，卷积核大小均为ir×ir；第一层输入为上一个迭代块的最终输出和γfu^td的并集，其中上标t表示转置，d表示子块d的计算结果，γ表示步长；第2～lr-1层之间以线性整流函数(rectifiedlinearunit,relu)连接，且每层输入均为之前网络层输出的并集；最后一层输出子块r的计算结果r；完整的子块r用如下非线性映射函数表示：

其中，θ^r表示该迭代块中子块r的内部参数；fr(·)表示子块r所训练的非线性映射；

c)子块g用于根据输入数据特征自动选择阈值θ，并通过软阈值算子更新当前迭代块的最终输出首先，对输入的子块r的计算结果r的所有元素取绝对值，然后经过全局均值池化(globalaveragepooling,gap)得到中间量a；同时另一条传输路径中，会将全局均值池化之后的结果，输入到一个小型的全连接网络之中，这个全连接网络以非线性激活函数sigmoid作为最后一步，得到一个尺度数值β∈(0,1)，子块g得到的阈值为θ＝a×β；

定义软阈值算子soft为：

soft(·,θ)＝sgn(·)×max(0,|·|-θ)，(5)

其中，sgn(·)为符号函数；子块g的最终输出结果为当前迭代块的最终输出完整的子块t用如下非线性映射函数表示：

其中，θ^t表示该迭代块中子块t的内部参数；ft(·)表示子块t所训练的非线性映射；

综上，将以上三个子块d、r和g级联，单个迭代块整体可以用如下非线性映射函数组表示：

其中，表示该迭代块中所有子块内部参数的集合；f(·)表示子块的级联；表示所训练的从到的非线性映射，是各子块非线性映射ft(r|θ^t)的组合；

设计的基于稀疏表示的深度学习网络模型可整体表示为：

其中，θ表示整体网络内部参数的集合；foverall(y|θ)表示所训练的从欠采样时域信号y到网络最终输出值的所有迭代块级联的非线性映射；

所述网络的反馈功能是网络求解目标值的重要过程，在网络模型构建过程中，通过将网络的输出值与全采样波谱x进行比较并反馈梯度来更新迭代模块的参数，使网络输出值更逼近全采样时域信号对应的全采样波谱；

所述损失函数可以定义为：

其中，θ表示整体网络内部参数的集合；||·||2表示向量的二范数；k表示第k个迭代块，k＝1,2,...,k，k表示迭代块的总数；σ表示求和运算。

在步骤4)中，所述求解基于稀疏表示的深度学习网络的最优参数可采用深度学习中表现较好的adam优化器，利用步骤2)中生成的训练集进行网络训练，通过最小化步骤3)中的损失函数来得到最优目标参数集合

在步骤5)中，所述将待重建的欠采样的磁共振时域信号输入网络重建磁共振波谱是将待重建的欠采样的磁共振信号作为已训练好的网络的输入，依照训练好的相对最优参数进行网络的前向传播，即可得到重建后的磁共振波谱

本发明提出了一种基于对磁共振波谱进行稀疏约束的深度学习磁共振波谱重建方法。本方法首先利用磁共振波谱的时域信号的指数函数特性，仿真产生完整的时域信号；接着对时域信号进行欠采样，训练集由欠采样的时域信号、对应的欠采样模板和全采样的时域信号对应的波谱组成，前二者作为网络的输入，后者作为标签；然后依照稀疏表示重建欠采样数据的形式设计深度学习神经网络结构，用上述的训练集求解网络最优参数构成重建模型；最后将欠采样时域数据输入已训练的网络得到重建后的磁共振波谱。

本发明结合了深度学习和传统迭代方法的优势，既延续了前者优异的重建时间表现，又基于后者具有相对可靠的理论支撑，可进行快速且高质量的磁共振波谱重建，且重建结果优于前沿的密集连接卷积神经网络(xiaoboqu,yihuihuang,hengfalu,tianyuqiu,diguo,tatianaagback,vladislavorekhov,zhongchen,"acceleratednuclearmagneticresonancespectroscopywithdeeplearning,"angewandtechemieinternationaledition,doi:10.1002/anie.201908162,2019.)。本发明通过约束磁共振频域信号的稀疏性，以传统最优化方法为指导设计深度神经网络，具有重建速度快、重建质量高和网络可解释性强的特点。

附图说明

图1是实施例中采用的15％采样率的欠采样模板。

图2是基于稀疏表示的深度学习网络结构。其中，(a)为整体网络结构与第k+1个迭代块结构的示意图；(b)为自动阈值选择小型全连接网络结构示意图。

图3是胞质cd79b蛋白质的全采样¹h-¹⁵n异核单量子关系(hsqc)标签波谱和15％采样率下的重建波谱。其中，(a)是全采样标签波谱，(b)是本发明的重建波谱，(c)是前沿的密集连接卷积神经网络的重建波谱(xiaoboqu,yihuihuang,hengfalu,tianyuqiu,diguo,tatianaagback,vladislavorekhov,zhongchen,"acceleratednuclearmagneticresonancespectroscopywithdeeplearning,"angewandtechemieinternationaledition,doi:10.1002/anie.201908162,2019.)。ppm是化学位移的单位，在数值上为10^-6。

具体实施方式

以下实施例将结合附图对本发明作进一步说明。本发明实施例用指数函数生成全采样的时域信号对应的波谱作为训练集标签，欠采样的时域信号、对应的欠采样模板作为训练集输入，通过若干次迭代训练得到最优网络参数，最后将需要重建的欠采样数据输入网络得到重建后的磁共振波谱。

以下给出具体实施例。

本发明实施例包括以下步骤：

第一步：利用指数函数生成磁共振波谱的时域信号

本实施例总计生成40000个自由感应衰减信号，根据指数函数生成磁共振的全采样时域信号它可以被表示为有限个指数函数的叠加：

其中，表示长度为n的复数向量；sn的下标n表示信号向量s的第n个元素；j表示第j个指数函数，j＝1,2,...,j，j表示指数函数的个数(也即谱峰个数)；aj、fj、τj和φj分别表示第j个指数函数的幅值、中心频率、衰减因子和相位；δt表示两次采样的时间间隔。实施例中，n＝256，谱峰个数为1～10，对固定的谱峰个数将生成4000个自由感应衰减信号，幅值aj的取值范围[0.05,1]，频率fj的取值范围[0.01,0.99]，衰减因子τj的取值范围[10.0,179.2]，相位值φj固定为0，两次采样的时间间隔δt＝1。

第二步：建立包含全采样时域信号对应的波谱、欠采样时域信号和对应的欠采样模板的训练集

为时域中的欠采样并在未采样点填零的算子，其对应的欠采样模板为u。对全采样时域信号s进行欠采样操作可以得到欠采样时域信号y：

其中，yn与un中的下标n表示欠采样时域信号y与欠采样模板u的第n个元素，算子的具体操作为：若un＝1，s的第n个元素保留到y，表示s的第n个元素有采样到；若un＝0，则对应的yn＝0，表示s的第n个元素没有采样到。这一过程可以表示为先对全采样时域信号进行傅里叶变换得到对应的全采样波谱x，再将欠采样时域信号y，全采样波谱x和欠采样模板u共同组成训练集其中欠采样时域信号y和全采样波谱x分别作为神经网络的输入数据和标签。

欠采样模板如图1所示。在图1中，欠采样模板u中白色的是采样点，表示该位置对应的数据被采样；黑色表示没有采样到的点，该位置对应的数据丢失。实施例中采样率为15％。

第三步：设计基于稀疏表示的深度学习网络模型、网络的反馈功能及损失函数

深度学习模型以迭代块为核心，通过叠加若干个迭代块作为整体的网络结构，每个迭代块包含三个子块。以第k+1个迭代块为例，网络结构的解释如下：

a)子块d用于更新网络中间量d。它由6层卷积神经网络密集连接构成，卷积核大小均为3×3。第1层输入为uf^hxk和欠采样时域信号y的并集，其中f表示傅里叶变换，上标h表示复共轭转置，xk表示上一个迭代块(即第k个迭代块)的最终输出。第2～5层之间以线性整流函数(rectifiedlinearunit,relu)连接，且每层输入均为之前网络层输出的并集。第6层输出子块d的计算结果dk+1。完整的子块d用如下非线性映射函数表示：

其中，表示第k+1个迭代块中子块d的内部参数；fd(·)表示子块d所训练的非线性映射。

b)子块r用于更新网络中间量r。它由6层卷积神经网络密集连接构成，卷积核大小均为3×3。第1层输入为第k个迭代块的输出xk和γfu^tdk+1的并集，其中上标t表示转置，dk+1表示子块d的输出，γ表示步长；第2～5层之间以线性整流函数(rectifiedlinearunit,relu)连接，且每层输入均为之前网络层输出的并集；第6层输出子块r的计算结果rk+1。完整的子块r用如下非线性映射函数表示：

c)子块g用于根据输入数据特征自动选择阈值θ，并通过软阈值算子更新当前迭代块的最终输出xk+1。首先，对输入的子块r的计算结果rk+1的所有元素取绝对值，然后经过全局均值池化(globalaveragepooling,gap)得到中间量ak+1；同时另一条传输路径中，会将全局均值池化之后的结果，输入到一个小型全连接网络之中，这个全连接网络以非线性激活函数sigmoid作为最后一步，得到一个尺度数值βk+1∈(0,1)，在实施例中，这个全连接网络共有两个隐藏层，分别由2个和1个神经元构成。子块g得到的阈值为θk+1＝ak+1×βk+1。

定义软阈值算子soft为：

soft(·,θ)＝sgn(·)×max(0,|·|-θ)，(5)

其中，sgn(·)为符号函数。子块g的最终输出结果为当前迭代块的最终输出xk+1。完整的子块t用如下非线性映射函数表示：

综上，将以上3个子块d、r和g级联，第k+1个迭代块整体可以用如下非线性映射函数组表示：

其中，θk+1表示第k+1个迭代块中所有子块内部参数的集合；f(·)表示子块的级联；foverall(xk|θk+1)表示所训练的从xk到xk+1的非线性映射，是各子块非线性映射的组合。

设计的基于稀疏表示的深度学习网络具体结构如图2所示。整个网络可表示为：

其中，θ表示整体网络内部参数的集合；foverall(y|θ)表示所训练的从y到网络最终输出值的所有迭代块级联的非线性映射。

所述网络的反馈功能是网络求解目标值的重要过程，在网络模型构建过程中，通过将网络的输出值与全采样波谱x进行比较并反馈梯度来更新迭代模块的参数，使网络输出值更逼近全采样时域信号对应的全采样波谱。

深度学习网络的损失函数定义如下：

其中，θ表示整体网络内部参数的集合；||·||2表示向量的二范数；k表示第k个迭代块，k＝1,2,...,k，k表示迭代块的总数，本方法中k＝10；σ表示求和运算。

第四步：训练基于稀疏表示的深度学习网络的最优参数

采用深度学习中表现较好的adam优化器(diederikkingmaandjimmyba,“adam:amethodforstochasticoptimization,”arxivpreprintarxiv:1412.6980,2014.)，通过最小化步骤三中的损失函数训练可得到最优目标参数集合

第五步：对欠采样的磁共振信号进行重建

将待重建的欠采样的磁共振信号作为已训练好的网络的输入，依照训练好的相对最优参数进行网络的前向传播，得到重建后的磁共振波谱

在实施例中，网络的输入为15％采样率的胞质cd79b蛋白质的¹h-¹⁵n异核单量子关系(hsqc)数据，数据维度为116×256。

在实施例中，胞质cd79b蛋白质的全采样¹h-¹⁵n异核单量子关系(hsqc)标签波谱和15％采样率下的重建波谱如图3。

可以看出，利用基于稀疏表示的深度学习磁共振波谱重建方法可以快速重建出高质量的磁共振波谱，且重建结果优于前沿的密集连接卷积神经网络(xiaoboqu,yihuihuang,hengfalu,tianyuqiu,diguo,tatianaagback,vladislavorekhov,zhongchen,"acceleratednuclearmagneticresonancespectroscopywithdeeplearning,"angewandtechemieinternationaledition,doi:10.1002/anie.201908162,2019.)。

本发明提出了一种利用磁共振波谱的稀疏性的深度学习磁共振波谱重建方法。这种基于对磁共振波谱进行稀疏约束的深度学习磁共振波谱重建方法具有重建速度快、重建质量高和网络可解释性强的特点。