SODL-IR-FISTA: sparse online dictionary learning with iterative reduction FISTA for cone-beam X-ray luminescence computed tomography

Cone beam X-ray luminescence computed tomography (CB-XLCT) is an emerging imaging technique with potential for early 3D tumor detection. However, the reconstruction challenge due to low light absorption and high scattering in tissues makes it a difficult inverse problem. In this study, the online dictionary learning (ODL) method, combined with iterative reduction FISTA (IR-FISTA), has been utilized to achieve high-quality reconstruction. Our method integrates IR-FISTA for efficient and accurate sparse coding, followed by an online stochastic approximation for dictionary updates, effectively capturing the sparse features inherent to the problem. Additionally, a re-sparse step is introduced to enhance the sparsity of the solution, making it better suited for CB-XLCT reconstruction. Numerical simulations and in vivo experiments were conducted to assess the performance of the method. The SODL-IR-FISTA achieved the smallest location error of 0.325 mm in in vivo experiments, which is 58% and 45% of the IVTCG-L1 (0.562 mm) and OMP-L0 (0.721 mm), respectively. Additionally, it has the highest DICE similarity coefficient, which is 0.748. The results demonstrate that our approach outperforms traditional methods in terms of localization precision, shape restoration, robustness, and practicality in live subjects.


Introduction
With the continuous development of X-ray excitable nanophosphors, X-ray luminescence computed tomography (XLCT) has garnered significant attention due to its promising performance [1].XLCT represents a hybrid imaging modality combining X-ray computed tomography (CT) with optical imaging [2].High-resolution X-ray structural imaging with high sensitivity and high specificity optical molecular imaging can be naturally integrated in XLCT [3].Compared to CT, XLCT is capable of acquiring molecular-level information, offering advantages in sensitivity and specificity.Compared with other optical molecular tomography (OMT) modalities, XLCT offers the following advantages.Firstly, in XLCT, the nanophosphors used can withstand multiple X-ray excitations and generate visible or near-infrared light (NIR) under X-ray irradiation.Its luminescence intensity can reach the single photon level, which is one order of magnitude higher than bioluminescence tomography (BLT) [4,5].Secondly, compared with Cerenkov luminescence tomography (CLT) [6,7], XLCT avoids the potential impact of radioactive substances.Moreover, compared to fluorescence molecular tomography (FMT) [8,9], the energy attenuation of X-rays absorbed by biological tissues is minimal, which can avoid the interference of spontaneous fluorescence and background fluorescence, achieving high resolution and high signal-to-noise ratio [10].Therefore, XLCT has become a promising imaging technology that can be applied in basic research, drug development and clinical studies [11,12].
Since the introduction of the first XLCT imaging system by Xing et al., XLCT systems have been widely studied [13].XLCT systems are categorized into pencil-beam/narrow-beam XLCT (PB/NB-XLCT) and cone-beam XLCT (CB-XLCT) based on the X-ray excitation and scanning methods [14].PB/NB-XLCT emits X-rays in a strip shape and can achieve higher spatial resolution through translational scanning, but its data collection time is long and is not suitable for rapid biomedical applications [15].CB-XLCT emits X-rays in a cone shape, eliminating the translational scanning process required for PB-XLCT and NB-XLCT, and significantly speeding up the imaging process at the expense of spatial resolution.Therefore, CB-XLCT is more suitable for fast biomedical applications, such as small animal whole-body imaging [16].
However, the CB-XLCT reconstruction is an ill-posed problem due to the complex excitation process and the high scattering of light during propagation in biological tissues [17].To address these issues and enhance reconstruction quality, several methods have been proposed by researchers in recent years.For instance, some scholars have attempted to employ regularization methods to address this problem.Tzoumas et al. used L 1 -norm as convex relaxation of L 0 -norm to promote the sparsity of the solution [18].Zhang et al. proposed a sparse nonconvex L p (0 < p < 1) regularization method for efficient reconstruction [19].Liu et al. proposed a method combining L 1 and total variation regularization (TV) to improve the imaging quality of CB-XLCT [15].Gao et al. proposed a primordial dual Newton conjugate layer method based on restart strategy on the basis of L 1 norm to alleviate ill-posedness [14].An elastic net-L 1 L 2 method based on the combination of L 1 and L 2 regular terms was proposed for CB-XLCT reconstruction by Zhao et al. [17].Beyond regularization, alternative approaches have been explored.Pu et al. demonstrated the use of compressed sensing (CS) techniques in CB-XLCT reconstruction to counteract the impact of sparse distribution of nanophosphors in vivo [20].Guo et al. proposed the utilization of a mixed proton propagation model and Lasso-LSQR algorithm to enhance reconstruction accuracy [21].Chen et al. proposed an algorithm based on dictionary learning and group structure strategy (DLGS) to effectively capture sparse features and enhance the interrelationship of elements in the same group [22].In addition to the above methods, the reconstruction performance can also be improved by using optical properties, multi-spectral information and permission source regions [8,23,24].
In this paper, inspired by the application of dictionary learning in compressed sensing (CS), an online dictionary learning (ODL) method leveraging iterative reduction FISTA (IR-FISTA) is introduced to enhance the reconstruction quality of CB-XLCT.In the field of sparse representations, there are the two well-known computational approaches: Orthogonal Matching Pursuit (OMP) as a representative of greedy algorithms and FISTA as a representative of convex relaxations.Compared to OMP, which is used in the sparse coding stage of DLGS, the convex relaxation algorithms represented by FISTA are recommended for many applications due to their generality and reliability.Typically, they offer better accuracy, but they are relatively slower in terms of speed [25].Therefore, we use IR-FISTA, which is much faster than FISTA, especially when the Lipschitz constant is large or the solution is highly sparse [26], aligning well with the demands of the CB-XLCT reconstruction problem.In the dictionary update stage, we adopt an online algorithm based on stochastic approximations, processing one sample at a time, updating dictionary by sequentially minimizing a quadratic local surrogate of the expected cost of the objective function, which makes the updating process more effective compared to global surrogate [27].
In order to evaluate the reconstruction performance of our method, a series of numerical simulations and simulation experiments were carried out.Our method was compared with L 0norm-based orthogonal matching pursuit (OMP-L 0 ) [28], and the incomplete variables truncated conjugate gradient with L 1 -norm (IVTCG-L 1 ) [29] in terms of location accuracy, target shape, robustness, dual-source resolution and in vivo practicability.
The rest of this paper is structured as follows.Section 2 provides a detailed description of the CB-XLCT imaging model and the proposed method.In Section 3, the numerical simulation and in vivo experiments are designed, and the evaluation indexes are proposed.Section 4 summarizes the experimental results obtained using different reconstruction methods.Finally, Section 5 discusses and summarizes the main findings of this work.

Methodology
CB-XLCT imaging mainly includes two parts: forward problem and reverse reconstruction.The forward problem involves calculating the optical signal distribution.In this process, the light emitted by the target light source and its transmission to the surface of the organism are considered.Knowledge of the location of the target light source, the anatomical structure of the biological tissue, and the optical parameters is required.A photon propagation model is utilized to describe this process.Based on the forward problem, the optical information of the biological surface is obtained, and combined with specific reconstruction algorithms, the position information and other quantitative information of the target are reconstructed as reverse reconstruction process of CB-XLCT.

Photon propagation model
According to the customized optical/CT dual modality imaging system, X-rays interact with the target to produce photons, and the excitable nanophosphors will produce visible light or NIR light after being stimulated by X-ray, the source energy density S(r) (W cm −3 ) can be calculated as follows [30]: where X(r) is the X-ray intensity (W cm −3 ) at position r, and ρ(r) is the density of the nanophosphor at position r, ε is the light yield of nanoparticles (mg/mL).In CB-XLCT, a conebeam incoherent X-ray source passing through biological tissues and exciting nanophosphors will produce photoelectric effect, which will make part of the X-ray electron flow absorbed or scattered.According to Beer-Lambert's law, when X-rays pass through tissues, the X-ray intensity distribution (energy attenuation) can be expressed as [31]: where X(r 0 ) is the intensity of the X-ray light source at the initial position r 0 , and µ t (τ) is the X-ray attenuation coefficient at position τ, which can be obtained with attenuation-based CT technology.Then we analyze the transmission of emitted photons, the transmission of visible light or NIR light in biological tissue is often modeled by radiative transfer equation (RTE) [5].However, due to the high computational cost of RTE, in general, the transmission of light in tissues with high scattering and low absorption can be simplified as the following steady-state diffusion equation (DE) model with Robin-type boundary condition, which is defined as follows [32,33]: where , which is the diffusion coefficient, and µ a (r) and µ ′ s (r) are the absorption and reduced scattering coefficients at the position r, respectively; Φ(r) is the photon flux, and Ω is the domain of the imaged object and ∂Ω is the boundary of Ω; v(r) is the outward unit normal vector on ∂Ω, and k(r, n, n ′ ) is the boundary mismatch factor, which depends on the refractive index n in Ω and n ′ in the surrounding medium.
Based on the finite element method (FEM) [34], Eq. ( 3) can be transformed to obtain the following matrix equation: with where ρ is the nanophosphors density to be reconstructed; M i,j and F i,j are the elements of matrix M and F, respectively; ψ i and ψ j denote the corresponding elements of the test function; X i,j (r) is the X-ray intensity for each vertex.
Since the matrix M is positive definite, Eq. ( 4) can be rewritten as: where and Φ is the photon measurements vector on the object surface, and its size is m × 1.In the practical application of CB-XLCT, the noise of the XLCT imaging system needs to be considered.Equation ( 6) can be written as: where b represents the photon measurements vector with noise, σ is the system noise, and x is the unknown nanophosphors distribution in the imaging object.

Reconstruction based on the SODL-IR-FISTA method
The objective of CB-XLCT reconstruction is to estimate the distribution of nanophosphor x from the measured photon density Φ.However, due to the severe scattering of visible and NIR light in biological tissue, the system matrix G becomes ill-conditioned, making it unrealistic to solve Eq. ( 6) directly.Considering the distribution of nanophosphors is small and sparse [35], to obtain a stable and robust solution, the L 0 regularization term can be used to find the solution by minimizing the following objective function [19]: where K is the sparsity level of the reconstruction result.The CB-XLCT reconstruction problem is transformed into an optimization problem of x in Eq. ( 8).
The SODL-IR-FISTA method uses the dictionary (matrix G can be regarded as a redundant dictionary) to continuously learn the sparse representation of the data, and by updating the dictionary, the updated dictionary can better represent the sparse data.Specifically, the solution is compartmentalized into three stages: sparse coding stage, dictionary update stage and re-sparse stage.The IR-FISTA is used in the sparse coding stage, which accelerates the solution under the premise of ensuring accuracy.In the dictionary update stage, an online algorithm based on stochastic approximation is used to further capture sparse features.Considering that the distribution of nanophosphor usually aggregates, the clustering algorithm in statistics is used for reference in the re-sparse stage to identify and eliminate outliers, as well as to control the sparsity of the solution, which can make the solution more accurate.The general flow of SODL-IR-FISTA is summarized in Algorithm 1.

Sparse coding stage
The sparse coding stage aims to address the objective function described in Eq. ( 8).However, due to its NP-hard nature, achieving a stable solution through conventional optimization methods is challenging.The convex relaxation algorithm represented by FISTA finds an approximation solution by transforming Eq. ( 8) into the following equation: where D is the learned dictionary (initially D = G) and λ is the regularization parameter, the above equation is also known as the LASSO problem.

Algorithm 1. SODL-IR-FISTA for CB-XLCT reconstruction
FISTA is celebrated for its simplicity and rapid convergence properties, with convergence speed decreasing quadratically with the number of the iterations.Although FISTA only needs to perform matrix-vector multiplications twice in each iteration, the computational cost of matrix-vector multiplications escalates significantly in large-scale problems involving an ill-conditioned matrix.This is because a large condition number of matrix would result in a large number of iterations.Consequently, FISTA may be inefficient for large-scale problems [36].Note that in the iterations of FISTA, a significant portion of the components of x typically remain zero following each iteration, particularly as the iterates approach the solution.This means that updating these components is "wasteful".Based on this observation, Wang et al. designed an iterative reduction algorithm [26], as shown below.
Using the current Karush-Kuhn-Tucker (KKT) information [37] to estimate the working set W k of the current iteration.
where τ is a small positive number.Moreover, let Even though we perform a proximal iteration solely on the working set W k , the computational expense of identifying this set equals that of performing a proximal iteration on a complete index set.Therefore, an effective approach involves executing multiple proximal iterations subsequent to identifying the working set.Due to this, after obtaining W k , we solve the following simplified subproblem: and let x k+1 . Hence, we have the procedures of iterative reduction FISTA as Algorithm 2.

Algorithm 2. Iterative reduction FISTA (IR-FISTA)
Moreover, let g(f , x) is used to represent the KKT residual of the LASSO problem, while g(f , x * ) can be used to measure the accuracy of the solution x * .For the FISTA algorithm, let where Ω(x) = λx 1 .And let The specific process of FISTA algorithm can be found in Algorithm 3.

Dictionary update stage
During the dictionary update stage, our method employs a block-coordinate descent approach coupled with warm restarts.One of its primary benefits is its lack of parameters, eliminating the need for any tuning of learning rates.
The new dictionary is obtained by the last updated dictionary and the obtained solution as warm restart, so that: Furthermore, the procedure necessitates storing only the matrices in R k×k and in R m×k , not all the vectors x k , where R is the set of real numbers, R N is the vector space composed of all N-dimensional column vectors and R M×N is the vector space composed of all matrices of size M × N. Concretely, in the dictionary update stage, this method sequentially updates each column of D. In practice, the vectors x k are sparse and the coefficients of the matrix A k are often concentrated on the diagonal, which makes the block-coordinate descent more efficient [27].Algorithm 4 describes the flow of the dictionary update algorithm.

Re-sparse stage
At this stage, the primary objective is to re-sparse the obtained solution and filter out the points that fail to meet the specified criteria.In most cases, the target data points should be closer to each other in space, and the points that do not meet the requirements may be far away from other points.Inspired by data mining, a re-sparse method is designed.We use the distance-based outlier detection algorithm for reference.By calculating the average distance from each point to other points and determining the deviation points according to a threshold, we can identify potential outliers.When there are multiple light sources, the non-zero data points in the solution are divided into several specific groups by distance-based clustering algorithm.This clustering ensures that data points within the same group exhibit similar spatial attributes, while those in different groups exhibit different attributes.Specifically, after outlier detection, for the remaining non-zero elements, the distance between adjacent elements is compared to determine whether they belong to the same group.If the distance between two adjacent points is less than the set value, they are considered to belong to the same group, otherwise they are not considered to belong to the same group.After the clustering process, several groups are selected as the basis for multi-target reconstruction according to the target number, and the elements of other groups are treated as outliers for subsequent elimination.Three possible scenarios of re-sparse stage are shown in Fig. 1.

Experiment design
In order to evaluate the effectiveness of the proposed method in CB-XLCT reconstruction, several numerical simulations and in vivo experiment were designed.Specifically, in terms of location accuracy, morphological recovery, robustness and in vivo practicability, two existing algorithms: OMP-L 0 and IVTCG-L 1 were compared.All experiments and procedures were conducted utilizing a laptop equipped with Inter Core i5-9300 H CPU (2.40 GHz) and 8GB RAM.

Numerical simulation setup
The heterogeneous phantom model was used for numerical simulation experiments.The model was a cylinder with a simulation height of 30 mm and a radius of 10 mm, including five organs and tissues: heart, liver, lung, muscle and bone.The 3D view of the model is depicted in Fig. 2(A).Table 1 displays the optical parameters of the tissues and organs at a wavelength of 650 nm within the model [22].In the inverse reconstruction, COMSOL Multiphysics 5.6 (COMSOL, Inc., Burlington, Massachusetts) software [38] is used to discretize the model into a finite element mesh of 4626 nodes and 25840 tetrahedral elements, as shown in Fig. 2 [39] implements the Monte Carlo (MC) method to obtain the forward simulation results of single-target and dual-target.The forward simulation result of single-target is shown in Fig. 2(C).Three groups of experiments were set up in the model.In the single-target simulation experiment, a uniform sphere with a radius of 1 mm and a center of (-4.0, 5.0, 5.0) mm was used to simulate early tumors.In the dual-target experiment, two identical spheres with a radius of 1 mm were used at (-6.0, 3.0, 8.0) mm and (-6.0, 3.0, 16.0) mm to simulate tumors.In addition, a series of anti-noise simulation experiments were carried out to evaluate the robustness of the method.Gaussian noise at level of 5%, 10%, 15%, 20%, 25%, 30% was added in the single-source simulation experiment, and the reconstruction results under different Gaussian noise ratios were observed.

Implanted experiment
In order to evaluate the feasibility of SODL-IR-FISTA in vivo CB-XLCT, we collected data from female BALB/C nude mouse using an XLCT imaging system.All experimental procedures were approved by the Animal Ethics Committee of Northwest University, China.Specifically, we implanted a flexible plastic tube (with a radius of 1 mm and a height of 2 mm) containing Gd 2 O 2 S:Tb nanophosphors (about 20 µl) into the mouse body and placed it at (15.2, 16.5, 20.0) mm to simulate an early tumor.The optical parameters of the tissues at a wavelength of 620 nm within our experiment are shown in Table 2.The imaging system included a controller, a CMOS X-ray detector panel, a cone beam X-ray source, an electron-multiplying CCD (EMCCD) camera and a rotation stage, as shown in Fig. 3.The cone beam X-ray source (Apogee, Oxford Instruments, USA) excited nanophosphors in biological tissues with a voltage of 50kVp and current of 1 mA.The EMCCD camera was placed perpendicular to the X-ray source, and the optical data with a wavelength of 620 nm were captured in a manner with an integration time and binning of 30 s and 2 × 2, respectively.To obtain 3D volume data, the mouse was fixed on a rotating stage and rotated 360°at 1°intervals to obtain micro-CT images and X-ray projection images.And the CMOS X-ray detector panel was used to detect transmitted X-rays, and then obtained high-resolution CT imaging.Subsequently, the Feldkamp-David-Kress (FDK) method [40] was used to perform 3D reconstruction of the XCT projection, and the forward measurement data of the nanophosphor was mapped onto the surface.The mouse model was discretized into 8935 nodes and 49545 tetrahedral elements for reconstruction, as Fig. 4 illustrates.

Evaluation indexes
In order to quantify and analyze the positioning accuracy and shape recovery performance of different methods in CB-XLCT reconstruction, two commonly used evaluation indexes are adopted: location error (LE) [41] and Dice similarity coefficient (DICE) [42].
LE is the Euclidean distance between the reconstructed region center (x r , y r , z r ) and the actual region center (x a , y a , z a ).It is a quantitative measure to evaluate the positioning accuracy of the reconstruction result.The LE value is always greater than 0, and the smaller the LE value, the   higher the positioning accuracy of the reconstruction.
DICE measures the overlap between the reconstructed target region X r and the actual target region X a .The DICE value ranges from 0 to 1.The larger the value, the greater the similarity between the reconstructed target region and the actual target region, indicating that the higher the morphological recovery ability of the algorithm.

Single-target experiment
The results of the single light source reconstruction experiment are shown in Fig. 5.The first row displayed the 3D view of the reconstruction results of each method, and the reconstructed source is visualized as a green area.In the section view, the yellow-green area represents a sectional view of the reconstruction results of the single-target at the Z = 5 mm plane, where the white circle represents the actual position and edge of the real target.Moreover, Table 3 displays the quantitative evaluation results of the three reconstruction methods.It can be seen that our method obtains a lower LE compared to other methods, which indicates that the reconstruction center of the SODL-IR-FISTA is closer to the true position.In addition, the DICE value of the SODL-IR-FISTA is higher, indicating that our method reconstructs the target region with the largest overlap with the true target region and has a better effect in morphological restoration.The reconstruction results of the dual-target experiment by the three methods are shown in Fig. 6.The representation of both the reconstructed and actual targets follows the same format as in the single-target experiment.It can be seen from 3D view and section view (Y = 3 mm) that the reconstruction results of SODL-IR-FISTA are the most accurate and the characterization of the actual target is better.The quantitative analysis of the reconstruction results of the three methods is shown in Table 4.The results indicate that, in comparison to other two methods, the LE of each target result reconstructed by the SODL-IR-FISTA method is the smallest and the value of DICE is the largest.This shows that in the case of multiple targets, SODL-IR-FISTA still has good positioning and morphological reconstruction capabilities, and can accurately distinguish different targets.The results of the anti-noise experiment are shown in Fig. 7. Based on the single-target experiment, 5%, 10%, 15%, 20%, 25%, and 30% Gaussian noise are sequentially added to the measurement data to observe their influence on the reconstruction results.It can be observed from the figure that the fluctuations of LE and DICE are small, indicating that our method can maintain good reconstruction performance under the condition of significant noise, which proves that the method has strong robustness.

In vivo experiment results
The results of the in vivo experiment are shown in Fig. 8, the 3D view is expressed in the same way as the above experiments.In the section view, the red-yellow area is the reconstruction result, and the gray-white circular area represents the section of the real source.Observations from 3D view and section view (Z = 20 mm) indicate that the SODL-IR-FISTA method provides the optimal reconstruction effect.Concurrently, in order to quantitatively analyze the reconstruction results of different methods, LE and DICE were calculated and summarized in Table 5.From Table 5, the LE obtained by OMP-L 0 is 0.721 mm, which is more than twice that of SODL-IR-FISTA.The LE obtained by IVTCG-L 1 is 0.562 mm, while the SODL-IR-FISTA is only 0.325 mm.This indicates that the reconstruction result center of our method is closer to the actual position, which can also be confirmed by observing Fig. 8. And, compared with OMP-L 0 and IVTCG-L 1 , SODL-IR-FISTA has better performance in morphological recovery.The DICE of the reconstructed results is the highest among the three methods, with a value of 0.748, while the DICE of OMP-L 0 and IVTCG-L 1 are 0.392 and 0.380, respectively.These results indicate that, the spatial overlap between the reconstruction results of SODL-IR-FISTA and the real source is the highest, which means that its shape recovery is more accurate.

Discussion and conclusion
CB-XLCT uses a cone-beam geometry to cover the entire imaging object, which can simplify the scanning process and shorten the imaging time.Therefore, it is more suitable for biological applications, especially fast biomedical imaging.Due to the diffusion propagation and high scattering of photons in biological tissues, the reconstruction process of CB-XLCT is ill-posed, resulting in inferior image quality.To solve this problem, in this paper, a reconstruction method called SODL-IR-FISTA is proposed.Since the distribution of nanophosphors in the body is small and sparse relative to the overall imaging object, this method utilizes the feature of dictionary learning in continuously learning sparse representations of data, effectively capturing sparse features of data and achieving sparse reconstruction of CB-XLCT.The SODL-IR-FISTA consists of three stages: sparse encoding, dictionary update, and re-sparse.The sparse encoding stage adopts the IR-FISTA algorithm, enhancing computational speed without compromising accuracy.The dictionary updating stage uses an online algorithm based on stochastic approximation to update the dictionary efficiently.In the re-sparse stage, an algorithm to remove outliers was designed to obtain more sparse and accurate solutions.
To verify the effectiveness of the proposed method, a series of numerical simulations and in vivo experiments were designed and carried out, and the results were compared with two common reconstruction methods, OMP-L 0 and IVTCG-L 1 .The experimental results are as follows.Firstly, the experimental results of single-target and dual-target show that SODL-IR-FISTA has good performance in reconstruction positioning accuracy and morphological restoration.Subsequently, the anti-noise experiments show that SODL-IR-FISTA has good robustness.Ultimately, the in vivo experiments and their results show the feasibility of SODL-IR-FISTA in biomedical research and the potential for CB-XLCT reconstruction in preclinical research.
However, despite the excellent performance of SODL-IR-FISTA in CB-XLCT reconstruction, there are still some challenges that need to be addressed.Firstly, during the dictionary update stage, the online update algorithm is not fast enough, which affects the overall efficiency of the method.Secondly, although most of the algorithm parameters do not require excessive adjustment, the target number and threshold values still need to be manually determined and cannot be automatically adjusted.Finally, more prior information can be used to obtain better results, such as group sparsity or smoothing prior.
In summary, this study proposes SODL-IR-FISTA method based on IR-FISTA and online dictionary learning to improve the reconstruction effect of CB-XLCT.The experimental results show that compared with other typical regularization-based methods, this method provides better reconstruction quality in terms of positioning accuracy and shape restoration.In addition, this method has strong robustness.Therefore, we believe that this method has great potential for preclinical application of CB-XLCT and can promote its in vivo application.

Fig. 1 .
Fig. 1.Schematic diagram of re-sparse stage.The straight arrow represents a vector with many points on it, the circle indicates the point where the element value is greater than zero, and the black circle indicates represents outliers.(A) When target number is 1, only outlier detection is performed; (B) When the target number is greater than 1, perform clustering operation after outlier detection.The blue and orange circles form Group 1 and Group 2 respectively.(C) When the target number is greater than 1, select the appropriate group according to the target number.Here, takes the target number equal to 2 as an example, and classify Group 3 as outliers.

Fig. 2 .
Fig. 2. Heterogeneous cylindrical phantom for numerical simulation studies: (A) The 3D view of the cylindrical phantom.(B) The 3D view of meshed cylindrical phantom.(C) The forward simulation result of single-target.

Fig. 4 .
Fig. 4. (A) The 3D view of the mouse.(B) The 3D view of meshed biological tissues.

Fig. 5 .
Fig. 5. Reconstruction results of three methods for a single spherical target.

Fig. 8 .
Fig. 8. Reconstruction results of the in vivo experiment.