Abstract
Diffuse optical tomography is a novel molecular imaging technology for small animal studies. Most known reconstruction methods use the diffusion equation (DA) as forward model, although the validation of DA breaks down in certain situations. In this work, we use the radiative transfer equation as forward model which provides an accurate description of the light propagation within biological media and investigate the potential of sparsity constraints in solving the diffuse optical tomography inverse problem. The feasibility of the sparsity reconstruction approach is evaluated by boundary angular-averaged measurement data and internal angular-averaged measurement data. Simulation results demonstrate that in most of the test cases the reconstructions with sparsity regularization are both qualitatively and quantitatively more reliable than those with standard regularization. Results also show the competitive performance of the split Bregman algorithm for the DOT image reconstruction with sparsity regularization compared with other existing algorithms.
1. Introduction
Diffuse optical tomography (DOT) is an emerging imaging modality that has attracted much attention in clinical diagnosis, for example, in breast cancer detection, monitoring of infant brain tissue oxygenation level, and functional brain activation studies, cerebral hemodynamic, and so forth; compare [1–3]. With the use of near-infrared (NIR) light, DOT probes the optical properties, mainly the absorption coefficient and the scattering coefficient of human tissues. In experimental systems, a set of optical fibres and optodes are attached to the boundary of the object as measurement detectors and sources. NIR light as the inflow current is emitted by laser and guided by some fibre optics into the object, one position at a time. The light is transmitted, and then the outflow current is measured from all the measurement positions using light sensitive detectors.
By employing a setup comprising a set of external light sources and detectors, the optical properties of a tissue can be recovered by applying the principles of tomography. The difference in absorption or scattering between the normal and abnormal tissues can provide the imaging contrast for tissue diagnostics. In this context, the aim of diffuse optical tomography is to provide the spatial distribution of the absorption and scattering coefficients of a tissue.
The forward problem in DOT describes the photon propagation in tissues and the inverse problem involves estimating the absorption and scattering coefficients of tissues from light measurements on the surface. There are many models describing the light propagation, for example, the Fokker-Planck equation, differential approximations of the radiative transfer equation (RTE) and so forth (cf. [4, 5]), and the diffusion equation (DA), which is accepted as a popular forward model describing light propagation in tissues. DA is a low-order approximation of RTE. However, the validation of DA breaks down in the following situations. First, since DA is only valid in the tissues with high scattering and low absorption, while in many tissues of human body with low scattering, such as skeleton, joint, DA is largely limited. Second, DA is not suitable to describe the light propagation in regions near to the light source, which will result in the large error when using DA as the forward model to image small animal or small tissue. For more on DA, we refer to [1, 6, 7]. Due to these limitations of DA, RTE is a more natural choice to model the light propagation. In this paper, we will present numerical evidence to show that DOT based on RTE can achieve satisfactory resolution.
Since DOT suffers from severe ill-posedness caused by noise and incomplete measurement data, its efficient, stable, and accurate numerical treatment is very challenging. Due to its immense range of prospective applications, there has been a vast amount of research work on mathematical as well as practical aspects of the inverse problem. In particular, the design of efficient and stable numerical algorithms has received considerable attention.
In the linearized Jacobian-based methods, a popular approach for solving the minimization problem is to treat it as a nonlinear least-square problem that can be solved with many standard optimization techniques [8], among which Levenberg-Marquardt method (LM) [9] is perhaps the most commonly used one, which is a Gauss-Newton method with regularization. Based on RTE, the authors in [10] discussed in detail the LM method for parameter estimation problems. We will describe the LM method briefly in Section 3. In the nonlinear gradient-based methods, we directly minimize the nonlinear functional; in the optimization process only the gradient of the nonlinear functional needs to be computed, for example, the BFGS (Broyden-Fletcher-Goldfarb-Shanno) method and the L-BFGS method [8, 11].
The appropriate choice of regularization depends on a priori knowledge of the domain and the inclusions. Although regularization is usually a natural choice for its simplicity, it is not the optimal strategy. For example, when the coefficient distribution is sparse or discontinuous, it is well known that the regularization method [12, 13] or the bounded total variation (TV) regularization method [14, 15] is more efficient than regularization.
Sparse reconstruction has attracted much researchers’ attention, especially since Donoho et al. [16, 17] established the theory of compressive sensing (CS). CS describes the sparse reconstruction problem as an quasinorm optimization [18]. However, Donoho in [19] proved that the regularization can also obtain the sparsest solution and proposed the equivalent condition between the regularization and the regularization. In [20], Candes et al. studied the stable signal recovery from incomplete and inaccurate measurements and reduced a computationally difficult problem to the basis pursuit problem. Because the sparsity regularization is nonsmooth, it is still a challenge to find efficient methods to solve this convex basis pursuit optimization problem, and the choice of techniques for solving it becomes crucial. Classical gradient-based methods usually bring high computational burden [21]. Motivated by the inverse problem in imaging [22–25], the authors used and developed the Bregman iteration technique for the regularization problems and proved that the Bregman iteration method is an effective way to solve the norm minimization problems.
In this paper, we adopt the split Bregman method for solving sparsity regularization problems. The split Bregman method is a simple and efficient algorithm which can split the minimization into and functionals and significantly reduce the computational burden [24]. We will introduce the split Bregman algorithm in detail in Section 3.
The organization of this paper is as follows. Section 2 introduces the DOT forward problem, briefly presents the radiative transfer equation (RTE), and describes the imaging modality of the forward problem. Section 3 presents in detail the reconstruction algorithm, including the Levenberg-Marquardt algorithm and the split Bregman algorithm. Section 4 provides several simulations to show the validity of the sparsity regularization reconstruction. Section 5 gives some conclusion statements.
2. Forward Problem
In this section, we formulate our problem of RTE based diffuse optical tomography.
2.1. Radiative Transfer Equation
Photon propagation in tissues can be described by the radiative transfer equation. Let , or , denote the physical domain of the medium with boundary , the unit sphere, the unit outer normal vector, and the outgoing and incoming boundaries defined by
The variables and denote the spatial position and the angular direction. Then we consider the following boundary value problem (BVP) of RTE: where is the total attenuation coefficient, describes the probability that a photon is absorbed in unit length, its reciprocal being the absorption mean free path, and is the scattering coefficient, describing the probability that a photon is scattered in unit length, its reciprocal being the scattering mean free path. Further, is the radiance and is the internal light source. In this paper, we consider the case with no light source inside ; . is the inflow current on , and the boundary condition (3) implies that once a photon escapes the domain , it does not reenter it; compare [26]. is the scattering operator. Denote by the infinitesimal area element on the unit sphere . Then is defined as with a nonnegative normalized phase function: In many applications, the function is independent of . However, in our general study, we allow to depend on . Indeed, we can even allow to be a general function of , , and , that is, in the form . The scattering phase function describes the probability that a photon with an initial direction will have a direction after a scattering event. In DOT, one typical example is the Henyey-Greenstein phase function; compare [27]. Consider where the parameter is the anisotropy factor of the scattering medium. Note that for isotropic scattering, for forward scattering, and for backward scattering. In biomedical imaging problems, the scattering is strongly forward peaked and is close to 1.
2.2. Forward Problem
Similar to the X-ray CT, DOT experiments acquire the current distribution of detectors on the boundary under the multi-incidents [28]. The experimental procedure of acquiring potential measurements is as follows. First, set a set of laser devices and detectors on the boundary of the object. Then launch an incident impulse from one laser device and record the resulting measurements from all the detectors. In order to gather enough data information, repeat this procedure on other laser devices.
We can model this procedure mathematically. Before doing this, let us introduce some assumptions on the domain and the coefficients, assumed to be valid throughout the rest of this paper.
We assume that the absorption and scattering coefficients are approximated piecewise-constant. It is known that the BVP (2) and (3) has a unique solution [29], and the solution depends continuously on the input current on the incoming boundary .
The forward problem of DOT is to determine the outgoing current on the detectors when the incident impulse and the absorption and scattering coefficients are known. Excite the domain with a sequence of incident impulses , , and get a sequence of measurements corresponding to each incident impulse, with its component , , being the measurement value on detectors. Then a mathematical description of such an experiment is provided by a sequence of forward operators: which maps prescribed optical parameters to the corresponding measurements. Here, denotes the th forward operator with respect to the th incident impulse and the resulting detected measurement data on detectors. Then, for , the domain of the operator is defined as is a column vector representing the measurement data on detectors.
We mention that, in DOT, the measured quantity is the excitance on the boundary of the domain. Due to the limit of measurement techniques of the optical devices, the angularly resolved measurement data is not practical. In this study, for , we use the boundary angularly averaged data as our measurements [30]: where is the solution of the BVP (2) and (3) corresponding to the optical coefficients and the th incident impulse and is the unit outer normal vector, whereas the set In [30], a detailed analysis is given on properties of the forward operator, including the Lipschitz continuity and the Fréchet differentiability. The BVP (2) and (3) needs to be solved by numerical methods, such as the finite difference method or the finite element method in spatial discretization, and the finite element method or approximation in angular discretization. For the simulation results in this paper, we use the RTE2DMATLAB codes [6], in which the finite element method in both spatial and angular spaces is used to discretize the forward mapping. We refer the reader to [4, 6, 31–34] for details about the finite element implementation in both spatial and angular space.
3. Inverse Problems
In practice, due to the limitations of the experimental environment and the laboratory equipment, the measurement we get usually contains noise. Here, we assume that, for , the actual measurements have the noise level ; that is, , where represents the actual measurement data and represents the true data corresponding to the true optical coefficients. Then our inverse problem of DOT is to determine such that the following nonlinear equations hold: for .
In this paper, we only reconstruct the scattering coefficient, while assuming that the distribution of the total attenuation coefficient is known.
As is typical for many inverse problems, the DOT inverse problem is ill-posed. In order to reconstruct the optical parameter stably, regularization is required. We minimize the following Tikhonov functional: over the admissible set for the coefficient . Here, is a regularization penalty functional that enforces a priori knowledge on the optical parameter to be reconstructed and is the regularization parameter used to trade off the discrepancy term (the first item of ) that incorporates the information contained in the data and [35]. Then we analyse the minimization problem:
We consider two popular regularization formulations: the norm penalty and the sparsity constraint regularization. We will describe these two approaches below.
3.1. Standard Reconstruction
First, we use the traditional norm squared penalty, which consists of minimizing the following functional: with a specified based on a priori information on the solution. In the statistical inversion framework, the corresponding prior constructions are known as smoothness priors [36]. To compute a minimizer of problem (14), many iterative regularization methods are available. We use the Levenberg-Marquardt method based on linearization. Specifically, for every , the forward operator is linearized around some initial guess ; that is, where is the Fréchet derivative of with respect to the coefficient at and denotes the Taylor remainder for the linearization around . Then substituting the above linearized expression into the functional and ignoring the higher-order remainder , we get a linearized problem The Euler equation of the discrete problem is that is, with the identity matrix. The system can be solved directly to get a new estimate for based on the initial guess . Then we iteratively update the reconstruction by taking the solution as an initial guess. In practice, the iterative procedure achieves required accuracy within a few iterations. The complete algorithm is given in Algorithm 1. The stopping criterion can be defined based on monitoring the relative change of consecutive iterations.
| ||||||||||||||||||||||||
3.2. Sparsity Reconstruction
In the sparsity reconstruction, the functional to be minimized is of the form We will apply the Bregman framework to solve (20). The key to our method is to “decouple” the and portions in (20). Rather than (20), we will consider the constrained problem [23]: To solve the above minimization problem, the corresponding unconstrained optimization problem is where is the split parameter. Then we can iteratively solve the following subproblems [37, 38]: The minimization of the above subproblems can be iteratively solved by splitting it into the minimizations of and separately. This suggests the following steps.
Step 1. Consider
Step 2. Consider
Step 3. Consider
For the solving of Step 1, we can use an iterative method, for example, the Landweber iteration method, the Levenberg-Marquardt method, and so forth. Since the Levenberg-Marquardt method has a higher convergence rate than the conventional Landweber iteration method, we use the Levenberg-Marquardt method on Step 1. Thus, we will solve a minimization problem as follows.
Step 1. Consider To solve Step 1, we solve the explicitly given variational equation as follows: Step 2 is an norm regularization problem and it can be solved efficiently through the shrinkage operator; that is, where the shrinkage operator
From the three steps, we can see that the speed of the split Bregman method is largely dependent on the speed of solving Step 1, where the computation of the Jacobian matrix is very time-consuming. We list the split Bregman method in Algorithm 2.
| ||||||||||||||||
3.3. Analysis of Standard Reconstruction and Sparsity Reconstruction
The appropriate choice of regularization depends on a priori knowledge of the solutions. Although -norm regularization is a natural choice for its simplicity, it is not always the optimal strategy.
For the DOT reconstruction, the sought-for optical coefficients distribution usually consists of an essentially uninteresting background plus some small inclusions. Thus we will require the solution of the reconstruction in a sparse coefficient vector form; that is, the coefficient vector of the optical parameter contains only a finite number of nonzero elements.
The standard reconstruction and sparsity reconstruction add -norm penalization and -norm penalization into the problem, respectively, which allow additional constraints or prior information towards the approximate solution. Figure 1(a) illustrates how -norm penalization leads to sparse solutions. Take a linear problem as an example; suppose we are looking for a solution of the linear equation in a model space with two degrees of freedom and the line consists of that satisfy . To find the solution with smallest -norm, we can imagine taking a small circle around the origin and increase its radius until it first touches the solution line : the tangent point is the minimum; see Figure 1(b). Obviously, the point of solution has both and components. In a similar way, to find the solution with smallest -norm, we take a small -ball around the origin and increase its radius until it first touches : the touching point is the minimum -norm solution, which is sparser than the solution achieved with smallest -norm, because it is on one axis only; only one component is nonzero.
(a)
(b)
There is also an extreme penalty term, say, -norm penalization, which is simply the number of nonzero coefficients and leads to the sparsest solution. Nevertheless, the optimization problem with this penalty is not computationally tractable. Hence -norm penalization is preferred in practical problems. In addition, it has been proven that, for some large matrices , if there exist sufficiently sparse solutions, the sparsest solution can be achieved by the -norm minimization [19, 20]. Under certain conditions, the penalty term can provide an accurate result even with limited observations [17]. These phenomena are further investigated in different fields such as electrical impedance tomography and tomographic inversion [39, 40].
4. Numerical Implementations
In this section, we report simulation results on numerical examples in two dimensions (2D). We perform the simulations on a 3.0 GHz PC with 8 GB RAM in MATLAB 2013b environment under Windows 7. In our simulations, the scattering coefficient is reconstructed both with the standard regularization reconstruction and with the sparsity regularization reconstruction. We use two kinds of measurements for reconstruction: boundary angular-averaged measurements and internal angular-averaged measurements. The boundary measurements are the excitance received by detectors attached to the boundary of tissue. Due to the limit of measurement techniques of the optical devices, we cannot accurately receive the excitance from all angles; instead we receive a boundary angular-averaged data which can be regarded as the integration of radiance on the boundary of the domain in all outgoing directions. The internal angular-averaged measurements are the same as boundary angular-averaged measurements except for the location of detectors which we assumed to be inside the domain.
We design different simulations to demonstrate the following points.
First, for sparse coefficient distribution, the sparsity regularization reconstruction localizes the location of the inclusion better than standard regularization reconstruction, especially when there is noise and the data are incomplete. See Figures 4, 9, and 10.
Second, the proposed sparsity regularization reconstruction works better when than . See Figures 6 and 7.
Third, when the internal angular-averaged measurements are available, the proposed sparsity regularization reconstruction works better with internal angular-averaged measurements than with boundary angular-averaged measurements. See Figures 13 and 14.
Last, the split Bregman algorithm can efficiently solve the DOT image reconstruction problem with sparsity regularization. The results reconstructed by the Bregman algorithm are more accurate than those achieved by three other algorithms. See Figures 16 and 17.
We mention here that, in the DOT reconstruction problems, the measurement data are usually synthesized from the numerical solution of the forward problem. Under this situation, the phenomenon of inverse crime will happen especially when the same discretization is used for the forward and inverse process, because it will make the ill-posedness of the inverse problem not evident [41, 42]. Hence, in order to avoid the inverse crime, we will use different discretization meshes in the forward and inverse problems.
Now let us state our numerical experiments. Before stating the details of every experiment, we first state some common experiment settings that will be used in all of our experiments. Our purpose is to reconstruct the scattering coefficient based on the BVP (2) and (3), where we assume that the following 2D simulations are all performed on a unit circular domain centered at (0 mm, 0 mm), corresponding to the the internal light source , the inflow current for every incident impulse that is settled as The boundary value is a piecewise linear function whose spatial support is and achieves the value at the center node of ; here denotes the finite element through which the th incident impulse passes. The direction of the th incident impulse is node of the angular mesh that points approximately from the center of to the center of .
If we denote the domains that contain the inclusions as , then the background domain is . Then the true value of absorption coefficient in BVP (2) and (3) is defined as Since the absorption coefficient is assumed to be known in our scattering reconstruction problems, we only give the mathematical formulation of the absorption coefficient. While in order to compare explicitly the true value and the reconstructed value of scattering coefficient, we demonstrate the true value and the reconstructed results of scattering coefficient in figures.
In order to make discretization in angle, we can divide the angular space uniformly into directions with equal interval length as shown in Figure 2. And we set in all examples; that is, there are 32 angular directions.
In the first example, Figure 3, a small inclusion with 0.1 mm diameter, is centered at . 12 sources and 12 detectors are equally spaced attached to the boundary of the circle domain. The simulated boundary angular-averaged measurements performing on the boundary are generated with a spatial mesh (a) of 1097 nodes and 2104 elements and reconstruction mesh (b) having 286 nodes and 526 elements.
(a)
(b)
(c)
(d)
(a)
(b)
(c)
(d)
The image from standard regularization reconstruction without noise is displayed in (d) in Figure 3. The images from sparsity regularization reconstruction without or with noise are displayed in Figure 4 with the noise defined as , where is the signal-to-noise ratio and is a Gaussian random variable with zero mean and unity variation. In Figure 4, we vary from 0 percent to 10 percent and plot the sparsity regularization reconstruction image. Comparing images from standard regularization and sparsity regularization, we conclude that sparsity regularization can efficiently localize the sparse inclusion with a few measurements, even in the presence of different degrees of noise, while standard regularization fails to reconstruct the sparse inclusion even without noise.
In the second example, we compare the sparsity regularization reconstruction and the standard regularization reconstruction on squares of three different sizes.
In Figure 5, (a), (c), and (e) show the forward mesh of three different sizes squares. (b), (d), and (f) show the true distribution of the scattering coefficient in three different experiments.
(a)
(b)
(c)
(d)
(e)
(f)
(a)
(b)
(c)
(d)
(e)
(f)
(a)
(b)
(c)
(d)
(e)
(f)
Comparing (a), (c), and (e) in Figures 6 and 7 and (b), (d), and (f) in Figures 6 and 7, clearly, the sparsity regularization can localize the position of the inclusion better than the standard regularization in all the three kinds of squares and has a clear contrast with the backgrounds. But the standard regularization also has its advantages. Comparing (a) and (b), (c) and (d), and (e) and (f) in Figure 6, when , there are many blurred dots in the standard regularization reconstructed images; we can not identify the distribution domain of the reconstructed scattering coefficient, but the value of the scattering coefficient in standard regularization reconstruction is much closer to the true value of the scattering coefficient than that in the sparsity regularization reconstruction. However, this advantage is broken when ; see Figure 7.
On the other hand, the proposed sparsity regularization reconstruction performs better when than when in localizing the position of the inclusions as well as identifying the value of the scattering coefficient.
In conclusion, when the anisotropic factor is small, the standard regularization can identify the value of the scattering coefficient better than the sparsity regularization. Although there are some blurred dots in standard regularization reconstructed images, one can remove them by using multilevel approach. On the other hand, the sparsity regularization can localize the position of the inclusion better than standard regularization and has a clear contrast, especially in the forward-peaking regime with big anisotropic factor .
In the third example (Figure 8), three-letter inclusions (a) are put in the domain to evaluate our reconstruction algorithms; the forward mesh (a) has 1077 nodes and 2072 elements; the reconstruction mesh (b) has 280 nodes and 518 elements.
(a)
(b)
(c)
(a)
(b)
(c)
(d)
(a)
(b)
(c)
(d)
In Figures 9 and 10, the reconstructions with and 0.9, respectively, are carried out to demonstrate that sparsity regularization ((a), (c) in Figures 9 and 10) in general reconstructs better images than standard regularization ((b), (d) in Figures 9 and 10).
Comparing (a) and (c) in Figures 9 and 10, respectively, we find that the proposed sparsity regularization reconstruction performed better when the inverse problem is more severely ill-posed due to fewer measurements.
In conclusion, efficient algorithm and proper regularization, for example, sparsity regularization reconstruction for sparse coefficient distribution, are essential to recover high-resolution images. On the other hand, we can notice that there are blurred dots in the middle of the reconstructed images in (a) and (c) in Figures 9 and 10, due to the ill-posedness of inverse problem. These blurred dots have small contrast and thus can be removed through sparsity regularization reconstruction by choosing proper and or, beyond our paper, through multilevel approach [43].
Next, using sparsity regularization reconstruction, we will show the superiority of the internal angular-averaged measurements over the boundary angular-averaged measurements with two relatively complicated cases. The detectors for internal angular-average measurements in two cases lie on a 0.95-radius circle equidistant from the center at (0 mm, 0 mm). Case 1 is shown in the first column in Figures 13 and 14; Case 2 is shown in the second column in Figures 13 and 14.
Figure 11 shows the meshes for Case 1. (a) shows the forward mesh; it has 1131 nodes and 2160 elements. (b) shows the inverse mesh; it has 296 nodes and 540 elements. (c) shows the true distribution of the scattering coefficient. Figure 12 shows the meshes for Case 2. (a) shows the forward mesh; it has 1347 nodes and 2592 elements. (b) shows the inverse mesh; it has 350 nodes and 648 elements.
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
(d)
(a)
(b)
(c)
(d)
Comparing (a) and (c) in Figures 13 and 14, we find that, for Case 1, the internal angular-averaged measurements ((c) in Figures 13 and 14) can localize the position of the outside inclusion, although the internal inclusion is localized blurred while the boundary angular-averaged measurements ((a) in Figures 13 and 14) can not completely localize any inclusion.
For the complicated Case 2 ((b) and (d) in Figures 13 and 14), we increase the measurement data by increasing the number of detectors to 16, respectively, in order to get more information from the boundary measurements. Both of the two kinds of measurements can not accurately identify the inclusion, but the internal angular-averaged measurements ((d) in Figures 13 and 14) performed relatively better than the boundary angular-averaged measurements ((b) in Figures 13 and 14). A justification of this phenomenon is that the energy of the incident current is decayed much due to the large and relatively frequently change of the scattering coefficient of the inclusion. One can alleviate this phenomenon by increasing the number of detectors or measurement data. In conclusion, the internal data can better-pose the inverse problem.
As the last simulation, we compute Case 3 in Figure 15 to investigate the performance of split Bregman algorithm for DOT image reconstruction with sparsity regularization. The results are compared with Gauss-Newton algorithm [44] and the state-of-the-art algorithms including GPSR [45] and YALL1 [46]. To deal with the nondifferentiability of the absolute value at in Gauss-Newton algorithm, we replace by , . We choose in Gauss-Newton algorithm.
(a)
(b)
(c)
(a)
(b)
(c)
(d)
(a)
(b)
(c)
(d)
In Case 3, two circle inclusions are embedded at the top and bottom of the circle domain with 1137 nodes and 2168 elements for the forward mesh and 298 nodes and 542 elements for the inverse mesh. We use 132 boundary angular-averaged measurements with 0.1% Gaussian noise for reconstruction. The reconstructed images are shown in Figures 16 and 17, and we summarize the results with in Table 1, containing the parameters for the four methods, computational time, data misfit , relative solution error norm (RE), and signal-to-noise ratio (SNR). The relatively optimal parameters are chosen empirically.
The RE is calculated as and SNR is calculated as
The reconstruction results in Figures 16 and 17 show that all of the four algorithms can locate the position of two inclusions. The split Bregman and GPSR algorithm can achieve more accurate shape of the inclusions, and the reconstructed scattering coefficients achieved by the two algorithms are closer to the true value. The Gauss-Newton algorithm can find the approximate location of the inclusions, but the reconstructed circle inclusion is smaller in size and the reconstructed scattering coefficient departs from actual ones.
The results displayed in Table 1 also demonstrate that the split Bregman algorithm performs better. It leads to the lowest RE and data misfit with less calculation. The GPSR algorithm achieves higher SNR at the expense of about 25% extra computational time. The Gauss-Newton algorithm spends least computational time, but the reconstruction error is not very satisfactory.
These results are justified by the fact that the split Bregman algorithm decouple the sparsity reconstruction problem into and portions leading to a better compromise in the efficiency and quality of the reconstructed optical parameters.
5. Conclusion
In this paper, by using the image modalities in DOT, we employ the sparsity regularization method on the RTE-based coefficient identification problems, which is proven to perform in general better than the standard regularization reconstruction, especially for sparse distribution coefficient and large noise, and in the forward-peaking regime with big anisotropic factor . On the other hand, we construct cases and compare the reconstruction results with boundary and internal measurement to test the validity of the proposed method; results show that the proposed method is practicable and feasible; it performs steadily with various measurement data; meanwhile, the internal measurement can better-pose the inverse problem and achieve more accurate results. However, we usually cannot obtain the internal measurement data in practice. Hence, our method cannot reconstruct inclusions with complicated internal structure accurately with small amount of boundary angular-averaged measurements. One can alleviate this phenomenon by increasing the number of detectors or measurement data. In the further work we will consider seeking for multi-imaging modality which can further improve the inversion quality with boundary angular-averaged measurements. Results show the competitive performance of the split Bregman algorithm for the DOT image reconstruction compared with other existing algorithms.
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgments
The authors would like to express their gratitude to the referees for their careful reading and detailed suggestions, which helped to improve the revised version of this paper. This work is supported by the National Nature Science Foundation of China under Grant no. 91230119, no. 81071207, and no. 81271622, and by “the Fundamental Research Funds for the Central Universities” (HIT.IBRSEM.B.201416).