Generalized conditional gradient method with adaptive regularization parameters for fluorescence molecular tomography

Yi Chen; Yi Chen; Mengfei Du; Mengfei Du; Jun Zhang; Jun Zhang; Gege Zhang; Gege Zhang; Linzhi Su; Linzhi Su; Linzhi Su; Kang Li; Kang Li; Fengjun Zhao; Huangjian Yi; Lin Wang; Xin Cao; Xin Cao; Xin Cao

doi:10.1364/OE.486339

1. Introduction

Fluorescence molecular imaging (FMI) is a non-invasive imaging technique that detects the surface fluorescence distribution emitted from fluorescent sources within biological tissue by using a high-sensitivity detectors such as charge-coupled device (CCD) or scientific complementary metal–oxide–semiconductor (sCMOS) camera, photomultiplier tubes (PMT), and other photosensitive technologies [1–3]. However, due to the absorption and scattering of light transmission [4,5], FMI can only obtain the photon distribution information on the surface of the object, but three-dimensional (3D) spatial information of the fluorescent target cannot be obtained [6]. Therefore, based on FMI, fluorescence molecular tomography (FMT) has been developed as an imaging modality to obtain the spatial distribution and quantitative analysis of interior fluorescent probes via the reconstruction method, which overcomes the difficulty of quantifying FMI [7,8]. Because of its high specificity, strong sensitivity, and low cost, FMT has been widely applied in preclinical study and diagnosis based on various small animal models [9–11].

Due to the severe scattering of photons propagation through heterogeneous tissues, the complexity of the photon propagation model, and the highly ill-posed inverse problem, FMT reconstruction is still a challenging problem [12–14]. To solve these critical problems, researchers have proposed many different model-based optimization methods, which need appropriate priors or penalties to promote reconstruction and restrict the search space to a specific solution set. One of the effective strategies is to utilize the anatomical information of different biological tissues, which is obtained by computed tomography (CT) or magnetic resonance imaging (MRI), and is used as the prior information of the photon propagation model to construct a more accurate forward model [15,16], so the spatial resolution can be promoted. Furthermore, the ${L_p}$-norm regularization ($p \in [{0,\,2} ]$) of the unknown fluorescent source is also used to constrain the FMT reconstruction [17–19]. Mathematically, ${L_0}$-norm is a sparsest constraint, but the reconstruction algorithm based on ${L_0}$-norm is a combinatorial optimization problem and an NP-hard problem [20,21]. When $p = 2$, the over-smoothness of the ${L_p}$-norm will cause reconstruction artifacts and exacerbates the noise effect [22,23]. When the value of p is small ($p \in ({0,\,1} ]$), the ${L_p}$-norm tend to produce an over-sparse and incomplete boundary of the reconstruction target [24]. To overcome the drawbacks of the classical sparsity regularization, Zou et al. proposed an elastic-net (EN) regularization [25], which combines ${L_1}$ and ${L_2}$ norms with different weights. It not only ensures the sparsity, but also improves the smoothness of the reconstruction source [26]. Liu et al. applied EN regularization to FMT reconstruction for the first time and achieved good reconstruction performance [27]. Later, Wang et al. proposed an adaptive parameter search elastic net (APSEN) method for FMT reconstruction [26], which addresses the issue of parameter selection in reconstruction. However, the coordinate descent method is adopted by APSEN to optimize the objective function, and the optimal solution is solved linearly along the coordinate axis, which is easy to fall into local optimization when solving large-scale sparse models.

It should be noted that the reconstruction performance of the regularization method will be affected by the regularization parameters [28]. The larger the regularization parameters, the greater the influence of the regularization term on the reconstruction results. However, the optimal selection of regularization parameters is usually unknown, which depends on the specific reconstruction problem and the properties of the sparse vector [29]. Recently, some methods for selecting the optimal regularization parameters have been proposed, such as L-curve, U-curve, and cross-validation [30–32], which are of great significance to improve the reconstruction performance of FMT.

In addition to numerical methods, deep learning has also been introduced into FMT reconstruction. Guo et al. proposed an end-to-end 3D depth encoder network [33], which can significantly improve image quality and reduce reconstruction time. Meng et al. proposed a local connection network based on the K-nearest neighbor to improve the morphological reconstruction performance of FMT [34]. Zhang et al. proposed a three-dimensional fusion double-sampling convolution neural network to achieve FMT ultra-high spatial resolution reconstruction [35]. The deep learning method reconstructs the fluorescent source by directly establishing an end-to-end mapping model on a large dataset, which can greatly eliminate the modeling error [36]. Nevertheless, the deep learning method has two common shortcomings. On the one hand, the trained neural network can only be used for specific imaging objects, which means its generalization ability is weak [37–39]. On the other hand, compared with traditional mathematical methods, the interpretability of neural network is poor.

In this work, a generalized conditional gradient method for adaptive regularization parameters (GCGM-ARP) is proposed. The ill-conditioned inverse problem is interpreted as a non-differentiable object function with EN regularization, in which the EN regularization is a combination of ${L_1}$ and ${L_2}$-norms to balance the sparsity and shape recovery of the reconstructed fluorescent source. GCGM-ARP adopts L-curve to effectively generate regularization parameters suitable for different fluorescence distributions. The objection function is solved using the generalized conditional gradient method (GCGM). The per-iteration cost of the GCGM is particularly cheap and globally convergent, making it widely applied in fast sparse approximation. The method consists of two minimization sub-problems: determining the descent direction and determining the step size. The minimization process of determining the descent direction is simplified to a shrinkage process by the iterative soft thresholding algorithm (ISTA) [40], and the step size sub-problem can be efficiently solved by the alternating direction method of multipliers (ADMM) [41].

To assess the performance of the GCGM-ARP method, a series of numerical simulations and in vivo experiments were carried out. Iterative shrinkage with ${L_1}$-norm (IS-${L_1}$) [13], incomplete variables truncated conjugate gradient with ${L_1}$-norm (IVTCG-${L_1}$) [42], Nesterov’s method with EN regularization (N-EN) [27], and APSEN method were used for comparison. The results of the numerical experiments showed that GCGM-ARP performed the lowest location error (LE) and relative intensity error (RIE), as well as the highest Dice coefficient (Dice), as compared to the other methods. These findings indicated the superiority of the GCGM-ARP methods in terms of localization, shape recovery, and dual-source resolution. Moreover, the in vivo experiments further verified the practical applicability of GCGM-ARP method.

The novelty of this paper lies in the proposal of a GCGM-ARP method based on the EN regularization model. Firstly, in comparison to other methods based on the EN regularization model applied to FMT, the principle of the L-curve in the GCGM-ARP method is easier to understand and use in parameter selection. Secondly, the GCGM algorithm avoids falling into local optimization by selecting the gradient descent direction.

The remainder of this paper is structured as follows: the Section 2 introduces the FMT reconstruction model and GCGM-ARP reconstruction algorithm. The Section 3 introduces the evaluation index, the design of experiments and the result of experiments. Finally, the Section 4 gives discussion and conclusion of our work.

2. Methodology

2.1 Photon propagation model

In highly scattering media, such as biological tissues, the photon propagation model can be approximated by the diffusion equation (DE). In steady-state FMT with point excitation sources, the photon propagation can be described by a coupled diffusion equation with Robin boundary condition, which is defined as [43,44]:

(1)$$\left\{ {\begin{array}{l} {\nabla \cdot [{{D_x}(r )\nabla {\varPhi _x}(r )} ]- {\mu_{ax}}(r ){\varPhi _x}(r )}\\ { ={-} \Theta \delta ({r - {r_l}} ),\,({r \in \Omega } )}\\ {\nabla \cdot [{{D_m}(r )\nabla {\varPhi _m}(r )} ]- {\mu_{am}}(r ){\varPhi _m}(r )}\\ { ={-} {\varPhi _x}(r )\eta {\mu_{af}}(r ),({r \in \Omega } )}\\ {2{D_{x,m}}(r )\nabla {\varPhi _{x,m}}(r )+ q{\varPhi _{x,m}}(r )= 0,({r \in \partial \Omega } )} \end{array}} \right.$$

where $\nabla $ is the gradient operator, $\varOmega $ denotes the entire domain of the problem, and r denotes the position of the nodes in the $\varOmega $; ${r_l}$ is the position of a point excitation source with an amplitude of $\Theta $ of the mean free path of photon transmission located below the surface of $\varOmega $; $\eta $ is the quantum efficiency, and q is the optical reflectivity. The subscripts x and m denote the excitation wavelength and emission wavelength, respectively. ${\varPhi _x}(r )$ and ${\varPhi _m}(r )$ denote the excitation and emission photon flux density at position r, respectively. ${D_{x,m}} = 1/3[{{\mu_{ax,am}} + ({1 - g} ){\mu_{sx,sm}}} ]$ denotes the diffusion coefficient, where g denotes the anisotropy parameter, ${\mu _{ax,am}}$ and ${\mu _{sx,sm}}$ denote the absorption coefficient and scattering coefficient respectively. $\eta {\mu _{af}}(r )$ is the unknown fluorophore distribution to be reconstruct, where ${\mu _{af}}$ is absorption coefficient of the fluorescence agent. $\Theta \delta ({r - {r_l}} )$ represents the excitation light which is considered as the point source, where $\delta ({r - {r_l}} )$ represents the Dirac function.

By using the finite element method (FEM) to discretize the photon propagation model, the partial differential Eq. (1) is linearized into the following linear equations:

(2)$$\left\{ {\begin{array}{l} {{M_x}{\varPhi _x} = {S_x}}\\ {{M_m}{\varPhi _m} = {G_x}X} \end{array}} \right.$$

where ${M_x}$ and ${M_m}$ is forward system matrices of the excitation and emission propagation process, respectively. The vector ${S_x}$ represents the distribution of discrete excitation point source. The matrix ${G_x}$ is obtained according to the discretization of fluorescent source, which can be calculated as follows:

(3)$${G_x}({i,j} )= \mathop \smallint \nolimits_\varOmega ^{} {\varPhi _x}(r ){{B}_i}(r ){{B}_i}(r )dr$$

where ${{B}_i}(r )$ and ${{B}_i}(r )$ denote the base functions of the node i and node j, respectively. The vector X represents the reconstructed fluorescent source. ${M_m}$ is symmetric positive matrix [45], so the relationship between ${\varPhi _m}$ and X in Eq. (2) can be expressed as a Matrix-form equation:

(4)$${\varPhi _m} = {M_m}^{ - 1}GX = {A_m}X$$

Because only partial nodes can be detected, the immeasurable nodes are removed from ${\varPhi _m}$ and the corresponding rows are removed from ${A_m}$. The final matrix equation can be expressed as:

(5)$$\varPhi = AX$$

where $\varPhi \in {\mathrm{\mathbb{R}}^{M \times 1}}$ represents the measured photon flux at the boundary of the biological tissue, $A \in {\mathrm{\mathbb{R}}^{M \times N}}$ denotes system matrix, and $X \in {\mathrm{\mathbb{R}}^{N \times 1}}$ denotes the distribution of internal fluorescent source.

In general, the inverse problem of Eq. (5) is ill-posed, so it is impractical to solve X directly. To alleviate the ill-posedness of the inverse problem, the EN regularization was adopted in FMT reconstruction, which is defined as:

(6)$$\mathop {\min }\limits_X E(X )= \frac{1}{2}\left\|AX - \varPhi\right\| _F^2 + \alpha {\left\|X\right\|_1} + \frac{\beta }{2}\left\|X\right\|_2^2$$

where $E(X )$ represents the objective function. $\left\| {\cdot}\right\|_{F}$ represents the Frobenius norm; $\alpha > 0$ and $\beta > 0$ represent the regularization parameters of ${L_1}$-norm and ${L_2}$-norm, respectively. The objective function in Eq. (6) can be equivalently expressed as:

(7)$$\mathop {\min }\limits_X E(X )= \frac{1}{2}\left\|AX - \varPhi\right\| _F^2 + R(X )$$

where $R(X )= \,\alpha {\left\|X\right\|_1} + \frac{\beta }{2}\left\|X\right\|_2^2$ is the EN regularization. Compared with the classical ${L_1}$-norm sparse regularization, there is an additional ${L_2}$-norm in $R(X )$. And EN regularization has the advantages of both ${L_1}$-norm regularization and ${L_2}$-norm regularization [46]. Moreover, statistically, EN regularization is more stable than the classical ${L_1}$-norm regularization [47]. So, it is more suitable for ill-conditioned inverse problem.

Obviously, EN regularization is non-differentiable, so the simple gradient descent minimization method cannot be used to solve the problem. For this reason, here, the GCGM-ARP algorithm was utilized to address the inverse problem.

2.2 Reconstruction based on GCGM-ARP scheme

The GCGM algorithm is proposed by Bredies K et al. to address the following form of minimization problem [48]:

(8)$$\mathop {\min }\limits_X F(X )+ G(X )$$

where F is assumed to be smooth, but G is assumed to be proper, convex, lower semi-continuous and coercive in Hilbert space. More detailed conditions about the functional G that need to be satisfied can be found in [48].

In this work, GCGM is considered for FMT reconstruction. The objective function in Eq. (6) is rewritten as follows:

(9) $$E(X )= F(X )+ G(X )$$

where:

(10)$$\begin{array}{l} {F(X )= \frac{1}{2}\left\|AX - \varPhi\right\| _F^2} - \Psi (X )\\ {G(X )= R(X )+ \Psi (X )}\\ {\Psi (X )= \frac{\lambda }{2}\left\|X\right\|_2^2 - \frac{\beta }{2}}\left\|X\right\|_2^2,\; \; \lambda > 0\; \end{array}$$

There are two reasons for $\Psi (X )= \frac{\lambda }{2}\left\|X\right\|_2^2 - \frac{\beta }{2}\left\|X\right\|_2^2$. Firstly, in this case, $\Psi (X )$ has possess of desirable properties, such as proper, convex, lower semi-continuous and coercive in Hilbert space. Secondly, it can be found that ISTA can be applied to the minimization of the objective function $E(X )$, as the both sub-problems of GCGM involve ${L_1}$-norm regularization.

Thus, Eq. (8) can be expressed as:

(11)$$\mathop {\min }\limits_X \{{(X )= F(X )+ G(X )} \}$$

Next, the GCGM from the literature [48] is showed in Algorithm 1.

It is worth noting that throughout this paper, the operator symbol $\langle\cdot \; ,\; \; \cdot\rangle $ represents the inner product.

In Algorithm 1, there are two sub-problem of minimization that need to be solved, which are the determining direction sub-problem and the determining step size sub-problem. To explain our proposed method more intuitively, these two sub-problems will be discussed respectively.

oe-31-11-18128-i001

2.2.1 Sub-problem of determining the descent direction

Since the Fréchet derivative of F is used in Algorithm 1, its mathematical formula is given:

(12)$$F^{\prime}(X )= {A^T}({AX - \varPhi } )- ({\lambda - \beta } )X.$$

Thus, the descent direction ${Y^k}$ in the Algorithm 1 can be rewritten as follows:

(13)$$\mathop {\min }\limits_Y \langle{A^T}({A{X^k} - \varPhi } )- ({\lambda - \beta } ){X^k},Y\rangle + \alpha {\left\|Y\right\|_1} + \frac{\lambda }{2}\left\|Y\right\|_2^2$$

The Eq. (13) can be solved by explicit calculation. Then, the explicit componentwise expression can be obtained:

(14)$$Y_i^k + \frac{\alpha }{\lambda }sign({Y_i^k} )= {\left( {\frac{1}{\lambda }{A^T}({\varPhi - A{X^k}} )+ \left( {1 - \frac{\beta }{\lambda }} \right){X^k}} \right)_i}$$

where the subscript i represents the $i$-th component. The solution of the Eq. (14) can be expressed by soft threshold function ${S_\gamma }$, which is defined as:

(15)$${S_\gamma }(t )= \left\{ {\begin{array}{l} {t - \frac{\alpha }{\lambda },\; \; t \ge \frac{\alpha }{\lambda }}\\ {0,\; \; - \frac{\alpha }{\lambda } \le t \le \frac{\alpha }{\lambda }}\\ {t + \frac{\alpha }{\lambda },\; \; t \le - \frac{\alpha }{\lambda }} \end{array}} \right.$$

where $\gamma = \frac{\alpha }{\lambda }$ and $t \in \mathrm{\mathbb{R}}$.

It is obvious that determining the descent direction sub-problem can be addressed by ISTA, which can be formulated as:

(16)$${Y^k} = {S_\gamma }\left( {\frac{1}{\lambda }{A^T}({\varPhi - A{X^k}} )+ \left( {1 - \frac{\beta }{\lambda }} \right){X^k}} \right)$$

2.2.2 Sub-problem of determining step size

In this sub-problem, our objective function can be simplified to the following form:

(17)$$\mathop {min }\limits_Z F({{X^k} + WZ} )+ G({{X^k} + WZ} )$$

where ${X^k}$ is a constant, and to avoid confusion, we let ${X^k} = \textrm{X}$. The symbol $W \in {\mathrm{\mathbb{R}}^{N \times N}}$ is a diagonal matrix, which can be defined as follows:

(18)$$W = \left( {\begin{array}{ccccc} {Y_1^k - X_1^k}&{}&{}&{}&{}\\ {}&{Y_2^k - X_2^k}&{}&{}&{}\\ {}&{}& \ddots &{}&{}\\ {}&{}&{}& \ddots &{}\\ {}&{}&{}&{}&{Y_N^k - X_N^k} \end{array}} \right)$$

Combining Eq. (6) and Eq. (10), then $G(X )$ can be formulated as:

(19)$$G(X )= \,\alpha {\left\|X\right\|_1} + \frac{\lambda }{2}\left\|X\right\|_2^2$$

ADMM is a method for solving convex optimization problems. It breaks the objective function of the original problem into smaller sub-problems, then solves each sub-problem in parallel, and finally balances the solution of each sub-problem to obtain the global solution of the original problem [49]. It can be noticed that the objective function of the Eq. (17) is non-differentiable but convex, so ADMM is applied to address the minimization problem in the Eq. (16).

An auxiliary variable is defined as $S = \textrm{X} + WZ$. According to ADMM, Eq. (17) can be reformulated as an optimization problem with equality constraints:

(20)$$\begin{array}{r} {\frac{1}{2}\left\|A({\textrm{X} + WZ} )- \varPhi _F\right\|^2 + \frac{\lambda }{2}\left\|\textrm{X} + WZ\right\|_2^2 - \frac{\beta }{2}\left\|\textrm{X} + WZ\right\|_2^2 + \alpha {\left\|S\right\|_1} + \frac{\beta }{2}\left\|S\right\|_2^2,}\\ {s.t\,\,S = \textrm{X} + WZ\,} \end{array}$$

The augmentation Lagrange function of the Eq. (20) can be expressed as:

(21)$${L_\rho }({S,Z,V} )= \begin{array}{c} {\frac{1}{2}\left\|A({\textrm{X} + WZ} )- \varPhi\right\| _F^2 + \frac{{\lambda - \beta }}{2}\left\|\textrm{X} + WZ\right\|_2^2 - \frac{\beta }{2}\left\|\textrm{X} + WZ\right\|_2^2 + \alpha {\left\|S\right\|_1}}\\ { + \frac{\beta }{2}\left\|S\right\|_2^2 + {V^T}({S - \textrm{X} - WZ} )+ \frac{\rho }{2}\left\|S - \textrm{X} - WZ\right\|_2^2} \end{array}$$

where $\rho > 0$ is a penalty parameter, V is a Lagrange multiplier. The Eq. (21) is solved by the following iterative scheme:

(22)$$\left\{ \begin{array}{l} {S^{k + 1}} = \mathop {\textrm{argmin}}\limits_S {L_\rho }({S,{Z^k},{V^k}} )= \\ \mathop {\textrm{argmin}}\limits_S \left\{ {\begin{array}{c} {\alpha {\left\|S\right\|_1} + \frac{\beta }{2}\left\|S\right\|_2^2 + {V^T}({\,S - \textrm{X} - WZ} )}\\ { + \frac{\rho }{2}\left\|S - \textrm{X} - WZ\right\|_2^2} \end{array}} \right\}\\ {Z^{k + 1}} = \mathop {\textrm{argmin}}\limits_Z {L_\rho }({{S^k},Z,{V^k}} )= \\ \mathop {\textrm{argmin}}\limits_Z \left\{ {\begin{array}{c} {\frac{1}{2}\left\|A({\textrm{X} + WZ} )- \varPhi\right\| _F^2 + \frac{{\lambda - \beta }}{2}\left\|\textrm{X} + WZ\right\|_2^2}\\ {{V^T}({\,S - \textrm{X} - WZ} )+ \frac{\rho }{2}\left\|S - \textrm{X} - WZ\right\|_2^2} \end{array}} \right\}\\ {V^{k + 1}} = {V^k} - \rho ({{S^{k + 1}} - \textrm{X} - W{Z^{k + 1}}} )\end{array} \right.$$

where ${Z^{k + 1}}$ can be solved directly by derivation, but ${S^{k + 1}}$ can not. However, the soft threshold function in Eq. (15) can be adopted to address this optimization problem. First, the explicit calculation result of the ${S^{k + 1}}$ minimization problem is as follows:

(23)$${S^{k + 1}} + \frac{\alpha }{{\beta + \rho }}sign({{S^{k + 1}}} )= \frac{1}{{\beta + \rho }}({\rho \textrm{X} + WZ - V} )$$

By comparing Eq. (13), the solution of Eq. (22) can be obtained, as follows:

(24)$${S^{k + 1}} = {S_{\frac{\alpha }{{\beta + \rho }}}}\left( {\frac{1}{{\beta + \rho }}({\rho \textrm{X} + WZ - V} )} \right).$$

So, the sub-problem of determining the step size is solved.

oe-31-11-18128-i002

2.2.3 Adaptive regularization parameters

As described above, it is difficult to select the optimal regularization parameters in FMT reconstruction. If $\beta $ in Eq. (6) is forced to equal 0, the adaptive EN regularization will reduce into the adaptive lasso. According to the opinion of Zou and Hastie [25], it can be easily proved that, no matter what the value of $\beta $ is, the adaptive EN regularization will be reduced to an adaptive lasso [50]. Thus, regularization parameter $\beta $ will be set empirically, and only regularization parameter $\alpha $ needs to be adjusted adaptively.

For the optimization problem based on Eq. (6), the L-curve will be applied to adaptively adjust the regularization parameter $\alpha $. It can show the tradeoff between the size of the regularization solution and its fitness to the given data, when the regularization parameters change. In general, the L-curve is a log-log plot, which is composed of the solution norm ${X_2}$ and the corresponding residual norm $\left\|AX - \varPhi\right\|_{2}$.

Based on above statement, the complete procedure of GCGM-ARP strategy is summarized in Algorithm 2.

3. Experiments and results

In this section, numerical simulations and in vivo experiment were performed to validate the reconstruction performance of the GCGM-ARP strategy. Furthermore, IS-${L_1}$, IVTCG-${L_1}$, N-EN, APSEN methods were selected for comparison in terms of accuracy, shape recovery, and in vivo practicability. All the programs for the reconstruction algorithm were implemented using MATLAB (2019b) on desktop computer with Intel Core i3-10100 CPU (3.60 GHz) and 8GB RAM.

3.1 Regularization parameters and evaluation index

Suppose ${\eta _1}$ and ${\eta _2}$ are the regularization parameter of ${L_1}$ and ${L_2}$ norms used in the algorithm for the above comparison. In order to ensure the convergence of all algorithms, the maximum number of iterations was set to 1000 and the error thresholds were set to 1e-6 through experience. Moreover, the optimal regularization parameters of GCGM-ARP method were determined by L-curve, and the best regularization parameters of IS-${L_1}$, IVTCG-${L_1}$, and N-EN methods were determined according to related work and experience [13,27,42]. APSEN method can adaptively adjust regularization [26]. The optimal regularization parameters of different algorithms in different experiments have been shown in Table 1.

Table 1. Regularization parameter $\boldsymbol{\alpha }$ for different methods in different experiments

View Table | View all tables in this article

To quantitatively assess the accuracy of FMT reconstruction using different methods, LE, Dice, and RIE were adopted as quantitative indexes in this work.

LE measures the Euclidean distance between the real source center and the reconstructed source center, which is defined as:

(25)$$LE = \left\|{L_{real}} - {L_{rcon}}\right\|_2$$

where ${L_{real}}$ and ${L_{rcon}}$ denote the center coordinates of the real fluorescent source and the reconstructed source region, respectively. The lower the LE, the more accurate the reconstruction location.

Dice is introduced to assess the spatial overlap performance of the real source region and the reconstructed source region, which is defined as:

(26)$$Dice = \; \frac{{2|{{R_{real}} \cap {R_{rcon}}} |}}{{|{{R_{real}}} |+ |{{R_{rcon}}} |}}$$

where ${R_{real}}$ and ${R_{rcon}}$ represent the reconstructed region and the real fluorescent region, respectively. The higher the Dice index, the better the morphological reconstruction. Specifically, the reconstructed source region is determined by the non-zero tetrahedron region based on X.

RIE is used to evaluate the intensity deviation between the real source intensity and the reconstructed source intensity, which is defined as:

(27)$$RIE = \; \frac{{|{{I_{real}} - {I_{rcon}}} |}}{{{I_{real}}}}$$

where ${I_{real}}$ and ${I_{rcon}}$ is the intensity of real fluorescent source and the reconstructed source, respectively. The smaller the RIE, the better the fluorescent intensity recovery of the reconstructed source.

3.2 Experimental setting

3.2.1 Numerical simulation

The heterogeneous model based on a cylinder with radius of 10 mm and height of 30 mm was adopted in the numerical simulations. It consisted of five biological tissues: muscle, heart, lung, liver, and bone, as shown in Fig. 1(a). In this work, all the experiments were carried out under the excitation wavelength of 680 nm and the emission wavelength of 750 nm. The optical parameters of each biological tissue have been listed in detail in Table 2 [51].

Fig. 1. (a) is the 3D view of the cylinder model with single-source. (b) is the mesh for inverse problem. (c) shows the position of 4 excitation sources. (d) is the photon distribution simulated on the surface by the single fluorescence source.

Download Full Size | PDF

Table 2. Optical parameters in numerical simulation

View Table | View all tables in this article

In the forward process, a mesh consisting of 29415 nodes and 168046 tetrahedral elements was used to simulate the photon propagation. The intensity of each source is 1 nw/mm³, and the fluorescence distribution of the surface was simulated using the molecular optical simulation environment (MOSE) [51] based on Monte Carlo method, as shown in Fig. 1(d). In the inverse process, the cylinder model was segmented into a mesh using Comsol Multiphysics software [52]. The mesh included 4626 nodes and 25840 tetrahedral elements, as shown in Fig. 1(b). The excitation sources were located at four positions along the $Z = 15mm$ plane, as illustrated in Fig. 1(c).

We designed single-source and dual-source simulation experiment to evaluate the performance of the GCGM-ARP reconstruction method. In the single-source simulation experiment, one spherical fluorescence source with a radius of 1 $mm$ was placed at (-1, -1, 15) $mm$. In the dual-source simulation experiment, two cylindrical sources with a radius of 1 $mm$ and a height of 2 $mm$, which were placed at the center coordinates of S1(2, 1, 10) $mm$ and S2(2, 1, 18) $mm$, respectively.

As noise is inevitable in FMT, the anti-noise experiment is also designed to evaluate the robustness of our method. Gaussian noise of 5%, 15%, 25% was added to the measurement data based on single-source numerical simulation, and then reconstructed using different methods to observe the robustness and accuracy of these methods under different Gaussian noise ratios.

3.2.2 In vivo experiment

With the guidelines of the Animal Ethics Committee of the Northwestern University of China, in vivo experiment was carried out to further investigate the practical performance of GCGM-ARP. In vivo experimental data set was collected from an adult BALB/c mouse by a dual-modality FMT/CT system, and the chief components of the dual-modality FMT/CT system was exhibited in Fig. 2(a). To minimize the suffering of mouse, all animal experiments were performed under isoflurane gas anesthesia (3% isoflurane-air mixture). The detailed collection process was introduced as follows.

Fig. 2. (a) The chief components of the dual-modality FMT/CT system. (b) Reconstruction workflow of FMT in vivo experiment, including data collection, data processing, and FMT reconstruction.

Download Full Size | PDF

First, a spherical fluorescent bead with a radius of 1$\,mm$ containing Cy5.5 solution was implanted into the abdominal cavity of mouse as a fluorescence target. The fluorescent bead is wrapped in a plastic material, which can be easily detected by micro-CT to locate the real fluorescent region. Six hours later, a 680 nm continuous wave semiconductor laser was used to provide excitation illumination, and the surface fluorescence image with a 120° field of view was collected by a thermoelectric cooled electron multiplying charge coupled device (EMCCD) camera (-80$\mathrm{\circ{C}}$, iXonEM + 888) with an exposure time of 1 s. The emission light was captured and restrained the noise by a 750 ± 10 nm bandpass filter. After the acquirement of fluorescence image, the structure information of mouse needs to be collected by Micro-CT system (tube voltage of 60kVp, x-ray power of 40W). The CT image was processed into 3D volume data, and Amria 5.2 (Amria, Visage Imaging, Australia) was used to segment the main organs, including muscle, lung, heart, stomach, liver, and kidney, and then integrate into the xenogeneic mouse model. The true central position of the fluorescent target is (17.5, 21.5, 9.0)$\,mm$, and the optical parameter for different organs were from the literature [53]. Thus, the structure of the in vivo experimental process can be summarized, as shown in Fig. 2(b).

3.3 Experimental results

3.3.1 Single-source simulation reconstruction

In the spherical single-source simulation experiment, the 3D and slice views of the reconstructed results with different methods were shown in Fig. 3. In the 3D view, different organs were depicted in different colors, and the reconstruction source was depicted in red. In the slice view, the white circle represents the shape and region of the real source, while the red region represents the source of the reconstruction. By observing the slice view, it can be found that although all the methods can approximately reconstruct the fluorescence source region, the overlapping between GCGM-ARP and real source is more, so a better morphological reconstruction is obtained. The quantitative analysis of the reconstruction results of the four methods was listed in Table 3. From these results, it is obvious that our method achieves the smallest LE and RIE, and the largest Dice coefficient. This indicated that GCGM-ARP had possessed of superior positioning ability and shape recovery ability.

Fig. 3. (a), (c), (e), (g) and (i) show the 3D reconstruction results of the IS-${L_1}$, IVTCG-${L_1}$, N-EN, APSEN, and GCGM-ARP, respectively. (b), (d), (f), (h), and (j) show the transverse slice views in the Z = 15 mm of the corresponding four methods.

Download Full Size | PDF

Table 3. Quantitative results in the single-source numerical simulation.

View Table | View all tables in this article

3.3.2 Dual-source simulation reconstruction

In Fig. 4, the 3D view of the reconstruction result in the first row, and the other rows correspond to transverse view, sagittal view, and coronal view, respectively. The enlarged view next to the slice view provides a close view of the distribution of the reconstructed area. The reconstruction sources obtained by IS-${L_1}$ and IVTCG-${L_1}$ were over-sparse, although N-EN and APSEN showed improved performance, the reconstruction result obtained by GCGM-ARP reconstruction was closest to the real sources. Obviously, the dual-source reconstruction performance of GCGM-ARP is better than other methods. The detailed reconstruction results of the four algorithms were shown in Table 4, which further confirms our observations. The LE of S1 and S2 obtained by GCGM-ARP method is much lower than that of other methods, which was 0.432 mm and 0.314 mm, and the Dice is up to 0.839 and 0.834, which demonstrated that GCGM-ARP achieved more accurate dual-source localization and morphological recovery.

Fig. 4. Reconstruction results of different method for dual-source simulation.

Download Full Size | PDF

Table 4. Quantitative results in the dual-source numerical simulation.

View Table | View all tables in this article

3.3.3 Robustness experiment

The influence of different intensity of Gaussian noise on the reconstruction results of the four methods is shown in the Fig. 5. Due to the influence of background fluorescent signal, the reconstruction results of the four methods are different from the single-source numerical simulation experiments without noise. However, under different ratios of Gaussian noise, compared with other methods, GCGM-ARP shows excellent reconstruction performance, with the lowest LE and RIE, and the largest Dice. Consequently, anti-noise experiments showed that GCGM-ARP is the most robust among these methods.

Fig. 5. The quantitative analysis of different methods under different Gaussian noise levels.

Download Full Size | PDF

3.3.4 In vivo experiment reconstruction

The 3D view and slice view were displayed in Fig. 6, which prove the practicality of our methods in vivo. According to Fig. 6, it was showed that there was a deviation between the reconstruction position of IS-${L_1}$ and IVTCG-${L_1}$, and there were many reconstruction artifacts produced by N-EN. While the performance of APSEN improved, it still cannot match the accuracy achieved by GCGM-ARP, which achieved more accurate reconstruction, including the most accurate spatial location, and introduced less reconstruction artifacts. Here, it should be noted that, because of the unknown real source intensity, all the evaluation indicators related to the source intensity (RIE) were unavailable in vivo experiments. As shown in Table 5, the GCGM-ARP method indicated the best accuracy with the least LE and largest Dice similarity. These quantitative results further indicated the superior performance of GCGM-ARP method in obtaining the morphology and localization of fluorescence probe distribution in mouse.

Fig. 6. The 3D view and slice view at the $Z = 9.0\; mm$ plane of reconstruction results obtained by different method for the in vivo experiment.s

Download Full Size | PDF

Table 5. Quantitative results of the in vivo experiment.

View Table | View all tables in this article

4. Discussion and conclusion

In this work, the EN regularization with adaptive parameters via GCGM was proposed to improve the recovery the 3D distribution of the fluorescent source. The regularization parameter was determined by the L-curve to select the optimal regularization parameter according to different data. Compared with the classical ${L_1}$-norm and ${L_2}$-norm, EN regularization offered greater flexibility to balance the sparsity and smoothness of the reconstruction source. Nevertheless, it also increased the computational complexity, so the GCGM-ARP with iteration was employed to enhance the accuracy and reduce the computational complexity. Specifically, GCGM-ARP splits the inverse problem into two sub-problems, namely gradient descent direction and step size, which were addressed by ISTA and ADMM. In addition, GCGM-ARP is simple and easy to implement and has convergence, which can guarantee the stability solution for FMT.

The effectiveness of GCGM-ARP was demonstrated through three groups of numerical simulation experiments and one in vivo experiment. In the numerical simulation experiments, the reconstruction results and quantitative analysis indicated that GCGM-ARP outperformed IS-${L_1}$, IVTCG-${L_1}$, N-EN, and APSEN method in terms of positioning accuracy, shape recovery, and dual-source positioning ability. Moreover, in the robustness test, GCGM-ARP showed high reconstruction accuracy and morphological recovery ability despite the increased ill-posedness of FMT reconstruction. The in vivo experiment further proved the superiority and practicability of the method. It was worth noting that because of the noise in the in vivo experiment and the error caused by organ segmentation, the results of reconstruction were worse than that of using the same method in numerical simulation, but GCGM-ARP still obtained satisfactory results. Overall, these experiments demonstrated the effectiveness of GCGM-ARP in improving reconstruction accuracy, spatial resolution, and dual-source resolution.

Although GCGM-ARP had achieved satisfactory reconstruction results, there are still some limitations. Firstly, the speed of ISTA algorithm used in GCGM-ARP can be further optimized by incorporating acceleration methods such as gradient projection to accelerate convergence. Secondly, the accuracy of optical parameters estimation can significantly affect the reconstruction results. Currently, the optical parameters used in this study were estimated based on the relevant tissues in the relevant literature [53], but using near infrared imaging method to measure the optical properties of tissue and background optical properties can improve the accuracy of reconstruction results. Additionally, in the in vivo experiment, the specific quantification of the effect of low-fluorescence in the non-target region on the captured fluorescence distribution has not been conducted. Further analysis of the effect of low- fluorescence can better simulate the complexity of the environment in vivo and enhance the reliability of in vivo experiments. Finally, The L-curve composed of a series of discrete points limits the accuracy of the optimal regularization parameters, thus affecting the reconstruction results. Therefore, developing a method for selecting more accurate optimal parameters under time constraint will be the focus of our work in the future.

In summary, among the four methods used in the experiments, the reconstruction results of IS-${L_1}$ and IVTCG-${L_1}$ are poor in location and morphology of the reconstructed source. The reconstruction results of N-EN method are show higher accuracy in location and shape recovery compared to the former two methods. The reconstruction performance of APSEN is similar to that of GCGM-ARP, and better than that of N-EN. However, GCGM-ARP achieves lower LE and RIE, as well as a higher Dice, which proves its better performance in localization, morphology and fluorescent intensity reconstruction. This work has the potential to promot the application of FMT in preclinical or clinical biology.

Funding

National Natural Science Foundation of China (61701403, 61806164); Key Research and Development Program of Shaanxi Province (2019GY215, 2021ZDLSF06-04); China Postdoctoral Science Foundation (2018M643719); Graduate Innovation Program of Northwest University (CX2023185).

Disclosures

The authors declare no conflicts of interest.

Data availability

Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon request.

References

1. V. Ntziachristos, “Fluorescence molecular imaging,” Annu. Rev. Biomed. Eng. 8(1), 1–33 (2006). [CrossRef]

2. M. Koch, P. Symvoulidis, and V. Ntziachristos, “Tackling standardization in fluorescence molecular imaging,” Nat. Photonics 12(9), 505–515 (2018). [CrossRef]

3. B. Zhu, J. C. Rasmussen, and E. M. Sevick-Muraca, “A matter of collection and detection for intraoperative and noninvasive near-infrared fluorescence molecular imaging: To see or not to see?” Med. Phys. 41(2), 022105 (2014). [CrossRef]

4. P. Mohajerani and V. Ntziachristos, “An inversion scheme for hybrid fluorescence molecular tomography using a fuzzy inference system,” IEEE Trans. Med. Imaging 35(2), 381–390 (2015). [CrossRef]

5. Z. Hu, C. Fang, B. Li, Z. Zhang, C. Cao, M. Cai, S. Su, X. Sun, X. Shi, and C. Li, “First-in-human liver-tumour surgery guided by multispectral fluorescence imaging in the visible and near-infrared-I/II windows,” Nat. Biomed. Eng. 4(3), 259–271 (2020). [CrossRef]

6. R. Weissleder and M. Nahrendorf, “Advancing biomedical imaging,” Proc. Natl. Acad. Sci. 112(47), 14424–14428 (2015). [CrossRef]

7. A. Ale, V. Ermolayev, E. Herzog, C. Cohrs, M. H. De Angelis, and V. Ntziachristos, “FMT-XCT: in vivo animal studies with hybrid fluorescence molecular tomography–X-ray computed tomography,” Nat. Methods 9(6), 615–620 (2012). [CrossRef]

8. Z. Hu, M. Zhao, Y. Qu, X. Zhang, M. Zhang, M. Liu, H. Guo, Z. Zhang, J. Wang, and W. Yang, “In vivo 3-dimensional radiopharmaceutical-excited fluorescence tomography,” J. Nucl. Med. 58(1), 169–174 (2017). [CrossRef]

9. V. Ntziachristos, “Going deeper than microscopy: the optical imaging frontier in biology,” Nat. Methods 7(8), 603–614 (2010). [CrossRef]

10. X. Shi, Z. Zhang, Z. Zhang, C. Cao, Z. Cheng, Z. Hu, J. Tian, and N. Ji, “Near-infrared window II fluorescence image-guided surgery of high-grade gliomas prolongs the progression-free survival of patients,” IEEE Trans. Biomed. Eng. 69(6), 1889–1900 (2021). [CrossRef]

11. C. Cao, Z. Jin, X. Shi, Z. Zhang, A. Xiao, J. Yang, N. Ji, J. Tian, and Z. Hu, “First clinical investigation of near-infrared window IIa/IIb fluorescence imaging for precise surgical resection of gliomas,” IEEE Trans. Biomed. Eng. 69(8), 2404–2413 (2022). [CrossRef]

12. X. Song, D. Wang, N. Chen, J. Bai, and H. Wang, “Reconstruction for free-space fluorescence tomography using a novel hybrid adaptive finite element algorithm,” Opt. Express 15(26), 18300–18317 (2007). [CrossRef]

13. D. Han, J. Tian, S. Zhu, J. Feng, C. Qin, B. Zhang, and X. Yang, “A fast reconstruction algorithm for fluorescence molecular tomography with sparsity regularization,” Opt. Express 18(8), 8630–8646 (2010). [CrossRef]

14. L. Kong, Y. An, Q. Liang, L. Yin, and J. Tian, “Reconstruction for Fluorescence Molecular Tomography via Adaptive Group Orthogonal Matching Pursuit,” IEEE Transactions on Biomedical Engineering 1, 1 (2020). [CrossRef]

15. F. Stuker, C. Baltes, K. Dikaiou, D. Vats, L. Carrara, E. Charbon, J. Ripoll, and M. Rudin, “Hybrid small animal imaging system combining magnetic resonance imaging with fluorescence tomography using single photon avalanche diode detectors,” IEEE Trans. Med. Imaging 30(6), 1265–1273 (2011). [CrossRef]

16. S. Jiang, J. Liu, Y. An, G. Zhang, J. Ye, Y. Mao, K. He, C. Chi, and J. Tian, “Novel l 2, 1-norm optimization method for fluorescence molecular tomography reconstruction,” Biomed. Opt. Express 7(6), 2342–2359 (2016). [CrossRef]

17. J. Shi, F. Liu, G. Zhang, J. Luo, and J. Bai, “Enhanced spatial resolution in fluorescence molecular tomography using restarted L1-regularized nonlinear conjugate gradient algorithm,” J. Biomed. Opt. 19(4), 046018 (2014). [CrossRef]

18. D. Zhu and C. Li, “Nonconvex regularizations in fluorescence molecular tomography for sparsity enhancement,” Phys. Med. Biol. 59(12), 2901–2912 (2014). [CrossRef]

19. J. Baritaux, K. Hassler, and M. Unser, “An efficient numerical method for general L (p) regularization in fluorescence molecular tomography,” IEEE Trans. Med. Imaging 29(4), 1075–1087 (2010). [CrossRef]

20. Y. Yuan, H. Guo, H. Yi, J. Yu, X. He, and X. He, “Correntropy-induced metric with Laplacian kernel for robust fluorescence molecular tomography,” Biomed. Opt. Express 12(10), 5991–6012 (2021). [CrossRef]

21. H. Zhang, X. He, J. Yu, X. He, H. Guo, and Y. Hou, “L1-L2 norm regularization via forward-backward splitting for fluorescence molecular tomography,” Biomed. Opt. Express 12(12), 7807–7825 (2021). [CrossRef]

22. Bangerth, and Wolfgang, “A Framework for the Adaptive Finite Element Solution of Large-Scale Inverse Problems,” Siam Journal on Scientific Computing 30, 2965–2989 (2011). [CrossRef]

23. Y. An, J. Liu, G. Zhang, J. Ye, Y. Du, Y. Mao, C. Chi, and J. Tian, “A novel region reconstruction method for fluorescence molecular tomography,” IEEE Trans. Biomed. Eng. 62(7), 1818–1826 (2015). [CrossRef]

24. J. Shi, F. Liu, J. Zhang, J. Luo, and J. Bai, “Fluorescence molecular tomography reconstruction via discrete cosine transform-based regularization,” J. Biomed. Opt. 20(5), 055004 (2015). [CrossRef]

25. H. Zou and T. Hastie, “Regularization and variable selection via the elastic net,” Journal of the royal statistical society: series B (statistical methodology) 67(2), 301–320 (2005). [CrossRef]

26. H. Wang, C. Bian, L. Kong, Y. An, Y. Du, and J. Tian, “A Novel Adaptive Parameter Search Elastic Net Method for Fluorescent Molecular Tomography,” IEEE Trans. Med. Imaging 40(5), 1484–1498 (2021). [CrossRef]

27. Y. Liu, J. Liu, Y. An, and S. Jiang, “Novel regularized sparse model for fluorescence molecular tomography reconstruction,” in International Conference on Innovative Optical Health Science (SPIE, 2017), pp. 39–44.

28. H. F. Schepker, C. Bockelmann, A. Dekorsy, A. Bartels, D. Trede, and K. S. Kazimierski, “C-Curve: A Finite Alphabet Based Parameter Choice Rule for Elastic-Net in Sporadic Communication,” IEEE Commun. Lett. 18(8), 1443–1446 (2014). [CrossRef]

29. P. C. Hansen, “The L-curve and its use in the numerical treatment of inverse problems,” (1999).

30. C. Huang, H. Meng, Y. Gao, S. X. Jiang, and J. Tian, “Fast and robust reconstruction method for fluorescence molecular tomography based on deep neural network,” in Imaging, Manipulation, and Analysis of Biomolecules, Cells, and Tissues XVII (2019).

31. L. Zhao, V. K. Lee, S.-S. Yoo, G. Dai, and X. Intes, “The integration of 3-D cell printing and mesoscopic fluorescence molecular tomography of vascular constructs within thick hydrogel scaffolds,” Biomaterials 33(21), 5325–5332 (2012). [CrossRef]

32. M. Li, X. Cao, F. Liu, B. Zhang, J. Luo, and J. Bai, “Reconstruction of fluorescence molecular tomography using a neighborhood regularization,” IEEE Trans. Biomed. Eng. 59(7), 1799–1803 (2012). [CrossRef]

33. L. Guo, F. Liu, C. Cai, J. Liu, and G. Zhang, “3D deep encoder–decoder network for fluorescence molecular tomography,” Opt. Lett. 44(8), 1892–1895 (2019). [CrossRef]

34. H. Meng, Y. Gao, X. Yang, K. Wang, and J. Tian, “K-nearest Neighbor Based Locally Connected Network for Fast Morphological Reconstruction in Fluorescence Molecular Tomography,” IEEE Trans. Med. Imaging 39(10), 3019–3028 (2020). [CrossRef]

35. P. Zhang, G. Fan, T. Xing, F. Song, and G. Zhang, “UHR-DeepFMT: ultra-high spatial resolution reconstruction of fluorescence molecular tomography based on 3-D fusion dual-sampling deep neural network,” IEEE Trans. Med. Imaging 40(11), 3217–3228 (2021). [CrossRef]

36. P. Zhang, C. Ma, F. Song, G. Fan, Y. Sun, Y. Feng, X. Ma, F. Liu, and G. Zhang, “A review of advances in imaging methodology in fluorescence molecular tomography,” Phys. Med. Biol. 67(10), 10TR01 (2022). [CrossRef]

37. Y. Chen, W. Li, M. Du, L. Su, H. Yi, F. Zhao, K. Li, L. Wang, and X. Cao, “Elastic net-based non-negative iterative three-operator splitting strategy for Cerenkov luminescence tomography,” Opt. Express 30(20), 35282–35299 (2022). [CrossRef]

38. X. Cao, X. Wei, F. Yan, L. Wang, and X. He, “A Novel Stacked Denoising Autoencoder-Based Reconstruction Framework for Cerenkov Luminescence Tomography,” IEEE Access 7, 1 (2019). [CrossRef]

39. Y. Chen, M. Du, W. Li, L. Su, H. Yi, F. Zhao, K. Li, L. Wang, and X. Cao, “ABPO-TVSCAD: alternating Bregman proximity operators approach based on TVSCAD regularization for bioluminescence tomography,” Phys. Med. Biol. 67(21), 215013 (2022). [CrossRef]

40. D. Zhu and C. Li, “Nonuniform update for sparse target recovery in fluorescence molecular tomography accelerated by ordered subsets,” Biomed. Opt. Express 5(12), 4249–4259 (2014). [CrossRef]

41. Y. Jiao, Q. Jin, X. Lu, and W. Wang, “Alternating direction method of multipliers for linear inverse problems,” SIAM J. Numer. Anal. 54(4), 2114–2137 (2016). [CrossRef]

42. X. He, J. Liang, X. Wang, J. Yu, X. Qu, X. Wang, Y. Hou, D. Chen, F. Liu, and J. Tian, “Sparse reconstruction for quantitative bioluminescence tomography based on the incomplete variables truncated conjugate gradient method,” Opt. Express 18(24), 24825–24841 (2010). [CrossRef]

43. J. H. Lee, A. Joshi, and E. M. Sevick-Muraca, “Fully adaptive finite element based tomography using tetrahedral dual-meshing for fluorescence enhanced optical imaging in tissue,” Opt. Express 15(11), 6955–6975 (2007). [CrossRef]

44. A. Joshi, W. Bangerth, and E. M. Sevick-Muraca, “Adaptive finite element based tomography for fluorescence optical imaging in tissue,” Opt. Express 12(22), 5402–5417 (2004). [CrossRef]

45. C. Qin, J. Tian, X. Yang, K. Liu, G. Yan, J. Feng, Y. Lv, and M. Xu, “Galerkin-based meshless methods for photon transport in the biological tissue,” Opt. Express 16(25), 20317–20333 (2008). [CrossRef]

46. H. Li and L. Ding, “Generalized conditional gradient method for elastic-net regularization,” Journal of Computational and Applied Mathematics 403, 113872 (2022). [CrossRef]

47. B. Jin, D. A. Lorenz, and S. Schiffler, “Elastic-net regularization: error estimates and active set methods,” Inverse Problems 25(11), 115022 (2009). [CrossRef]

48. K. Bredies, D. A. Lorenz, and P. Maass, “A generalized conditional gradient method and its connection to an iterative shrinkage method,” Comput. Optim. Appl. 42(2), 173–193 (2009). [CrossRef]

49. S. Boyd, N. Parikh, E. Chu, B. Peleato, and J. Eckstein, “Distributed optimization and statistical learning via the alternating direction method of multipliers,” FNT in Machine Learning 3(1), 1–122 (2011). [CrossRef]

50. H. Zou and H. H. Zhang, “On the adaptive elastic-net with a diverging number of parameters,” Ann. Statist. 37(4), 1733 (2009). [CrossRef]

51. S. L. Jacques, “Optical Properties of Biological Tissues: A Review,” Phys. Med. Biol. 58(11), R37–R61 (2013). [CrossRef]

52. B. Parvitte, C. Risser, R. Vallon, and V. Zéninari, “Quantitative simulation of photoacoustic signals using finite element modelling software,” Appl. Phys. B 111(3), 383–389 (2013). [CrossRef]

53. H. Yi, D. Chen, W. Li, S. Zhu, X. Wang, J. Liang, and J. Tian, “Reconstruction algorithms based on l1-norm and l2-norm for two imaging models of fluorescence molecular tomography: a comparative study,” J. Biomed. Opt. 18(5), 056013 (2013). [CrossRef]

Methods		IS- $L_{1}$	IVTCG- $L_{1}$	N-EN	APSEN	GCGM-ARP
Single-source		$η_{1} :$ 2e-10	$η_{1} :$ 1e0	$η_{1} :$ 1e-5 $η_{2} :$ 1e-6	$η_{1} :$ 0.5 $η_{2} : 2$ e-9	$α :$ 8e-6 $β :$ 1e-5
Dual-source		$η_{1} :$ 2e-10	$η_{1} :$ 1e0	$η_{1} :$ 1e-5 $η_{2} :$ 1e-6	$η_{1} :$ 0.9 $η_{2} :$ 1.1e-9	$α :$ 1e-5 $β :$ 1e-5
Dual-source	5%	$η_{1} :$ 2e-10	$η_{1} :$ 1e0	$η_{1} :$ 1e-5 $η_{2} :$ 1e-6	$η_{1} :$ 0.5 $η_{2} : 2$ e-9	$α :$ 2e-6 $β :$ 1e-5
Robustness	15%	$η_{1} :$ 2e-10	$η_{1} :$ 1e0	$η_{1} :$ 1e-5 $η_{2} :$ 1e-6	$η_{1} :$ 0.6 $η_{2} :$ 1.67e-9	$α :$ 1e-5 $β :$ 1e-5
Robustness	25%	$η_{1} :$ 2e-10	$η_{1} :$ 1e0	$η_{1} :$ 1e-5 $η_{2} :$ 1e-6	$η_{1} :$ 0.8 $η_{2} :$ 1.25e-9	$α :$ 1e-4 $β :$ 1e-5
In vivo experiment		$η_{1} :$ 2e-10	$η_{1} :$ 1e0	$η_{1} :$ 1e-5 $η_{2} :$ 1e-6	$η_{1} :$ 0.5 $η_{2} : 2$ e-9	$α :$ 1e-5 $β :$ 1e-5

Tissue	$μ_{α x} (m m^{- 1})$	$μ_{s x}^{'} (m m^{- 1})$	$μ_{a m} (m m^{- 1})$	$μ_{s m}^{'} (m m^{- 1})$	g
Muscle	0.0745	0.4115	0.0474	0.3122	0.97
Bone	0.0521	2.4415	0.0326	2.114	0.93
Liver	0.3016	0.6676	0.1921	0.6023	0.93
Lung	0.1681	2.1569	0.1045	2.0477	0.93
Heart	0.0504	0.9437	0.0331	0.8203	0.90

Method	True center (mm)	Reconstructed center (mm)	LE (mm)	Dice	RIE
IS- $L_{1}$	(-1, -1, 15)	(-0.198, -0.721, 15.310)	0.904	0.496	0.302
IVTCG- $L_{1}$	(-1, -1, 15)	(-0.829, -0.643, 14.311)	0.794	0.625	0.824
N-EN	(-1, -1, 15)	(-0.983, -1.193, 14.541)	0.498	0.693	0.553
APSEN	(-1, -1, 15)	(-8.145, -1.077, 14.737)	0.331	0.795	0.337
GCGM-ARP	(-1, -1, 15)	(-0.835, -0.884, 14.820)	0.270	0.952	0.181

Method	True center (mm)	Reconstructed center (mm)	LE (mm)	Dice	RIE
IS- $L_{1}$	(2, 1, 10)	(1.010, 1.075, 10.620)	1.170	0.462	0.833
IS- $L_{1}$	(2, 1, 18)	(1.313, 0.349, 17.485)	1.077	0.446	0.833
IVTCG- $L_{1}$	(2, 1, 10)	(1.654, 0.732, 9.262)	0.858	0.548	0.749
IVTCG- $L_{1}$	(2, 1, 18)	(2.973, 1.407, 17.718)	1.092	0.298	0.749
N-EN	(2, 1, 10)	(2.524, 0.648, 9.925)	0.636	0.545	0.731
N-EN	(2, 1, 18)	(1.851, 1.265, 17.510)	0.576	0.538	0.731
APSEN	(2, 1, 10)	(1.766, 0.571, 9.975)	0.489	0.783	0.674
APSEN	(2, 1, 18)	(1.850, 1.264, 17.517)	0.571	0.698	0.674
GCGM-ARP	(2, 1, 10)	(1.848, 1.339, 10.220)	0.432	0.839	0.502
GCGM-ARP	(2, 1, 18)	(1.832, 0.802, 18.177)	0.314	0.834	0.502

Method	True center (mm)	Reconstructed center (mm)	LE (mm)	Dice
IS- $L_{1}$	(17.5, 21.5, 9.0)	(17.998, 20.641, 8.616)	1.065	0.254
IVTCG- $L_{1}$	(17.5, 21.5, 9.0)	(16.652, 21.356, 8.296)	1.112	0.254
N-EN	(17.5, 21.5, 9.0)	(17.625, 20.999, 8.456)	0.760	0.325
APSEN	(17.5, 21.5, 9.0)	(17.518, 20.846, 9.178)	0.678	0.625
GCGM-ARP	(17.5, 21.5, 9.0)	(17.893, 21.477, 8.818)	0.434	0.659

Generalized conditional gradient method with adaptive regularization parameters for fluorescence molecular tomography

Abstract

1. Introduction

2. Methodology

2.1 Photon propagation model

2.2 Reconstruction based on GCGM-ARP scheme

2.2.1 Sub-problem of determining the descent direction

2.2.2 Sub-problem of determining step size

2.2.3 Adaptive regularization parameters

3. Experiments and results

3.1 Regularization parameters and evaluation index

3.2 Experimental setting

3.2.1 Numerical simulation

3.2.2 In vivo experiment

3.3 Experimental results

3.3.1 Single-source simulation reconstruction

3.3.2 Dual-source simulation reconstruction

3.3.3 Robustness experiment

3.3.4 In vivo experiment reconstruction

4. Discussion and conclusion

Funding

Disclosures

Data availability

References

Data availability

Cited By

Figures (6)

Tables (5)

Equations (27)

Optics Express