Self-calibrated algorithms for diffuse optical tomography and bioluminescence tomography using relative transmission images

doi:10.1364/BOE.3.002794

« Show journal navigation

Self-calibrated algorithms for diffuse optical tomography and bioluminescence tomography using relative transmission images

Mohamed A. Naser, Michael S. Patterson, and John W. Wong »View Author Affiliations

Biomedical Optics Express, Vol. 3, Issue 11, pp. 2794-2808 (2012)
http://dx.doi.org/10.1364/BOE.3.002794

View Full Text Article

Acrobat PDF (1727 KB)

Journal Search
Article Lookup

Article Tools

Citations

Alert me when this article is cited
Export Citation/Save

Article Navigation

▸ Article Outline
– Abstract
1. Introduction
2. DOT algorithm
3. BLT algorithm
4. Results and discussion
5. Conclusion
– References and links

Article Tools
Full Text
Article Info
References
Cited By
Figures
Metrics
Download PDF

Viewing Equations in the
HTML Article

Optimal Viewing Recommendations

Full Text
Article Info
References (38)
Cited By
Figures (5)
Metrics

Abstract

Reconstruction algorithms for diffuse optical tomography (DOT) and bioluminescence tomography (BLT) have been developed based on diffusion theory. The algorithms numerically solve the diffusion equation using the finite element method. The direct measurements of the uncalibrated light fluence rates by a camera are used for the reconstructions. The DOT is self-calibrated by using all possible pairs of transmission images obtained with external sources along with the relative values of the simulated data and the calculated Jacobian. The reconstruction is done in the relative domain with the cancelation of any geometrical or optical factors. The transmission measurements for the DOT are used for calibrating the bioluminescence measurements at each wavelength and then a normalized system of equations is built up which is self-calibrated for the BLT. The algorithms have been applied to a three dimensional model of the mouse (MOBY) segmented into tissue regions which are assumed to have uniform optical properties. The DOT uses the direct method for calculating the Jacobian. The BLT uses a reduced space of eigenvectors of the Green's function with iterative shrinking of the permissible source region. The reconstruction results of the DOT and BLT algorithms show good agreement with the actual values when using either absolute or relative data. Even a small calibration error causes significant degradation of the reconstructions based on absolute data.

1. Introduction

Bioluminescence imaging (BLI) is a noninvasive technique that can monitor biological processes at a molecular level and visualize disease progression and response to treatment [1

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

–4

J. Tian, J. Bai, X. P. Yan, S. Bao, Y. Li, W. Liang, and X. Yang, “Multimodality molecular imaging,” IEEE Eng. Med. Biol. Mag. 27(5), 48–57 (2008). [CrossRef] [PubMed]

]. In its simplest form, a sensitive camera is used to image light that reaches the surface from the internal bioluminescence source. These images can provide useful information but they are not quantitative. Bioluminescence tomography (BLT) attempts to reconstruct the three dimensional distribution of the bioluminescent source (i.e. optical power per unit volume) and to provide true quantitative information about its magnitude and location [5

A. J. Chaudhari, F. Darvas, J. R. Bading, R. A. Moats, P. S. Conti, D. J. Smith, S. R. Cherry, and R. M. Leahy, “Hyperspectral and multispectral bioluminescence optical tomography for small animal imaging,” Phys. Med. Biol. 50(23), 5421–5441 (2005). [CrossRef] [PubMed]

–9

H. Dehghani, S. C. Davis, S. Jiang, B. W. Pogue, K. D. Paulsen, and M. S. Patterson, “Spectrally resolved bioluminescence optical tomography,” Opt. Lett. 31(3), 365–367 (2006). [CrossRef] [PubMed]

]. Many different strategies have been proposed for BLT [10

W. Cong, G. Wang, D. Kumar, Y. Liu, M. Jiang, L. V. Wang, E. A. Hoffman, G. McLennan, P. B. McCray, J. Zabner, and A. Cong, “Practical reconstruction method for bioluminescence tomography,” Opt. Express 13(18), 6756–6771 (2005). [CrossRef] [PubMed]

–18

N. V. Slavine, M. A. Lewis, E. Richer, and P. P. Antich, “Iterative reconstruction method for light emitting sources based on the diffusion equation,” Med. Phys. 33(1), 61–68 (2006). [CrossRef] [PubMed]

] but, in general, three components are necessary: a forward model of light propagation from the internal source to the animal surface, a calibration of the camera that links the image (e.g. counts per pixel) to the emittance (optical power per unit area), and a reconstruction algorithm, usually iterative, that estimates the source distribution that minimizes the difference between the observed images and those calculated by application of the forward model.

The second component, camera calibration, typically requires three sub-components: an assumption that the emission from the animal surface is Lambertian [19

X. Chen, X. Gao, D. Chen, X. Ma, X. Zhao, M. Shen, X. Li, X. Qu, J. Liang, J. Ripoll, and J. Tian, “3D reconstruction of light flux distribution on arbitrary surfaces from 2D multi-photographic images,” Opt. Express 18(19), 19876–19893 (2010). [CrossRef] [PubMed]

], modeling of the camera optics that transfers the light to the image plane, and measurement of the absolute camera response on a pixel-by-pixel basis. Recent work by Kienle and Foschum showed that the Lambertian approximation is a poor description of the true surface emission [20

A. Kienle and F. Foschum, “250 years Lambert surface: does it really exist?” Opt. Express 19(5), 3881–3889 (2011). [CrossRef] [PubMed]

]. In addition, measurement of the camera response and modeling of the optics are challenging technical problems. In this paper we propose an alternative strategy for BLT. In addition to the bioluminescence image, a transmission image (or images) is acquired using an external fiber optic source (or sources) of known optical power. The reconstruction algorithm uses the ratio of the internal source image to the external source image as input data. Because the angular dependence of the emission at any point on the surface is very likely the same for both sources, and because the same camera is used to capture both images, the camera calibration problem is reduced to a simple mapping of the animal surface to the camera image plane. In our system [21

J. Pekar, “Multispectral bioluminescence tomography with x-ray CT spatial priors” Ph.D. thesis (McMaster University, Hamilton, 2011), Open Access Dissertations and Theses, Paper 4329, http://digitalcommons.mcmaster.ca/opendissertations/4329.

] the surface anatomy is acquired by an integral x-ray CT scanner but other methods, such as structured light imaging, have been applied [22

M. A. J. Chaudhari, A. A. Joshi, F. Darvas, and R. M. Leahy, “A method for atlas-based volumetric registration with surface constraints for optical bioluminescence tomography in small animal imaging,” Proc. SPIE 6510, 651024, 651024–10 (2007). [CrossRef]

Several investigators have concluded that better BLT reconstructions can be achieved if forward calculations take into account the true three dimensional distribution of scattering and absorption coefficients [23

H. Dehghani, B. W. Pogue, S. C. Davis, and M. S. Patterson, “Modeling and image reconstruction in spectrally resolved bioluminescence tomography,” Proc. SPIE 6434, 64340V, 64340V–9 (2007). [CrossRef]

,24

G. Alexandrakis, F. R. Rannou, and A. F. Chatziioannou, “Tomographic bioluminescence imaging by use of a combined optical-PET (OPET) system: a computer simulation feasibility study,” Phys. Med. Biol. 50(17), 4225–4241 (2005). [CrossRef] [PubMed]

]. In previous papers [25

M. A. Naser and M. S. Patterson, “Algorithms for bioluminescence tomography incorporating anatomical information and reconstruction of tissue optical properties,” Biomed. Opt. Express 1(2), 512–526 (2010). [CrossRef] [PubMed]

,26

M. A. Naser and M. S. Patterson, “Improved bioluminescence and fluorescence reconstruction algorithms using diffuse optical tomography, normalized data, and optimized selection of the permissible source region,” Biomed. Opt. Express 2(1), 169–184 (2011). [CrossRef] [PubMed]

] we have described a strategy whereby this information can be obtained by diffuse optical tomography (DOT). Transmission images acquired with external fiber optics sources and the BLI camera can be used to reconstruct the optical properties which, in turn, are used in the BLT algorithm. Our previous algorithms assumed that a perfect camera calibration (as described in the previous paragraph) was available to link the transmission images to the actual emittance produced by each source. In this paper we propose an alternative strategy for DOT: instead of “absolute data”, ratio images or “relative data” obtained with different source positions are used. If, for example, nine external sources are employed and nine transmission images are acquired, 36 ratio images can be formed and used as input to a DOT algorithm. The relative power of each of the fiber optic sources must also be known. In the simulations presented in this paper we demonstrate that the optical property reconstructions based on these ratios are generally not as accurate as those based on absolute calibration. However, the BLT solutions based on the relative data are comparable in accuracy to those based on absolute data. In addition, we show by example that the relative method outperforms the absolute method if the calibration process is imperfect.

2. DOT algorithm

Although we will use only the continuous wave (CW) solution, for generality we begin with the frequency domain diffusion equation and the Robin-type boundary condition at a specific wavelength given by [27

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

,28

A. D. Klose, “Transport-theory-based stochastic image reconstruction of bioluminescence sources,” J. Opt. Soc. Am. A 24(6), 1601–1608 (2007). [CrossRef]

]

[- \nabla \cdot κ (r, λ) \nabla + (μ_{a} (r, λ) + i \frac{ω}{c})] ϕ (r, ω, λ) = s (r, ω, λ) (r \in Ω),

(1)

[1 + 2 κ (r, λ) ξ \hat{n} (r) μ_{a} (r, λ)] ϕ (r, λ, ω) = 0 (r \in \partial Ω),

(2)

where

φ (r, λ, ω)

is the light fluence rate at position r, wavelength λ, and light modulation frequency f given by

ω = 2 π f

; Ω and

\partial Ω

are the domain and its boundary, respectively, and the source distribution is given by

s (r, ω, λ)

. The spatial distribution of the tissue optical properties at wavelength λ is given by the absorption coefficient

μ_{a} (r, λ)

and the diffusion coefficient

κ (r, λ)

, where the diffusion coefficient is defined by κ(r, λ) =

{[3 (μ_{a} (r, λ) + {μ^{'}}_{s} (r, λ))]}^{- 1}

, where

{μ^{'}}_{s} (r, λ)

is the reduced scattering coefficient

{μ^{'}}_{s} = μ_{s} (1 - g),

μ_{s}

is the scattering coefficient, and g is the anisotropy factor.

\hat{n} (r)

is a unit vector pointed outwardly normal to

\partial Ω

, and ξ is derived from Fresnel’s law as [29

H. Dehghani, M. E. Eames, P. K. Yalavarthy, S. C. Davis, S. Srinivasan, C. M. Carpenter, B. W. Pogue, and K. D. Paulsen, “Near infrared optical tomography using NIRFAST: Algorithm for numerical model and image reconstruction,” Commun. Numer. Methods Eng. 25(6), 711–732 (2009). [CrossRef] [PubMed]

]

ξ = ((2 / (1 - R_{0})) - 1 + {| \cos θ_{c} |}^{3}) / (1 - {| \cos θ_{c} |}^{2})

, where the critical angle

θ_{c} = \sin^{- 1} (1 / n)

and

R_{0} = {(n - 1)}^{2} / {(n + 1)}^{2}

, n is the tissue refractive index, and the refractive index of the surrounding medium is assumed to be 1.

Equation (1) can be solved numerically using the method of finite elements (FE) [30

S. R. Arridge, M. Schweiger, M. Hiraoka, and D. T. Delpy, “A finite element approach for modeling photon transport in tissue,” Med. Phys. 20(2), 299–309 (1993). [CrossRef] [PubMed]

–32

H. Jiang, “Frequency-domain fluorescent diffusion tomography: a finite-element-based algorithm and simulations,” Appl. Opt. 37(22), 5337–5343 (1998). [CrossRef] [PubMed]

], and the discrete version of Eq. (1) is given by

A (λ) ϕ (λ) = s (λ),

(3)

where A is the forward matrix obtained by discretizing the diffusion operator in Eq. (1). Assuming the object can be segmented into N regions where each element within the region has the same scattering and absorption coefficients, the number of unknown optical coefficients is 2N and the optical properties vector can be written in its transpose form as

μ^{T} = [{μ^{'}}_{s 1}, \dots, {μ^{'}}_{s N}, μ_{a 1}, \dots, μ_{a N}],

where

{μ^{'}}_{s n}

and

μ_{a n}

are the scattering and absorption coefficients of region n. Differentiating both sides of Eq. (3) with respect to

μ_{n}

gives

\frac{\partial A}{\partial μ_{n}} ϕ + A \frac{\partial ϕ}{\partial μ_{n}} = \frac{\partial s}{\partial μ_{n}} = 0,

(4)

where the differentiation of the source with respect to the optical properties gives zero, since the source is not a function of the optical properties. The derivative of the forward matrix A with respect to

μ_{n}

is obtained by changing

μ_{n}

\pm Δ μ_{n}

and calculating the difference in the forward matrix

Δ A

. The first order derivative of the light fluence rate with respect to the optical properties is obtained from Eq. (4):

\frac{\partial ϕ}{\partial μ_{n}} = - A^{- 1} \frac{\partial A}{\partial μ_{n}} ϕ .

(5)

The first order derivative of the light fluence rate with respect to the optical properties at the detectors around the object is given by

\frac{\partial ϕ_{d}}{\partial μ_{n}} = \frac{\partial ϕ (d)}{\partial μ_{n}} = - A^{- 1} (d, :) \frac{\partial A}{\partial μ_{n}} ϕ,

(6)

where d is the detector index and

ϕ_{d}

is a vector of length

N_{d}

, where

N_{d}

is the number of detectors and the matrix

A^{- 1} (d, :)

consists of the d rows and all columns of the matrix

A^{- 1}

. If the detectors do not coincide with the FE nodes at the boundary, interpolation can be used to get the value from the surrounding nodes. It is assumed that the values of the light fluence rate at the detector locations can be obtained from the measurements of the corresponding pixels of the CCD camera. For a phasor representation of the light fluence rate where

ϕ_{d} = φ_{d} e^{i θ_{d}}

, Eq. (6) can be written as

\begin{array}{l} \frac{\partial ϕ_{d}}{\partial μ_{n}} = \frac{\partial φ_{d}}{\partial μ_{n}} e^{i θ_{d}} + i φ_{d} e^{i θ_{d}} \frac{\partial θ_{d}}{\partial μ_{n}} \\ \Rightarrow \frac{\partial ϕ_{d}}{\partial μ_{n}} / ϕ_{d} = \frac{\partial φ_{d}}{\partial μ_{n}} / φ_{d} + i \frac{\partial θ_{d}}{\partial μ_{n}} = \frac{\partial \log (φ_{d})}{\partial μ_{n}} + i \frac{\partial θ_{d}}{\partial μ_{n}} . \end{array}

(7)

Therefore

\begin{array}{l} \frac{\partial \log (φ_{d})}{\partial μ_{n}} = real (\frac{\partial ϕ_{d}}{\partial μ_{n}} / ϕ_{d}), \\ \frac{\partial θ_{d}}{\partial μ_{n}} = imag (\frac{\partial ϕ_{d}}{\partial μ_{n}} / ϕ_{d}) . \end{array}

(8)

The light fluence rate data captured by the detectors around the object can be approximated by a first order Taylor expansion such that

\begin{array}{l} \log (φ_{d}) = \log (φ_{d 0}) + \sum_{n} \frac{\partial \log (φ_{d})}{\partial μ_{n}} δ μ_{n}, \\ θ_{d} = θ_{d 0} + \sum_{n} \frac{\partial θ_{d}}{\partial μ_{n}} δ μ_{n}, \end{array}

(9)

where

φ_{d 0}

and

θ_{d 0}

are the magnitude and phase of the light fluence rate corresponding to the initial guess

μ_{0}

. In matrix form, Eq. (9) is given by

(\begin{matrix} \log (φ_{d}) \\ θ_{d} \end{matrix}) = (\begin{matrix} \log (φ_{d 0}) \\ θ_{d 0} \end{matrix}) + (\begin{matrix} \frac{\partial \log (φ_{d})}{\partial μ_{1}} & \dots & \frac{\partial \log (φ_{d})}{\partial μ_{2 N}} \\ \frac{\partial θ_{d}}{\partial μ_{1}} & \dots & \frac{\partial θ_{d}}{\partial μ_{2 N}} \end{matrix}) (\begin{matrix} \partial μ_{1} \\ \begin{array}{l} ⋮ \\ ⋮ \end{array} \\ \partial μ_{2 N} \end{matrix})

(10)

which can be concisely written as

Φ_{d} = Φ_{d 0} + J d μ .

(11)

For multiple external sources, measured data are obtained for each source and so Eq. (11) for

N_{s}

sources is modified to

(\begin{matrix} Φ_{d}^{1} \\ ⋮ \\ Φ_{d}^{N_{s}} \end{matrix}) = (\begin{matrix} Φ_{d 0}^{1} \\ ⋮ \\ Φ_{d 0}^{N_{s}} \end{matrix}) + (\begin{matrix} J_{1} \\ ⋮ \\ J_{N_{s}} \end{matrix}) d μ,

(12)

or, concisely,

{\bar{Φ}}_{d} = {\bar{Φ}}_{d 0} + \bar{J} d μ

. This equation can be solved iteratively for

d μ

such that

d μ = {({\bar{J}}^{T} \bar{J})}^{- 1} ({\bar{J}}^{T} ({\bar{Φ}}_{d} - {\bar{Φ}}_{d 0})) .

(13)

Equation (13) can be used to reconstruct the optical properties for absolute (or completely calibrated) data, where the magnitude and phase of the light fluence rate at the object boundary can be obtained by heterodyne detection. In the case of non-contact camera detection, usually only the magnitude is measured (i.e., f = 0).

If the measured light fluence rate is not calibrated, the magnitude of the light fluence rate at the detectors is multiplied by an unknown vector of calibration coefficients α and the phase is shifted by an unknown vector of phase delay β. Therefore the measured light fluence rate

ϕ_{m}

is given by

\begin{array}{l} ϕ_{m} = φ_{m} e^{i θ_{m}} = α φ_{d} e^{i (θ_{d} + β)} \\ \Rightarrow (\begin{matrix} \log (φ_{m}) \\ θ_{m} \end{matrix}) = (\begin{matrix} \log (φ_{d}) \\ θ_{d} \end{matrix}) + (\begin{matrix} \log (α) \\ β \end{matrix}), \\ or Φ_{m} = Φ_{d} + (\begin{matrix} \log (α) \\ β \end{matrix}), \end{array}

(14)

and for multiple source measurements

(\begin{matrix} Φ_{m}^{1} \\ \begin{array}{l} ⋮ \end{array} \\ Φ_{m}^{N_{s}} \end{matrix}) = (\begin{matrix} Φ_{d}^{1} + (\begin{matrix} \log (α) \\ β \end{matrix}) \\ ⋮ \\ Φ_{d}^{N_{s}} + (\begin{matrix} \log (α) \\ β \end{matrix}) \end{matrix}) .

(15)

A system of equations to reconstruct the optical properties as in Eq. (12) can be obtained by using the relative uncalibrated measured data to remove the magnification factor vector α and the phase delay vector β such that

(\begin{matrix} Φ_{m}^{2} - Φ_{m}^{1} \\ Φ_{m}^{3} - Φ_{m}^{1} \\ ⋮ \\ Φ_{m}^{N_{s}} - Φ_{m}^{1} \\ Φ_{m}^{3} - Φ_{m}^{2} \\ Φ_{m}^{4} - Φ_{m}^{2} \\ ⋮ \\ Φ_{m}^{N_{s}} - Φ_{m}^{2} \\ ⋮ \\ Φ_{m}^{N_{s}} - Φ_{m}^{N_{s} - 1} \end{matrix}) = (\begin{matrix} Φ_{d}^{2} - Φ_{d}^{1} \\ Φ_{d}^{3} - Φ_{d}^{1} \\ ⋮ \\ Φ_{d}^{N_{s}} - Φ_{d}^{1} \\ Φ_{d}^{3} - Φ_{d}^{2} \\ Φ_{d}^{4} - Φ_{d}^{2} \\ ⋮ \\ Φ_{d}^{N_{s}} - Φ_{d}^{2} \\ ⋮ \\ Φ_{d}^{N_{s}} - Φ_{d}^{N_{s} - 1} \end{matrix}) = (\begin{matrix} Φ_{d 0}^{2} - Φ_{d 0}^{1} \\ Φ_{d 0}^{3} - Φ_{d 0}^{1} \\ ⋮ \\ Φ_{d 0}^{N_{s}} - Φ_{d 0}^{1} \\ Φ_{d 0}^{3} - Φ_{d 0}^{2} \\ Φ_{d 0}^{4} - Φ_{d 0}^{2} \\ ⋮ \\ Φ_{d 0}^{N_{s}} - Φ_{d 0}^{2} \\ ⋮ \\ Φ_{d 0}^{N_{s}} - Φ_{d 0}^{N_{s} - 1} \end{matrix}) + (\begin{matrix} J_{2} - J_{1} \\ J_{3} - J_{1} \\ ⋮ \\ J_{N_{s}} - J_{1} \\ J_{3} - J_{2} \\ J_{4} - J_{2} \\ ⋮ \\ J_{N_{s}} - J_{2} \\ ⋮ \\ J_{N_{s}} - J_{N_{s} - 1} \end{matrix}) d μ,

(16)

or, concisely,

{\bar{\bar{Φ}}}_{m} = {\bar{\bar{Φ}}}_{d 0} + \bar{\bar{J}} d μ

Therefore, the optical properties can be reconstructed iteratively from the relative measured data

{\bar{\bar{Φ}}}_{m}

and the relative simulated data

{\bar{\bar{Φ}}}_{d 0}

at the initial guess

μ_{0}

by using the following equation:

d μ = {({\bar{\bar{J}}}^{T} \bar{\bar{J}})}^{- 1} ({\bar{\bar{J}}}^{T} ({\bar{\bar{Φ}}}_{m} - {\bar{\bar{Φ}}}_{d 0})) .

(17)

Equation (17) can be used to reconstruct the optical properties by using the relative measurement, and the reconstructed optical properties at iteration i are given by

μ_{i} = μ_{i - 1} + d μ_{i} .

(18)

The relative simulated data

{\bar{\bar{Φ}}}_{d 0}

are recalculated at every iteration, using the updated optical properties

μ_{i}

. The relative Jacobian

\bar{\bar{J}}

can be updated at every iteration if there is high nonlinearity or less frequently if there are no large changes in the optical properties during these iterations. The solution is obtained when the maximum relative error between the relative measured data

{\bar{\bar{Φ}}}_{m}

and the relative simulated data

{\bar{\bar{Φ}}}_{d 0}

goes below a certain limit.

3. BLT algorithm

The diffusion equation in Eq. (3) can be written as [26

]

φ (λ) = G (λ) s (λ),

(19)

where

G (λ)

is the Green’s function matrix at wavelength λ and s(λ) is the bioluminescence source distribution. Assuming that the source magnitude is not a function of wavelength or, equivalently, the source spectrum is known in the range of measurement, the light fluence rate at the detectors on the object surface is given by

φ (d; λ) = G (d, R; λ) s (R) or φ_{d} (λ) = G_{d R} (λ) s_{R},

(20)

where d is the detector index at the boundary

\partial Ω

, and R is the permissible source region index which contains the index of nodes of the finite element mesh inside the domain Ω.

The light fluence rate measured by the CCD camera is related to the light fluence rate at the “detectors” on

\partial Ω

φ_{m} = α \cdot φ_{d} \Rightarrow φ_{m} = α \cdot G_{d R} s_{R},

(21)

where α is an unknown scaling vector that represents all calibrations required to obtain the right values of the light fluence rate at the detectors. Using Eq. (21), the logarithm of the measured light fluence rates corresponding to the bioluminescence source

s_{R}

and the external sources

s_{1}, s_{2}, \dots, s_{N_{s}}

are given by

\begin{array}{l} \log φ_{m}^{s_{R}} = \log α + \log (G_{d R} s_{R}) \\ \log φ_{m}^{s_{1}} = \log α + \log (G_{d} s_{1}) \\ \log φ_{m}^{s_{2}} = \log α + \log (G_{d} s_{2}) \\ ⋮ \\ \log φ_{m}^{s_{N_{s}}} = \log α + \log (G_{d} s_{N_{s}}) \end{array}

(22)

where

G_{d} = G (d, :)

is the Green's function where all columns of the matrix are included, which means the permissible region is the whole object. The external sources

s_{1}, s_{2}, \dots, s_{N_{s}}

used for the DOT reconstruction and the measurements resulted from these sources are also used for calibrating the bioluminescence measurements. The calibrating vector α can be removed from Eq. (21), and a calibrated equation can be obtained as

G_{d R} s_{R} = \exp (\log φ_{m}^{s_{R}} - \frac{1}{N_{s}} \sum_{i = 1}^{N_{s}} (\log φ_{m}^{s_{i}} - \log (G_{d} s_{i}))) \Rightarrow G_{d R} s_{R} = {\bar{φ}}_{m}^{s_{R}},

(23)

where

{\bar{φ}}_{m}^{s_{R}}

are the calibrated bioluminescence measurements. For multi-spectral measurements, Eq. (23) for N wavelengths in a normalized form can be given by

(\begin{matrix} G_{d R} (λ_{1}) / \max ({\bar{φ}}_{m}^{s_{R}} (λ_{1})) \\ G_{d R} (λ_{2}) / \max ({\bar{φ}}_{m}^{s_{R}} (λ_{2})) \\ ⋮ \\ G_{d R} (λ_{N}) / \max ({\bar{φ}}_{m}^{s_{R}} (λ_{N})) \end{matrix}) s_{R} = (\begin{matrix} {\bar{φ}}_{m}^{s_{R}} (λ_{1}) / \max ({\bar{φ}}_{m}^{s_{R}} (λ_{1})) \\ {\bar{φ}}_{m}^{s_{R}} (λ_{2}) / \max ({\bar{φ}}_{m}^{s_{R}} (λ_{2})) \\ ⋮ \\ {\bar{φ}}_{m}^{s_{R}} (λ_{N}) / \max ({\bar{φ}}_{m}^{s_{R}} (λ_{N})) \end{matrix}),

(24)

or, concisely,

{\bar{G}}_{d R} s_{R} = {\bar{φ}}_{m}

. A symmetric system of equations can be obtained from Eq. (24):

({\bar{G}}_{d R}^{T} {\bar{G}}_{d R}) s_{R} = {\bar{G}}_{d R}^{T} {\bar{φ}}_{m} .

(25)

As described in detail in a previous paper [33

M. A. Naser and M. S. Patterson, “Bioluminescence tomography using eigenvectors expansion and iterative solution for the optimized permissible source region,” Biomed. Opt. Express 2(11), 3179–3193 (2011). [CrossRef] [PubMed]

], Eq. (25) can be solved by assuming that the bioluminescence source

s_{R}

is written as a superposition of the eigenvectors of the matrix

[{\bar{G}}_{d R}^{T} {\bar{G}}_{d R}]

corresponding to non-zero eigenvalues such that

s_{R} = \sum_{i = 1}^{M} a_{i} v_{i}

, where

v = [v_{1}, v_{2}, \dots, v_{M}]

is a set of eigenvectors of the matrix

[{\bar{G}}_{d R}^{T} {\bar{G}}_{d R}]

corresponding to non-zero eigenvalues and M is the number of chosen eigenvectors. The vector

a = {[v^{T} ({\bar{G}}_{d R}^{T} {\bar{G}}_{d R}) v]}^{- 1} [v^{T} ({\bar{G}}_{d R}^{T}) {\bar{φ}}_{m}]

consists of the coefficients of the expansion [33

]. The equation described above is solved several times for different permissible source regions starting from the whole object size and ending with a few nodes. For every solution, the objective function

f = {‖ {\bar{G}}_{d R} (:, R) s_{R} - {\bar{φ}}_{m} ‖}^{1}

, which gives the difference between the calculated and measured light fluence rate, is calculated. After all iterations, the solution that minimizes the objective function is considered the best solution. The shrinking of the permissible region after every iteration is done by removing source nodes from the domain R corresponding to low source magnitude as described in [33

]

4. Results and discussion

To demonstrate these algorithms, a highly detailed anatomy for a laboratory mouse (MOBY) [34

W. P. Segars, B. M. Tsui, E. C. Frey, G. A. Johnson, and S. S. Berr, “Development of a 4-D digital mouse phantom for molecular imaging research,” Mol. Imaging Biol. 6(3), 149–159 (2004). [CrossRef] [PubMed]

] is used to create the FE mesh. The 3D object (Fig. 1(a) ) represents a finite element mesh of a section of a mouse abdomen where adipose, liver and spleen, bowel and intestine, kidney, and bone tissues are identified and assigned realistic wavelength-dependent optical properties. It is assumed that all elements within each region have the same optical properties. There are 9 optical fibers used as external sources located on the ventral surface of the mouse as shown in Fig. 1(b). It is assumed that these fibers produce Gaussian sources located one transport length inside the object (~1 mm) with a maximum value of 1 nW/mm³ and standard deviation of 1 mm. The measurements corresponding to these sources are used for the DOT reconstruction and BLT calibration at all wavelengths. The detectors shown in Fig. 1(b) are 2183 triangular surface elements and form a field of view of 225 degrees on the dorsal surface. It is assumed that the light fluence rate at these detectors can be obtained from CCD camera images by appropriate geometrical mapping [35

A. Da Silva, M. Leabad, C. Driol, T. Bordy, M. Debourdeau, J. M. Dinten, P. Peltié, and P. Rizo, “Optical calibration protocol for an x-ray and optical multimodality tomography system dedicated to small-animal examination,” Appl. Opt. 48(10), D151–D162 (2009). [CrossRef] [PubMed]

Fig. 1 (a) A 3D finite element mesh for the object used in the simulation. The mesh is a regular Delaunay tetrahedral with 17966 nodes and 99140 tetrahedrons generated by iso2mesh [36

Q. Fang and D. Boas, “Tetrahedral mesh generation from volumetric binary and gray-scale images,” in IEEE International Symposium on Biomedical Imaging: from Nano to Macro, 2009. ISBI '09 (IEEE, 2009), pp. 1142–1145.

]. The object represents a section of a mouse abdomen which illustrates 5 different tissue regions: (1) adipose, (2) liver and spleen, (3) bowel and intestine, (4) kidneys, and (5) bone. (b) A schematic of the object showing the sources and detectors set up and the coordinate system.

In the first step, the optical properties of all tissue regions at the wavelengths 580, 590, 600, 610, 620, and 630 nm were reconstructed using the DOT algorithm described above. The light fluence rate at each detector is determined for excitation of each source fiber. However these measurements are not calibrated and are shifted in both magnitude and phase from the correct values

Φ_{d}

as shown in Eq. (15). Equations (16) and (17) show that the optical properties can be reconstructed by using relative measured data, relative simulated data, and the relative Jacobian as shown in Eq. (16). The use of relative values in Eq. (16) cancels out the shift in magnitude and phase of the measured data from the correct calibrated data shown in Eq. (15). The DOT algorithm described above can be used for both frequency domain

f > 0

or continuous wave (CW)

f = 0

. In case of CW, the phase θ and the derivatives

\frac{\partial θ}{\partial μ}

in Eq. (10) will be zeros, which can be removed from the equations or left, as there is no effect on the results. The initial estimates used for scattering and absorption for all tissues at all wavelengths are 1 mm⁻¹ and 0.01 mm⁻¹ respectively. The DOT algorithm stops when the maximum relative error of the measured and simulated light fluence rates is less than 1% or the number of iterations exceeds 60.

Figure 2 shows a comparison of the relative percentage error in reconstructing the optical properties for CW at wavelength 590 nm using absolute data, relative data, and imperfect calibrated data. The forward model, Eq. (3), is used to simulate the measured data using the actual optical properties and 2% Gaussian noise is added to the calculated data to simulate realistic conditions. The true scattering and absorption coefficients were calculated using the data reported in [24

,37

S. A. Prahl, “Optical properties spectra,” (Oregon Medical Laser Clinic, 2001), http://omlc.ogi.edu/spectra/index.html.

] and are shown in Tables 1 and 2 respectively. The corresponding relative errors of reconstructing scattering, absorption, and effective attenuation

μ_{e f f} = \sqrt{3 ({μ^{'}}_{s} + μ_{a}) μ_{a}}

are shown in the figure. The errors using relative data are larger than when using absolute data. The maximum relative errors obtained in the range 580 nm to 630 nm are 2% for scattering, 3% for absorption, and 1% for effective attenuation when using absolute data while the maximum errors are 9% for scattering, 12% for absorption, and 4% for effective attenuation when using relative data. The simulation took 32 minutes for 60 iterations on a pc with Intel processor with a clock of 3.33 GHz and 24 GB RAM. The Jacobian was recalculated using the new reconstructed optical properties every 10 iterations.

Fig. 2 The relative percent error in the reconstructions of scattering

{μ^{'}}_{s}

, absorption

μ_{a}

,and effective attenuation

μ_{e f f}

, for the 5 tissue regions (1) adipose, (2) liver, (3) bowel, (4) kidney, and (5) bone at 590 nm using absolute, relative, and imperfect calibrations. Note the split vertical scale.

Table 1 True [24

,37

S. A. Prahl, “Optical properties spectra,” (Oregon Medical Laser Clinic, 2001), http://omlc.ogi.edu/spectra/index.html.

] reduced scattering coefficients (mm⁻¹) for each region in the heterogeneous object at six different wavelengths

	580 nm	590 nm	600 nm	610 nm	620 nm	630 nm
Adipose	1.30	1.29	1.28	1.27	1.26	1.25
Liver	0.79	0.77	0.76	0.75	0.74	0.72
Bowel	1.37	1.35	1.32	1.29	1.27	1.24
Kidney	2.80	2.73	2.66	2.60	2.53	2.47
Bone	3.08	3.01	2.93	2.86	2.80	2.73

Table 2 True [24

,37

S. A. Prahl, “Optical properties spectra,” (Oregon Medical Laser Clinic, 2001), http://omlc.ogi.edu/spectra/index.html.

] absorption coefficients (mm⁻¹) for each region in the heterogeneous object at six different wavelengths

	580 nm	590 nm	600 nm	610 nm	620 nm	630 nm
Adipose	0.0073	0.0043	0.0021	0.0014	0.0010	0.0008
Liver	0.6474	0.3993	0.1899	0.1200	0.0824	0.0647
Bowel	0.0198	0.0128	0.0063	0.0040	0.0028	0.0023
Kidney	0.1209	0.0746	0.0356	0.0226	0.0156	0.0123
Bone	0.1040	0.0670	0.0325	0.0207	0.0142	0.0112

In Fig. 2 we also show an example of the results obtained when the calibration process is imperfect. There are many possible sources of error in this multi-step process but here we use a simple illustration. One way to calibrate the camera is to use a reference object with known optical properties [38

L. Hervé, A. Koenig, A. Da Silva, M. Berger, J. Boutet, J. M. Dinten, P. Peltié, and P. Rizo, “Noncontact fluorescence diffuse optical tomography of heterogeneous media,” Appl. Opt. 46(22), 4896–4906 (2007). [CrossRef] [PubMed]

]. We have assumed that this reference object has the same external shape as the Moby phantom and uniform internal optical properties. We further assume that the reduced scattering coefficient is nominally 1.00 mm⁻¹ but actually 1.05 mm⁻¹ and the absorption coefficient is nominally 0.0100 mm⁻¹ but actually 0.0105 mm⁻¹. We performed forward calculations of the transmission for each source-detector pair under nominal and actual conditions. The ratio is the imperfect calibration factor used in all subsequent simulations including the BLT calculations. This 5% error in scattering and absorption coefficients is relatively small and likely underestimates the cumulative error in a real calibration process. With this imperfect calibration the maximum relative errors are 8% for scattering, 55% for absorption, and 35% for effective attenuation.

For BLT, to simulate the measured data

φ_{m}^{s_{R}}

using Eq. (23), the true source distribution

s_{R}

and the Green's function

G_{d R}

calculated using the true optical properties at each wavelength given in Tables 1 and 2 were used. 2% Gaussian noise was added to simulate realistic data:

φ_{m}^{s_{R}}

in Eq. (23). Two cases are simulated for

φ_{m}^{s_{R}}

in Eq. (23): the first case is to set each element in α equal to 1 in Eq. (22) when the data are calibrated and match the correct values at the detectors around the object. In this case there is no need to use external sources for calibrations and

{\bar{φ}}_{m}^{s_{R}} = φ_{m}^{s_{R}}

in Eq. (23). The second case is to set each element in α equal to 1000 in Eq. (22) to simulate uncalibrated data and normalization with external sources is needed as described in Eq. (23). To simulate the measured data due to external sources

φ_{m}^{s_{i}}

, the Green's functions

G_{d}

calculated using the true optical properties at each wavelength were used in Eq. (22). The external sources are simulated by Gaussian sources as described before and the elements of α equal 1000. In Eq. (23), the Green’s functions

G_{d R}

and

G_{d}

were calculated using the reconstructed optical properties, by the DOT algorithm described above, using absolute data, relative data, and imperfect calibrated data. The permissible region R was initially chosen to be a cylinder of 20 mm diameter and height 20 mm with its center located at (0, 3, 3). At each iteration, the permissible region is slightly shrunk by an arbitrary factor γ The permissible region initially containing 5010 mesh nodes is finally shrunk to 10 nodes in 60 iterations, γ = (5010/10)^(1/60) = 1.1092. The 60 iterations took approximately 30 minutes on a PC with Intel processor with a clock of 3.33 GHz and 24 GB RAM. The execution time varies greatly according to the choice of mesh size, the volume of the initial permissible region, and number of iterations. It is useful to start with a small size mesh and fewer iterations to get a quick initial solution before solving the complete problem. This can provide a more focused permissible region and hence greatly reduce the simulation time and give more accurate results.

Figure 3 shows transverse and coronal cross sections for the actual and reconstructed source distributions. The transverse section passes through the two kidneys at

z = 3 mm

and the coronal section passes through the two kidneys at

y = 3 mm

. The actual sources are spheres of 4 mm diameters located in the two kidneys at (−5, 3, 3) and (6, 3, 3) with magnitudes of 2 nW/ mm³ and 1 nW/ mm³ respectively and total spatially integrated powers of 67.2 nW and 39 nW respectively.

Fig. 3 BLT reconstruction of two uniform bioluminescence sources localized in the left and right kidneys (4 mm diameters). The “anatomy” is shown in Fig. 1. The actual sources in transverse section (a) and coronal section (b) have magnitudes 2 nW/ mm³ and 1 nW/ mm³ for the left and right source respectively. The reconstructed sources are shown in (c) and (d) using absolute data, (e) and (f) using relative data, and (g) and (h) using imperfect calibrated data. The total powers of the actual and reconstructed sources and the normalized magnitude errors of the reconstructed sources using different reconstruction methods are shown in the figure.

Two objective metrics have been calculated to assess the accuracy of the BLT reconstruction of each source. The first is simply the total power integrated over three dimensions,

\sum_{i} s (i) V (i)

, where

V (i)

is the volume of the element i in the FE mesh and the summation includes all elements that include the source. The second we define as the “normalized magnitude error” calculated according to

\sum_{i} | s^{actual} (i) - s^{reconst} (i) | V (i) / \sum_{i} s^{actual} (i) V (i) .

This element-by-element summation is sensitive to the shape and position of the reconstructed source and is a more stringent metric than the total power.

The reconstructed sources in Figs. 3(c) and 3(d) using absolute data and in 3(e) and 3(f) using relative data show good agreement with the actual sources in Figs. 3(a) and 3(b) in terms of location, total power reconstructed (within 0.75%) and magnitudes with normalized magnitude errors of 0.27. There are no significant differences in the reconstruction results when using absolute or relative data. The reconstructed sources using imperfect calibrated data are shown in Figs. 3(g) and 3(h) with total power reconstructed (within 6%) and normalized magnitude error of 0.54.

Figure 4 shows the result of the reconstruction of large uniform sources that fill the entire two kidneys; the actual sources in Figs. 4(a) and 4(b) are of uniform magnitude 1 nW/mm³. The total powers of the sources are 308.3 nW in the left kidney and 316.3 nW in the right kidney. The reconstructed sources in Figs. 4(c) and 4(d) using absolute data and in 4(e) and 4(f) using relative data show comparable results in terms of the total power reconstructed which is within 0.63% of the actual values and the normalized magnitude errors which are 0.14 and 0.17 respectively. The reconstructed sources using imperfect calibration in Figs. 4(g) and 4(h) have total power reconstructed within 6.3% and normalized magnitude error of 0.54.

Fig. 4 BLT reconstruction of large distributed sources that fill the entire two kidneys. The magnitude of the actual sources in (a) and (b) is 1 nW/mm³, and the reconstructed sources are shown in (c) and (d) using absolute data, (e) and (f) using relative data, and (g) and (h) using imperfect calibrated data. The total powers of the actual and reconstructed sources and the normalized magnitude errors of the reconstructed sources using different reconstruction methods are shown in the figure.

Figure 5 shows the BLT reconstruction for a large non-uniform source that fills the entire right kidney. There is a 4 mm diameter “hot spot” of 2 nW/mm³ located at (6, 3, 3) inside the kidney while the rest of the organ “background” of 1 nW/mm³ as shown in Figs. 5(a) and 5(b). The total power of the actual source is 355.3 nW which is reconstructed within 1% when using absolute data in Figs. 5(c) and 5(d) or relative data in 5(e) and 5(f). The hot spot can be recognized in both cases but the location is slightly more accurate when absolute data are used. The normalized magnitude error is 0.2 using absolute data and 0.25 while using relative data. The total power is reconstructed within 5.3% using imperfect calibration and the normalized magnitude error is 0.39. The hot spot location is less accurate than that using absolute or relative data.

Fig. 5 BLT reconstruction of a large non-uniform source distribution that fills the entire right kidney. The actual source in (a) and (b) has a hot spot of 4 mm diameter with magnitude 2 nW/mm³ which is twice the magnitude of the background. The reconstructed source in (c) and (d) using the absolute data, (e) and (f) using the relative data, and (g) and (g) using the imperfect calibrated data. The total powers of the actual and reconstructed sources and the normalized magnitude errors of the reconstructed sources using different reconstruction methods are shown in the figure.

5. Conclusion

Self-calibrated algorithms for three dimensional DOT and BLT based on the finite element solution of the diffusion equation have been developed. Both algorithms avoid absolute calibration of the detectors as by using ratios, camera images can be used for reconstructing the optical properties and bioluminescence sources. The algorithms have been applied to simulations based on a highly detailed anatomy for a laboratory mouse (MOBY) and show good agreement between the actual values of optical properties and bioluminescence sources and the reconstructions. The DOT algorithm used the direct method to calculate the Jacobian. Using relative values for the measured data, simulated data, and Jacobian showed good reconstruction results compared to using absolute values. The same measurements due to external sources used for DOT at each wavelength were also used to calibrate the bioluminescence measurements. The reconstruction results for BLT using relative calibration show good agreements with the results using absolute data. The BLT algorithms used the eigenvectors expansion of the Green's function as a reduced space for the reconstruction, and iterative shrinking of the permissible source region. Small multiple sources and large homogeneous and inhomogeneous source distributions were reconstructed using both absolute and relative data. The results showed good agreement between the reconstructed and actual sources in terms of location and total power and magnitude. The BLT reconstructions based on relative data are slightly inferior to those based on absolute data obtained with perfect calibration. However, any calibration procedure will be imperfect, and we have demonstrated that even small errors in the calibration can lead to performance that is worse than that based on relative data.

Acknowledgment

This work was supported in part by National Institutes of Health grant RO1 CA158100.

References and links

1.	V. Ntziachristos, J. Ripoll, L. V. Wang, and R. Weissleder, “Looking and listening to light: the evolution of whole-body photonic imaging,” Nat. Biotechnol. 23(3), 313–320 (2005). [CrossRef] [PubMed]
2.	R. Weissleder and U. Mahmood, “Molecular imaging,” Radiology 219(2), 316–333 (2001). [PubMed]
3.	J. K. Willmann, N. van Bruggen, L. M. Dinkelborg, and S. S. Gambhir, “Molecular imaging in drug development,” Nat. Rev. Drug Discov. 7(7), 591–607 (2008). [CrossRef] [PubMed]
4.	J. Tian, J. Bai, X. P. Yan, S. Bao, Y. Li, W. Liang, and X. Yang, “Multimodality molecular imaging,” IEEE Eng. Med. Biol. Mag. 27(5), 48–57 (2008). [CrossRef] [PubMed]
5.	A. J. Chaudhari, F. Darvas, J. R. Bading, R. A. Moats, P. S. Conti, D. J. Smith, S. R. Cherry, and R. M. Leahy, “Hyperspectral and multispectral bioluminescence optical tomography for small animal imaging,” Phys. Med. Biol. 50(23), 5421–5441 (2005). [CrossRef] [PubMed]
6.	S. Ahn, A. J. Chaudhari, F. Darvas, C. A. Bouman, and R. M. Leahy, “Fast iterative image reconstruction methods for fully 3D multispectral bioluminescence tomography,” Phys. Med. Biol. 53(14), 3921–3942 (2008). [CrossRef] [PubMed]
7.	D. C. Comsa, T. J. Farrell, and M. S. Patterson, “Quantification of bioluminescence images of point source objects using diffusion theory models,” Phys. Med. Biol. 51(15), 3733–3746 (2006). [CrossRef] [PubMed]
8.	C. Kuo, O. Coquoz, T. L. Troy, H. Xu, and B. W. Rice, “Three-dimensional reconstruction of in vivo bioluminescent sources based on multispectral imaging,” J. Biomed. Opt. 12(2), 024007 (2007). [CrossRef] [PubMed]
9.	H. Dehghani, S. C. Davis, S. Jiang, B. W. Pogue, K. D. Paulsen, and M. S. Patterson, “Spectrally resolved bioluminescence optical tomography,” Opt. Lett. 31(3), 365–367 (2006). [CrossRef] [PubMed]
10.	W. Cong, G. Wang, D. Kumar, Y. Liu, M. Jiang, L. V. Wang, E. A. Hoffman, G. McLennan, P. B. McCray, J. Zabner, and A. Cong, “Practical reconstruction method for bioluminescence tomography,” Opt. Express 13(18), 6756–6771 (2005). [CrossRef] [PubMed]
11.	X. He, J. Liang, X. Qu, H. Huang, Y. Hou, and J. Tian, “Truncated total least squares method with a practical truncation parameter choice scheme for bioluminescence tomography inverse problem,” Int. J. Biomed. Imaging 2010, 12 (2010). [CrossRef] [PubMed]
12.	H. Huang, X. Qu, J. Liang, X. He, X. Chen, D. Yang, and J. Tian, “A multi-phase level set framework for source reconstruction in bioluminescence tomography,” J. Comput. Phys. 229(13), 5246–5256 (2010). [CrossRef]
13.	J. Feng, K. Jia, G. Yan, S. Zhu, C. Qin, Y. Lv, and J. Tian, “An optimal permissible source region strategy for multispectral bioluminescence tomography,” Opt. Express 16(20), 15640–15654 (2008). [CrossRef] [PubMed]
14.	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 sparse a priori information,” Opt. Express 17(10), 8062–8080 (2009). [CrossRef] [PubMed]
15.	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]
16.	N. Cao, A. Nehorai, and M. Jacobs, “Image reconstruction for diffuse optical tomography using sparsity regularization and expectation-maximization algorithm,” Opt. Express 15(21), 13695–13708 (2007). [CrossRef] [PubMed]
17.	X. Gu, Q. Zhang, L. Larcom, and H. Jiang, “Three-dimensional bioluminescence tomography with model-based reconstruction,” Opt. Express 12(17), 3996–4000 (2004). [CrossRef] [PubMed]
18.	N. V. Slavine, M. A. Lewis, E. Richer, and P. P. Antich, “Iterative reconstruction method for light emitting sources based on the diffusion equation,” Med. Phys. 33(1), 61–68 (2006). [CrossRef] [PubMed]
19.	X. Chen, X. Gao, D. Chen, X. Ma, X. Zhao, M. Shen, X. Li, X. Qu, J. Liang, J. Ripoll, and J. Tian, “3D reconstruction of light flux distribution on arbitrary surfaces from 2D multi-photographic images,” Opt. Express 18(19), 19876–19893 (2010). [CrossRef] [PubMed]
20.	A. Kienle and F. Foschum, “250 years Lambert surface: does it really exist?” Opt. Express 19(5), 3881–3889 (2011). [CrossRef] [PubMed]
21.	J. Pekar, “Multispectral bioluminescence tomography with x-ray CT spatial priors” Ph.D. thesis (McMaster University, Hamilton, 2011), Open Access Dissertations and Theses, Paper 4329, http://digitalcommons.mcmaster.ca/opendissertations/4329.
22.	M. A. J. Chaudhari, A. A. Joshi, F. Darvas, and R. M. Leahy, “A method for atlas-based volumetric registration with surface constraints for optical bioluminescence tomography in small animal imaging,” Proc. SPIE 6510, 651024, 651024–10 (2007). [CrossRef]
23.	H. Dehghani, B. W. Pogue, S. C. Davis, and M. S. Patterson, “Modeling and image reconstruction in spectrally resolved bioluminescence tomography,” Proc. SPIE 6434, 64340V, 64340V–9 (2007). [CrossRef]
24.	G. Alexandrakis, F. R. Rannou, and A. F. Chatziioannou, “Tomographic bioluminescence imaging by use of a combined optical-PET (OPET) system: a computer simulation feasibility study,” Phys. Med. Biol. 50(17), 4225–4241 (2005). [CrossRef] [PubMed]
25.	M. A. Naser and M. S. Patterson, “Algorithms for bioluminescence tomography incorporating anatomical information and reconstruction of tissue optical properties,” Biomed. Opt. Express 1(2), 512–526 (2010). [CrossRef] [PubMed]
26.	M. A. Naser and M. S. Patterson, “Improved bioluminescence and fluorescence reconstruction algorithms using diffuse optical tomography, normalized data, and optimized selection of the permissible source region,” Biomed. Opt. Express 2(1), 169–184 (2011). [CrossRef] [PubMed]
27.	S. R. Arridge, “Optical tomography in medical imaging,” Inverse Probl. 15(2), R41–R93 (1999). [CrossRef]
28.	A. D. Klose, “Transport-theory-based stochastic image reconstruction of bioluminescence sources,” J. Opt. Soc. Am. A 24(6), 1601–1608 (2007). [CrossRef]
29.	H. Dehghani, M. E. Eames, P. K. Yalavarthy, S. C. Davis, S. Srinivasan, C. M. Carpenter, B. W. Pogue, and K. D. Paulsen, “Near infrared optical tomography using NIRFAST: Algorithm for numerical model and image reconstruction,” Commun. Numer. Methods Eng. 25(6), 711–732 (2009). [CrossRef] [PubMed]
30.	S. R. Arridge, M. Schweiger, M. Hiraoka, and D. T. Delpy, “A finite element approach for modeling photon transport in tissue,” Med. Phys. 20(2), 299–309 (1993). [CrossRef] [PubMed]
31.	M. Schweiger, S. R. Arridge, M. Hiraoka, and D. T. Delpy, “The finite element method for the propagation of light in scattering media: boundary and source conditions,” Med. Phys. 22(11), 1779–1792 (1995). [CrossRef] [PubMed]
32.	H. Jiang, “Frequency-domain fluorescent diffusion tomography: a finite-element-based algorithm and simulations,” Appl. Opt. 37(22), 5337–5343 (1998). [CrossRef] [PubMed]
33.	M. A. Naser and M. S. Patterson, “Bioluminescence tomography using eigenvectors expansion and iterative solution for the optimized permissible source region,” Biomed. Opt. Express 2(11), 3179–3193 (2011). [CrossRef] [PubMed]
34.	W. P. Segars, B. M. Tsui, E. C. Frey, G. A. Johnson, and S. S. Berr, “Development of a 4-D digital mouse phantom for molecular imaging research,” Mol. Imaging Biol. 6(3), 149–159 (2004). [CrossRef] [PubMed]
35.	A. Da Silva, M. Leabad, C. Driol, T. Bordy, M. Debourdeau, J. M. Dinten, P. Peltié, and P. Rizo, “Optical calibration protocol for an x-ray and optical multimodality tomography system dedicated to small-animal examination,” Appl. Opt. 48(10), D151–D162 (2009). [CrossRef] [PubMed]
36.	Q. Fang and D. Boas, “Tetrahedral mesh generation from volumetric binary and gray-scale images,” in IEEE International Symposium on Biomedical Imaging: from Nano to Macro, 2009. ISBI '09 (IEEE, 2009), pp. 1142–1145.
37.	S. A. Prahl, “Optical properties spectra,” (Oregon Medical Laser Clinic, 2001), http://omlc.ogi.edu/spectra/index.html.
38.	L. Hervé, A. Koenig, A. Da Silva, M. Berger, J. Boutet, J. M. Dinten, P. Peltié, and P. Rizo, “Noncontact fluorescence diffuse optical tomography of heterogeneous media,” Appl. Opt. 46(22), 4896–4906 (2007). [CrossRef] [PubMed]

OCIS Codes
(170.3880) Medical optics and biotechnology : Medical and biological imaging
(170.6960) Medical optics and biotechnology : Tomography

ToC Category:
Image Reconstruction and Inverse Problems

History
Original Manuscript: July 13, 2012
Revised Manuscript: October 3, 2012
Manuscript Accepted: October 9, 2012
Published: October 11, 2012

Citation
Mohamed A. Naser, Michael S. Patterson, and John W. Wong, "Self-calibrated algorithms for diffuse optical tomography and bioluminescence tomography using relative transmission images," Biomed. Opt. Express 3, 2794-2808 (2012)
http://www.opticsinfobase.org/boe/abstract.cfm?URI=boe-3-11-2794

Sort: Author | Year | Journal | Reset

References

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]
R. Weissleder and U. Mahmood, “Molecular imaging,” Radiology219(2), 316–333 (2001). [PubMed]
J. K. Willmann, N. van Bruggen, L. M. Dinkelborg, and S. S. Gambhir, “Molecular imaging in drug development,” Nat. Rev. Drug Discov.7(7), 591–607 (2008). [CrossRef] [PubMed]
J. Tian, J. Bai, X. P. Yan, S. Bao, Y. Li, W. Liang, and X. Yang, “Multimodality molecular imaging,” IEEE Eng. Med. Biol. Mag.27(5), 48–57 (2008). [CrossRef] [PubMed]
A. J. Chaudhari, F. Darvas, J. R. Bading, R. A. Moats, P. S. Conti, D. J. Smith, S. R. Cherry, and R. M. Leahy, “Hyperspectral and multispectral bioluminescence optical tomography for small animal imaging,” Phys. Med. Biol.50(23), 5421–5441 (2005). [CrossRef] [PubMed]
S. Ahn, A. J. Chaudhari, F. Darvas, C. A. Bouman, and R. M. Leahy, “Fast iterative image reconstruction methods for fully 3D multispectral bioluminescence tomography,” Phys. Med. Biol.53(14), 3921–3942 (2008). [CrossRef] [PubMed]
D. C. Comsa, T. J. Farrell, and M. S. Patterson, “Quantification of bioluminescence images of point source objects using diffusion theory models,” Phys. Med. Biol.51(15), 3733–3746 (2006). [CrossRef] [PubMed]
C. Kuo, O. Coquoz, T. L. Troy, H. Xu, and B. W. Rice, “Three-dimensional reconstruction of in vivo bioluminescent sources based on multispectral imaging,” J. Biomed. Opt.12(2), 024007 (2007). [CrossRef] [PubMed]
H. Dehghani, S. C. Davis, S. Jiang, B. W. Pogue, K. D. Paulsen, and M. S. Patterson, “Spectrally resolved bioluminescence optical tomography,” Opt. Lett.31(3), 365–367 (2006). [CrossRef] [PubMed]
W. Cong, G. Wang, D. Kumar, Y. Liu, M. Jiang, L. V. Wang, E. A. Hoffman, G. McLennan, P. B. McCray, J. Zabner, and A. Cong, “Practical reconstruction method for bioluminescence tomography,” Opt. Express13(18), 6756–6771 (2005). [CrossRef] [PubMed]
X. He, J. Liang, X. Qu, H. Huang, Y. Hou, and J. Tian, “Truncated total least squares method with a practical truncation parameter choice scheme for bioluminescence tomography inverse problem,” Int. J. Biomed. Imaging2010, 12 (2010). [CrossRef] [PubMed]
H. Huang, X. Qu, J. Liang, X. He, X. Chen, D. Yang, and J. Tian, “A multi-phase level set framework for source reconstruction in bioluminescence tomography,” J. Comput. Phys.229(13), 5246–5256 (2010). [CrossRef]
J. Feng, K. Jia, G. Yan, S. Zhu, C. Qin, Y. Lv, and J. Tian, “An optimal permissible source region strategy for multispectral bioluminescence tomography,” Opt. Express16(20), 15640–15654 (2008). [CrossRef] [PubMed]
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 sparse a priori information,” Opt. Express17(10), 8062–8080 (2009). [CrossRef] [PubMed]
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]
N. Cao, A. Nehorai, and M. Jacobs, “Image reconstruction for diffuse optical tomography using sparsity regularization and expectation-maximization algorithm,” Opt. Express15(21), 13695–13708 (2007). [CrossRef] [PubMed]
X. Gu, Q. Zhang, L. Larcom, and H. Jiang, “Three-dimensional bioluminescence tomography with model-based reconstruction,” Opt. Express12(17), 3996–4000 (2004). [CrossRef] [PubMed]
N. V. Slavine, M. A. Lewis, E. Richer, and P. P. Antich, “Iterative reconstruction method for light emitting sources based on the diffusion equation,” Med. Phys.33(1), 61–68 (2006). [CrossRef] [PubMed]
X. Chen, X. Gao, D. Chen, X. Ma, X. Zhao, M. Shen, X. Li, X. Qu, J. Liang, J. Ripoll, and J. Tian, “3D reconstruction of light flux distribution on arbitrary surfaces from 2D multi-photographic images,” Opt. Express18(19), 19876–19893 (2010). [CrossRef] [PubMed]
A. Kienle and F. Foschum, “250 years Lambert surface: does it really exist?” Opt. Express19(5), 3881–3889 (2011). [CrossRef] [PubMed]
J. Pekar, “Multispectral bioluminescence tomography with x-ray CT spatial priors” Ph.D. thesis (McMaster University, Hamilton, 2011), Open Access Dissertations and Theses, Paper 4329, http://digitalcommons.mcmaster.ca/opendissertations/4329 .
M. A. J. Chaudhari, A. A. Joshi, F. Darvas, and R. M. Leahy, “A method for atlas-based volumetric registration with surface constraints for optical bioluminescence tomography in small animal imaging,” Proc. SPIE6510, 651024, 651024–10 (2007). [CrossRef]
H. Dehghani, B. W. Pogue, S. C. Davis, and M. S. Patterson, “Modeling and image reconstruction in spectrally resolved bioluminescence tomography,” Proc. SPIE6434, 64340V, 64340V–9 (2007). [CrossRef]
G. Alexandrakis, F. R. Rannou, and A. F. Chatziioannou, “Tomographic bioluminescence imaging by use of a combined optical-PET (OPET) system: a computer simulation feasibility study,” Phys. Med. Biol.50(17), 4225–4241 (2005). [CrossRef] [PubMed]
M. A. Naser and M. S. Patterson, “Algorithms for bioluminescence tomography incorporating anatomical information and reconstruction of tissue optical properties,” Biomed. Opt. Express1(2), 512–526 (2010). [CrossRef] [PubMed]
M. A. Naser and M. S. Patterson, “Improved bioluminescence and fluorescence reconstruction algorithms using diffuse optical tomography, normalized data, and optimized selection of the permissible source region,” Biomed. Opt. Express2(1), 169–184 (2011). [CrossRef] [PubMed]
S. R. Arridge, “Optical tomography in medical imaging,” Inverse Probl.15(2), R41–R93 (1999). [CrossRef]
A. D. Klose, “Transport-theory-based stochastic image reconstruction of bioluminescence sources,” J. Opt. Soc. Am. A24(6), 1601–1608 (2007). [CrossRef]
H. Dehghani, M. E. Eames, P. K. Yalavarthy, S. C. Davis, S. Srinivasan, C. M. Carpenter, B. W. Pogue, and K. D. Paulsen, “Near infrared optical tomography using NIRFAST: Algorithm for numerical model and image reconstruction,” Commun. Numer. Methods Eng.25(6), 711–732 (2009). [CrossRef] [PubMed]
S. R. Arridge, M. Schweiger, M. Hiraoka, and D. T. Delpy, “A finite element approach for modeling photon transport in tissue,” Med. Phys.20(2), 299–309 (1993). [CrossRef] [PubMed]
M. Schweiger, S. R. Arridge, M. Hiraoka, and D. T. Delpy, “The finite element method for the propagation of light in scattering media: boundary and source conditions,” Med. Phys.22(11), 1779–1792 (1995). [CrossRef] [PubMed]
H. Jiang, “Frequency-domain fluorescent diffusion tomography: a finite-element-based algorithm and simulations,” Appl. Opt.37(22), 5337–5343 (1998). [CrossRef] [PubMed]
M. A. Naser and M. S. Patterson, “Bioluminescence tomography using eigenvectors expansion and iterative solution for the optimized permissible source region,” Biomed. Opt. Express2(11), 3179–3193 (2011). [CrossRef] [PubMed]
W. P. Segars, B. M. Tsui, E. C. Frey, G. A. Johnson, and S. S. Berr, “Development of a 4-D digital mouse phantom for molecular imaging research,” Mol. Imaging Biol.6(3), 149–159 (2004). [CrossRef] [PubMed]
A. Da Silva, M. Leabad, C. Driol, T. Bordy, M. Debourdeau, J. M. Dinten, P. Peltié, and P. Rizo, “Optical calibration protocol for an x-ray and optical multimodality tomography system dedicated to small-animal examination,” Appl. Opt.48(10), D151–D162 (2009). [CrossRef] [PubMed]
Q. Fang and D. Boas, “Tetrahedral mesh generation from volumetric binary and gray-scale images,” in IEEE International Symposium on Biomedical Imaging: from Nano to Macro,2009. ISBI '09 (IEEE, 2009), pp. 1142–1145.
S. A. Prahl, “Optical properties spectra,” (Oregon Medical Laser Clinic, 2001), http://omlc.ogi.edu/spectra/index.html .
L. Hervé, A. Koenig, A. Da Silva, M. Berger, J. Boutet, J. M. Dinten, P. Peltié, and P. Rizo, “Noncontact fluorescence diffuse optical tomography of heterogeneous media,” Appl. Opt.46(22), 4896–4906 (2007). [CrossRef] [PubMed]

Cited By

Alert me when this paper is cited

OSA is able to provide readers links to articles that cite this paper by participating in CrossRef's Cited-By Linking service. CrossRef includes content from more than 3000 publishers and societies. In addition to listing OSA journal articles that cite this paper, citing articles from other participating publishers will also be listed.

Click here to see a list of articles that cite this paper