A Novel Adaptive Parameter Search Elastic Net Method for Fluorescent Molecular Tomography

Fluorescence molecular tomography (FMT) is a new type of medical imaging technology that can quantitatively reconstruct the three-dimensional distribution of fluorescent probes in vivo. Traditional Lp norm regularization techniques used in FMT reconstruction often face problems such as over-sparseness, over-smoothness, spatial discontinuity, and poor robustness. To address these problems, this paper proposes an adaptive parameter search elastic net (APSEN) method that is based on elastic net regularization, using weight parameters to combine the L1 and L2 norms. For the selection of elastic net weight parameters, this approach introduces the L0 norm of valid reconstruction results and the L2 norm of the residual vector, which are used to adjust the weight parameters adaptively. To verify the proposed method, a series of numerical simulation experiments were performed using digital mice with tumors as experimental subjects, and in vivo experiments of liver tumors were also conducted. The results showed that, compared with the state-of-the-art methods with different light source sizes or distances, Gaussian noise of 5%–25%, and the brute-force parameter search method, the APSEN method has better location accuracy, spatial resolution, fluorescence yield recovery ability, morphological characteristics, and robustness. Furthermore, the in vivo experiments demonstrated the applicability of APSEN for FMT.

tumors were also conducted. The results showed that, compared with the state-of-the-art methods with different light source sizes or distances, Gaussian noise of 5%-25%, and the brute-force parameter search method, the APSEN method has better location accuracy, spatial resolution, fluorescence yield recovery ability, morphological characteristics, and robustness. Furthermore, the in vivo experiments demonstrated the applicability of APSEN for FMT. Index Terms-Fluorescence molecular tomography, adaptive parameter search, elastic net.

I. INTRODUCTION
F LUORESCENCE molecular tomography (FMT) is a new type of medical imaging technology that can quantitatively reconstruct the three-dimensional (3-D) distributions of fluorescent probes in vivo. It also addresses the problem of insufficient depth resolution, which is caused by the lack of depth information in traditional fluorescence imaging methods [1]- [4]. FMT has been applied to preclinical and preliminarily clinical applications in cancer diagnosis and treatment [5], [6].
In recent years, research on FMT has mainly been focused on improving the 3-D reconstruction quality, especially the location accuracy and morphological recovery ability [7], [8]. Generally, FMT reconstruction requires the establishment of a correspondence between the known surface fluorescence distribution of the imaged object and the unknown distribution of the fluorescent probes inside it; then, the unknown information is obtained using a specific solution method. However, FMT reconstruction is highly ill-posed due to the high scattering effect of fluorescence in biological tissues [9]- [11]. To solve this problem, researchers have developed many optimization methods, such as Lp norm regularization (p = (0, 2]), in which the Lp norm of the unknown fluorescent probe is used to constrain the FMT reconstruction [12]- [14]. However, different p-values suit different reconstruction situations. A large p-value (for example, p = 2) is more appropriate for morphological reconstruction; in this case, the reconstructed region includes most of the ground truth but also introduces many reconstructed artifacts, causing over-smoothness [15]. Small p-values (p ≤ 1) are usable for an accurate reconstruction of part of the ground truth without the reconstructed artifacts. However, there are some difficulties, such as oversparseness, incomplete reconstruction of the fluorescent probe, and the lack of detailed information on the boundary [16]. Other available regularization methods are total variation norm regularization [17], group sparse norm regularization [18], etc. These methods add other prior information to the original Lp This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ norm regularization method, such as the aggregation of fluorescent probe distribution and connectivity within subregions. These methods have achieved relatively satisfied results, but they have high computational complexity and are not very accurate.
On the other hand, the reconstruction effect of the regularization method is usually related to the choice of the regularization parameter λ. Generally, when the regularization parameter is small, there are many artifacts in the reconstructed image owing to the noise amplification effect [18]. When the regularization parameter is large, the reconstructed image will be over-smoothed or over-sparse. Researchers have proposed several methods for selecting the appropriate regularization parameters, such as the L-curve, U-curve, cross-validation, and discrepancy principle methods [19]- [22]. CR. Vogel noted out that the choice of regularization parameters should be related to the specific inverse problem [23].
Recently, elastic net (EN) regularization has been proposed to address the problem of over-sparseness or over-smoothness of traditional regularization methods [24]. The method uses the weight parameters to combine the L1 norm and L2 norm, which improves the reconstruction effect by balancing the weight parameters while also alleviating the discontinuity of the reconstructed tumor region. Recent studies have shown that EN possesses certain advantages when used to solve the problem of Lp norm regularization, and it has attracted wide research attention. For example, T. Rymarczyk et al. presented an effective algorithm based on multiple ENs to solve the inverse problem of the electrical impedance tomography [25]; P. Causin et al. applied the EN regularization theory of entropy to diffuse optical tomography applications [26]; T. Nguyen-Duc et al. used EN regularization for better approximation of the sparse input of high-undersamplingdynamic magnetic resonance imaging (MRI) data [27]. Some researchers have also explored the application of EN in FMT reconstruction [28]. However, choosing the right balance of weight parameters is a major concern in EN.
In this work, an adaptive parameter search elastic net (APSEN) method is proposed to optimize the FMT reconstruction problem. APSEN is based on the L1 norm method and optimizes the over-sparseness results of the L1 norm reconstruction by adding the L2 norm to improve the quality of tumor region morphological reconstruction. APSEN uses the L2 norm of the residual vector to discriminate, which can effectively yield parameters suitable for the distribution of various types of fluorescent probes. By introducing the L0 norm of valid reconstruction results and using the proportional search and coordinate descent algorithm, APSEN can rapidly and accurately reconstruct fluorescent probe distributions in vivo. Based on the results of the previous step during the iteration of the solution, the APSEN method adaptively constrains the selection of the parameters for the next step.
To evaluate the performance of the APSEN method, we performed numerical simulation experiments and in vivo experiments for liver tumors. Comparisons were made to the L1-based iterative shrinkage (IS-L1) [29], the L2-based Tikhonov (Tikhonov-L2) [15], the fused least absolute shrinkage and selection operator (LASSO) method (FLM) [18], and Nesterov's method-based elastic net (N-EN) [28]. The results demonstrate that APSEN significantly improves the reconstruction accuracy and shape recovery of the fluorescent probe distribution compared with the above methods.
The structure of this paper is organized as follows. Section II introduces the forward model of FMT and the reconstruction algorithm based on APSEN. Section III presents the process and results of numerical simulations and in vivo experiments. Section IV discusses the proposed APSEN method and summarizes the conclusions.

A. Forward Model of FMT
For a steady-state FMT that has the point excitation sources, photon propagation in biological tissues can be described by a pair of coupled diffusion equations with the Robin-type boundary condition [30]. This can be expressed as follows where r denotes the nodes inside the problem domain , and r l denotes the positions of point excitation sources with the amplitude ; r l is placed on one mean free-path of photon transport beneath the surface of . Subscripts x and m denote the excitation and emission wavelengths, respectively.
x,m (r) denotes the photon flux density at the position r inside the region . ημ a f (r) is the fluorescent source to be reconstructed, where η is the quantum efficiency. D x,m = 1/3(μ ax,am + (1-g)μ sx,sm ) denotes the diffusion coefficient, where g is the anisotropy parameter; μ ax,am and μ sx,sm are the absorption and scattering coefficients, respectively. q denotes the optical reflective index of the biological tissues.
Based on the finite element method [31], we can discretize the light propagation model and linearize the partial differential equation (1) into the following linear equation where Y = [y 1 , . . . , y N ] ∈ R N denotes the photon segment on the surface of the measured object. A = f T 1 , . . . , f T N ∈ R N× p denotes the system weight matrix, where f l = a 1 , . . . , a p ∈ R p is the feature vector, and X = x 1 , . . . , x p ∈ R p is the fluorescence intensity distribution in biological tissues. Therefore, restoring the fluorescence intensity distribution X in the above linear matrix equation is the purpose of solving the FMT inverse problem. The details of the forward model can be found in [31].

B. Inverse Problem of FMT Reconstruction
As described above, the linear relationship between the fluorescent sources inside the tissue and the surface photon density has been established. Due to the inverse problem being ill-posed, the unique solution is absent. To solve the inverse problem based on (2), we propose an APSEN method.
The objective function in traditional FMT reconstruction can be expressed as the sum of the least squares and regularization terms [12], [29]. This is defined as where E(X) denotes the objective function. α and β denote the regularization parameters. This combination of different norms, which is called EN regularization, can combine the advantage of each regularization and find a balance between the sparseness and smoothness of the image. If β = 1, equation (3) becomes the L1 norm regular expression (lasso regularization); if β = 0, equation (3) will become the L2 norm regular expression (ridge regression). We use the coordinate descent algorithm [32] to solve the objective equation (3), which uses the EN penalty term. A detailed procedure is introduced in Section II.B.1. Obviously, multiple regularization parameters need to be determined to obtain the optimized solution from equation (3). In this work, we introduce the L0 norm of the valid reconstruction result and the L2 norm of the residual vector as the standard and propose a novel APSEN method that adaptively adjusts the regularization parameters. Detailed information regarding the method is introduced in Section II.B.2. The flow chart of the APSEN method is shown in Fig. 1. 1) EN Based on the Coordinate Descent Algorithm: Because of the complexity of the regularization term in EN regularization, we used a fast and accurate iterative algorithm based on the coordinate descent algorithm to solve equation (3). The EN regularization is summarized in Algorithm 1, and the details of the method are explained below.
Because of the complex regularization term, we used the coordinate descent method to solve equation (3) and simplify the operation process. We have N values of y i , where y i denotes the value of the i -th term in the measured emission light distribution Y . The derivative of equation (3) is where x l denotes the value of the l-th term in the fluorescence intensity distribution X, a il denotes the i -th term of the lth feature vector in the system weight matrix A, and x k l denotes x l in the k-th iteration. To simplify (4), we let Through w s , we find that whereỹ a i j x k j is the fitted value excluding the contribution from x k l , andŷ i is the current fit of the model for the i -th observation; hence, r i = y i −ŷ i is the current residue. Therefore, we can rewrite w s as where the first term on the right-hand side in equation (6) is the gradient of the loss with respect to x k l . From equation (6), it can be found that the coordinate reduction is computationally effective. For more efficient calculation, the first term on the right-hand side in equation (6) can be written as where a l , Y = N i=1 a il y i . Thus, the ideal way to perform the coordinate descent method is to use the sparseness of the solution to improve the calculation efficiency by summing the non-zero terms. Furthermore, in this case, scaling the variable will not change the sparsity.

Algorithm 1 EN
Input: matrix A, measured emission light distribution Y , regularization parameters α and β Initialization: fluorescence intensity distribution X 0 = ∅, iterator number k = 0, tolerance for the optimization tol, maximum number of iterations I for the stopping criterion.
a l , a t x k t End while End while Output: X = X * To solve equation (4), we can use the soft-thresholding operator proposed by Friedman et al. [33] 1 Therefore, when equation (8) equals zero, the updated value x k+1 l on the corresponding gradient direction is achieved. According to equation (8), the iterative result x k+1 l is as follows Therefore, we finally acquire the value of the iterative result x k+1 l . If the change of another complete cycle in X is less than tol or the number of iterations exceeds I , the operation is completed, otherwise, the process repeats.
2) Adaptive Parameter Search EN: The FMT reconstruction based on (3) creates a problem in that the parameters are difficult to select. To optimize the parameter selection problem based on (3), we proposed APSEN which can adaptively adjust the regularization parameters to improve the parameter selection and final results of FMT reconstruction in the region of interest (ROI). APSEN is summarized in Algorithm 2.
It can be seen from (3) that EN has two regularization parameters, and the two-dimensional (2-D) parameter search process is complicated and computationally intensive. After determining that the parameter αβ of the L1 norm term in (3) is unchanged, we increase the proportion of the L2 norm term of the objective function by changing β. In the iterative search, we introduce the L0 norm of the valid reconstruction results. The L0 norm of valid reconstruction results reflects the degree of reconstruction of the corresponding parameters and changes with the parameters. Equation (3) shows that when β decreases, the L2 norm proportion increases, the reconstruction result X tends to be smooth, and the number of valid reconstruction results increases accordingly. Then, we use the residual vector as the criterion for judging the reconstruction result and select the solution with the smallest L2 norm of the residual vector in the search path as the optimal solution with the current path. The core idea of APSEN is to add the L2 norm term based on the L1 norm term solution, use the L0 norm of valid reconstruction results as a new parameter search criterion to constrain the direction of parameter search and perform a fast search, and simultaneously use the residual vector as a judgment standard to select the optimal parameter in each epoch search path. The expression for the L0 norm of valid reconstruction results is as follows where EN(•) is the EN regularization (Algorithm 1) which is mentioned in Section II.B.1 with a matrix A and a measured emission light distribution Y as inputs. Further, N † denotes the number of positive items obtained by running the EN regularization when the parameters are α † and β † .
Referring to the previous work, the FMT is reconstructed using L1 norm regularization to determine parameter α * , which ensures the accuracy of the reconstructed spatial information of APSEN and obtains the target value N * as the starting point of the search by (10). APSEN converts the 2-D parameter search into a one-dimensional parameter search by determining parameter αβ = α * and then adjusting β. We want to control the accuracy of the search by setting the search step size ratio R, which controls the end of the search path for each epoch. The parameter acquisition method used in APSEN is introduced in Section III.B.1. The L0 norm of valid reconstruction results at the end of each search epoch can be computed by using where N s denotes the search target number of each search epoch and N * denotes the number of valid reconstruction results of the corresponding parameter α * and β * . APSEN uses an adaptive objective function parameter search ratio N s /N m−1 in the iteration. N s /N m−1 is the ratio of the target value of the L0 norm needed to search the L0 norm of the valid reconstruction result obtained in m − 1 iterations, and parameters can be adaptively adjusted based on the result of each iteration. If N s /N m−1 > 1, β value will decrease and the proportion of the L2 norm term in the objective function will increase. Meanwhile N s /N m−1 < 1, β value will increase and the proportion of the L2 norm term in the objective function will decrease. This ratio simultaneously constrains the search direction to reach N s quickly. Hence, the search method is as follows where s m denotes the search step size of β at step m in each search epoch. At the beginning of the search, we set the initial search step size s 0 . Using the search method described above, we determine the weight of the L1 norm term and increase the weight of the L2 norm term to improve the final image smoothing effect. To eliminate the possibility of fluctuations at the end of the ratio search algorithm, which will cause the search results to not converge, we have added interruption conditions I * .
In each epoch of the search, APSEN searches for the target value N s according to a predetermined step size and calculates the L2 norm size of the corresponding residual vector in the determined ROI. The ROI is selected as follows where d denotes the ROI parameter and j denotes the maximum number of eligible xi in X. We define the ROI to distinguish the reconstruction region from the background region, and we calculate and judge whether the corresponding light intensity value of each node is greater than d times (d < 1) the maximum light intensity value of the corresponding solution.
If it is greater than d times the maximum light intensity value, we consider the point in the ROI. The corresponding residual vector is given by where r denotes the residual vector, ROI(•) denotes the corresponding calculation result of equation (13), and A R O I denotes the A matrix corresponding to the ROI. APSEN selects the solution corresponding to the smallest value with r 2 in each epoch as the optimal solution in the current search path. If the r 2 in the next epoch of search does not change or reaches an acceptable range r 2 < ε, the algorithm is terminated and the optimal solution X * is output.
In this study, we set the search step size ratio R = 3 and s * = 0.0001, and the ROI parameter d = 0.03. In most instances of this study, the optimal reconstruction results can be obtained within 10 iterations. When I * > 10, the subsequent iterative calculation only produces a slight influence on the reconstruction result, which is time-consuming. Finally, we utilized I * = 10 as the maximum number of iterations. For the selection of the parameter n, we basically adopt the same parameter selection as the comparison method, as shown in Section III.B.3. As described in the algorithm 2, different n can lead to different L0 norm of the initial vector. However, according to the algorithm debugging, under the condition of controlling other parameters unchanged, the influence of different n on reconstruction accuracy is within 10%.

III. EXPERIMENTS AND RESULTS
This section describes the numerical simulations and in vivo experiments are designed to verify the reconstruction performance of APSEN. Specifically, the location accuracy, spatial resolution, fluorescence yield recovery ability, morphological recovery ability, robustness, and practicality of APSEN were evaluated. This section shows the parameter acquisition in the APSEN, the results of comparison with the brute-force parameter search method and the results of four experiments, which includes dual-sources reconstruction, multi-size single-source reconstruction, anti-noise ability experiments, and in vivo reconstruction experiments. All the MATLAB and Python programs were run on a desktop computer with an Intel(R) Core(TM) i7-6700 CPU (3.40 GHz) and 16GB RAM.   the mouse model are the heart, liver, lungs, kidneys, and muscles. We used spheres with different diameters and positions to simulate the fluorescence light source in a digital mouse [34]. The single-source and dual-sources are shown in Fig. 2.
The fluorescence yield was set to 0.5 mm −1 , and the power of the excitation light source was set to 0.02 W. The excitation light sources were on the same plane as the fluorescence sources. We captured the fluorescence distribution of the phantom surface opposite to the excitation light source within a 160 • field of view (FOV). To facilitate the calculations involved, Amira 5.2 (Thermo Fisher Scientific, USA) was used to scatter the selected torso section of the digital mouse model into a uniform tetrahedral mesh. The specific numbers of nodes, numbers of tetrahedral elements, and tumor center coordinates are shown in Table I below.
To verify the performance of APSEN more practically, we performed three simulation experiments and compared APSEN with IS-L1, Tikhonov-L2, FLM, and N-EN. We first designed a dual-sources simulation experiment with different edge-to-edge distances (EEDs). The radius of the dual-sources tumors was 0.5 mm, and the EED was set to 2 mm, 3 mm and 3.7 mm. The EED simulation experiment was used to evaluate Fig. 3. Configuration of multi-modality imaging tomography system for FMI/CT data acquisition. The system includes 1: EMCCD, 2: X-ray tube, 3: X-ray detector, 4: small animal gas anesthesia machine, 5: small animal bed, and 6: rotating gantry.

TABLE II OPTICAL ABSORPTION AND SCATTERING COEFFICIENTS OF BIOLOGICAL TISSUES IN NUMERICAL AND in vivo EXPERIMENTS
the positioning accuracy and spatial resolution of APSEN. We then designed simulation experiments with single sources of different sizes. We placed a single-source with a radius of 0.5 mm, 1.0 mm, or 1.5 mm at the same tumor center coordinates. Therefore, the morphological recovery ability of APSEN was evaluated by simulating single sources of different sizes. To eliminate the influence of the mesh size on the final result, a different mesh size was used in the dual-sources experiment in which EED was set to 3.7 mm compared to that used in the single-source simulation experiment and the other dual-sources experiments. Finally, we designed a simulation experiment to test the anti-noise ability of APSEN. In particular, 5%, 15%, and 25% Gaussian noise was added in the 1.0 mm single-source simulation experiment. Anti-noise ability simulation experiments were performed to evaluate the robustness of APSEN. To better calculate the fluorescence yield recovery ability more accurately, we did not normalize the data before the reconstruction procedure, and for better comparison between methods, the same ROI was used in all reconstruction methods.
2) In Vivo Experiment: In vivo experiments were performed by preparing 4-week-old BALB/c nude mice (Vital River Laboratory Animal Technology Co. Ltd., Beijing, China) to establish the tumor model. All animal protocols were approved by the Institutional Animal Care and Use Committee at Peking University (Permit Number 2011-0039), and all procedures were in accordance with approved guidelines.  The mice were intraperitoneally injected with 2% pentobarbital sodium at a dosage of 40 mg/kg. After the mice were anesthetized, we implanted Huh7-GFP-fLuc cells into their livers to establish orthotopic liver cancer models. After 12 days, three tumor-bearing mice models were successfully obtained. To evaluate the feasibility of this method in vivo, we injected 2 nmol MMPSense 750 FAST (Perkin Elmer, MA, USA) into these three mice via the tail vein 6 h before imaging. Studies have shown that matrix metalloproteinases (MMPs) can be over-expressed in various tumor tissues, so we used MMPSense to activate the MMPs to visualize the tumors. We utilized the multi-modality imaging tomography system, which was developed by the key laboratory of molecular imaging of Chinese Academy of Sciences, to image these mice [35]. As shown in Fig. 3, the multi-modality imaging tomography system mainly consists of two subsystems: a fluorescence acquisition system and a CT system. The fluorescence acquisition system consists of a 750 nm continuous wave laser and thermoelectric cooling electron multiplying charge coupled device (EMCCD) cameras (iXonem + 888, Andor, UK). The CT system is composed of an X-ray tube (UltraBright, Oxford Instrument, USA) and an X-ray detector (CMOS Flat-panel Detector, Hamamatsu, Japan).
Data acquisition for the in vivo experiments included fluorescence acquisition and structural data acquisition. First, fluorescence data were acquired using a fluorescence acquisition system. Then, CT imaging was performed on each mouse by using the CT system to obtain the structural data of the mouse. After multi-modality raw data acquisition, the mouse was imaged by MRI (1.5T, M3TM, Aspect Imaging, Israel) to provide the location of the tumor, as MRI is the gold standard for locating tumors. After the data collection was complete, we proceeded with data processing. We segmented the major organs, including the muscle, bones, lungs, heart, and liver, and integrated these organs into the mouse model. Finally, we mapped and registered the fluorescence image to the surface of the mouse torso. The structure of the in vivo experimental process is shown in Fig. 4. The optical parameters of the different organs used in the numerical simulations and in vivo experiments were obtained from [36], as shown in Table II.
3) The Evaluation Indexes: To evaluate the accuracy of FMT reconstruction quantitatively in terms of the source location as well as the fluorescence yield and shape recovery in the ROI, the location error (LE), relative intensity error (RIE) and dice similarity (Dice) were used as the quantitative indexes. LE measures the distance variation between the centers of the actual region and the reconstructed region, respectively, and is defined as (15) where L t denotes the center coordinates of the true region and L r denotes the center coordinates of the reconstructed region in the ROI. RIE measures the intensity variation between the true region and reconstructed regions and is given by where I t denotes the true fluorescence yield of sources and I r denotes the mean fluorescence yield of reconstructed sources in the ROI. If RIE is close to 0, the reconstructed result reveals a better fluorescence yield recovery. Dice is the similarity between the reconstructed and true fluorescence regions [18]. Thus, Dice is computed as where ROI (X) denotes the reconstructed result of ROI and S denotes the true fluorescence region. Note that Dice measures the similarity of the shapes of two objects and ranges from 0 to 1. If Dice is close to 1, the reconstructed result coincides well with the true region; otherwise, if Dice is close to 0, the reconstructed result reveals poor shape recovery.

1) Parameter Acquisition:
We used all the single-source phantoms to select the search step size ratio parameter in each search epoch. Using single sources with radii of 0.5 mm, 1.0 mm and 1.5 mm, respectively, we tested to select the search step size ratio R from 2 to 5. In this work, to reduce the artifact error, LE was calculated using only the reconstruction results within 3 mm of the center of the light source. The quantitative analysis of the reconstruction results is presented in Fig. 5.
The quantitative analysis results showed that all search step size ratio settings could more effectively search for better target solutions. In particular, for the quantitative analysis results of RIE, when R increases or decreases, APSEN's search path will be changed correspondingly, resulting in different β obtained from the search. When radius = 1.0 mm and R = 3, we get β around 0.4, and L0 norm obtained in APSEN iteration is closer to the number of the nodes of the true source, and the total intensity of the reconstructed light source is also close to the true light source intensity. According to the calculation formula of RIE, RIE is the smallest in this case. Similarly, when radius = 1.5 mm and R = 5, we get β around 0.6, and the RIE is the smallest in this case. But in the evaluation of the reconstruction results, we believe that the reconstruction results are better if two results which in LE, RIE and Dice are better than other parameter results in quantitative analysis. When R was set to 3, the searched optimal solution had smaller LE and RIE and larger Dice coefficient. Therefore, when R was 3, more accurate reconstruction results could be obtained. In summary, R was set to 3 in the following experiments.

2) Comparison With the Brute-Force Parameter Search
Method: To prove the validity of the adaptive parameter search method, a comparison between the adaptive and brute-force parameter search methods was performed. In the comparison experiment using the brute-force parameter search method, the same initial L1 norm parameter with APSEN was adopted, and β was divided into 10 groups from 0.1 to 1 for the reconstruction experiments. The quantitative analysis of the reconstruction results is shown in Fig. 6 (a). According to the evaluation method in the previous section, the quantitative analysis results showed that the adaptive parameter search method could effectively search for the optimal solution. For example, when β = 0.5, the reconstruction results of the 1.5 mm single-source experiment are similar to the adaptive parameter search reconstruction results: however, in the 0.5 mm and 1.0 mm single-source experiments as well as the 0.5 mm dual-sources experiment, the adaptive parameter search reconstruction results are better than the β = 0.5 reconstruction results. When β = 1, the expression of the EN regularization term changes to LASSO regularization, and the reconstruction results of the 0.5 mm single-source experiment are similar to the adaptive parameter search reconstruction results. However, in the 1.0 mm and 1.5 mm single-source experiments as well as the 0.5 mm dual-sources experiment, the adaptive parameter search reconstruction results are significantly better than the reconstruction results with β = 1.
All the comparison experiments demonstrate that the optimal parameters searched vary with the data set.
3) Dual-Sources Reconstruction: To verify the location accuracy and spatial resolution of APSEN, we set up a dual-sources reconstruction experiment. We set two sources S1 and S2 with different EEDs in the liver of the simulated mouse. The specific values of the simulated mouse and coordinates of the two source centers are shown in Table I. We utilized IS-L1, Tikhonov-L2, FLM, and N-EN for comparison. All the comparison methods and the initial L1 norm parameter of APSEN were tested according to the parameters given in the corresponding references [15], [18], [28], [29]. The initial parameter selections were: IS-L1: α = 2e −10 ; Tikhonov-L2: α = 1e −9 ; FLM: The regularization parameters of all methods in the experiment were manually optimized by brute-force parameter searches of different orders of magnitude based on the initial parameter selection [29]. Fig. 7 shows the reconstruction results obtained using different methods on the dual-sources simulated mouse. The shape on the axial plane proves the validity and accuracy of the reconstruction method. The quantitative analysis of the dual-sources simulated mouse reconstruction results is shown in Table III. Under different EEDs, APSEN has the best results in terms of the total location accuracy (S1 + S2) and Dice coefficient among the five methods, as well as better RIE results. By comparing the reconstruction results of EED = 3.7 mm and EED = 3 mm and 2 mm, we find that the accuracy of APSEN is not affected by grid density. At the same time, analyzing the reconstruction time, we find that the EED = 3.7 mm is less sparse than the grids with EED = 3 mm and 2 mm, which means that there are less unknowns to be calculated, so the calculation time for EED = 3.7 mm is the least. For EED = 3 mm and 2 mm, the grid density is close. When EED = 2 mm, the measurable surface fluorescence distribution generated by two relatively close light sources during the simulation have a considerable overlapping, so the surface fluorescence distribution area is relatively less than the case of EED = 3 mm, which indicate the number of effectively solved equations is also less than EED = 3 mm. Therefore, it takes less time for the reconstruction of EED = 2 mm than EED = 3 mm. Thus, APSEN could successfully improve the location accuracy and spatial resolution.
4) Multi-Sized Single-Source Reconstruction: To measure the morphological recovery ability of APSEN, we set up single sources of different sizes in the reconstruction experiments.
We compared the proposed APSEN to the four methods described previously in terms of shape restoration ability (Fig. 8). The 3-D reconstruction results and shape on the axial plane proved the effectiveness of the method in morphological restoration. The quantitative analysis of the single-source simulation reconstruction results is shown in Table IV. IS-L1 and FLM generated overly convergent source boundaries, and Tikhonov-L2 produced over-smoothing artifacts around the single source. N-EN and APSEN could obtain satisfactory results under different single-source conditions; nevertheless, APSEN produced more morphologically accurate results and fewer image artifacts. The 3-D mesh reconstructed by APSEN exhibited a good overlap with the actual mesh. APSEN had a smaller LE, and the reconstructed light source center was closer to the true value. RIE was also smaller for APSEN than for the other methods, suggesting that the fluorescence yield of the reconstructed light source was closer to the true value. More importantly, according to the Dice coefficient, the reconstruction area of APSEN was closer to that of the actual source area. Therefore, APSEN has good morphological recovery ability.

5) Anti-Noise Ability Experiment:
To examine the robustness of the APSEN, we set up an anti-noise ability experiment. Using the 1.0mm single-source simulation as an example, the effects of adding 5%, 15%, and 25% Gaussian noise on the reconstruction performance of the five methods were compared, as shown in Fig. 9. We performed independent simulation reconstruction for the five reconstruction methods. The simulated reconstruction source was displayed in 3-D and axial views.
Owing to the influence of the Gaussian noise, the reconstructed light source was further away from the center of the actual light source. The quantitative analysis of anti-noise ability verification results is shown in Table V. ASPEN exhibits a superior reconstruction effect, with a superior RIE value, the lowest LE value and the highest Dice value among the investigated methods. The anti-noise ability verification results demonstrate that APSEN is the most robust of these approaches.

6) In Vivo Reconstruction:
To evaluate the utility of the proposed APSEN for in vivo imaging, we performed in vivo FMT on three orthotopic liver tumor-bearing mice. According to the registration of the mouse body boundaries, the FMT images reconstructed by the five methods were compared to MRI data. As shown in Fig. 10, the 3-D reconstruction results and the shape of the axial plane prove the effectiveness of the different methods. APSEN has the smallest LE and a tumor area matching similarity that is significantly better than those of the other methods. The quantitative analysis also verified 7) Reconstruction Efficiency Analysis: To evaluate the time efficiency of the proposed APSEN, we recorded the actual running times of all methods. The running time of the EN algorithm in the adaptive and brute-force parameter search methods was analyzed as shown in Fig. 6 (b). These results prove that the adaptive parameter search method requires less calculation time than the brute-force parameter search method and has a better reconstruction effect. According to the experimental results in Tables III-VI, APSEN spent at most twice the running time to reach the optimal parameter selection compared to N-EN. However, compared to the manual optimization of the brute-force parameter search and entire FMT reconstruction process, this running time was acceptable. In particular, the APSEN could use an acceptable account of time to obtain better results in terms of reconstruction accuracy, robustness, and practicality.

IV. DISCUSSION AND CONCLUSION
FMT can quantitatively reconstruct the 3-D fluorescence distribution inside an organism and is proposed to solve the problem of insufficient depth resolution in traditional FMI methods. To improve the quality of 3-D reconstruction further, we proposed a novel APSEN method to achieve morphological FMT reconstruction with higher accuracy of tumor regions. The EN regularization combines the conventional L1 and L2 norms and reconciles the smoothness and sparsity of the reconstruction result. Owing to the complex regularization term, we used the coordinate descent method to simplify the operation process. The EN regularization simultaneously increased the problem of parameter selection. We further proposed an adaptive parameter search method that used the L0 norm of valid reconstruction results as a new parameter search criterion to constrain the direction of the parameter search and to perform a fast search and simultaneously used the residual vector as a judgment standard to select the optimal parameter in each epoch search path. This method reduced the number of parameters that needed to be adjusted.
For comparison, we used the IS-L1, Tikhonov-L2, FLM, and N-EN approaches to verify the performance of APSEN and conducted numerical simulation and in vivo experiments. The experimental results showed that 1) APSEN could guarantee a high reconstruction accuracy and spatial resolution. 2) Compared with the traditional methods, APSEN had advantages in tumor morphological reconstruction. 3) In the antinoise ability experiment, APSEN had a higher reconstruction accuracy and morphological recovery ability, than the other approaches, demonstrating its robustness. 4) The simulation results revealed that APSEN has a better fluorescence yield recovery ability than the other methods and acceptable solution time, indicating its practicality. 5) The in vivo experiments verified the significant advantages of APSEN in tumor detection, demonstrating the practicality and potential of APSEN in biomedical research. 6) APSEN has reasonable time efficiency. Although APSEN achieved better results than the other approaches considered, some concerns remain to be addressed. The starting point of the parameter search is based on the results of L1 norm regularization; therefore, the accuracy of the reconstruction is based on the selection of the starting point. Moreover, our adaptive parameter search method is currently only applicable to the EN regularization reconstruction algorithm based on the coordinate descent method. In the future, it would be beneficial to develop a more general adaptive parameter search.
In summary, we proposed an APSEN method that can provide precise morphological FMT reconstruction of tumor areas. The EN method simultaneously reconciled the smoothness and sparseness of the reconstruction results using the coordinate descent method to simplify the calculation process. The L0 norm of valid reconstruction results was used as a search criterion to constrain the search direction, and the residual vector was used as a judgment standard. The method was subsequently evaluated by performing numerical simulations and in vivo experiments. The reconstruction results showed that compared to the traditional methods, the proposed method had better results in terms of localization, spatial resolution, fluorescence yield recovery, and morphological reconstruction; it also simultaneously had better robustness and acceptable time efficiency. Future work will focus on developing a more general adaptive parameter search method, as well as the pre-clinical or clinical biomedical application of the proposed method.