Automatic selection of regularization parameters for dynamic fluorescence molecular tomography: a comparison of L-curve and U-curve methods

Maomao Chen; Han Su; Yuan Zhou; Chuangjian Cai; Dong Zhang; Jianwen Luo

doi:10.1364/BOE.7.005021

1. Introduction

Fluorescence molecular tomography (FMT) is an emerging imaging technique that allows for noninvasive and quantitative reconstruction of three-dimensional (3D) distribution of fluorescent agents in the biological body in vivo [1]. By adding time as a new dimension, dynamic FMT (DFMT) can be used for pharmacokinetic studies [2, 3], which reflect the absorption, distribution and excretion characteristics of fluorescent agents inside the biological tissues. Hence, DFMT is a potential method for tumor detection, treatment monitoring and drug delivery studies [4–6].

In order to monitor the metabolic process of fluorescent agents in the body, the imaged object is continuously monitored for tens of minutes. During each minute, a set of fluorescent measurements is acquired [2]. The reconstruction of DFMT images consists of two parts [3]. Firstly, conventional FMT reconstruction methods are used to obtain a time series of 3D FMT images. Secondly, the FMT images are fitted into a two-compartment model to provide the pharmacokinetic parameters of each voxel in the imaged volume. Therefore, high-quality reconstruction of each FMT image is crucial for DFMT problems. In order to obtain a meaningful solution, the FMT reconstruction problem must be regularized, for example, using Tikhonov regularization, and a regularization parameter λ is required. In a typical DMFT implementation, the imaged object is monitored for more than 50 minutes, and the fluorescent concentration varies with time. Thus, each set of the fluorescent measurements should be considered as an independent FMT reconstruction problem, in which a particular regularization parameter should be used. As a consequence, tens of different regularization parameters are required and are very difficult to choose manually.

A number of methods for automatic selection of regularization parameters have been proposed over years. The unbiased predictive risk estimator (UPRE) method [7] is based on minimizing the statistical estimator of the mean squared norm of predictive error. The discrepancy principle (DP) method [8, 9] is to seek a regularization parameter $λ$ , with which the regularized residual equals to the data variance. These two methods require accurate estimation of the system noise, which is very difficult for DFMT problems. As an alternative to UPRE, the method of generalized cross validation [10, 11] (GCV) does not require prior knowledge of the variance of the white noise. However, previous studies have reported that GCV may fail at high noise levels, especially when the noises of the measurements are correlated [12, 13]. The f-slope method [14] plots the solution norm (norm of the regularized solution) against $\ln (1 / λ)$ and the regularization parameter is selected at the flattest part of the curve. The f-slope method considers only the solution norm but not the residual norm (norm of the differences between the collected measurements and the predicted measurements using the regularized solution). It has been reported that this method may overestimate the regularization parameter for implementations with small noise levels [15].

For Tikhonov regularization, the L-curve method [12, 13] is commonly used for automatic selection of the regularization parameter. It balances the residual norm and the solution norm by selecting the regularization parameter at the point corresponding to the maximum curvature of the L-curve. The L-curve method has been applied to optical tomography problems [15], and shows good repeatability and reconstruction quality. However, it may fail to find a pronounced corner for highly ill-posed problems [15, 16]. Moreover, for diffuse optical tomography (DOT), the L-curve method yields over-smooth solution [17].

As an alternative to L-curve, the U-curve method has been proposed [18, 19], and its feasibility for FMT problems has been studied [20]. It has been proved that the U-curve always has a local minimum [18]. In addition, the U-curve method is relatively computationally efficient [20], because the interval where the local minimum exists can be determined.

The parametric images are obtained by curve fitting the fluorescent concentrations at a series of time points to a two-compartment model for the DFMT problems. The quality of the parametric images relies heavily on the accuracy of each FMT image. Therefore, in this work, we quantitatively evaluated the parametric images obtained with the L-curve and U-curve methods in numerical simulations, phantom experiments and in vivo experiments. We found that the U-curve method resulted in higher accuracy, noise-resistance and robustness in parametric imaging. It could be a proper choice for automatic selection of the regularization parameters in DFMT problems.

This paper is organized as follows. Section 2.1 describes the noncontact, full-angle FMT and X-ray computed tomography (XCT) imaging system. Section 2.2 introduces the reconstruction method for DFMT problems. Sections 2.3 and 2.4 present the selection of regularization parameters with the L-curve and U-curve methods, respectively. Section 2.5 introduces the discrete Picard condition (DPC), which is employed to estimate the accuracy of the selected parameters. The setups and results of the numerical simulations, phantom experiments and in vivo experiments are shown in Sections 3.1-3.3. Finally, Section 4 discusses the results and draws the conclusion.

2. Methods

2.1 Data acquisition system

A hybrid FMT/XCT system [21] is used to obtain the DFMT boundary measurements. As shown in Fig. 1, a free-space and full-angle FMT system is used to acquire the fluorescence data sets, while an XCT system is used in the in vivo experiments to obtain the anatomical information of the small animal, which provides structural priors to constrain the reconstruction results. For DFMT problems, the small animal is fixed on the stage and continuously rotated for K minutes to monitor the metabolic process of the fluorescent agents. During each minute, in order to collect the fluorescent measurements, the small animal is rotated by 360° and the fluorescent data sets of S projections (uniformly distributed within 360°) are acquired by a charge coupled device (CCD) camera. Therefore, totally P = KS projections are obtained during the entire DFMT acquisition process.

Fig. 1 Schematic of the hybrid FMT/XCT system. The FMT system is used to acquire the DFMT measurements, while the XCT system obtains the anatomical information of the small animal.

Download Full Size | PDF

2.2 Reconstruction method for DFMT problems

2.2.1 FMT images reconstruction

For FMT problems, a coupled diffusion equation [22] is used to describe the propagation of excitation and emission light. Using Green's function theory, the fluorescence measurements $Φ_{m} (r_{d}, r_{s})$ detected at a point $r_{d}$ due to a excitation source at $r_{s}$ can be written as [23]

Φ_{m} (r_{d}, r_{s}) = Θ \int_{Ω} G_{m} (r_{d,} r) x (r) G_{x} (r, r_{s}) d^{3} r

where x(r) denotes the fluorescent yield to be reconstructed. The Green's function

G_{m} (r_{d}, r)

describes the light propagation from an arbitrary position r inside the medium to the detector point r_d at the emission wavelength.

G_{x} (r, r_{s})

stands for the light propagation from the source point r_s to position r at the excitation wavelength;

Θ

accounts for the unknown gain and attenuation factors of the system.

In order to reduce the influence of the heterogeneity, the normalized Born approximation is employed as follows [24]

Η (r_{d}, r_{s}) = \frac{Φ_{m} (r_{d}, r_{s})}{Φ_{x} (r_{d}, r_{s})} = Θ \int_{Ω} \frac{G_{m} (r_{d}, r) x (r) G_{x} (r, r_{s})}{G_{x} (r_{d}, r_{s})} d^{3} r

where

Φ_{x} (r_{d}, r_{s})

denotes the excitation measurements corresponding to

Φ_{m} (r_{d}, r_{s})

.

By discretizing the imaged geometry Ω into N voxels, the forward model can be transformed into a linear system as [25]

b_{k} = W x_{k}

where

b_{k} = {[b_{k} (r_{1}), b_{k} (r_{2}), \dots, b_{k} (r_{M})]}^{T}

is the vector of boundary measurement predicted by the Born ratio model. The index

k (k = 1, 2, \dots K)

denotes the kth set of the fluorescent measurements and FMT images, and M is the number of the measurements.

x_{k} = {[x_{k} (r_{1}), x_{k} (r_{2}), \dots, x_{k} (r_{N})]}^{T}

is the vector of the fluorescence distribution. W is the weight matrix with a size of

M \times N

and its entries are defined by [25]

W_{i j} = Δ V Θ \frac{G_{m} (r_{d}^{i}, r^{j}) G_{x} (r^{j}, r_{s}^{i})}{G_{x} (r_{d}^{i}, r_{s}^{i})}

here,

Δ V

is the volume of each voxel. The index

i (i = 1, 2, \dots, M)

denotes the particular source-detector pair corresponding to the ith element of predicted boundary measurements, and the index

j (j = 1, 2, \dots, N)

corresponds to a particular voxel in the discretized geometry

Ω

.

Finally, the FMT image $x_{k}$ can be solved using the Tikhonov regularization method. The objective function is given as follows [3]

Ψ (x_{k}) = {‖ b_{k} - W x_{k} ‖}_{2}^{2} + λ_{k} {‖ L x_{k} ‖}_{2}^{2}

where

λ_{k}

and

b_{k} = {[b_{k} (r_{1}), b_{k} (r_{2}), \dots, b_{k} (r_{M})]}^{T}

are the regularization parameter and the boundary measurements during the kth minute, respectively. The regularization matrix L is used to constrain the reconstruction results. In the numerical simulations and phantom experiments, the structural priors are not used and L is set to be an identity matrix. In the in vivo experiments, L is constructed as a Laplacian-type matrix [26] based on the structural information obtained from the XCT system described in Section 2.1.

2.2.2 Two-compartment model and curve fitting

Indocyanine green (ICG) is used as the fluorescent agent in this paper. Its metabolic process can be approximated as a two-compartment model described by a double-exponential, four-parameter function [27, 28]

I (t) = - A \exp (- α t) + B \exp (- β t)

where I(t) denotes the fluorescence concentration at time t. A and B are the zero-time intercepts, which constitute the initial fluorescence concentration.

α

and

β

are the uptake and excretion rates respectively, which describe the wash-in and wash-out processes of ICG in organs and tissues [4].

In DFMT problems, after the FMT images are obtained, the pharmacokinetic parameters for each voxel of the geometry are calculated using a nonlinear least-square curve fitting method. Let $φ_{j} = {[A_{j}, B_{j}, α_{j}, β_{j}]}^{T}$ denote the pharmacokinetic parameter vector for voxel j, the DFMT result can be obtained as follows [2,3]

x_{k} (φ_{j}) = - A_{j} \exp (- α_{j} t_{k}) + B_{j} \exp (- β_{j} t_{k})

where

t_{k}

denotes the time point k, and

x_{k} (φ_{j})

is the fluorescence concentration for voxel j at

t_{k}

.

2.3 L-curve method

For Tikhonov regularization, the L-curve method is generally employed to select the regularization parameter. With a log-log plot of the solution norm $η (λ)$ against the corresponding residual norm $ρ (λ)$ , $λ_{L}$ is determined corresponding to the point of the maximum curvature on the L-curve.

In order to calculate the solution norm and the residual norm, singular value decomposition (SVD) is employed. By considering that the regularization matrix L is square and invertible, the objective function Eq. (5) can be transformed into the standard form as

Ψ ({\bar{x}}_{k}) = {‖ b_{k} - \bar{W} {\bar{x}}_{k} ‖}_{2}^{2} + λ_{k} {‖ {\bar{x}}_{k} ‖}_{2}^{2}

where

\bar{W} = W L^{- 1}

, and

{\bar{x}}_{k} = L x_{k}

. The SVD of the transformed weight matrix

\bar{W}

is described as follows [12]

\bar{W} = U S V^{T} = \sum_{i = 1}^{n} u_{i} σ_{i} v_{i}^{T}

where S is a diagonal matrix that contains the singular values

σ_{i}

, which are nonnegative and nonincreasing numbers, i.e.,

σ_{1} \geq σ_{2} \geq \dots \geq σ_{n} \geq 0

.

u_{i}

and

v_{i}

are the left and right singular vectors of the corresponding ith columns of U and V, respectively.

Using the SVD theory, the solution norm $η (λ)$ and the residual norm $ρ (λ)$ can be analytically obtained as [12]

η (λ) = {‖ \bar{x} (λ) ‖}_{2}^{2} = \sum_{i = 1}^{n} \frac{σ_{i}^{2} f_{i}^{2}}{{(λ^{2} + σ_{i}^{2})}^{2}}

ρ (λ) = {‖ b - \bar{W} \bar{x} (λ) ‖}_{2}^{2} = \sum_{i = 1}^{n} \frac{λ^{4} f_{i}^{2}}{{(λ^{2} + σ_{i}^{2})}^{2}} + {‖ r_{⊥} ‖}_{2}^{2}

where

f_{i} = u_{i}^{T} b_{k}

is the Fourier coefficients, and

r_{⊥}

is the least squares residual, which stands for the component of the boundary measurements

b_{k}

orthogonal to the left singular vectors

u_{1}, \dots, u_{n}

[12, 18]. The Matlab Regularization Toolbox is used to obtain the L-curve parameters [29].

2.4 U-curve method

The U-curve is the plot of the sum of the reciprocals of the solution norm $η (λ)$ and the residual norm $ρ (λ)$ , i.e., [20]

U (λ) = \frac{1}{η (λ)} + \frac{1}{ρ (λ)}

A typical U-curve mainly consists of three parts: the left and right vertical parts, which are dominated by the residual norm and the solution norm, respectively, and the horizontal part in middle, which reflects a relative balance between $ρ (λ)$ and $η (λ)$ . It has been proved [18] that the function $U (λ)$ strictly decreases in the interval of $(0, σ_{n}^{2 / 3})$ and strictly increases in the interval of $(σ_{1}^{2 / 3}, \infty)$ , where $σ$ is the singular value of the transformed weight matrix $\bar{W}$ as described in Eq. (9). Therefore, a local minimum always exists in the range $λ \in (σ_{n}^{2 / 3}, σ_{1}^{2 / 3})$ , and the minimum point is selected to be $λ_{U}$ . In this work, all the singular values in the range $(σ_{n}^{2 / 3}, σ_{1}^{2 / 3})$ are tested to locate the optimal parameter $λ_{U}$ .

2.5 Discrete Picard condition

To estimate the reliability of the selected regularization parameters, the DPC is employed. It is defined as follows [30]: if the Fourier coefficients $| U_{i}^{T} b |$ averagely decay to zero faster than the singular values $σ_{i}$ , the regularized solutions $x_{λ}$ is guaranteed to have approximately the same properties as the exact solution $x_{0}$ . For DFMT problems, the right-hand-side b is contaminated with measurement noise, approximation errors and rounding errors, and the boundary measurements rarely satisfy the DPC. Therefore, it is essential to select a proper regularization parameter, which can fairly balance the residual norm and the solution norm in the regularized solution. The regularization parameter should fulfill the DPC, and locates before the point where the data becomes dominated by errors. Generally, a smaller regularization parameter improves the resolution of the reconstructed images while increases the noises of the images.

In this paper, the Fourier coefficients $| U_{i}^{T} b |$ , the singular values $σ_{i}$ , the corresponding quotients $| U_{i}^{T} b | / σ_{i}$ , and the regularization parameters $λ_{L}$ and $λ_{U}$ are plotted on the same graph to evaluate how well the selected parameters satisfy the DPC [20].

3. Simulations, experiments and results

3.1 Simulations

3.1.1 Setups

The simulations are carried out on a cylinder with a height of 3 cm and a radius of 1.5 cm, as depicted in Fig. 2(a). Two tubes filled with ICG with an edge-to-edge distance (EED) of 0.5 cm are placed inside the cylinder, as shown in Fig. 2(b). The height and the radius of the tubes are 1 cm and 0.2 cm, respectively. FMT images $x_{k} (k = 1, \dots, 60)$ in 60 minutes (L = 60) are generated according to the metabolic curves shown in Fig. 2(c). The curves of tubes 1 and 2 are set to mimic the metabolic process of the fluorescent agents within 60 minutes inside the liver and lungs [31], respectively. The corresponding pharmacokinetic parameters are listed in Table 1 [31]. A CCD camera is simulated to generate the fluorescent data sets. The boundary measurements $b_{k} (k = 1, \dots, 60)$ in each minute are calculated using Eq. (3). Five levels of zero-mean white Gaussian noise with signal-to-noise ratios (SNRs) of 40 dB, 35 dB, 30 dB, 25 dB and 20 dB respectively are added to generate the noisy measurements.

Fig. 2 Setups of the simulations and phantom experiments. (a) 3D geometric configurations of the simulations and phantom experiments. The height and the radius of the cylinder are 3 cm and 1.5 cm, respectively. Two cylindrical fluorescent tubes with an EED of 0.5 cm are included. The height and the radius of the tubes are 1 cm and 0.2 cm, respectively. (b) Cross-sectional image corresponding to the red line depicted in (a). (c) Metabolic curves of the simulations and phantom experiments. The fluorescent concentrations for the tubes in the simulations are depicted as red solid and blue dashed lines, respectively. The 10 groups of fluorescent concentrations for the tubes in the phantom experiments are shown as red stars and blue circles, respectively.

Download Full Size | PDF

Table 1. Pharmacokinetic parameters in the simulations and phantom experiments

View Table | View all tables in this article

In the reconstruction, the cylinder is discretized into $N = 17, 856$ voxels, and the number of the measurements obtained during each minute is $M = 16, 424$ . The absorption coefficient $μ_{a}$ and reduced scattering coefficient ${μ^{'}}_{s}$ are set to be 0.02 cm⁻¹ and 10 cm⁻¹, respectively.

3.1.2 Results

60 regularization parameters are obtained with the two methods at each noise level, and the means and the standard deviations of the parameters are listed in Table 2. The mean parameter of the L-curve method increases from 1.15 × 10² to 2.44 × 10³ when the SNR drops from 40 dB to 20 dB, while the mean parameters of the U-curve method are very stable for all the noise levels.

Table 2. Means and standard deviations of the regularization parameters obtained with the two methods at different noise levels

View Table | View all tables in this article

The L-curve and U-curve at the time point of $t_{k} = 10$ min at the SNR of 40 dB are depicted in Figs. 3(a) and 3(b), respectively. Both the maximum curvature of the L-curve and the minimum point of the U-curve can be located. Figure 3(c) shows the DPC plot in the simulations. The regularization parameters $λ_{L}$ and $λ_{U}$ are shown as black dashed and red solid horizontal lines, respectively. It can be observed that the Fourier coefficients $| u_{i}^{T} b |$ (green stars) averagely decay to zero faster than the corresponding singular values (blue dots) above about 10⁻¹, which is closer to $λ_{U}$ . It indicates that the U-curve method may obtains higher imaging resolution. Similar results can be obtained from the regularization parameters at other time points.

Fig. 3 (a) L-curve plots on log-log scale. (corner at λ_L = 1.40 × 10²). (b) U-curve plots on log-log scale (corner at λ_U = 2.09 × 10⁻²). (c) The DPC plot for the simulations, where i denotes the number of the singular values. The Fourier coefficients, the singular values and their corresponding quotients are depicted as green crosses, blue dots and red circles, respectively. The L-curve and U-curve coefficients are shown as black dashed and red solid horizontal lines, respectively.

Download Full Size | PDF

In order to compare the imaging resolution of the two methods, the FMT images at the SNR of 40 dB at the time point of $t_{k} = 10$ min are presented in Fig. 4. 40 dB is selected because both methods obtain good parametric images. The 3D geometric configuration is shown in Fig. 4(a), and the 3D FMT images reconstructed using the regularization parameters obtained with the L-curve and U-curve methods are shown in Figs. 4(b) and 4(c), respectively. Figures 4(d)-4(f) are the cross-sectional images corresponding to the slices depicted by the red lines in Fig. 4(a)-4(c). Although the images obtained with the U-curve method are noisier than those with the L-curve method (see the noise in the background in Fig. 4(c)), the two fluorescent tubes can be better distinguished by the U-curve method. It indicates that the U-curve method obtains higher imaging resolution. The FMT images in Fig. 4(b) and 4(c) are reconstructed with the parameters of $λ_{L} = 1.40 \times 10^{2}$ and $λ_{U} = 2.09 \times 10^{- 2}$ , respectively. Similar results can be obtained from the FMT images at other time points (not shown).

Fig. 4 FMT images reconstructed using the regularization parameters obtained with the L-curve and U-curve methods at the time point of $t_{k} = 10$ min in the simulations. (a) The 3D geometric configuration of the simulations. (b)-(c) The 3D FMT images reconstructed with the L-curve and U-curve methods, respectively. (d)-(f) The cross-sectional images corresponding to the slices depicted by the red lines shown in (a)-(c), respectively. All the images are normalized to the maximum of the results.

Download Full Size | PDF

The intensity profiles along the x axis (i.e., y = 0) in the cross-sectional images shown in Figs. 4(d)-4(f) is depicted in Fig. 5. The true profile is shown as black dashed lines. The profiles obtained with the L-curve and U-curve methods are plotted as blue and red solid lines, respectively. The profiles further demonstrate that the U-curve method obtains higher imaging resolution.

Fig. 5 Intensity profiles along the x axis (i.e., y = 0) in the cross-sectional images shown in Figs. 4(d)-4(f). The true profile is shown as black dashed lines. The profiles obtained with the L-curve and U-curve methods are plotted as blue and red solid lines, respectively.

Download Full Size | PDF

The parametric images reconstructed with the two methods at the SNRs of 40 dB, 30 dB and 20 dB respectively are depicted in Fig. 6. Both methods can obtain acceptable results at the noise level of 40 dB, as shown in Figs. 6(b1)-6(b4) and 6(e1)-6(e4). When the SNR decreases to 30 dB, the images obtained with the L-curve method obviously deteriorate. The two tubes cannot be distinguished in Fig. 6(c3). By contrast, the U-curve method can still obtain reasonable parametric images, as shown in Figs. 6(f1)-6(f4). Actually, the U-curve method obtains acceptable results even when the SNR decreases to 25 dB (not shown). Finally, when the SNR drops to 20 dB, the $α$ image obtained with the U-curve method is too blurred to recognize the two tubes, as depicted in Fig. 6(g3). The results illustrate that the U-curve method performs better than the L-curve method in high-noise situations.

Fig. 6 Cross-sectional parametric images of the simulations at different noise levels. (a1)-(a4) The true parametric images. (b1)-(d4) The images reconstructed with the L-curve method at different noise levels. (e1)-(g4) The images reconstructed with the U-curve method at different noise levels.

Download Full Size | PDF

The normalized root mean square error ( $N R M S E = {‖ \hat{e} - e ‖}_{2} / {‖ e ‖}_{2}$ , where $e$ and $\hat{e}$ stands for the true and reconstructed parametric images, respectively) is used to quantitatively estimate the relative errors between the true and the reconstructed images. A smaller NRMSE indicates better accuracy of the reconstructed image [31]. As listed in Table 3, at all the noise levels, the U-curve method obtains smaller NRMSEs than the L-curve method.

Table 3. NRMSEs of the parametric images obtained with the two methods at different noise levels

View Table | View all tables in this article

Table 4 lists the average pharmacokinetic parameters for the two tubes obtained with the two methods at the noise level of 40 dB. Because the pharmacokinetic parameters A and B have arbitrary units (a.u.) without knowing the gain of the FMT system, the two parameters are normalized to the mean value of A in tube 1 [31]. Compared with those of the L-curve method, the pharmacokinetic parameters obtained with the U-curve method are closer to the true parameters. It indicates higher reconstruction accuracy of the U-curve method.

Table 4. Pharmacokinetic parameters and correlation coefficients in the simulations

View Table | View all tables in this article

Figure 7 depicts the reconstructed metabolic curves. In order to evaluated how well the reconstructed curves match with the true curves, the correlation coefficients (CC) are calculated as

C C = \frac{\sum_{i = 1}^{n} (X_{r, i} - {\bar{X}}_{r}) (X_{t, i} - {\bar{X}}_{t})}{\sqrt{\sum_{i = 1}^{n} {(X_{r, i} - {\bar{X}}_{r})}^{2} \cdot \sum_{i = 1}^{n} {(X_{t, i} - {\bar{X}}_{t})}^{2}}}

where

X_{t}

and

X_{r}

denote the true and reconstructed metabolic curves, respectively. A larger CC indicates better reconstruction accuracy. The CC results are listed in Table 4. The CCs obtained with the U-curve method (0.999 and 0.999) are larger than those with the L-curve method (0.998 and 0.998). It means the reconstructed curves obtained with the U-curve method match better with the true curves.

Fig. 7 Metabolic curves in the simulations. (a) The metabolic curves of tube 1. The true curve and the curves obtained with the L-curve and U-curve methods are depicted as green solid, blue dotted and red dashed lines, respectively. (b) The metabolic curves of tube 2. The true curve and the curves obtained with the L-curve and U-curve methods are depicted as green solid, blue dotted and red dashed lines, respectively.

Download Full Size | PDF

3.2 Phantom experiments

3.2.1 Setups

The phantom experiments are performed with a noncontact, full-angle FMT imaging system shown in Fig. 1 (XCT imaging is not used in the phantom experiments). The geometric setups of the phantom experiments are identical to those of the simulations, as shown in Figs. 2(a) and 2(b). In order to obtain dynamic measurements, ten sets of the fluorescent measurements are acquired by filling the two tubes with 10 groups of ICG concentrations, as indicated by the red stars and blue circles in Fig. 2(c). Considering that the ICG concentrations vary much faster during the absorption stage than during the excretion stage, the time intervals are set to be 2 minutes for the first 5 groups in the increasing stage of the metabolic curves, and 10 minutes for the last 5 groups in the decreasing stage of the curves.

The excitation is applied by a 300 W Xenon lamp (MAX-302, Asahi Spectra, Torrance, CA, USA) with a power density of 0.03 mW/cm². The excitation light is filtered through a 770 ± 6 nm band-pass excitation filter (XBPA770, Asahi Spectra, Torrance, CA, USA). On the opposite side of the excitation source, the emitted ICG fluorescence is filtered with an 840 ± 6 nm band-pass emission filter (FF01-840/12-25, Semrock, Rochester, NY, USA), and detected by a 512 × 512 pixel, −70°C cooled CCD camera (iXon DU-897, Andor Technologies, Belfast, Northern Ireland, UK). The exposure time of the CCD camera is 1 s, and the CCD binning is set to be 1 × 1.

The imaged phantom cylinder is filled with 1% intralipid; therefore the optical coefficients are identical to those of the simulations ( ${μ^{'}}_{s} = 10.0 c m^{- 1}$ , $μ_{a} = 0.02 c m^{- 1}$ ). The reconstructed region is discretized into $N = 18, 692$ voxels, and the number of the measurements obtained during each minute is $M = 16, 200$ .

3.2.2 Results

10 regularization parameters are obtained with the L-curve ( $λ_{L}$ ) and U-curve ( $λ_{U}$ ) methods in the phantom experiments. The means and standard deviations of the parameters are listed in Table 5. The mean of $λ_{L}$ (3.62 × 10¹) is significantly larger than that of $λ_{U}$ (9.83 × 10⁻²).

Table 5. Means and standard deviations of the regularization parameters in the phantom experiments

View Table | View all tables in this article

The log-log plots of the L-curve and U-curve at the time point of $t_{k} = 10$ min are shown in Figs. 8(a) and 8(b), respectively. It can be observed that, the L-curve dose not exhibit a neat corner in Fig. 8(a), while the minimum point of the U-curve can be easily detected in Fig. 8(b). The DPC plot is shown in Fig. 8(c), and $λ_{L}$ and $λ_{U}$ are shown as black dashed and red solid horizontal lines, respectively. The Fourier coefficients $| u_{i}^{T} b |$ (green stars) averagely decay to zero faster than the corresponding singular values (blue dots) above about 10⁻², which is close to $λ_{U}$ . It indicates that $λ_{U}$ satisfies the DPC and the U-curve method can obtain good reconstruction images. By comparison, $λ_{L}$ is significantly larger than 10⁻², which means that the L-curve method may not obtain reasonable reconstruction results. Similar results can be obtained from the regularization parameters at other time points (not shown).

Fig. 8 (a) L-curve plots on log-log scale (corner at λ_L = 3.53 × 10¹). (b) U-curve plots on log-log scale (corner at λ_U = 9.35 × 10⁻²). (c) The DPC plot in the phantom experiments, where i denotes the number of the singular values. The Fourier coefficients, the singular values and their corresponding quotients are depicted as green crosses, blue dots and red circles, respectively. The L-curve and U-curve coefficients are shown as black dashed and red solid horizontal lines, respectively.

Download Full Size | PDF

Figures 9(b) and 9(c) depict the FMT images reconstructed with the parameters of $λ_{L} = 3.53 \times 10^{1}$ and $λ_{U} = 9.35 \times 10^{- 2}$ at the time point of $t_{k} = 10$ min, respectively. The L-curve method fails to reconstruct a reasonable FMT images as shown in Fig. 9(b). It may result from the wrong detection of the maximum curvature point of the L-curve, as depicted in Fig. 8(a). On the contrary, the two fluorescent tubes can be clearly distinguished in the FMT images reconstructed with the U-curve method, as shown in Fig. 9(c). Figures 9(d)-9(f) are the cross-sectional images corresponding to the slices depicted by the red lines in Fig. 9(a)-9(c). According to the results, the U-curve method achieves higher imaging resolution. Similar results can be obtained from the FMT images at other time points (not shown).

Fig. 9 FMT images reconstructed using the regularization parameters obtained with the L-curve and U-curve methods at the time point of $t_{k} = 10$ min in the phantom experiments. (a) The 3D geometric configuration of the phantom experiments. (b) and (c) The 3D FMT images reconstructed with the L-curve and U-curve methods, respectively. (d)-(f) The cross-sectional images corresponding to the slices depicted by the red lines shown in (a)-(c), respectively. All the images are normalized to the maximum of the results.

Download Full Size | PDF

Figure 10 shows the cross-sectional parametric images, which correspond to the slices indicated by the red lines in Figs. 9(b) and 9(c). The fluorescent tubes cannot be observed in the DFMT images obtained with the L-curve method (Figs. 10(b1)-10(b4)), while the two tubes can be clearly recognized in the DFMT images obtained with the U-curve method (Figs. 10(c1)-10(c4)). As depicted in Table 6, the U-curve method obtains a smaller mean value of NRMSEs, which indicates higher accuracy of the parametric images.

Fig. 10 Cross-sectional parametric images of the phantom experiments. (a1)-(a4) True parametric images. (b1)-(b4) Parametric images reconstructed with the L-curve method. (c1)-(c4) Parametric images reconstructed with the U-curve method.

Download Full Size | PDF

Table 6. NRMSEs of the parametric images obtained with the two methods in the phantom experiments

View Table | View all tables in this article

Table 7 lists the average pharmacokinetic parameters and the CCs of the two tubes obtained with the L-curve and U-curve methods in the phantom experiments, and Fig. 11 depicts the corresponding metabolic curves. The average parameters A and B are normalized to the mean value of A in tube 1 [31]. The U-curve method shows smaller differences of the pharmacokinetic parameters and larger CCs of the metabolic curves. It indicates that the U-curve method achieves higher reconstruction accuracy.

Table 7. Pharmacokinetic parameters and correlation coefficients in the phantom experiments

View Table | View all tables in this article

Fig. 11 Metabolic curves in the phantom experiments. (a) The metabolic curves of tube 1. The true curve and the curves obtained with the L-curve and U-curve methods are depicted as green solid, blue dotted and red dashed lines, respectively. (b) The metabolic curves of tube 2. The true curve and the curves obtained with the L-curve and U-curve methods are depicted as green solid, blue dotted and red dashed lines, respectively.

Download Full Size | PDF

3.3 In vivo experiments

3.3.1 Setups

In vivo experiments are conducted based on the hybrid FMT/XCT system described in Sec. 2.1. A healthy BALB/c nude mouse with an age of about 8 weeks is anesthetized and fixed on the rotation stage. A bolus of ICG (0.1mL, 50μg/mL) is injected via the tail vain. During the acquisition of DFMT measurements, the mouse is monitored for 50 minutes ( $K = 50$ ). In each minute, the mouse is rotated by 360° with an angular increment of 15°, and obtains S = 24 projections of the fluorescent measurements. Therefore, 1,200 projections (P = KS = 50 × 244) are acquired for the entire DFMT measurement. The system configurations of the DFMT measurements for the in vivo experiments are identical to those for the phantom experiments, as described in section 3.2.1. All procedures involving animals were in accordance with the ethical standards of the Ethics Committee of Tsinghua University.

The XCT images are collected after the fluorescence data acquisition is finished. An XCT contrast agent for hepatobiliary, Fenestra LC (Advanced Research Technologies, Montreal, Canada) is injected at a dose of 15 mL/kg body weight through the tail vein. One hour after the injection, the XCT images are collected to provide structural prior information. The X-ray tube works at 45 kVp and 1 mA during the scan, and the XCT images are collected by a complementary metal oxide semiconductor (CMOS) flat-panel detector (C7921-02, Hamamatsu, Japan).

Figure 12(a) shows the coronal XCT image of the chest region of the mouse with a height of 2.5 cm. The XCT images are manually segmented into five regions: heart, liver, lungs, kidneys and the background. The cross-sectional images in the kidney and lung regions indicated by the red solid and yellow dashed lines in Fig. 12(a) are shown in Figs. 12(b) and 12(c), and the corresponding segmentation results are depicted in Figs. 12(d)-12(e). A heterogeneous mouse model is created based on the segmentation results by assigning different optical properties to relevant regions, and a Laplacian-type prior matrix [26] is constructed to constrain the reconstruction results. The optical coefficients used for different organs are listed in Table 8 [32, 33]. In the mouse experiments, totally $M = 12, 968$ measurements are acquired, and the reconstructed 3D volume is discretized into $N = 6, 062$ voxels.

Fig. 12 XCT results of the mouse experiments. (a) Coronal XCT image in the chest region of the mouse. (b) and (c) Transversal XCT images indicated by the red solid and yellow dashed lines in (a), respectively. (d) and (e) Manually segmented organs corresponding to (b) and (c).

Download Full Size | PDF

Table 8. Setups of the optical properties in different regions

View Table | View all tables in this article

3.3.2 Results

50 regularization parameters are obtained with the L-curve ( $λ_{L}$ ) and U-curve ( $λ_{U}$ ) methods in the in vivo experiments. The means and standard deviations of the parameters are listed in Table 9. The mean of $λ_{L}$ (7.83 × 10⁶) is about 10⁸ times larger than that of $λ_{U}$ (9.02 × 10⁻²).

Table 9. Means and standard deviations of the regularization parameters in the in vivo experiments

View Table | View all tables in this article

Figures 13(a) and 13(b) show the L-curve and U-curve at the time point of $t_{k} = 10$ min in the in vivo experiments, and both methods show the corners. The DPC plot is illustrated in Fig. 13(c), where $λ_{L}$ and $λ_{U}$ are shown as black dashed and red solid horizontal lines, respectively. Fourier coefficients $| u_{i}^{T} b |$ (green crosses) averagely decay to zero faster than the corresponding singular values (blue dots) above 10⁻¹, which is very close to $λ_{U}$ . By contrast, $λ_{L}$ is about 10⁸ times larger than 10⁻¹.

Fig. 13 (a) L-curve plots on log-log scale (corner at λ_L = 7.73 × 10⁶). (b) U-curve plots on log-log scale (corner at λ_U = 7.88 × 10⁻²). (c) The DPC plot in the in vivo experiments, where i denotes the number of the singular values. The Fourier coefficients, the singular values and their corresponding quotients are depicted as green crosses, blue dots and red circles, respectively. The L-curve and U-curve coefficients are shown as black dashed and red solid horizontal lines, respectively.

Download Full Size | PDF

Figures 14 shows the cross-sectional images of the pharmacokinetic parameters (A, B, α, and β) in the liver region, which is indicated by the red solid line in Fig. 12(a). The shapes of the organs are well reconstructed in all parametric images, because the results are constrained by the Laplacian-type prior matrix. The images obtained with the L-curve method are over smooth inside the organs, due to the large regularization parameters, as shown in Figs. 14(a)-14(d). However, differences can be observed inside the organs in the parametric images obtained with the U-curve method, as shown in Figs. 14(e)-14(h). This may be beneficial for the detection of pathological tissues inside the organs, for example, hepatic hemangiomas inside the liver.

Fig. 14 Cross-sectional images of the parametric images in the liver region of the mouse corresponding to the slice shown in Fig. 11(b). (a)-(d) Parametric images reconstructed using the regularization parameters obtained with the L-curve method. (e)-(h) Parametric images reconstructed using the regularization parameters obtained with the U-curve method.

Download Full Size | PDF

Table 10 lists the average pharmacokinetic parameters and the CCs of the liver region in the in vivo experiments, and Fig. 15 depicts the corresponding metabolic curves. Because the true pharmacokinetic parameters of the liver are difficult to obtain in the in vivo experiments, the parameters estimated using a fiber system [2] are used as the reference. The reconstructed parameters A and B are normalized to the mean value of A [31]. The pharmacokinetic parameters and the metabolic curves obtained with the U-curve method match better with those of the reference. It indicates higher reconstruction accuracy of the U-curve method.

Table 10. Pharmacokinetic parameters and correlation coefficients in the in vivo experiments

View Table | View all tables in this article

Fig. 15 Metabolic curves in the in vivo experiments. The true curve and the curves obtained with the L-curve and U-curve methods are depicted as green solid, blue dotted and red dashed lines, respectively.

Download Full Size | PDF

4. Discussion and conclusion

DFMT is a promising technique to study the pharmacokinetic property of the fluorescent agents in biological body, and the accuracy of the parametric images rely on the high-quality reconstruction of each FMT image. For typical DFMT implementations, the imaged object is monitored for more than 50 minute. During each minute, a set of fluorescent measurements is acquired, and the corresponding FMT image is reconstructed. Therefore, it is essential to find a proper method for automatic selection of the regularization parameters. For Tikhonov regularization, the L-curve method is generally employed to determine the regularization parameters. However, it may fail to find an appropriate λ when the point of the maximum curvature is difficult to detect [15, 16], and the reconstructed images may tend to be over-smooth in some cases [17, 19].

In previous studies, the U-curve method has been proposed and tested in some numerical applications [18, 19] and FDOT problems [20]. In this study, the feasibility of the U-curve criterion for DFMT problems is investigated, and the performances are compared with those of the L-curve method. The U-curve method has four major advantages.

Firstly, the U-curve method can obtain higher imaging resolution. Generally, a larger regularization parameter reduces the resolution of the reconstructed images, and a smaller parameter increases the noise of the images. Therefore, a proper regularization parameter should neither be too large to affect the imaging resolution, nor be too small to fail in satisfying the DPC. The DPC plot is useful for objective evaluation of the qualities of the regularization parameters. As depicted in Figs. 3, 8 and 13, the L-curve method tends to obtain large regularization parameters. It explains the poor resolution and over-smoothness of the images obtained with the L-curve method. By comparison, the U-curve regularization parameters are always located in a suitable range, and obtain better imaging resolution (Figs. 4(c) and 9(c)). In the in vivo experiments, differences can be observed inside the organs of the DFMT images reconstructed with the U-curve method, as shown in Figs. 14(e)-14(h), which may be beneficial for the detection of pathological tissues inside the organs.

Secondly, the U-curve method achieves higher accuracy. The errors of the parametric images obtained with the U-curve method are smaller than those with the L-curve method, as demonstrated by the NRMSEs in the simulations and phantom experiments (Tables 2 and 6). Moreover, The pharmacokinetic parameters obtained with the U-curve method are closer to the true (or reference) parameters (Tables 4, 7 and 10), and the corresponding metabolic curves match better with the true (or reference) curves (Figs. 7, 11 and 15) in all the simulations, phantom experiments and in vivo experiments. It is worth mentioning that, because the pharmacokinetic parameters A and B have arbitrary units (a.u.) without knowing the gain of the FMT system, the two parameters are normalized to the mean value of A. Therefore, the differences between these two methods are less significant for the parameter A.

Thirdly, the U-curve method shows better noise-resistance. As depicted by the simulation results, the L-curve method cannot distinguish the two tubes at the SNR of 30 dB, while the images obtained with the U-curve method are still acceptable. The U-curve method performs better than the L-curve method at high noise levels.

Finally, the U-curve method shows better robustness. The L-curve method searches for the point where the curvature reaches its maximum, or the point where the curve is the closest to the origin. As a consequence, it may fail to determine the corner when the log-log plot is over-flat. This phenomena is observed in the results of the phantom experiments, where the L-curve method cannot find appropriate regularization parameters (Fig. 8(a)) or reconstruct reasonable DFMT images (Fig. 9(b)). By contrast, it has been proved [18] that the U-curve method always has a local minimum in the interval $λ \in (σ_{n}^{2 / 3}, σ_{1}^{2 / 3})$ . This has been verified in this paper. The corner of the U-curve can be easily determined in all the simulations, phantom experiments and in vivo experiments (Figs. 3(b), 8(b) and 13(b)), and the corresponding solutions are satisfactory.

In conclusion, compared with the L-curve method, the U-curve method obtains higher resolution of the reconstruction images, better accuracy of the pharmacokinetic parameters and better robustness. It could be a proper choice for automatic selection of the regularization parameters for the DFMT problems.

Funding

This work is supported by the National Natural Science Foundation of China under Grant Nos. 81227901, 81271617, 61322101 and 61361160418; and the National Major Scientific Instrument and Equipment Development Project under Grant No. 2011YQ030114.

References and links

1. V. Ntziachristos, J. Ripoll, L. V. Wang, and R. Weissleder, “Looking and listening to light: the evolution of whole-body photonic imaging,” Nat. Biotechnol. 23(3), 313–320 (2005). [CrossRef] [PubMed]

2. G. Zhang, F. Liu, B. Zhang, Y. He, J. Luo, and J. Bai, “Imaging of pharmacokinetic rates of indocyanine green in mouse liver with a hybrid fluorescence molecular tomography/x-ray computed tomography system,” J. Biomed. Opt. 18(4), 040505 (2013). [CrossRef] [PubMed]

3. G. Zhang, F. Liu, H. Pu, W. He, J. Luo, and J. Bai, “A Direct Method With Structural Priors for Imaging Pharmacokinetic Parameters in Dynamic Fluorescence Molecular Tomography,” IEEE Trans. Biomed. Eng. 61(3), 986–990 (2014). [CrossRef] [PubMed]

4. A. B. Milstein, K. J. Webb, and C. A. Bouman, “Estimation of kinetic model parameters in fluorescence optical diffusion tomography,” J. Opt. Soc. Am. A 22(7), 1357–1368 (2005). [CrossRef] [PubMed]

5. B. Alacam, B. Yazici, X. Intes, S. Nioka, and B. Chance, “Pharmacokinetic-rate images of indocyanine green for breast tumors using near-infrared optical methods,” Phys. Med. Biol. 53(4), 837–859 (2008). [CrossRef] [PubMed]

6. B. Alacam and B. Yazici, “Direct reconstruction of pharmacokinetic-rate images of optical fluorophores from NIR measurements,” IEEE Trans. Med. Imaging 28(9), 1337–1353 (2009). [CrossRef] [PubMed]

7. S. Vatankhah, V. E. Ardestani, and R. A. Renaut, “Application of the chi(2) principle and unbiased predictive risk estimator for determining the regularization parameter in 3-D focusing gravity inversion,” Geophys. J. Int. 200(1), 265–277 (2014). [CrossRef]

8. C. D. Brooks and P. K. Lamm, “A discrepancy principle for generalized local regularization of linear inverse problems,” J. Inverse Ill-Posed Probl. 22(1), 95–119 (2014). [CrossRef]

9. D. Peykov and R. W. Paynter, “On the choice of tuning parameters for use with Robust GCV, Modified GCV and the Discrepancy Principle in the inversion of ARXPS data,” J. Electron. Spectrosc. 197, 93–101 (2014). [CrossRef]

10. G. H. Golub, M. Heath, and G. Wahba, “Generalized Cross-Validation as a Method for Choosing a Good Ridge Parameter,” Technometrics 21(2), 215–223 (1979). [CrossRef]

11. M. A. Lukas, “Strong robust generalized cross-validation for choosing the regularization parameter,” Inverse Probl. 24(3), 034006 (2008). [CrossRef]

12. P. C. Hansen and D. P. Oleary, “The Use of the L-Curve in the Regularization of Discrete Iii-Posed Problems,” SIAM J. Sci. Comput. 14(6), 1487–1503 (1993). [CrossRef]

13. P. C. Hansen, “Analysis of Discrete Ill-Posed Problems by Means of the L-Curve,” SIAM Rev. 34(4), 561–580 (1992). [CrossRef]

14. L. M. Wu, “A parameter choice method for Tikhonov regularization,” Electron. Trans. Numer. Anal. 16, 107–128 (2003).

15. T. Correia, A. Gibson, M. Schweiger, and J. Hebden, “Selection of regularization parameter for optical topography,” J. Biomed. Opt. 14(3), 034044 (2009). [CrossRef] [PubMed]

16. M. Hanke, “Limitations of the L-curve method in ill-posed problems,” BIT 36(2), 287–301 (1996). [CrossRef]

17. J. P. Culver, R. Choe, M. J. Holboke, L. Zubkov, T. Durduran, A. Slemp, V. Ntziachristos, B. Chance, and A. G. Yodh, “Three-dimensional diffuse optical tomography in the parallel plane transmission geometry: evaluation of a hybrid frequency domain/continuous wave clinical system for breast imaging,” Med. Phys. 30(2), 235–247 (2003). [CrossRef] [PubMed]

18. D. Krawczyk-Stando and M. Rudnicki, “Regularization parameter selection in discrete ill-posed problems - The use of the U-curve,” Int. J. Appl. Math Comp. 17(2), 157–164 (2007).

19. D. Krawczyk-Stando and M. Rudnicki, “The Use of L-Curve and U-Curve in Inverse Electromagnetic Modelling,” Stud. Comp. Intell. 119, 73–82 (2008).

20. J. Chamorro-Servent, J. Aguirre, J. Ripoll, J. J. Vaquero, and M. Desco, “Feasibility of U-curve method to select the regularization parameter for fluorescence diffuse optical tomography in phantom and small animal studies,” Opt. Express 19(12), 11490–11506 (2011). [CrossRef] [PubMed]

21. X. Guo, X. Liu, X. Wang, F. Tian, F. Liu, B. Zhang, G. Hu, and J. Bai, “A combined fluorescence and microcomputed tomography system for small animal imaging,” IEEE Trans. Biomed. Eng. 57(12), 2876–2883 (2010). [CrossRef] [PubMed]

22. S. R. Arridge, “Optical tomography in medical imaging,” Inverse Probl. 15(2), R41–R93 (1999). [CrossRef]

23. G. Zhang, X. Cao, B. Zhang, F. Liu, J. Luo, and J. Bai, “MAP estimation with structural priors for fluorescence molecular tomography,” Phys. Med. Biol. 58(2), 351–372 (2013). [CrossRef] [PubMed]

24. A. Soubret, J. Ripoll, and V. Ntziachristos, “Accuracy of fluorescent tomography in the presence of heterogeneities: study of the normalized Born ratio,” IEEE Trans. Med. Imaging 24(10), 1377–1386 (2005). [CrossRef] [PubMed]

25. D. Hyde, E. Miller, D. H. Brooks, and V. Ntziachristos, “A statistical approach to inverting the Born ratio,” IEEE Trans. Med. Imaging 26(7), 893–905 (2007). [CrossRef] [PubMed]

26. S. C. Davis, H. Dehghani, J. Wang, S. Jiang, B. W. Pogue, and K. D. Paulsen, “Image-guided diffuse optical fluorescence tomography implemented with Laplacian-type regularization,” Opt. Express 15(7), 4066–4082 (2007). [CrossRef] [PubMed]

27. H. Shinohara, A. Tanaka, T. Kitai, N. Yanabu, T. Inomoto, S. Satoh, E. Hatano, Y. Yamaoka, and K. Hirao, “Direct measurement of hepatic indocyanine green clearance with near-infrared spectroscopy: separate evaluation of uptake and removal,” Hepatology 23(1), 137–144 (1996). [CrossRef] [PubMed]

28. M. Gurfinkel, A. B. Thompson, W. Ralston, T. L. Troy, A. L. Moore, T. A. Moore, J. D. Gust, D. Tatman, J. S. Reynolds, B. Muggenburg, K. Nikula, R. Pandey, R. H. Mayer, D. J. Hawrysz, and E. M. Sevick-Muraca, “Pharmacokinetics of ICG and HPPH-car for the detection of normal and tumor tissue using fluorescence, near-infrared reflectance imaging: a case study,” Photochem. Photobiol. 72(1), 94–102 (2000). [CrossRef] [PubMed]

29. P. C. Hansen, “Regularization tools version 4.0 for matlab 7.3,” Numer. Algorithms 46(2), 189–194 (2007). [CrossRef]

30. P. C. Hansen, “The Discrete Picard Condition for Discrete Ill-Posed Problems,” BIT 30(4), 658–672 (1990). [CrossRef]

31. G. Zhang, H. Pu, W. He, F. Liu, J. Luo, and J. Bai, “Full-direct method for imaging pharmacokinetic parameters in dynamic fluorescence molecular tomography,” Appl. Phys. Lett. 106(8), 081110 (2015). [CrossRef]

32. M. J. Niedre, G. M. Turner, and V. Ntziachristos, “Time-resolved imaging of optical coefficients through murine chest cavities,” J. Biomed. Opt. 11(6), 064017 (2006). [CrossRef] [PubMed]

33. D. Hyde, R. Schulz, D. Brooks, E. Miller, and V. Ntziachristos, “Performance dependence of hybrid x-ray computed tomography/fluorescence molecular tomography on the optical forward problem,” J. Opt. Soc. Am. A 26(4), 919–923 (2009). [CrossRef] [PubMed]

Regions	A(a.u)	B(a.u)	α(min⁻¹)	β(min⁻¹)
Tube 1	1.0	1.0	0.435	0.011
Tube 2	0.8	0.8	0.296	0.020

	L-curve		U-curve
	Mean	Standard deviation	Mean	Standard deviation
40 dB	1.15 × 10²	2.66 × 10¹	2.13 × 10⁻²	7.67 × 10⁻⁴
30 dB	4.23 × 10²	3.61 × 10²	2.13 × 10⁻²	7.00 × 10⁻⁴
20 dB	2.44 × 10³	4.16 × 10²	2.13 × 10⁻²	7.20 × 10⁻⁴

	L-curve					U-curve
	A	B	α	β	Mean	A	B	α	β	Mean
40dB	0.98	0.99	1.96	2.03	1.49	0.77	0.77	1.67	1.65	1.22
30dB	1.01	1.05	1.92	2.55	1.63	0.75	0.80	1.91	1.75	1.30
20dB	1.73	2.13	2.83	3.48	2.54	0.86	1.03	2.49	2.10	1.62

	Tube 1					Tube 2
	A (a.u)	B (a.u)	α (min-1)	β (min-1)	CC	A (a.u)	B (a.u)	α (min-1)	β (min-1)	CC
True	1.000	1.000	0.435	0.011	/	0.800	0.800	0.296	0.020	/
L-curve	1.000	0.997	0.394	0.013	0.998	0.777	0.772	0.261	0.024	0.998
U-curve	1.000	1.006	0.447	0.011	0.999	0.813	0.815	0.291	0.020	0.999

L-curve						U-curve
A	B	α	β	Mean		A	B	α	β	Mean
1.980	2.058	4.488	1.866	2.598		1.244	1.252	2.392	1.318	1.552

Automatic selection of regularization parameters for dynamic fluorescence molecular tomography: a comparison of L-curve and U-curve methods

Abstract

1. Introduction

2. Methods

2.1 Data acquisition system

2.2 Reconstruction method for DFMT problems

2.2.1 FMT images reconstruction

2.2.2 Two-compartment model and curve fitting

2.3 L-curve method

2.4 U-curve method

2.5 Discrete Picard condition

3. Simulations, experiments and results

3.1 Simulations

3.1.1 Setups

3.1.2 Results

3.2 Phantom experiments

3.2.1 Setups

3.2.2 Results

3.3 In vivo experiments

3.3.1 Setups

3.3.2 Results

4. Discussion and conclusion

Funding

References and links

Cited By

Figures (15)

Tables (10)

Equations (13)

Biomedical Optics Express

	Mean	Standard deviation
L-curve	3.62 × 10¹	7.83 × 10⁻¹
U-curve	9.83 × 10⁻²	3.22 × 10⁻³

Regions	Heart	Liver	Lungs	Kidneys	Background
$μ_{a} (c m^{- 1})$	0.350	0.500	0.250	0.175	0.300
$μ_{s}^{'} (c m^{- 1})$	23	13	30	20	10

	Mean	Standard deviation
L-curve	7.83 × 10⁶	2.09 × 10⁶
U-curve	9.02 × 10⁻²	8.09 × 10⁻³

Regions	A(a.u)	B(a.u)	α(min-1)	β(min-1)	CC
Reference	0.926	1.171	0.260	0.012	/
L-curve	1.000	1.282	0.303	0.019	0.902
U-curve	1.000	1.175	0.226	0.019	0.976