Reconstruction of localized fluorescent target from multi-view continuous-wave surface images of small animal with lp sparsity regularization

Shinpei Okawa; Tatsuya Ikehara; Ichiro Oda; Yukio Yamada

doi:10.1364/BOE.5.001839

1. Introduction

Fluorescence diffuse optical tomography (FDOT) is a non-invasive biomedical imaging modality [1,2]. FDOT uses near-infrared (NIR) light propagating diffusely in biological media which strongly scatter and weakly absorb NIR light [3, 4]. FDOT reconstructs the distributions of the concentration and properties of fluorophores, such as lifetime and quantum yield in biological media by solving an inverse problem with the input of the fluorescent light detected at the surfaces of the media. One of the applications of FDOT is molecular imaging of small animals for drug developments. Delivery of drugs labeled with fluorophores in small animals can be monitored by FDOT [5]. Screening of breast cancer by FDOT is also expected. Some fluorophores targeting tumors enhance the image contrast of the tumors without relying on the endogenous components such as hemoglobin, water or fat [6]. Quantum yield and lifetime reconstructed by FDOT reflect conditions of diseased tissues and provide information useful for diagnoses [7,8].

FDOT has many advantages over other imaging modalities such as x-ray CT, MRI, PET, and so on. Measuring instruments for FDOT are smaller in size and easier to use than the other modalities. However, FDOT has a disadvantage of a low spatial resolution. Fluorophore concentrations reconstructed by FDOT usually have broader distributions than the true ones. Measurement noises also affect the reconstructed image quality. Especially for continuous-wave (CW) mode FDOT, it is difficult to reconstruct fluorophores embedded deeply in an object at the correct positions and with the correct concentrations. It is desired to obtain high quality images with high spatial resolutions and good quantitativeness using CW instruments which are the easiest to be handled and cheapest to be built among biomedical optical imaging techniques.

Various FDOT systems have been developed in recent years. Graves et al. developed a system which measured fluorescent light by detectors aligned at an object surface opposite to that illuminated by raster scanned light source, and reconstructed 3-D fluorophore distributions [9]. Delioanis et al. reported a system having excitation light sources and detectors which rotated around an object to obtain projections over 360° angles [10]. These systems can obtain large amount of measurement data for reconstruction of high quality images with the costs of long measurement time and heavy computational load. To shorten the measurement time, we have developed a CW measurement system that acquires fluorescence data at multiple angles of view by use of reflecting mirrors and wide illuminations of excitation light from multiple angles. This system can obtain measurement data from a whole body surface of a small animal in a short period of time and has been built with a low cost. The computational load can be reduced due to a small number of independent measurement data obtained with this system. However, images reconstructed from the measurement data might have a low spatial resolution and quantitativeness due to the ill-posedness which is partly attributed to a small number of independent measurement data. In this study we try to incorporate a regularization method into the image reconstruction process to obtain high quality images from the small number of independent measurement data.

Regularization techniques have been studied for image reconstruction of diffuse optical tomography (DOT) [11–13] which is an NIR imaging technique related to FDOT. It is demonstrated that efficient use of structural or anatomical prior information is effective to improve the quality of the DOT images. Spatially variant Tikhonov regularization was used to reduce high frequency noises in the reconstructed images by Pogue et al [14]. Boverman et al. evaluated the effect of prior segmentation of breast into grandular and adipose tissues [15]. Yalavarthy et al. regularized the inverse solution by use of MRI-derived breast geometry [16]. Anisotropic diffusion regularization with a priori edge information of the object was proposed by Douri et al. to preserve the edge of the inner structure [17]. An algorithm assuming that the targets have blocky structures for regularization was introduced by Hiltunen et al [18]. Mutual information and joint entropy were used by Panagiotou et al. to reflect the structure obtained from an alternative high resolution modality in the DOT images [19].

Sparsity regularization improves spatial resolutions of the DOT images. Cao et al. applied L₁ norm minimization by use of an expectation maximization algorithm for a linearized DOT inverse problem, and showed that reconstructed regions with abnormal optical properties were localized more sharply than other methods they used [20]. Shimokawa et al. proposed a method to reconstruct a sparse image based on hierarchical Bayesian estimation [21]. Sparsity regularization has been used for various inverse problems [22, 23].

Also in FDOT image reconstruction, regularization methods have been used. Freiberger et al. used regularization to obtain inverse solutions by a total variation method [24]. An efficient reconstruction with L₁ norm minimization was numerically investigated by Han et al. [25]. They also implemented a sparsity regularization with a shrinkage method [26]. Yi et al. reconstructed sparse images by L₁ norm minimization and Landweber iterative method with adaptive meshing of the finite element method (FEM) [27]. The L₁ sparsity constraint was also applied to fluorescence/bioluminescence diffuse optical tomography (F/BDOT) [28, 29]. Okawa et al. used recursive spatial filtering for F/BDOT to obtain localized fluorophore distributions [30].

In this paper, we try to implement the l_p (0 < p ≤ 1) sparsity regularization in the FDOT reconstruction algorithm by partial use of the focal underdetermined system solver (FOCUSS) algorithm [31, 32]. FOCUSS algorithm is proposed for solving the underdetermined linear inverse problems by use of recursive linear projection to obtain a sparse solution, which is used as a weight to improve the solution in the next iteration. The updating process is equivalent to minimizing the p-norm of the solution which is strictly constrained by the linear forward model.

The effects of the regularization method on the quality of the FDOT images reconstructed from multi-view CW surface fluorescence images of objects including a small animal are investigated. FDOT images are reconstructed by minimizing the residual error between the measurement and predicted data and the l_p norm of the concentration of the fluorophore simultaneously. Numerical simulations, phantom and mouse experiments demonstrate that the l_p (0 < p ≤ 1) sparsity regularization improves the localization of the fluorophore in the reconstructed images.

2. Method

2.1. Measurement system

An in vivo optical imaging system, Clairvivo OPT (Shimadzu Co., Kyoto, Japan), was used for acquisition of the fluorescence surface images of objects in this study. Figure 1 shows the schematics for acquiring the fluorescence surface images. The system consisted of a horizontal animal table, excitation light sources, a CCD camera, reflecting mirrors, a camera lens, an optical filter and eight white LED light sources illuminating the animal from four different angles. The table was made of plexiglass and was transparent to visible and NIR light. A small animal either anesthetized or sacrificed was laid down and fixed by adhesive tapes at a right position on the table without using an animal holder. The maximum sizes of an object were 40 mm, 30 mm and 120 mm in width, height and length, respectively, which were limited by the view of the CCD camera. A gas port for anesthesia was installed at an end of the table, and would be connected to the head of a small animal when in vivo measurements would be made.

Fig. 1 (a) Schematics of the CW measurement system to acquire the surface fluorescence images of the measured object, and the angular positions of (b) excitation and (c) emission light.

Download Full Size | PDF

For excitation, five CW diode lasers emitting light at the wavelength of 785 nm and a power of 0.07 mW/cm² illuminated the whole object from five different polar angles of 52, 120, 180, 240 and 308 degrees, i.e., Excitation 1 to Excitation 5 as shown in Fig. 1(b). Each laser light illuminated the object with a width of 40 mm and a length of more than 120 mm.

Fluorescence emission light was detected by a single CCD camera. An optical filter was placed in front of the CCD camera to cut the excitation light when the emission light was measured. The optical filter was a band-pass interference filter with a central wavelength of 845 nm and a bandwidth of 55 nm. Autofluorescence was observed even by use of the optical filter, but its effects on the fluorescence surface and reconstructed images were negligible. Two flat mirrors and two concave mirrors were placed around the object at the polar angles of 70, 164, 196 and 290 degrees as shown in Fig. 1(c), and five surface views of the fluorescence emission from the object, Emission 1 to Emission 5 including the polar angle of zero degree, were acquired simultaneously by the CCD camera. The five views were imaged in parallel at different positions on the CCD sensor without overlapping to achieve simultaneous acquisition. Totally, 25 fluorescence surface images were acquired by switching the position of the excitation light one by one from Excitation 1 to Excitation 5, and were used to reconstruct the distribution of the fluorophore concentration in the object. The CCD sensor had 1024×1024 pixels, and 4×4 binning of the pixels provided 256×256 data for one CCD image. For example, one surface view of the fluorescence emission from a cylindrical phantom used in the phantom experiments occupied 25×51 data, and the number of the data for five views for one excitation was 25×51×5. By summing up the data for the five angles of excitation, the total number of measurement data for the phantom experiments was 25×51×5×5= 31875 which was the input for one reconstructed image. This number of the input data is about 10 % of the maximum number of 256×256 provided by the CCD sensor, meaning a low efficiency of the usage of the sensor, which was caused by the simultaneous acquisition of five surface images. But, by limiting the number of the input data, the measurement time was shortened, the measurement noises were reduced and the computational load of image reconstruction was also reduced, resultantly.

The five emission images were aligned in the CCD sensor in the order of Emissions 2, 3, 1, 4 and 5 from the left to right as shown in Fig. 1(c), and corrected for distortion by the concave mirrors. When the reconstructed FDOT images were superimposed on the visible images the order was changed to Emissions 4, 5, 1, 2 and 3 for easy recognition as shown later in the sections of phantom and mouse experiments. We assumed that the emission light diffusively emitted from the surface of the measured object so that the directional dependence of the emission light from the curved surface of the object was ignored.

Eight white LEDs simultaneously illuminated the object by visible light from the four angles of Excitations 1, 2, 4 and 5. Visible light images of the object, to which the fluorescence surface images were superimposed, were also acquired by the CCD camera from the five angles of view same as the emission images. The optical filter was removed for visible light imaging.

2.2. Forward modeling

Light propagation in biological media is a phenomenon of the transport of radiative energy. NIR light is strongly scattered and weakly absorbed by biological tissues while propagating in the tissues. In this study, the photon diffusion equation (PDE), which is the approximation of the radiative transport equation, is used to describe light propagation in biological media [3].

The following PDEs are applied for CW excitation and emission light,

{- \nabla \cdot D_{x} (r) \nabla + μ_{a_{x}} (r) + ε N (r)} Φ_{x (r)} = q_{0} (r),

{- \nabla \cdot D_{m} (r) \nabla + μ_{a_{m}} (r)} Φ_{m} (r) = ε γ N (r) Φ_{x} (r)

where D = 1/(3μ′_s) is the diffusion coefficient with the reduced scattering coefficient μ′_s, μ_a the absorption coefficient of the background medium, and Φ the fluence rate at the position r with the subscripts x and m indicating excitation and emission light, respectively. ε and γ are the molar extinction coefficient and the quantum yield, respectively, and N is the fluorophore concentration. The light source q₀ is distributed uniformly in the area of the surface illuminated by the excitation light. The boundary conditions are given as −n · D_y∇Φ_y = Φ_y/(2A) where n is the vector outward normal to the surface of the medium, subscript y refers to x or m, and A is the parameter depending on the ratio of the refractive indices, n_r. In this study, n_r is give as 1.36, then A = 3. The forward solutions of light propagation are obtained by solving Eqs. (1) and (2) using the finite element method (FEM) [33].

Measured fluorescence intensities are given by the fluxes of the fluence rate of the fluorescence emission measured at the surface, Γ(r) = −n · D_m∇Φ_m = Φ_m/(2A), and are calculated from the solutions of Eqs. (1) and (2).

Now the calculated fluorescence intensity, Γ_i,j,k,l, is given for the i-th excitation, j-th angle of view, k-th pixel in an surface image, and the l-th voxel in a medium containing a unit fluorophore concentration where i and j = 1, 2,···,5, k = 1, 2,···,K, l = 1, 2,···,L. Here, K is the total number of the pixels in the surface image, and L is the total number of the voxels in the medium. Because Γ_i,j,k,l is the fluorescence intensity when the fluorophore of unit concentration is placed at the l-th voxel, it corresponds to the sensitivity of the l-th voxel to the measured fluorescence intensity. Then the measured fluorescence intensities represented by a vector, F, are calculated by the product of the sensitivity matrix, H, and the fluorophore concentration vector, f, consisting of the fluorophore concentrations at L voxels as Eq. (3),

F = H f,

where F and f are vectors of 5 × 5 × K and L components, respectively, and H is a matrix of 5 × 5 × K rows and L columns consisting of Γ_i,j,k,l, as Eq. (4),

H = (\begin{matrix} Γ_{1, 1, 1, 1} & \dots & Γ_{1, 1, 1, l} & \dots & Γ_{1, 1, 1, L} \\ Γ_{1, 1, 2, 1} & \dots & Γ_{1, 1, 2, l} & \dots & Γ_{1, 1, 2, L} \\ ⋮ & \dots & ⋮ & \dots & ⋮ \\ Γ_{1, 1, K, 1} & \dots & Γ_{1, 1, K, l} & \dots & Γ_{1, 1, K, L} \\ Γ_{1, 2, 1, 1} & \dots & Γ_{1, 2, 1, l} & \dots & Γ_{1, 2, 1, L} \\ ⋮ & \dots & ⋮ & \dots & ⋮ \\ Γ_{1, 5, K, 1} & \dots & Γ_{1, 5, K, l} & \dots & Γ_{1, 5, 1, L} \\ ⋮ & \dots & ⋮ & \dots & ⋮ \\ Γ_{5, 5, K, 1} & \dots & Γ_{5, 5, K, l} & \dots & Γ_{5, 5, K, L} \end{matrix})

2.3. Image reconstruction with l_p sparsity regularization

Reconstruction of the fluorophore distribution, f, is carried out by minimizing the residual error between the measured data represented by a vector M and the calculated data of the vector F = Hf, and is formulated by the following minimization problem with application of the l_p sparsity regularization,

min_{f} [{‖ M - H f ‖}^{2} + λ \cdot \sum_{l = 1}^{L} {| f_{l} |}^{p}],

where λ is a regularization parameter, and f_l is a component of the vector f. The first and second terms in the bracket of Eq. (5) are the residual error term and the sparsity regularization term, respectively, and the second term is introduced to minimize the l_p norm (0 < p ≤ 1) of the reconstructed f. By this regularization f becomes sparse as the exponential factor, p, approaches zero. Changing the magnitude of p adjusts the sparseness of f so that images can be modified by taking some prior information about the fluorescent target into account.

However, under the condition of 0 < p ≤ 1, there exists a difficulty in calculating gradients of the function in the bracket of Eq. (5) when a gradient-based optimization is used to solve Eq. (5), because |f_l|^p−1 goes to infinity at singular points where the values of f_l are close to zero. To overcome this difficulty, f_l is replaced by the following formulation with another variable z_l used in FOCUSS algorithm [31],

f_{l} = {| z_{l} |}^{2 / p} \cdot sgn (z_{l}) .

Then the minimization problem of Eq. (5) is rewritten as follows,

min_{z} [{‖ M - H f (z) ‖}^{2} + λ \cdot \sum_{l = 1}^{L} {| z_{l} |}^{2}],

where z represents a vector consisting of z_l. Now the gradient of the function in the bracket of Eq. (7) can be calculated without difficulty. A nonlinear conjugate gradient method [34, 35] is employed to solve Eq. (7) with the initial guess of f given by the solution using the Tikhonov regularization which is the case of p = 2 in Eq. (5) [21]. We have confirmed that reconstructed images with an initial guess of a uniform zero distribution are the same as those with the initial guess obtained by use of the Tikhonov regularization. Use of the solution with the Tikhonov regularization as the initial guess accelerates convergence of the image reconstruction process using the l_p sparsity regularization.

In image reconstruction, Intel^® Xenon processor W5580 (8M Cache, 3.20 GHz, 6.40 GT/s Intel^® QPI), RAM with 24 GB and Matlab with Parallel Computing Toolbox were used. Preprocessing including FEM meshing, generation of the sensitivity matrix, H, of Eq. (4) and singular value decomposition of H necessary for reconstruction using the Tikhonov regularization required about 24 hours. Once preprocessing was carried out image reconstruction using the Tikhonov regularization took about 40 s and fifty iterations for reconstruction using the l_p sparsity regularization took about 114 s, resulting in the total time of about 154 s.

3. Simulations and experiments

3.1. Numerical simulations

3.1.1. Single target

Simulation data of the measurement data, M, are generated by solving the fundamental equations (1) and (2) with 3D-FEM by use of COMSOL Multiphysics ver. 3.5a (COMSOL Inc., MA, USA). The object is a cylinder having a diameter of 25 mm and a height of 50 mm with the x-y-z coordinate as shown in Fig. 2. The cylinder is discretized into 24,485 tetrahedral elements. μ_a and μ′_s of the cylinder are given as 0.022 mm⁻¹ and 0.6 mm⁻¹, respectively. The optical properties are given by referring those of the rat liver [36], but μ_a is made smaller than that of rat liver by considering the existence of other tissues and organs. A single fluorescent target of indocyanine green (ICG) is assumed to be placed in the plane (x–y plane) perpendicular to the cylinder axis (z-axis) at a height of 25 mm (z = 0). ICG is a popular fluorescent diagnostic dye approved by Food and Drug Administration (FDA). The peak excitation and emission wavelengths of ICG are 780 and 830 nm, respectively [37–39]. As stated before, one surface view of the fluorescence emission from a cylindrical phantom has K = 25 × 51 = 1275 data, and the total number of the measurement data is 5 × 5 × 25 × 51 = 31875. In the image reconstruction, the cylinder is discretized into L = 25347 voxels with volumes of 1 mm³. Therefore, the matrix H has 31875 rows and 25347 columns. In this case, the number of the measurement data, 31875, is larger than that of the unknowns, 25347, and the ill-posedness of the inversion process might be relaxed.

Fig. 2 The geometry of a cylindrical object, and true and reconstructed x–y sectional images of the ICG concentration normalized by their maxima for Case (i). True images (first column), reconstructed images using the Tikhonov regularization (second column) and using the l_p sparsity regularization with p = 1 (third column) and p = 0.5 (fourth column) for the depths of the true target being 4 mm (top row), 6 mm (middle row) and 9 mm (bottom row).

Download Full Size | PDF

The size, position and ICG quantities (or concentrations) of the target are varied as the following three cases.

Case (i): The target has a volume of 1 mm³ and 100 pmol of ICG, and the depth, which means the distance from the center of the target to the nearest surface point in a radial direction, is varied as 4, 6, and 9 mm. When the coordinate of the target center is expressed by (x_c, y_c, z_c), these varied depths are described by (x_c, y_c, z_c) = (0.0, 8.5, 0.0), (0.0, 6.5, 0.0) and (0.0, 3.5, 0.0), respectively.
Case (ii): The target has a volume of 1.0 mm³, and the quantity of ICG is changed as 100, 10, and 1 pmol with its depth fixed as 6 mm.
Case (iii): The volume of the target is changed as 1, 8, and 27 mm³ with the quantity of ICG and depth of the target fixed as 100 pmol and 6 mm, respectively.

Shot noises with a standard deviation of $\sqrt{c}$ are added to the simulation data consisting of all the pixel data in the simulated surface image where c represents photon counts detected by a CCD camera. Image reconstruction is carried out using the Tikhonov regularization or the l_p sparsity regularization with p = 0.5 or 1. The size of the voxel for reconstruction is 1 mm × 1 mm × 1 mm. The reconstructed quantity of ICG in the target is evaluated with a volume-of-interest (VOI) analysis. The VOI which surrounds the true target position is a cube with a side of 5 mm, and the quantity of ICG in the VOI, R_VOI, is expressed by the following equation,

R_{VOI} = \int_{x_{c} - 2.5}^{x_{c} + 2.5} \int_{y_{c} - 2.5}^{y_{c} + 2.5} \int_{z_{c} - 2.5}^{z_{c} + 2.5} N (x, y, z) d z d y d x .

Sizes of the reconstructed targets are evaluated by the full width at half maximum (FWHM) of the profiles of the ICG quantity at each y-position, Q(y), which is calculated by integrating the ICG concentrations (quantity in a unit volume of 1.0 mm³), N(x, y, z), over a rectangular parallelepiped volume of 1 mm × 1 mm × 50 mm along the z-axis from z = −25 mm to z = 25 mm as expressed by the following,

Q (y) = \int_{- 0.5}^{0.5} \int_{- 25}^{25} N (x, y, z) d z d x .

3.1.2. Multiple targets

Spatial resolution and contrast of the reconstructed targets are evaluated for the cases of two targets which are positioned symmetrically to the y-axis with a distance of d at the fixed y-position in the x–y plane at z = 0. Both targets have the same volume of 1 mm³. The distance, d, and the ICG quantities in the targets are varied in the following cases (iv) and (v).

Case (iv): d is varied as 4.0, 6.0, and 8.0 mm with y_c = 6.5 mm, resulting in the coordinates of the target centers of (x_c, y_c, z_c) = (±2.0, 6.5, 0.0), (±3.0, 6.5, 0.0), and (±4.0, 6.5, 0.0). The quantities of the ICG in the both targets are 100 pmol equally.
Case (v): d is fixed as 6.0 mm with y_c = 6.5 mm, but the quantities of ICG in the targets are different. One target on the left hand side located at (x_c, y_c, z_c) = (−3.0, 6.5, 0.0) contains 100 pmol, and the other on the right hand side located at (x_c, y_c, z_c) = (3.0, 6.5, 0.0) contains 75, 50, or 25 pmol.

For the cases of multiple targets, the profile of the ICG quantity along the line connecting the two target centers, Q(x), is calculated in the manner similar to Eq. (9),

Q (x) = \int_{6}^{7} \int_{- 25}^{25} N (x, y, z) d z d y .

3.2. Phantom experiments

A phantom made of silicone resin with mixed pigment ink and kaolin particles was used for the experiment. The phantom was a cylinder with a diameter of 25 mm and a height of 150 mm, and had the optical properties of μ_a = 0.022 mm⁻¹ and μ′_s = 0.6 mm⁻¹, which were achieved by adjusting the concentrations of the pigment ink (0.31 wt%) and kaolin particle (1.95 wt%). A single fluorescent target with a diameter of 1 mm and a height of 3 mm was located in the plane at the height of 110 mm from the bottom of the phantom. The target contained 13 pmol ICG, and the depth of the target was changed as 4, 6, and 9 mm. The fluorescence surface images of the phantom were acquired using Clairvivo OPT. Image reconstruction was carried out for a partial region of the phantom with a diameter of 25 mm and a height of 50 mm, which included the target. The reconstructed target images were evaluated by a VOI analysis with the size of the VOI of 5 mm × 5 mm × 10 mm.

3.3. Mouse experiment

A hairless mouse (BALB/c, 34 weeks old) was used as an object, which was sacrificed with cervical dislocation, and fixed on the plexiglass table using pieces of adhesive tape at the neck and tail. A target made of a glass tube containing ICG and MRI-contrast agent (Magnevist^®) was embedded in the abdomen of the mouse. Magnevist^® (gadopentetate dimeglumine) is an FDA-approved agent and used clinically for MRI to visualize lesions with abnormal vascularity [40]. The target had inner and outer diameters of 1.0 mm and 1.3 mm, respectively, and a length of 3.8 mm with the ICG quantity of 15.6 pmol, and its position was checked from the MR image of the mouse. In addition, another tube with an inner diameter of 1.0 mm containing the MRI-contrast agent was attached on the skin of the left abdomen of the mouse as a marker for co-registration of FDOT and MR images. The fluorescence surface images of the mouse were acquired using Clairvivo OPT.

A 4.7 T MRI (Varian, Inc., CA, USA) instrument which was located about 50 m apart from Clairvivo OPT was used to obtain T1-weighted MR image of the mouse. The mouse fixed on the plexiglass table was transferred from Clairvivo OPT to the MRI instrument as soon as possible and the MR images were taken about one hour after the measurement of fluorescence surface images with keeping the posture of the mouse on the table.

For the mouse experiment, one fluorescence surface image had K = 532 data (about half that for the phantom experiments), and the total number of the measurement data used for image reconstruction was 5×5×532=13,300. Image reconstruction was carried out for the abdominal region of the mouse which was modeled with a total of L = 6979 cubic voxels with each voxel having a side of 1.0 mm. For forward calculation to provide the matrix H, an FEM mash having 125,613 tetrahedral elements was used to solve the photon diffusion equations (1) and (2). The FEM model of the mouse abdomen was generated from the MR image using a software, ZedView (LEXI Co., Ltd., Tokyo, Japan). The number of the measurement data of 13,300 was again larger than that of the unknowns of 6979. The optical properties of the liver of the mouse, μ′_s = 0.56 mm⁻¹ and μ_a = 0.072 mm⁻¹ [36], were used for the background optical properties, because the target tube was placed near the liver of the mouse.

The coordinates of the FEM model generated from the MR image and the fluorescence imaging system, Clairvivo Opt, were co-registered by a conventional affine transformation using the markers on the mouse skin which were observable in both the MR and fluorescence surface images. The reconstructed fluorescent tomographic images based on the FEM model were co-registered with and superimposed on the MR image.

This experiment was conducted under the approval of the committee of Shimadzu Corporation in accordance with the guideline for animal experiments.

The conditions of the numerical simulations, phantom experiments and mouse experiment are summarized in Tables 1 and 2.

Table 1. Summary of the conditions of simulations with single target and experiments

View Table | View all tables in this article

Table 2. Summary of the conditions of simulations with multiple targets

View Table | View all tables in this article

4. Results and discussions

4.1. Numerical simulations

4.1.1. Single target

Figure 2 shows the true and reconstructed images for Case (i) with the ICG concentrations normalized by their maxima. When the Tikhonov regularization is used, the target is reconstructed at almost the correct position, although the positions of the maximum concentration shift about 1 mm from the true position toward the cylinder axis for the cases of the target depths of 6 mm and 9 mm. The target is reconstructed much broader than the true image with FWHMs calculated from the Q(y) profiles (not shown) ranging from 4 mm to 6 mm. The deeper the target position is, the broader the target is reconstructed.

On the other hand, when the l_p sparsity regularization is used, the reconstructed targets are more sharply localized than using the Tikhonov regularization. The targets were localized at the correct positions with FWHMs of 3 mm and 1 mm when p = 1 and 0.5, respectively. A smaller value of p localizes the target more sharply, and the target reconstructed with p = 0.5 is as sharp as the true image.

Figure 3 shows R_VOI for Case (i). R_VOI using the Tikhonov regularization (blue dashed line) are about 30 % to 45 % of the true value of 100 pmol while they are 88 % to 99 % of the true value when the l_p sparsity regularization (red solid and green chained lines) is used. R_VOI decreases as the target depth increases, probably due to diffusive light propagation and small signal-to-noise ratios of the detected fluorescence emission light. This effect of the target depth is often observed in FDOT and DOT. The results for Case (i) indicate that the l_p sparsity regularization reconstructs the quantity of the fluorophore robustly even when the target is embedded deeply inside the medium.

Fig. 3 Reconstructed quantities of ICG in the VOI using the Tikhonov regularization (blue dashed line) and the l_p sparsity regularization with p = 1 (green chained line) and p = 0.5 (red solid line) as a function of the depth of the target.

Download Full Size | PDF

Figure 4 shows the reconstructed images of the normalized ICG concentration for Case (ii). Similarly for Case (i), the reconstructed targets are localized very well using the l_p sparsity regularization even when the true quantity of ICG reduces from 100 pmol to 1 pmol. The reconstructed target is sharper as p becomes smaller, and is slightly broader as the true quantity becomes smaller. FWHMs are approximately from 5 mm to 7 mm when the Tikhonov regularization is used, while they are approximately 1 mm to 2 mm when p =0.5, and approximately 3 mm when p =1.

Fig. 4 True and reconstructed x–y sectional images of the ICG concentration normalized by their maxima for Case (ii). True images (first column), reconstructed images using the Tikhonov regularization (second column) and the l_p sparsity regularization with p = 1 (third column) and p = 0.5 (fourth column) for the true quantity of ICG in the targets of 100 pmol (top row), 10 pmol (middle row) and 1 pmol (bottom row).

Download Full Size | PDF

The results of the VOI analysis, R_VOI, shown in Fig. 5 indicate that the reconstructed ICG quantity changes almost linearly to the change in the true ICG quantity. When using the l_p sparsity regularization, the value of R_VOI ranges from 84 % to 95 % of the true values which changes from 100 pmol to 1 pmol while it ranges from 31 % to 36 % of the true values when using the Tikhonov regularization. As the true value decreases, the percentage decreases slightly due to lower signal-to-noise ratios of the measurement noise.

Fig. 5 Reconstructed quantities of ICG in the VOI using the Tikhonov regularization (blue dashed line) and the l_p sparsity regularization with p = 1 (green chained line) and p = 0.5 (red solid line) as a function of the true fluorophore quantity.

Download Full Size | PDF

The true and reconstructed images of the normalized ICG concentration for Case (iii) are shown in Fig. 6. The reconstructed targets are well localized when using the l_p sparsity regularization with successful recovery of the ICG quantity. The reconstructed values of R_VOI are 90 % to 95 % and 83 % to 94 % of the true value for p =1 and 0.5, respectively. On the other hand, the values of R_VOI using the Tikhonov regularization are 36 % to 39 % of the true value. When the volume of the true target is 27 mm³, the reconstructed target image using the Tikhonov regularization looks the closest to the correct image, while those using the sparsity regularization with p = 1 and 0.5 look closest to the true images when the volumes of the true targets are 8 mm³ and 1 mm³, respectively. From these results, it seems that the appropriate p value depends on the target size, and prior information obtained by other imaging modalities will be useful for selecting an appropriate p value to reconstruct the target size faithfully.

Fig. 6 True and reconstructed x–y sectional images of the normalized ICG concentration for Case (iii). True images (first column), reconstructed images using the Tikhonov regularization (second column) and the l_p sparsity regularization with p = 1 (third column) and p = 0.5 (fourth column) for the true targets of 1 mm³ (top row), 8 mm³ (middle row) and 27 mm³ (bottom row).

Download Full Size | PDF

4.1.2. Multiple targets

Figure 7 shows the reconstructed images for Case (iv), and Fig. 8 plots the profiles of Q(x) along the line connecting the two target centers. When the distance between the targets is 4 mm, the targets are reconstructed separately only when p = 0.5. The maxima of Q(x) are 3.5, 16.0, and 28.4 pmol for the Tikhonov regularization, p = 1, and p = 0.5, respectively. The targets are reconstructed larger in volume and smaller in concentration than the true ones. The reconstructed quantities of the two targets are almost the same. When the distance between the two targets is 6 mm, the two targets are recognized in all of the reconstructed images. Especially, the two targets are clearly separated when p = 1 and 0.5. The maxima of Q(x) are 3.1, 12.0, and 25.3 pmol for the Tikhonov regularization, p = 1, and 0.5, respectively. The two targets are separated well regardless of the p value when the distance is 8 mm. The maxima of Q(x) are 3.7, 16.6, and 22.3 pmol for the Tikhonov regularization, p = 1, and 0.5, respectively. In all cases, the reconstructed quantities are almost the same for the two targets.

Fig. 7 True and reconstructed x–y sectional images of the normalized ICG concentration for Case (iv). True images (first column), reconstructed images using the Tikhonov regularization (second column) and the l_p sparsity regularization with p = 1 (third column) and p = 0.5 (fourth column) for the distances between the targets of 4 mm (top row), 6 mm (middle row) and 8 mm (bottom row).

Download Full Size | PDF

Fig. 8 The profiles of the reconstructed quantities of ICG along the line connecting the two target centers, Q(x), using the Tikhonov regularization (blue dashed line) and the l_p sparsity regularization with p = 1 (green chained line) and p = 0.5 (red solid line) for Case (iv) when the distance between the two targets is (a) 4mm, (b) 6 mm and (c) 8 mm. The true quantities are 100 pmol for both targets.

Download Full Size | PDF

Figure 9 shows the reconstructed images for Case (v), and Fig. 10 plots the profiles of Q(x) revealing two targets at the fixed positions with different quantities of ICG. The reconstructed ICG quantities of the targets on the right hand side are smaller than those on the left hand side corresponding to the true quantities listed in Table 2. When the Tikhonov regularization is used, the reconstructed targets are not clearly separated similarly to the result for Case (iv). As p decreases, the spatial resolution of the reconstructed targets is improved, and the difference in the reconstructed quantities of ICG between the two targets increases. When the true ICG quantity of the target on the right hand side is 25 pmol, it is difficult to find it in the reconstructed images when p = 1 and 0.5 while it is observable when the Tikhonov regularization is used. This is because the ICG concentrations in the reconstructed images are normalized by the maximum in each image. When the Tikhonov regularization is used, the maximum in the target on the left hand side is very small and close to that on the right hand side as seen in Fig. 10(c).

Fig. 9 True and reconstructed x–y sectional images of the ICG concentration normalized by the maxima of the reconstructed ICG concentration for Case (v) for Case (v). True images (first column), reconstructed images using the Tikhonov regularization (second column) and the l_p sparsity regularization with p = 1 (third column) and p = 0.5 (fourth column). The two targets are at the fixed positions but with different ICG quantities; 100 pmol and 75 pmol (top row), 100 pmol and 50 pmol (middle row), and 100 pmol and 25 pmol (bottom row).

Download Full Size | PDF

Fig. 10 Profiles of Q(x) reconstructed using the Tikhonov regularization (blue dashed line) and the l_p sparsity regularization with p = 1 (green chained line) and p = 0.5 (red solid line) for Case (v). The target on the left hand side contain 100 pmol of ICG while the targets on the right hand side contain (a) 75 pmol, (b) 50 pmol and (c) 25 pmol of ICG.

Download Full Size | PDF

In the cases of multiple targets, the maxima of the reconstructed concentration are smaller than those in the cases of single target. This can be caused by the nature of the l_p sparsity regularization. The l_p sparsity regularization term in Eq. (5) is minimized when the reconstructed fluorophore distribution of the target is localized in a small number of voxels (especially a single voxel) and/or when the reconstructed fluorophore concentrations in the target are close to zero. In the cases of single target, both the residual error and sparsity regularization terms in Eq. (5) can be simultaneously minimized by localizing the fluorophore distribution in a single voxel.

But, in the cases of multiple targets, if the fluorophore distribution is localized in a single voxel, the residual error term increases, and resultantly the sum of the residual error and the l_p sparsity regularization terms is not minimized. In order to minimize the sum, it is necessary to allow the fluorophores to exist at multiple voxels with lower fluorophore concentrations than those in the case of single target. So, the l_p sparsity regularization term becomes smaller by taking smaller fluorophore concentrations. Appropriate adjustment of the FEM mesh by use of sophisticated meshing such as dual-meshing or adaptive-meshing may solve this problem. [41, 42].

4.2. Phantom experiments

The fluorescence surface images superimposed on the visible light images of the phantom are shown in Fig. 11 for the target depths of (a) 4 mm, (b) 6 mm and (c) 9 mm. The color bars of the images show the photon counts of the fluorescence intensity normalized by the maximum for each target depth, and naturally the photon counts sharply decreased with the target depth. Also, the fluorescence surface emission from the target at the deep position was detected in a broader area than that from the target at the shallow position. The main noise of the signals acquired by Clairvivo Opt was the shot noise of the CCD sensor, and the signal-to-noise ratios (SNR) were estimated approximately as 100, 50 and 20 for the cases of the target depths of 4, 6, and 9 mm, respectively from the equation, $SNR = c_{\max} / \sqrt{c_{\max}}$ , where c_max was the maximum photon counts in the fluorescence surface image.

Fig. 11 Fluorescence surface images superimposed on the visible light images of the phantom acquired by Clairvivo OPT with 5 angular positions of the excitation light, Ex 1 to 5 (numbered in Fig. 1(b)), for the targets at the depths of (a) 4 mm, (b) 6 mm and (c) 9 mm. Five views of Emissions 4, 5, 1, 2, and 3 (numbered in Fig. 1(c)) are aligned from left to right in each images. Color bars indicate the photon counts per second.

Download Full Size | PDF

Figure 12 shows the true and reconstructed cross-sectional images of the ICG concentration normalized by their maxima, and Fig. 13 plots R_VOI as a function of the target depth. Similarly to the numerical simulations, the l_p sparsity regularization localized the targets more sharply than the Tikhonov regularization. Smaller p localized targets more sharply, and the deeper targets were reconstructed more broadly. The VOI analysis revealed that the recovered ICG quantities using the l_p sparsity regularization were approximately 18 pmol, 11 pmol and 5 pmol against the true value of 13 pmol for the target depths of 4 mm, 6 mm and 9 mm, respectively, while those using the Thikonov regularization were approximately 3.5 pmol, 4 pmol and 1.5 pmol, respectively. The reasons why the recovered ICG concentrations using the l_p sparsity regularization for the target depth of 4 mm were larger than the true value are not clear, but unintentional heterogeneity of the optical properties in the phantom and possible failure of the diffusion approximation for the case of the target depth of 4 mm may be the reasons.

Fig. 12 True and reconstructed x–y sectional images of the ICG concentration in the phantom experiments for the true targets at the depths of (a) 4 mm (top row), (b) 6 mm (middle row) and (c) 9 mm (bottom row); true images (first column), reconstructed images using the Tikhonov regularization (second column) and the l_p sparsity regularization with p = 1 (third column) and p = 0.5 (fourth column).

Download Full Size | PDF

Fig. 13 Reconstructed quantities of ICG in the VOI using the Tikhonov regularization (blue dashed line) and the l_p sparsity regularization with p = 1 (green chained line) and p = 0.5 (red solid line) as a function of the target depth of the target in the phantom experiment.

Download Full Size | PDF

4.3. Mouse experiment

The fluorescence surface images superimposed on the visible light images of the mouse are shown in Fig. 14. Strong fluorescence emissions were detected around the right ventral of the mouse. The SNR of the surface images in the mouse experiment was approximately 50. Figure 15 shows the reconstructed ICG concentration images superimposed on the MR images. The fluorophore concentrations are normalized by the maximum in each image as before, and the correct positions of the target are indicated in the MR images. The quality of the reconstructed images in the mouse experiment is similar to those in the numerical simulations and phantom experiments. The reconstructed target was localized well by use of the l_p sparsity regularization while it was reconstructed in a broader region by use of the Tikhonov regularization. Figure 16 shows the profiles of Q(y) along the y-axis passing through the reconstructed targets. The maximum of Q(y), 10 pmol, reconstructed by use of the l_p sparsity regularization with p = 0.5 was three times that with p = 1 and ten times that by use of the Tikhonov regularization.

Fig. 14 Fluorescence surface images of the mouse acquired by Clairvivo OPT with 5 angular positions of excitation light, Excitation 1 to 5 (numbered in Fig. 1(b)). Five views of Emission 4, 5, 1, 2, and 3 (numbered in Fig. 1(c)) are aligned from left to right in each image. Color bar indicates the photon counts per second.

Download Full Size | PDF

Fig. 15 MR images (first column) and reconstructed fluorescence tomographic images using the Tikhonov regularization (second column), and the l_p sparsity regularization with p = 1 (third column) and p = 0.5 (fourth column). The sagittal (top row), coronal (middle row) and axial (bottom row) images are shown. Red circles in the MR images indicate the correct position of the target.

Download Full Size | PDF

Fig. 16 Reconstructed quantities of ICG, Q(y), (right) along the y-axis indicated in the reconstructed axial MR image (left). The reconstructed values using the Tikhonov regularization (blue dashed line), and the l_p sparsity regularization with p = 1 (green chained line) and p = 0.5 (red solid line) are shown together with the true values (black dotted line).

Download Full Size | PDF

The spatial resolution was improved by use of the l_p sparsity regularization, although the position of the reconstructed target was about 2 mm away from the true position. The error in the reconstructed target position may have been caused by the difference between the background optical properties used in the reconstruction process and those of the actual mouse, i.e., the background optical properties were assumed homogeneous in the reconstruction process while they were heterogeneous in the actual mouse. Various factors such as movement, temperature and metabolism of mouse will affect measurements when a mouse is alive. In this study, a sacrificed mouse with an embedded glass tube containing ICG was used for the purpose of focusing discussions on the validation of the algorithm and the performances of the imaging system. The present imaging system required 3 minutes to obtain the fluorescent surface images for image reconstruction. Real-time tomographic imaging with shorter acquisition and reconstruction times than in the present study will be a next challenge.

4.4. Selection of regularization parameter

The reconstruction method with the l_p sparsity regularization in this study adjusts the degree of sparsity of fluorophore distributions by varying the value of p appearing in Eq. (6). As p becomes smaller, the power of z, 2/p, becomes larger. Then the variations in the values of z_l at different positions are enlarged more by smaller p when z_l are transformed to f_l. When z_l = 0.1 and 1.0 for example, p = 0.5 provides f_l = 0.0001 and 1.0, respectively, while p = 1 provides f_l = 0.01 and 1.0, respectively. Therefore, localization of the reconstructed target is improved with smaller p through the transformation of z_l to f_l.

The regularization term in Eq. (7) functions to reduce the effects of noises and mismatches between the forward model and actual situation, and smoothens z_l by the manner of the Tikhonov regularization. Throughout this paper, λ in Eq. (7) is given as 5 × 10⁻⁷ which was empirically determined prior to this study by numerical simulations where noisier measurement data were used. Reconstructed images using λ smaller than 5 × 10⁻⁷ (not shown) were not so much different from the images obtained here.

When λ was very small, the reconstructed images were found to be less affected by the noises contained in the measurement data and by the mismatches between the forward model and the actual geometry. This fact indicated that the noises and mismatches in this study were small enough to reconstruct good images without regularization. However, if noises in the data measured by other measuring systems are larger than those by Clairvivo Opt used in this study, larger λ must be used to suppress the effect of the noises, and resultantly the reconstructed target images will be blurred more. When the l_p sparsity regularization is used for noisier measurement data, λ should be selected by a certain method such as the L-curve method.

4.5. Comparison with other sparsity regularization methods

Various methods to reconstruct the sparse distribution of the fluorophores have been proposed. However, there is no standard sparsity regularization method currently, and we have found no report on FDOT reconstruction algorithm using the p-norm (0 < p < 1) as discussed in this paper. Usually, the sparse reconstruction is carried out with minimizing the 1-norm of the solution in the inverse problems for FDOT or DOT, and the non-differentiability of the 1-norm at its singular points requires various techniques proposed in the recent literatures as the following.

By use of barrier penalty, L1 sparsity constraint can be achieved [43]. To obtain the solution converging to an optimal point, the objective function needed a parameter t which gradually changed from zero to infinity to keep the changes in μ_a within a certain range. Although the rate of the change in the parameter t may affect the convergence process and quality of the reconstructed image, the effect of the sparsity constraint explained in the method was well demonstrated by phantom experiments and in vivo functional brain imaging.

A method, which used the L¹ sparsity constraint and total variation method simultaneously, was proposed [44]. The method achieved the sparsity constraint by use of the positivity constraint on the solution. The gradient of the L¹ norm of the fluorophore distribution was equal to the vector which had unities as the components when the reconstructed variables were larger than zero. Practically, the non-differentiability of the L¹ norm penalty near zero was overcome by the method. The method must need a threshold value for the positivity constraint, although the method to constrain the solution was not mentioned in the literature. The L¹ sparsity constraint worked well and improved the sparseness of the reconstructed fluorophore distribution.

A framework of the EM algorithm can be used to obtain the sparse solution which minimizes the L₁ norm. In the literature [20], a linearized inverse problem of DOT was reformulated as a statistic model by assuming that the absorption coefficient was a random variable. The image was reconstructed by maximizing the log-likelihood estimator with the EM algorithm. The maximization (M-step) was done with a soft-threshold method which required a parameter besides the regularization parameter to formulate the statistic model. Numerical simulation showed that the proposed L₁ norm minimization with the EM algorithm reconstructed more sharply localized distributions of the absorption coefficient than other methods such as the Tikhonov regularization and simultaneous iterative reconstruction technique.

The differences among the above methods can be seen in the formulation of the optimized cost functions and in the additional parameters. However, solutions which minimize the L₁ (L1, L¹ or 1-) norm by use of the various methods including the l_p sparsity regularization must be equivalent each other if they work appropriately under the same conditions of the measurement system. Every literature has reported similarly that the L₁ norm minimization has localized the target more sharply than the Tikhonov regularization. And this is consistent with the results in this paper. The major feature of the method proposed in this paper is that the method can be applied to the p-norm with 0 < p < 1 even at singular points where other methods fail to work. The merit of the proposed method is that the method can obtain the fluorophore distribution localized more sharply than other methods which minimize the L₁ norm. This can be a breakthrough to improve the spatial resolution of the reconstructed optical images in biomedicine, although the selection of an appropriate p value is a problem remained for future work.

5. Conclusions

The multi-view CW surface images of the fluorescence emission from the objects were efficiently acquired by the Clairvivo OPT system. The objects were easily installed in the system, the measurements were performed in a short period of time, and the obtained multi-view surface images were used as the input measurement data for image reconstruction of FDOT. The image reconstruction algorithm incorporating the l_p sparsity regularization was developed for efficient reconstruction of the fluorescence targets. Regularization was done by minimizing the l_p norm of the reconstructed distribution of the fluorophore concentration. The reconstructed distribution was transformed by introducing another variable for calculation of the gradient of the p-norm (0 < p ≤ 1) to overcome the non-differentiability at the singular points.

Numerical experiments have proved that the l_p sparsity regularization reconstructs the targets with improved localization and with much narrower distributions of the fluorophore than the Tikhonov regularization. Also the l_p sparsity regularization improves the correctness of the reconstructed fluorophore quantities.

Phantom and mouse experiments have demonstrated that the combination of the multi-view surface image acquisition system and the developed image reconstruction algorithm reconstructed the localized targets successfully. Co-registration of the reconstructed FDOT and MR images showed that FDOT reconstructed the fluorescence targets at the correct positions and with the correct sizes. Even though the heterogeneity of the background optical properties was not taken into account in the forward model for image reconstruction, the fluorescence targets were reconstructed almost at the correct positions. The l_p sparsity regularization alleviated the mismatches between the forward model and the real object. Through the numerical simulations, phantom and mouse experiments, the reconstructed distributions of the fluorophore with p < 1 localized the target more sharply than those with p = 1 as reported in some literatures previously.

The combination of the multi-view CW fluorescence measurements and the reconstruction algorithm with the l_p sparsity regularization will be a useful tool for FDOT in preclinical animal tests for developments of novel drugs.

Acknowledgments

This work was partly supported by JSPS KAKENHI Grant Number 24760314.

References and links

1. V. Ntziachristos, C.-H. Yung, C. Bremerand, and R. Weissleder, “Fluorescence molecular tomography resolves protease activity in vivo,” Nat. Med. 8(7), 757–760 (2002). [CrossRef] [PubMed]

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

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

4. A. P. Gibson, J. C. Hebden, and S. R. Arridge, “Recent advances in diffuse optical imaging,” Phys. Med. Biol. 50, R1–R43 (2005). [CrossRef] [PubMed]

5. R. Weissleder, “Molecular imaging in cancer,” Science 321, 1168–1171 (2006). [CrossRef]

6. D. J. Hawrysz and E. M. Sevick-Muraca, “Developments toward diagnostic breast cancer imaging using near-infrared optical measurements and fluorescent contrast agents,” Neoplasia 2(5), 388–417 (2000). [CrossRef]

7. K. Vishwanath, B. Pogue, and M.-A. Mycek, “Quantitative fluorescence lifetime spectroscopy in turbid media: comparison of theoretical, experimental and computational method,” Phys. Med. Biol. 47, 3387–3405 (2002). [CrossRef] [PubMed]

8. D. Y. Paithankar, U. A. Chen, B. W. Pogue, M. S. Patterson, and E. M. Sevick-Muraca, “Imaging of fluorescent yield and lifetime from multiply scattered light reemitted from random medium,” Appl. Opt. 36(10), 2260–2272 (1997). [CrossRef] [PubMed]

9. E. Graves, J. Ripoll, R. Weissleder, and V. Ntziachristos, “A submillimeter resolution fluorescence molecular imaging system for small animal imaging,” Med. Phys. 30, 901–911 (2003). [CrossRef] [PubMed]

10. N. C. Deliolanis, T. Lasser, D. Hyde, A. Soubret, J. Ripoll, and V. Ntziachristos, “Free-space fluorescence molecular tomography utilizing 360° geometry projections,” Opt. Lett. 32(4), 382–384 (2007). [CrossRef] [PubMed]

11. T. Yates, C. Hebdan, A. Gibson, N. Everdell, S. R. Arridge, and M. Douek, “Optical tomography of the breast using a multi-channel time-resolved imager,” Phys. Med. Biol. 50, 2503–2517 (2005). [CrossRef] [PubMed]

12. A. P. Gibson, T. Austin, N. L. Everdell, M. Schweiger, S. R. Arridge, J. H. Meek, J. S. Wyatt, D. T. Delpy, and J. C. Hebden, “Three-dimensional whole-head optical tomography for passive motor evoked responses in the neonate,” NueroImage 30, 521–528 (2006). [CrossRef]

13. J. C. Hebden, A. Gibson, R. M. Yusof, N. Everdell, E. M. C. Hillman, D. T. Delpy, S. R. Arridge, T. Austin, J. H. Meek, and J. S. Wyatt, “Three-dimensional optical tomography of the premature infant brain,” Phys. Med. Biol. 47, 4155–4166 (2002). [CrossRef] [PubMed]

14. B. W. Pogue, T. O. McBride, J. Prewitt, U. Lösterberg, and K. D. Paulsen, “Spatially variant regularization improves diffuse optical tomography,” Appl. Opt. 38(13), 2950–2961 (1999). [CrossRef]

15. G. Boverman, E. L. Miller, A. Li, Q. Zhang, T. Chaves, D. H. Brooks, and D. A. Boas, “Quantitative spectroscopic optical tomography of the breast guided by imperfect a priori structural information,” Phys. Med. Biol. 50, 3941–3956 (2005). [CrossRef] [PubMed]

16. P. K. Yalavarthy, B. W. Pogue, H. Dehghani, C. M. Carpenter, S. Jiang, and K. D. Paulsen, “Structural information within regularization matrices improves near infrared diffuse optical tomography,” Opt. Express 15(13), 8043–8058 (2007). [CrossRef] [PubMed]

17. A. Douiri, M. Schweiger, J. Riley, and S. R. Arridge, “Anisotropic diffusion regularization methods for diffuse optical tomography using edge prior information,” Meas. Sci. Tech. 18, 87–95 (2007). [CrossRef]

18. P. Hiltunen, D. Calvetti, and E. Somersalo, “An adaptive smoothness regularization algorithm for optical tomography,” Opt. Express 16(24), 19957–19977 (2008). [CrossRef] [PubMed]

19. C. Panagiotou, S. Somayajula, A. P. Gibson, M. Schweiger, R. M. Leahy, and S. R. Arridge, “Information theoretic regularization in diffuse optical tomography,” J. Opt. Soc. Am. A 26(5), 1277–1290 (2009). [CrossRef]

20. N. Cao, A. Nehorai, and M. Jacob, “Image reconstruction for diffuse optical tomography using sparsity regularization and expectation-maximization algorithm,” Opt. Express , 15(21), 13695–13708 (2007). [CrossRef] [PubMed]

21. T. Shimokawa, T. Kosaka, O. Yamashita, N. Hiroe, T. Amita, Y. Inoue, and M. Sato, “Hierarchical Bayesian estimation improves depth accuracy and spatial resolution of diffuse optical tomography,” Opt. Express 20(18), 20427–20446 (2012). [CrossRef] [PubMed]

22. P. Xu, Y. Tian, H. Chen, and D. Yao, “L_p Norm Iterative Sparse Solution for EEG Source Localization,” IEEE Trans. Biomed. Eng. 54(3), 400–409 (2007). [CrossRef] [PubMed]

23. P. M. Shankar and M. A. Neifeld, “Sparsity constrained regularization for multiframe image restoration,” J. Opt. Soc. Am. A 25(5), 1199–1214 (2008). [CrossRef]

24. M. Freiberger, C. Clason, and H. Scharfetter, “Total variation regularization for nonlinear fluorescence tomography with an augmented Lagrangian splitting approach,” Appl. Opt. 49(19), 3741–3747 (2010). [CrossRef] [PubMed]

25. D. Han, X. Yang, K. Liu, C. Qin, B. Zhang, X. Ma, and J. Tian, “Efficient reconstruction method for L1 regularization in fluorescence molecular tomography,” Appl. Opt 49(36), 6930–6937 (2010). [CrossRef] [PubMed]

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

27. H. Yi, D. Chen, X. Qu, K. Peng, X. Chen, Y. Zhou, J. Tian, and J. Liang, “Multilevel, hybrid regularization method for reconstruction of florescent molecular tomography,” Appl. Opt. 51(7), 975–986 (2012). [CrossRef] [PubMed]

28. P. Mohajerani, A. A. Eftekhar, J. Huang, and A. Adibi, “Optimal sparse solution for fluorescent diffuse optical tomography: theory and phantom experimental results,” Appl. Opt. 46(10), 1679–1685 (2007). [CrossRef] [PubMed]

29. Y. Lu, X. Zhang, A. Douraghy, D. Stout, J. Tian, T. F. Chan, and A. F. Chatziioannou, “Source reconstruction for spectrally-resolved bioluminescence tomography with aparse A priori information,” Opt. Express 17(10), 8062–8088 (2009). [CrossRef] [PubMed]

30. S. Okawa and Y. Yamada, “Reconstruction of fluorescence/bioluminescence sources in biological medium with spatial filter,” Opt. Express 18(12), 13151–13172 (2010). [CrossRef] [PubMed]

31. Z. He, A. Cichocki, R. Zdunek, and S. Xie, “Improved FOCCUS method with conjugate gradient iterations,” IEEE Trans. Signal Process. 57(1), 399–404 (2009). [CrossRef]

32. S. Okawa, Y. Hoshi, and Y. Yamada, “Improvement of image quality of time-domain diffuse optical tomography with l_p sparsity regularization,” Biomed. Opt. Express 2(12), 3334–3348 (2011). [CrossRef] [PubMed]

33. M. Schweiger, S. R. Arridge, and D. T. Delpy, “Application of the finite-element method for the forward and inverse model in optical tomography,” J. Math. Imaging Vis. 3, 263–283 (1993). [CrossRef]

34. C. R. Vogel, Computational Methods for Inverse Problems (Frontiers in Applied Mathematics) (SIAM, Philadelphia, 2002). [CrossRef]

35. S. R. Arridge, “A gradient-based optimization scheme for optical tomography,” Opt. Express 12(6), 213–226 (1998). [CrossRef]

36. W. F. Cheong, S. A. Prahl, and A. J. Welch, “A review of the optical properties of biological tissues,” IEEE J. Quantum Electron. 26, 2166–2185 (1990). [CrossRef]

37. A. B. Milstein, S. Oh, K. Webb, C. A. Bouman, Q. Zhang, D. A. Boas, and R. P. Millane, “Fluorescence optical diffusion tomography,” Appl. Opt. 42(16), 3081–3094 (2003). [CrossRef] [PubMed]

38. A. Marjono, A. Yano, S. Okawa, F. Gao, and Y. Yamada, “Total light approach of time-domain fluorescence diffuse optical tomography,” Opt. Express , 16(19), 15268–15285 (2008). [CrossRef] [PubMed]

39. B. Toczylowska, E. Zieminska, G. Goch, D. Milej, A. Gerega, and A. Liebert, “Neurotoxic effects of indocyanine green-cerebellar granule cell culture viability study,” Biomed. Opt. Express , 5(3), 800–816 (2014). [CrossRef] [PubMed]

40. J. Barkhausen, W. Ebert, J. F. Debatin, and H.-J. Weinmann, “Imaging of myocardial infarction: comparison of magnevist and gadophrin-3 in rabbits,” J. Am. Coll. Cardiol. 39(8), 1392–1398 (2002). [CrossRef] [PubMed]

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

42. M. Huang and Q. Zhu, “Dual-mesh optical tomography reconstruction method with a depth correction that uses a priori ultrasound information,” Appl. Opt. 43(8), 1654–1662 (2006). [CrossRef]

43. V. C. Kavuri, Z.-J. Lin, F. Tian, and H. Liu, “Sparsity enhanced spatial resolution and depth localization n diffuse optical tomography,” Biomed. Opt. Express , 3(5), 943–957 (2012). [CrossRef] [PubMed]

44. J. Dutta, S. Ahn, C. Li, S. R. Cherry, and R. M. Leahy, “Joint L¹ and total variation regularization for fluorescence molecular tomography,” Phys. Med. Biol. 57, 1459–1476 (2012). [CrossRef] [PubMed]

Reconstruction of localized fluorescent target from multi-view continuous-wave surface images of small animal with l_p sparsity regularization

Abstract

1. Introduction