CN112037300B - 基于交替方向乘子网络的光学重建方法及装置 - Google Patents

Abstract

本发明公开了一种基于交替方向乘子网络的光学重建方法，本发明方法基于传统光学成像的前向过程，以不同位置的光源进行仿真，以光学特征参数和解剖结构信息作为先验信息，得到成像目标表面光子能量分布与成像目标内部能量分布的线性关系，并将其转换为凸优化问题，将凸优化问题转化为求解L1正则化的极小化问题，对极小化问题进行迭代求解后，通过参数与层的转换构建交替方向乘子网络，对得到的交替方向乘子网络的结构层进行迭代求解，得到最优解，获得重建目标的三维分布。本发明避免了传统重建方法在重建过程中的不适定性问题，有效提高了光学分子成像的重建结果，在分子影像、重建算法等领域有重要的应用价值。

Description

基于交替方向乘子网络的光学重建方法及装置

技术领域

本发明属于分子影像技术领域，具体涉及一种基于交替方向乘子网络的光学重建方法及装置。

背景技术

光学分子成像(optical molecular tomography,OMT)是在利用光学成像技术获取表面二维光学信号分布的情况下，融合组织结构信息、生物组织特异性光学参数等信息，然后通过有效的光传输数学模型和重建算法，最终实现体内特异性探针三维分布的成像技术。光学分子断层成像根据表面二维光学信号的生成方式可以分为两种：一种是无需激发的生物自发光断层成像，主要包含：生物发光断层成像(Bioluminescence tomography,BLT)和契伦科夫光学成像(Cerenkov Luminescence tomography,CLT)；另一种是需要激发的激发荧光断层成像，主要包含：荧光分子断层成像技术(Fluorescence MolecularTomography,FMT)和X射线发光断层成像(X-Ray Luminescence Computed Tomography,XLCT)。相比于传统的医学成像设备，分子影像学技术有效地利用特异性的分子探针在分子水平上实时监测肿瘤的发生、发展过程，为疾病早期发现、精确诊断、准确治疗奠定技术基础。光学分子成像主要涉及两大问题：前向问题和逆向问题。前向问题主要是通过数学模型来描述光在生物组织中的传输过程，以得到表面节点的能量分布信息。逆向问题则是基于前向问题建立的模型，选择合适的重建算法，求解出荧光目标在生物体内的位置及浓度等一系列信息。但由于表面节点数量有限，远远低于整个成像物体的全部节点数，从而导致在逆向问题中，需求解的未知数远大于方程数，使其具有严重的不适定性。如何精确的重建出荧光目标的真实位置，是光学成像的核心问题。传统的重建算法多依赖于正则化方法以及多种参数的有效选择，这是影响重建精度与速度的一个关键因素。

综上所述，现有技术存在的问题是：光学分子成像中获得的表面节点数据与整体节点数量相差甚大，增大其逆向问题的不适定性，现有的重建方法无法完全解决这一问题。

发明内容

为解决现有技术中存在的不足，本发明提供了一种基于交替方向乘子网络的光学重建方法及装置，解决现有技术中光学分子成像获得的表面节点数据与整体节点数量相差大，导致其逆向问题的不适定性增加的技术问题。

为了解决上述技术问题，本发明采用如下技术方案予以实现：

一种基于交替方向乘子网络的光学重建模型构建方法，包括以下步骤：

步骤1、成像目标内包含重建目标，对成像目标进行CT扫描，获取成像目标的计算机成像数据；

步骤2、使用有限元剖分方法对得到的成像目标的计算机成像数据进行离散化处理，得到多个网格节点和多个四面体，得到成像目标的解剖结构信息；

步骤3、依据扩散光学层析成像算法，由成像目标的解剖结构信息得到成像目标的光学特性参数；

步骤4、根据成像目标的光学特性参数，基于光传输模型和有限元理论，以步骤2得到的网格节点为目标光源，构建成像目标表面光子能量分布与成像目标内部能量分布的线性关系；

步骤5、结合交替方向乘子算法，将步骤4得到的成像目标表面光子能量分布与成像目标内部能量分布的线性关系转换为凸优化问题，将得到的凸优化问题转化为求解L1正则化的极小化问题；

步骤6、对极小化问题进行迭代求解后，通过参数与层的转换构建交替方向乘子网络，对得到的交替方向乘子网络的结构层进行迭代求解，得到最优解，获得重建目标的三维分布；

所述步骤6包括以下子步骤：

步骤S61、结合曾广拉格朗日函数，对极小化问题以式(1)进行迭代求解：

其中，x_l ⁽ⁿ⁾为重建目标的光源分布信息，z_l ⁽ⁿ⁾为非线性变换参数、β_l ⁽ⁿ⁾为乘子更新参数，n表示迭代次数，D_l表示滤波操作的变换矩阵，ρ_l为惩罚项参数，S(·)是带参数λ_l/ρ_l,l∈{1,2,…,L}的非线性收缩函数，L为滤波器的个数，η_l为初始更新率且η_l>0，T为矩阵的转置，A为光传输的系统权重矩阵，Φ为表面能量分布；

步骤S62，将步骤S61求解得到的x_l ⁽ⁿ⁾、z_l ⁽ⁿ⁾、β_l ⁽ⁿ⁾转换为第n层的重建目标的光源分布信息X⁽ⁿ⁾、非线性变换参数Z⁽ⁿ⁾和乘子更新参数β⁽ⁿ⁾；

β⁽ⁿ⁾＝β^(n-1)+η⁽ⁿ⁾×z⁽ⁿ⁾

其中，为第n层的惩罚参数，S_PLF(·)为分段线性函数，p_i为预定义值且-1≤p_i≤1，i为乘子更新参数β⁽ⁿ⁾组成的数组中包含的数值个数，是滤波操作中第l个滤波器在第n层与p_i对应的值，N是分段线性函数取值范围中的最大值，η⁽ⁿ⁾是第n层的更新率且η⁽ⁿ⁾>0；

步骤S63、对步骤S62得到的X⁽ⁿ⁾、Z⁽ⁿ⁾和β⁽ⁿ⁾进行迭代求解，得到最优解，获得重建目标的三维分布。

本发明还包括如下技术特征：

具体的，所述步骤4中构建成像目标表面光子能量分布与成像目标内部能量分布的线性关系通过以下公式求得：

A·x＝Φ (3)

其中，A为光传输的系统权重矩阵，Φ为表面能量分布，x是成像目标表面光子能量分布。

具体的，所述步骤5中求解L1正则化的极小化问题的公式如下：

其中，λ是正则化参数，x_l是x的子集。

一种基于交替方向乘子网络的光学重建装置，包括图像采集模块及重建模块；

所述图像采集模块用于获得待重建的图像；

所述重建模块用于将待重建的图像输入至本发明所述的基于交替方向乘子网络的光学重建模型构建方法获得的基于交替方向乘子网络的光学重建模型中，获得重建目标的三维分布

与现有技术相比，本发明的有益效果是：

1.本发明适用于常见的光学分子成像模式，如：生物发光断层成像、契伦科夫光学成像、荧光分子断层成像、X射线发光断层成像等。

2.本发明采用基于交替方向乘子网络的深度学习方法，可直接根据采集到的成像目标的表面光子分布来重建真实光源的信息，通过训练正则化参数避免因参数选取不当带来的误差。既避免了XLCT重建的病态性问题，同时也提高了重建结果的精度和速度。

3.与现有的端到端的重建网络相比，本发明能够看到网络内部的操作过程，可解释性更强。

附图说明

图1为本发明的基于交替方向乘子网络的光学重建方法流程图；

图2为本发明的交替方向乘子网络的迭代求解机制图；

图3为经交替方向乘子网络进行重建后获得的重建结果示意图，其中，(a)为真实发光目标与重建的发光目标的二维展示；(b)为真实发光目标与重建的发光目标的三维展示。

以下结合附图和实施例对本发明的具体内容作进一步详细解释说明。

具体实施方式

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。

实施例1

本实施例给出一种基于交替方向乘子网络的光学重建模型构建方法，如图1所示，按照以下步骤进行：

步骤1、成像目标内包含重建目标，对成像目标进行CT扫描，获取成像目标的计算机成像数据；

在本实施例中，以简单圆柱仿体作为成像目标，对简单圆柱仿体进行CT扫描，获取简单圆柱仿体的计算机成像数据；

步骤2、使用有限元剖分方法对得到的成像目标的计算机成像数据进行离散化处理，得到多个网格节点和多个四面体，得到成像目标的解剖结构信息；

在本实施例中，使用有限元剖分方法对得到的成像目标的计算机成像数据进行离散化处理，得到9803个节点和54112个四面体，然后得到成像目标的解剖结构信息，如，形状、尺寸、位置等。

步骤3、依据扩散光学层析成像算法，由成像目标的解剖结构信息得到成像目标的光学特性参数；

利用步骤2获得的成像目标的解剖结构信息，结合成像目标的材质及含氧量等信息，通过扩散光学层析算法计算出成像目标的光学特性参数，具体包括吸收系数和散射系数。

步骤4、根据成像目标的光学特性参数，基于光传输模型和有限元理论，以步骤2得到的网格节点为目标光源，构建成像目标表面光子能量分布与成像目标内部能量分布的线性关系；

具体的，从步骤2得到的9803个节点中随机选取1个节点作为光源目标，由光传输模型结合光学特性参数得到的系统权重矩阵，通过光源目标分布信息和系统权重矩阵求得表面光子能量分布信息，得到成像目标表面光子能量分布与成像目标内部能量分布的线性关系：

A·x＝Φ

其中，A为光传输的系统权重矩阵，Φ为表面能量分布，x是成像目标表面光子能量分布。

步骤5、结合交替方向乘子算法，将步骤4得到的成像目标表面光子能量分布与成像目标内部能量分布的线性关系转化为求解L1正则化的极小化问题；求解L1正则化的极小化问题的公式如下：

其中，λl是正则化参数，

步骤6、对极小化问题进行迭代求解后，通过参数与层的转换构建交替方向乘子网络，对得到的交替方向乘子网络的结构层进行迭代求解，得到最优解，获得重建目标的三维分布，具体包括以下子步骤：

步骤S61、结合曾广拉格朗日函数，对极小化问题以下式进行迭代求解：

其中，x_l ⁽ⁿ⁾为重建目标的光源分布信息，z_l ⁽ⁿ⁾为非线性变换参数、β_l ⁽ⁿ⁾为乘子更新参数，n表示迭代次数，D_l表示滤波操作的变换矩阵，ρ_l为惩罚项参数，S(·)是带参数λ_l/ρ_l,l∈{1,2,…,L}的非线性收缩函数，L为滤波器的个数，η_l为初始更新率且η_l>0，T为矩阵的转置，A为光传输的系统权重矩阵，Φ为表面能量分布；在本实施例中，η_l取值为1。

步骤S62，将步骤S61求解得到的x_l ⁽ⁿ⁾、z_l ⁽ⁿ⁾、β_l ⁽ⁿ⁾转换为第n层的重建目标的光源分布信息X⁽ⁿ⁾、非线性变换参数Z⁽ⁿ⁾和乘子更新参数β⁽ⁿ⁾；

在经过重建操作后，重建出成像目标内部的光源分布信息。

在输入z^(n-1)和β^(n-1)的情况下，通过下式得出第n层的重建目标的光源分布信息X⁽ⁿ⁾：

在给定输入β^(n-1)的情况下，通过下式计算得到第n层的非线性变换参数Z⁽ⁿ⁾：

在给定输入z⁽ⁿ⁾和β^(n-1)的情况下，通过下式计算得到第n层的乘子更新参数β⁽ⁿ⁾：

β⁽ⁿ⁾＝β^(n-1)+η⁽ⁿ⁾×z⁽ⁿ⁾

在上述三个公式中，为第n层的惩罚参数，S_PLF(·)为分段线性函数，p_i为预定义值且-1≤p_i≤1，i为乘子更新参数β⁽ⁿ⁾组成的数组中包含的数值个数，是滤波操作中第l个滤波器在第n层与p_i对应的值，N是分段线性函数取值范围中的最大值，η⁽ⁿ⁾是第n层的更新率且η⁽ⁿ⁾>0；

其中，S_PLF(·)由一组控制点决定，用于将乘子更新参数β_l ⁽ⁿ⁾变为分段线性值，本实施例中η⁽ⁿ⁾取值为1。

步骤S63、对步骤S62得到的X⁽ⁿ⁾、Z⁽ⁿ⁾和β⁽ⁿ⁾进行迭代求解，得到最优解，获得重建目标的三维分布。

转换过程如图2所示：图中，X⁽³⁾，Z⁽³⁾，β⁽³⁾分别代表交替乘子网络第三层中的重建层、非线性变换层和乘子更新层。将第二层输出的Z⁽²⁾、X⁽²⁾和β⁽²⁾输入第三层，得到X⁽³⁾；以β⁽²⁾作为非线性变化层的输入，输出为Z⁽³⁾；以β⁽²⁾、X⁽³⁾和Z⁽³⁾作为乘子更新的输入，输出为β⁽³⁾。将求得的X⁽³⁾，Z⁽³⁾，β⁽³⁾继续作为下一层的输入，不断迭代求解至得到最优解，获得重建结果。

仿真结果分析：

经过仿真实验，按照上述发明步骤得到的重建结果如图3所示，其中，(a)显示为成像物体重建结果的二维展示效果，其中黑色圆圈代表真实发光目标，圆圈处的不规则形状为重建得到的发光目标；(b)为对应的重建发光目标的三维展示。

待重建目标的真实中心位置为(0,6,15)mm，重建光源的中心位置为(-0.3，6.38,15.5)mm，位置误差为：从而说明，基于本发明的光学重建算法，能够很好的重建出光源，为使用深度学习网络方法进行光学分子成像提供了新思路。

实施例2

一种基于交替方向乘子网络的光学重建装置，其特征在于，包括图像采集模块及重建模块；

图像采集模块用于获得待重建的图像；

重建模块用于将待重建的图像输入至本发明基于交替方向乘子网络的光学重建模型构建方法获得的基于交替方向乘子网络的光学重建模型中，获得重建目标的三维分布。

需要说明的是，本发明并不局限于上述实施例，在本发明公开的技术方案的基础上，本领域的技术人员根据所公开的技术内容，不需要创造性的劳动就可以对其中的一些技术特征作出一些替换和变形，这些替换和变形均在本发明的保护范围内。

Claims (4)

1.一种基于交替方向乘子网络的光学重建方法，其特征在于，包括以下步骤：

步骤1、成像目标内包含重建目标，对成像目标进行CT扫描，获取成像目标的计算机成像数据；

步骤2、使用有限元剖分方法对得到的成像目标的计算机成像数据进行离散化处理，得到多个网格节点和多个四面体，得到成像目标的解剖结构信息；

步骤3、依据扩散光学层析成像算法，由成像目标的解剖结构信息得到成像目标的光学特性参数；

步骤4、根据成像目标的光学特性参数，基于光传输模型和有限元理论，以步骤2得到的网格节点为目标光源，构建成像目标表面光子能量分布与成像目标内部能量分布的线性关系；

步骤5、结合交替方向乘子算法，将步骤4得到的成像目标表面光子能量分布与成像目标内部能量分布的线性关系转换为凸优化问题，将得到的凸优化问题转化为求解L1正则化的极小化问题；

步骤6、对极小化问题进行迭代求解后，通过参数与层的转换构建交替方向乘子网络，对得到的交替方向乘子网络的结构层进行迭代求解，得到最优解，获得重建目标的三维分布；

所述步骤6包括以下子步骤：

步骤S61、结合曾广拉格朗日函数，对极小化问题以式(1)进行迭代求解：

其中，x_l ⁽ⁿ⁾为重建目标的光源分布信息，z_l ⁽ⁿ⁾为非线性变换参数、β_l ⁽ⁿ⁾为乘子更新参数，n表示迭代次数，D_l表示滤波操作的变换矩阵，ρ_l为惩罚项参数，S(·)是带参数λ_l/ρ_l,l∈{1,2,…,L}的非线性收缩函数，L为滤波器的个数，η_l为初始更新率且η_l>0，T为矩阵的转置，A为光传输的系统权重矩阵，Φ为表面能量分布；

步骤S62，将步骤S61求解得到的x_l ⁽ⁿ⁾、z_l ⁽ⁿ⁾、β_l ⁽ⁿ⁾转换为第n层的重建目标的光源分布信息X⁽ⁿ⁾、非线性变换参数Z⁽ⁿ⁾和乘子更新参数β⁽ⁿ⁾；

其中，为第n层的惩罚参数，S_PLF(·)为分段线性函数，p_i为预定义值且-1≤p_i≤1，i为乘子更新参数β⁽ⁿ⁾组成的数组中包含的数值个数，是滤波操作中第l个滤波器在第n层与p_i对应的值，N是分段线性函数取值范围中的最大值，η⁽ⁿ⁾是第n层的更新率且η⁽ⁿ⁾>0；

步骤S63、对步骤S62得到的X⁽ⁿ⁾、Z⁽ⁿ⁾和β⁽ⁿ⁾进行迭代求解，得到最优解，获得重建目标的三维分布。

2.如权利要求1所述的基于交替方向乘子网络的光学重建方法，其特征在于，所述步骤4中构建成像目标表面光子能量分布与成像目标内部能量分布的线性关系通过以下公式求得：

A·x＝Φ (3)

其中，A为光传输的系统权重矩阵，Φ为表面能量分布，x是成像目标表面光子能量分布。

3.如权利要求1所述的基于交替方向乘子网络的光学重建方法，其特征在于，所述步骤5中求解的L1正则化的极小化问题的公式如下：

其中，λ是正则化参数，z_l为引入的辅助变量。

4.一种基于交替方向乘子网络的光学重建装置，其特征在于，包括图像采集模块及重建模块；

所述图像采集模块用于获得待重建的图像；

所述重建模块用于将待重建的图像输入至权利要求1-3任一项权利要求所述的基于交替方向乘子网络的光学重建模型构建方法获得的基于交替方向乘子网络的光学重建模型中，获得重建目标的三维分布。