Structural information within regularization matrices improves near infrared diffuse optical tomography

doi:10.1364/OE.15.008043

Structural information within regularization matrices improves near infrared diffuse optical tomography

Phaneendra K. Yalavarthy, Brian W. Pogue, Hamid Dehghani, Colin M. Carpenter, Shudong Jiang, and Keith D. Paulsen »View Author Affiliations

Optics Express, Vol. 15, Issue 13, pp. 8043-8058 (2007)
http://dx.doi.org/10.1364/OE.15.008043

View Full Text Article

Acrobat PDF (822 KB)

Journal Search
Article Lookup

Article Tools

Citations

Alert me when this article is cited
Export Citation/Save

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

Abstract

Near-Infrared (NIR) tomographic image reconstruction is a non-linear, ill-posed and ill-conditioned problem, and so in this study, different ways of penalizing the objective function with structural information were investigated. A simple framework to incorporate structural priors is presented, using simple weight matrices that have either Laplacian or Helmholtz-type structures. Using both MRI-derived breast geometry and phantom data, a systematic and quantitative comparison was performed with and without spatial priors. The Helmholtz-type structure can be seen as a more generalized approach for incorporating spatial priors into the reconstruction scheme. Moreover, parameter reduction (i.e. hard prior information) in the imaging field through the enforcement of spatially explicit regions may lead to erroneous results with imperfect spatial priors.

1. Introduction

Near Infrared (NIR) optical tomography is a method that uses light in the 650-900nm wave-length range to recover images of the internal spatial distribution of tissue optical properties, absorption (or chromophore concentrations) and scattering parameters[1–3

D. A. Boas, D. H. Brooks, E. L. Miller, C. A. DiMarzio, M. Kilmer, R. J. Gaudette, and Q. Zhang, “Imaging the body with diffuse optical tomography,” IEEE Sig. Proc. Mag. 18, 57–75 (2001). [CrossRef]

]. When imaged at multiple wavelengths, these quantities can be used to estimate tissue hemoglobin and water concentration[3

S. Srinivasan, B. W. Pogue, S. Jiang, H. Dehghani, C. Kogel, S. Soho, J. J. Gibson, T. D. Tosteson, S. P. Poplack, and K. D. Paulsen, “Interpreting Hemoglobin andWater Concentration, Oxygen Saturation and Scattering Measured in Vivo by Near-Infrared Breast Tomography,” Proc. Nat. Acad. Sci. U.S.A. 100, 12349–12354 (2003). [CrossRef]

] and have the advantage of being acquired non-invasively and without ionizing radiation. The imaging procedure can be rapidly or repeatedly applied to investigate physiological state, and systems can be integrated into conventional imaging platforms such as X-ray mammography, Ultrasound and MRI[4–9

Q. Zhu, N. G. Chen, and S. C. Kurtzman, “Imaging tumor angiogenesis by use of combined near-infrared diffusive light and ultrasound,” Opt. Lett. 28, 337–339 (2003). [CrossRef] [PubMed]

]. These hybrid systems have been shown to achieve superior performance in terms of resolution and quantitative accuracy which should provide more accurate physiological data from the tissue under investigation[4–16

Q. Zhu, N. G. Chen, and S. C. Kurtzman, “Imaging tumor angiogenesis by use of combined near-infrared diffusive light and ultrasound,” Opt. Lett. 28, 337–339 (2003). [CrossRef] [PubMed]

]. However, a fundamental question is how to utilize the spatial information from the clinical system to maximize the accuracy of NIR tomography. In this study, the ability to improve the quantitative accuracy of regions imaged with NIR tomography was investigated, in the setting where prior spatial information is readily available.

This work explores image reconstruction strategies that take advantage of multimodality image data, in particular, the combination of MRI with NIR optical tomography for breast cancer imaging. MRI provides structural information at high spatial resolution (near 1 mm), whereas NIR imaging has relatively poor spatial resolution (near 4-7 mm). Yet MR imaging would benefit from the molecular-specific signatures available through NIR[3

, 16–18

B. Brooksby, B. W. Pogue, S. Jiang, H. Dehghani, S. Srinivasan, C. Kogel, T. D. Tosteson, J. Weaver, S. P. Poplack, and K. D. Paulsen, “Imaging breast adipose and fibroglandular tissue molecular signatures using Hybrid MRI-Guided Near-Infrared Spectral Tomography,” Proc. Nat. Acad. Sci. U.S.A. 103, 8828–8833 (2006). [CrossRef]

], specifically tissue hemoglobin content, oxygenation level, and water, as well as scattering particle size and number density [19

X. Wang, B.W. Pogue, S. Jiang, X. Song, K.D. Paulsen, C. Kogel, S.P. Poplack, and W.A. Wells, , “Approximation of Mie scattering parameters in near-infrared tomography of normal breast tissue in vivo,” J. Biomed. Opt. 10, 051704:1-8 (2005). [CrossRef] [PubMed]

]. The inverse problem (image reconstruction procedure) in NIR imaging is known to be a non-linear, ill-posed and ill-conditioned problem [20

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

]. Use of structural information in NIR reconstruction schemes has been explored by several research groups. For example, Li et. al [7

A. Li, E. L. Miller, M. E. Kilmer, T. J. Brukilacchio, T. Chaves, J. Stott, Q. Zhang, T. Wu, M. Chorlton, R. H. Moore, D. B. Kopans, and D. A. Boas, “Tomographic optical breast imaging guided by three-dimensional mammography,” Appl. Opt. 42, 5181–5190 (2003). [CrossRef] [PubMed]

] have used the data derived from X-ray mammography for choosing different regularization parameters for the region of interest (ROI) and surrounding tissue, and have shown that the contrast and resolution of the reconstructed images can be improved. Srinivasan et al. [21

S. Srinivasan, B. W. Pogue, H. Dehghani, S. Jiang, X. Song, and K. D. Paulsen, “Improved quantification of small objects in near-infrared diffuse optical tomography,” J. Biomed. Opt. 9, 1161–1171 (2004). [CrossRef] [PubMed]

] have developed a three-step reconstruction process for improving the quantification accuracy of small-objects in NIR tomography, where they use the conventional NIR reconstructed images (first step) as a structural prior for the last two steps. Earlier papers have shown that optical contrast can be correlated to MR contrast [6

B. Brooksby, S. Jiang, C. Kogel, M. Doyley, H. Dehghani, J. B. Weaver, S. P. Poplack, B. W Pogue, and K. D Paulsen, “Magnetic resonance guided near infrared tomography of the breast,” Rev. Sci. Instrum. 75, 5262–5270 (2004). [CrossRef]

, 9

V. Ntziachristos, A. G. Yodh, M. Schnall, and B. Chance, “Concurrent MRI and diffuse optical tomography of breast after indocyanine green enhancement,” Proc. Nat. Acad. Sci. U.S.A. 97, 2767–2772 (2000). [CrossRef]

, 13

V. Ntziachristos, X. H. Ma, and B. Chance, “Time-correlated single photon counting imager for simultaneous magnetic resonance and near-infrared mammography,” Rev. Sci. Instrum. 69, 4221–4233 (1998). [CrossRef]

] and that structural MRI images can be used to reduce the number of unknown parameters to be estimated [14

B. Brooksby, H. Dehghani, B. W. Pogue, and K. D. Paulsen, “Near infrared (NIR) tomography breast image reconstruction with a priori structural information from MRI: algorithm development for reconstructing heterogeneities,” IEEE J. Sel. Top. Quantum Electron. 9, 199–209 (2003). [CrossRef]

]. The difficulty with parameter reduction approaches (referred to as hard-priors [22

H. Dehghani, B.W. Pogue, J. Shudong, B. Brooksby, and K. D. Paulsen, “Three-dimensional optical tomography: resolution in small-object imaging,” Appl. Opt. 42, 3117–3126 (2003). [CrossRef] [PubMed]

]) is the potential of introducing error by imposing incorrect model assumptions and introducing variation due to uncertainty in the prior information (even when the underlying model is appropriate). For example, the features which lead to contrast in one imaging system may not be spatially coregistered with those that produce contrast in another imaging system. Further, segmentation of congruent features always includes classification errors resulting from digitization. Recently, Boverman et. al [10

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

] showed that even imperfect priors which encode breast background structure improves anomaly localization, but at the expense of biasing the spectroscopic dimension of the image reconstruction. Another type of approach for constraining the problem, often described as soft-priors, penalizes the variation within regions which are assumed to have the same properties by controlling regularization. Brooksby et. al [15

B. Brooksby, S. Jiang, H. Dehghani, B. W. Pogue, K. D. Paulsen, J. Weaver, C. Kogel, and S. P. Poplack, “Combining near infrared tomography and magnetic resonance imaging to study in vivo breast tissue: implementation of a Laplacian-type regularization to incorporate magnetic resonance structure,” J. Biomed. Opt. 10, 051504:1-10 (2005). [CrossRef] [PubMed]

, 16

, 23

B. Brooksby, S. Srinivasan, S. Jiang, H. Dehghani, B. W. Pogue, K. D. Paulsen, J. Weaver, C. Kogel, and S. P. Poplack, “Spectral priors improve near-infrared diffuse tomography more than spatial priors,” Opt. Lett. 30, 1968–1970 (2005). [CrossRef] [PubMed]

] have developed a Laplacian type of regularization that allows intra-region variability. This is a method which works well even if the confidence in the prior structural information is low.

This paper develops a more generalized framework, known as soft-priors, for incorporating the structural priors into the NIR image reconstruction process, and explores a covariance-based constraint scheme adopted from finite differencing of the Helmholtz equation in addition to the soft and hard prior approaches noted above. The soft-prior approach allows optical property variation with a given region, reducing biases caused by the use of imperfect prior information. The results indicate that imperfect structural information can generate errors in the hard-priors case, whereas the soft-priors are able to quantify regions more appropriately. Simulation and experimental studies are performed to demonstrate the superior reconstruction image quality. These types of procedures are needed to improve NIR imaging both in terms of high spatial resolution available from MRI and high contrast inherent in the NIR signal.

2. Methods: mathematical formulation

2.1. Diffusion-based Light Transport Model

Light transport in breast tissue can be described accurately by the Diffusion Equation (DE), which is an approximation to the Radiative Transport Equation (RTE) [24

H. Jiang, K. D. Paulsen, U. Osterberg, B. W. Pogue, and M. S. Patterson, “Optical image reconstruction using frequency domain data: simulations and experiments,” J. Opt. Soc. Am. A 13, 253–266 (1996). [CrossRef]

]. In the frequencydomain, the DE is given by

- \nabla . D (r) \nablaΦ (r, ω) + (μ_{a} (r) + \frac{iω}{c}) Φ (r, ω) = Q_{o} (r, ω)

(1)

where Φ(r,ω) is the photon density at position r and light modulation frequency is given by ω (in this work,ω = 100MHz). The isotropic source term is represented by Q_o (r,ω) and speed of light in tissue by c. μ_a(r) is the optical absorption coefficient and D(r) is the optical diffusion coefficient, which is defined as

D (r) = \frac{1}{3 [μ_{a} (r) + {μ'}_{s} (r)]}

(2)

where μ´_s (r) is the reduced scattering coefficient, which is defined as μ´_s = μ_s(1-g). μ_s is the scattering coefficient and g is the anisotropy factor. A Robin-type (Type-III) boundary condition is applied to exactly model the refractive-index mismatch at the boundary [25

M. Schweiger, S. R. Arridge, M. Hiroaka, and D. T. Delpy, “The Finite Element Model for the Propagation of Light in Scattering Media: Boundary and Source Conditions,” Med. Phys. 22, 1779–1792 (1995). [CrossRef] [PubMed]

]. The boundary data for a frequency domain system are the amplitude and phase of the measured signal, which is used with a Finite Element Method (FEM) based reconstruction procedure to obtain the internal spatial distributions of μ_a and μ´_s.

2.2. Standard image reconstruction

The objective function (Ω) for this procedure can be written as

Ω = min_{D, μ_{a}} {{∥ y - F (D, μ_{a}) ∥}^{2} + λ {∥ (D, μ_{a}) - (D_{0}, μ_{a 0}) ∥}^{2}}

(3)

Where, F is the forward operator that generates the model response and y is the experimental measured data. ∥.∥² represents L2-Norm of the vector. This is also known as the Tikhonov approach [27

A. N. Tikhonov, “Regularization of mathematically incorrectly posed problems,” Soviet Math. Dokl. 4, 1624–1627 (1963).

], where λ is the regularization parameter that balances the current estimate of optical properties and the initial values, and the data-model misfit. More specifically, it is the ratio of the variances of the data noise and the parameter (λ = σ² _y/σ² _(D,μa))[28

P. K. Yalavarthy, B. W. Pogue, H. Dehghani, and K. D. Paulsen, “Weight-matrix structured regularization provides optimal generalized least-squares estimate in diffuse optical tomography,” Med. Phys. 34(6), 2085–2098 (2007). [CrossRef] [PubMed]

]. Minimizing Eq. (3) (by setting first derivatives with respect to μ_a and D to zero) leads to an update equation

(J^{T} J + λI) (δ D, {δμ}_{a}) = J^{T} (y - F (D, μ_{a})) - λ [(D, μ_{a}) - (D_{0}, μ_{a 0})]

(4)

Note that deriving this update equation is not so straight forward, the full derivation is given in Ref. [28

]. In Eq. 4, J is the Jacobian matrix and I is the identity matrix. Note that J^T J is ill-conditioned; λI stabilizes the matrix. However, a slight deviation from this update equation is generally employed, which is also known as the Levenberg-Marquardt (LM) regularization method[29

K. Levenberg, “A method for the solution of certain nonlinear problems in least squares,” Q. Appl. Math. 2, 164–168 (1944).

, 30

D. W. Marquardt, “An algorithm for least squares estimation of nonlinear parameters,” J. Soc. Ind. Appl. Math. 11, 431–441 (1963). [CrossRef]

], assuming [(δD, δμ_a) = (D,μ_a)-(D ₀,μ_a0)] leading to [31

H. Jiang, K. D. Paulsen, U. Osterberg, B.W. Pogue, and M. S. Patterson, “Simultaneous reconstruction of optical absorption and scattering maps in turbid media from near-infrared frequency-domain data,” Opt. Lett. 20, 2128–2130 (1995). [CrossRef] [PubMed]

]

(J^{T} J + 2 λI) (δ D, {δμ}_{a}) = J^{T} (y - F (D, μ_{a}))

(5)

Most of the literature reports λ^{^*} ≡ 2λ [20

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

, 31

], which is true only for the Levenberg-Marquardt minimization which does not involve the parameter field in the objective function (Eq. 3)[29

K. Levenberg, “A method for the solution of certain nonlinear problems in least squares,” Q. Appl. Math. 2, 164–168 (1944).

, 30

D. W. Marquardt, “An algorithm for least squares estimation of nonlinear parameters,” J. Soc. Ind. Appl. Math. 11, 431–441 (1963). [CrossRef]

]. A detailed discussion about different least-squares minimization methods is presented in Ref. [28

]. In this LMapproach, λ^{^*} typically starts being the ratio of the variances and is reduced at each of the iterations by a small factor (in here, it is √10 and also multiplied by the maximum of the diagonal values of J^T J[21

]). The iterative procedure is repeated until experimental data matches with modeled data within a preset value ε (≈ data noise level). In general, the initial values, (D ₀,μ_a0), are obtained from a pre-reconstruction step where the data is calibrated by a homogeneous fitting procedure[32

B. W. Pogue, K. D. Paulsen, H. Kaufman, and C. Abele, “Calibration of near-infrared frequency-domain tissue spectroscopy for absolute absorption coefficient quantitation in neonatal head-simulating phantoms,” J. Biomed. Opt. 5, 185–193 (2000). [CrossRef] [PubMed]

, 33

S. Jiang, B. W. Pogue, T. O. McBride, M. M. Doyley, S. P. Poplack, and K. D. Paulsen, “Near-infrared breast tomography calibration with optoelastic tissue simulating phantoms,” J. Electron. Imaging 12, 613–620 (2003). [CrossRef]

2.3. Inclusion of a priori information: Soft-Priors

The objective function with inclusion of prior information is given as[15

]

Ω = min_{D, μ_{a}} {{∥ y - F (D, μ_{a}) ∥}^{2} + λ {∥ L [(D, μ_{a}) - (D_{0}, μ_{a 0})] ∥}^{2}}

(6)

Here also λ is the ratio of the variance of the data noise to parameter field and L is a penalty matrix (dimensionless in all the cases considered in this work) which can be derived from MRI structural priors as indicated below. The update equation resulting from this procedure is:

(J^{T} J + λ L^{T} L) (δ D, {δμ}_{a}) = J^{T} (y - F (D, μ_{a})) - λ L^{T} L [(D, μ_{a}) - (D_{0}, μ_{a 0})]

(7)

In this work, each location in the computational discretized model is labeled according to tissue type (fatty, fibroglandular or tumor) determined from MRI T1-weighted images[15

, 16

, 23

]. It was also assumed that there is no covariance between the different regions of the imaging domain. Since the domain model does not itself change throughout the iterative reconstruction algorithm, the L-matrix is calculated before the reconstruction procedure and it is used through out the process to penalize the solution. In implementation of Eq. 7 for simultaneous reconstruction of D and μ_a, using the Jacobian form described in Ref. [18

H. Dehghani, B. W. Pogue, S. P. Poplack, and K. D. Paulsen, “Multiwavelength three-dimensional near-infrared tomography of the breast: initial simulation, phantom, and clinical results,” Appl. Opt. 42, 135–145 (2003). [CrossRef] [PubMed]

] leads to the block matrices

{[\frac{H_{D} 2}{H_{D μ_{a}}}| \frac{H_{D μ_{a}}}{H_{{μ_{a}}^{2}}}] + [\frac{(λ_{D}) L^{T} L}{0}| \frac{0}{(λ_{μ_{a}}) L^{T} L}]} [\begin{matrix} δk \\ δ μ_{a} \end{matrix}] = [\begin{matrix} {(J^{T})}_{D} (y - F (D, μ_{a})) \\ {(J^{T})}_{μ_{a}} (y - F (D, μ_{a})) \end{matrix}] - [\begin{matrix} (λ_{D}) L^{T} L (D_{i - 1} - D_{0}) \\ (λ_{μ_{a}}) L^{T} L (μ_{ai - 1} - μ_{a 0}) \end{matrix}]

(8)

where H represents Hessian matrix (J^T J). Since D and μ_a are considered to be independent of each other in the estimation procedure, the cross-terms in the combined regularization matrix (λL^T L) become zero. Two forms for the L-matrix are discussed in the subsections below, including the Laplacian and Helmholtz structures.

2.3.1. a Laplacian-type structured regularization matrix

The Laplace equation in parameter u(r) can be written as

\nabla^{2} u (r) = 0

(9)

A finite difference approximation to this equation for ‘N’ number of equi-space (step size = h) nodes can be written as [34

D. R. Lynch, Numerical partial differential equations for environmental scientists and engineers - A first practical course , (Springer, Edition 1, 2005). [PubMed]

]

\nabla^{2} u (r) h^{2} \approx u_{1} + u_{2} + \dots - {Nu}_{\frac{N}{2}} + \dots + u_{N - 1} + u_{N} = 0

(10)

Dividing the whole equation by ‘-N’ leads to

\frac{- u_{1}}{N} + \frac{- u_{2}}{N} + \dots + u_{\frac{N}{2}} + \dots + \frac{- u_{N - 1}}{N} + \frac{- u_{N}}{N} = 0

(11)

The L matrix is a matrix that relates each nodal property of the numerical model to all other nodes. Therefore given a node i within the mesh, its relationship to another node j having Laplacian structure (Eq. 11) within the same mesh can be given as [15

, 16

]

L_{ij} = {\begin{matrix} \begin{matrix} 0 & if i and j are not in the same region \end{matrix} \\ \begin{matrix} - \frac{1}{N} & if i and j are in the same region \end{matrix} \\ \begin{matrix} 1 & if i = j \end{matrix} \end{matrix}

(12)

where N is the number of finite element mesh nodes comprising a given region. In this case, L^T L approximates a second-order Laplacian smoothing operator within each region, and functionally works to average the update within a region, while allowing discontinuity between different regions.

2.3.2. a Helmholtz-type structured regularization matrix

The Helmholtz equation in parameter u(r) for a damped wave can be written as

\nabla^{2} u (r) - k^{2} u (r) = 0

(13)

where κ is the wave number, specifically κ = 1/l, where l covariance length[34

D. R. Lynch, Numerical partial differential equations for environmental scientists and engineers - A first practical course , (Springer, Edition 1, 2005). [PubMed]

]. l also represents the decay length scale over which the parameter u(r) has correlation. Making κ = 0, will give the Laplace equation (Eq. 9). A finite difference approximation to this equation for ‘N’ number of equi-space (step size = h) nodes can be written as [34

D. R. Lynch, Numerical partial differential equations for environmental scientists and engineers - A first practical course , (Springer, Edition 1, 2005). [PubMed]

]

(\nabla^{2} - κ^{2}) u (r) h^{2} \approx u_{1} + u_{2} + \dots + [- (N + {(κh)}^{2})] u_{\frac{N}{2}} + \dots + u_{N - 1} + u_{N} = 0

(14)

Dividing the whole equation by ‘-(N +(κh)²)’ gives

\frac{u_{1}}{- (N + {(κh)}^{2})} + \frac{u_{2}}{- (N + {(κh)}^{2})} + \dots + u_{\frac{N}{2}} + \dots + \frac{u_{N - 1}}{- (N + {(κh)}^{2})} + \frac{u_{N}}{- (N + {(κh)}^{2})} = 0

(15)

Writing this as a generalized L-matrix form similar to Eq. 12

L_{ij} = {\begin{matrix} \begin{matrix} 0 & if i and j are not in the same region \end{matrix} \\ \begin{matrix} - \frac{1}{N + {(κ h)}^{2}} & if i and j are in the same region \end{matrix} \\ \begin{matrix} 1 & if i = j \end{matrix} \end{matrix}

(16)

For the FEM nodes case, h is chosen to be the distance between the nodes. Moreover, κ = 1/l is generally chosen to be the inverse of the size of the feature (tumor in this case) in the imaging domain. This can be seen as using the best priors for the estimation of optical properties[28

], termed as best priors estimate (BPE). As the prior structural information is available through MRI, l is chosen to be the diameter of the target (tumor) in the BPE case. In this case, L^T L (L given by Eq. 16) approximates a second-order Helmholtz smoothing operator. To determine the effect of κ on the parameter reconstruction, different values for κ are chosen. It is shown that for small values of κ, which corresponds to a large correlation length (l), both Laplacian and Helmholtz structures recover the same optical property distribution.

2.4. Inclusion of a priori information: hard-priors

Reduction of parameter space to the number of regions segmented from a high resolution imaging modality is known as hard-priors[28

]. In this section, a detailed formulation for this approach is given. The estimation optical properties (D and μ_a) in this procedure is performed through LM minimization, which was discussed in detail in Ref. [28

]. The update equation in this case is given by Eq. 5, with Jacobian having dimension of (2^{^*}NM)×(2^{^*}NR) instead of (2^{^*}NM)×(2^{^*}NN). In here, NM, NN and NR represent the number of measurements, number of FEM nodes, and number of regions, respectively. The multiplication factor 2 in the dimensions results from the treatment of amplitude and phase separately for the frequency-domain signal, and the estimated parameters being D and μ_a. This transformation requires multiplication of original Jacobian, J, by region mapper matrix (R)[22

\tilde{J} = JR

(17)

where the R has dimension of (2^{^*}NN)×(2^{^*}NR), having a block form

R = [\begin{matrix} \sum \\ \sum \end{matrix}]

. J can be also seen as a block matrix with J = [J_D J_{μ_a} ]. If RL_j represented a segmented region (j) in the FEM mesh for j = 1, 2, …NR, then ∑ has the from

\sum_{ij} = {\begin{matrix} 1 if i \in R L_{j} \\ 0 otherwise \end{matrix}

(18)

Overall, the new Jacobian J͂matrix elements were produced by adding the sensitivity of nodes belonging to the same region[22

]. Using J͂ instead of J, in Eq. 5 gives the update of optical properties of the segmented regions (δ͂D,δμ͂_a) having dimension of (2^*NR)x1. To map back this solution vector to the original dimensions (2^*NN)×1, requires

(δ D, {δμ}_{a}) = R (\tilde{δD}, \tilde{δ μ_{a}})

(19)

It is also important to relaize that even though, NR«NM, solving this might still be illposed[20

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

], leading to LM update equation[28

3. Methods: simulations and experiments

3.1. Breast geometry - Effect of imperfect a priori information

The techniques described in Sec. 2 were used to reconstruct images from synthetic data with 1% random noise added. Numerical experiments using synthetic data generated on MRI T1-weighted breast images with incorrect priors to show the effectiveness of soft-priors. Figure 1(a) (first column) shows the original distribution of three tissue layers, namely fatty (μ_a = 0.006 mm^-1 and μ´_s = 0.6 mm^-1), fibro glandular (μ_a = 0.012 mm^-1 and μ´_s = 1.2 mm^-1), and tumor (μ_a = 0.018 mm^-1 and μ´_s = 1.8 mm^-1) for the breast geometry (labeled as 0, 1 and 2 respectively, in Fig. 1(a) first column). Sixteen light collection/delivery fibers were arranged equally spaced on a circle (indentions in Fig. 1(a)). In succession, one fiber was used as the source while all other fibers served as detectors which provided a total of 240 measurements. In these studies, the source was modeled as a Gaussian profile with a Full Width Half Maximum (FWHM) of 3 mm to most accurately represent the light applied in the clinical system[35

T. O. Mcbride, B. W. Pogue, S. Jiang, U. L. Osterberg, and K. D. Paulsen, “A parallel-detection frequencydomain near-infrared tomography system for hemoglobin imaging of the breast in vivo,” Rev. Sci. Instrum. 72, 1817–1824 (2001). [CrossRef]

], and was placed at a depth of one transport scattering distance from the tissue boundary. A mesh of 1785 nodes (corresponding to 3418 linear triangular elements) with uniform nodal density was used for the diffusion model and reconstruction calculations[36

P. K. Yalavarthy, H. Dehghani, B. W. Pogue, and K. D. Paulsen, “Critical computational aspects of near infrared circular tomographic imaging: Analysis of measurement number, mesh resolution and reconstruction basis,” Opt. Express 14, 6113–6127 (2006). [CrossRef] [PubMed]

]. A total of 7% of the glandular layer (label-1) FEM nodes were labeled (relative to the original glandular layer nodes) as fat (label-0) to introduce imperfect structural priors.

Fig. 1. (a) Simulated μ_a and μ´_s distributions from a breast (obtained from a volunteer) are shown in the first column. Optical properties for the region labeled ‘0’ (fat) are: μ_a = 0.006 mm^-1 and μ_s = 0.6 mm^-1. Region ‘1’ (fibroglandular) values are: μ_a = 0.012 mm^-1 and μ´_s = 1.2 mm^-1. Region ‘2’ (tumor) values are: μ_a = 0.018 mm^-1 and μ´_s = 1.8 mm^-1. Reconstructed μ_a and μ´_s images from different techniques with simulated data having 1% random noise and imperfect structural information in defining region ‘1’ (7% reduction compared to the original segmentation) are shown in the rest of the columns. The middle two columns use soft prior structural information while the last column shows the result with hard prior information. In the Helmholtz case, κ = 1/8 mm^-1 (BPE) was used. (b) Cross-sectional plots, along the dotted line in the actual image (see first column of (a)), of true and reconstructed μ_a and μ´_s distributions.

The same initial estimates (optical properties of region ‘0’) were used as homogeneous starting conditions. The iterative procedure was stopped once the data-model misfit (residual) did not improve by more than 2% when compared with the previous iteration. The starting value for λ is chosen to be 25000 and 75 for μ_a and D respectively, derived from the noise characteristics[28

], for Eq. 5. For the case of Eq. 7, a starting value of 10 was chosen and was decreased as the iterations progressed, this procedure is outlined in Ref. [21

Table 1. Mean and standard deviation of the reconstructed (a) μ_a and (b) μ´_s values in different regions (labeled in first column of Fig. 1(a)) recovered with simulated data having 1% random noise and imperfect structural information defining region ‘1’ (7% reduction compared to the original segmentation). The corresponding reconstructed images are shown in Fig. 1(a)

Methods	Region-0	Region-1	Region-2
Actual	0.006	0.012	0.018
Laplacian	0.0064±0.0010	0.0117±0.0018	0.0156±0.0018
Helmholtz (κ = 1/8)	0.0062±0.0011	0.0120±0.0020	0.0156±0.0017
Hard Priors	0.006	0.0118	0.0843

Methods	Region-0	Region-1	Region-2
Actual	0.6	1.0	1.8
Laplacian	0.63±0.10	1.13±0.18	1.67±0.23
Helmholtz (κ = 1/8)	0.64±0.09	1.09±0.16	1.64±0.22
Hard Priors	0.64	1.13	0.23

3.2. Breast geometry - effect of data noise level

The techniques described in Sec. 2 were used to reconstruct images from synthetic data with 1, 3, 5 and 10%Gaussian distributed noise to see the effect of data noise level on the reconstruction techniques. The breast geometry (and the optical properties) were equivalent to the previous section (Fig. 1(a)), however, perfect spatial priors were used. The same FEMmesh as described above was employed in the forward and reconstruction problems. Optical properties of region ‘0’ were used as initial guess for the iterative procedure. The regularization parameter (λ) and stopping criterion was chosen according to each data noise level.

3.3. Phantom studies

A multi-layered gelatin phantom (86 mm diameter, 25 mm height) was fabricated with different optical properties using heated mixtures of water (80%), gelatin (20%) (G2625, Sigma Inc.), India ink for absorption, and TiO₂ (titanium oxide powder, Sigma Inc.) for scattering[15

] (see Fig. 3). Different layers of gelatin were constructed by successively hardening gel solutions containing different amounts of ink and TiO₂. A cylindrical hole (diameter of 16 mm and height of 24 mm) was filled with liquid to mimic the tumor. The first column in Fig. 4 shows the axial cross-section of three-layers of the phantom (Fig. 3) where the region labeled ‘0’ has the optical properties (μ_a = 0.0065 mm^-1 and μ´_s = 0. 65 mm^-1), similar to the typical fatty layer in the breast[6

] and a thickness of 10 mm. The fibroglandular layer (diameter 76mm) has optical properties (region labeled ‘1’) of μ_a = 0.01 mm^-1 and μ´_s = 1.0 mm^-1. The tumor (represented by the region labeled ‘2’) has a diameter of 16 mm with optical properties of μ_a = 0.02 mm^-1 and μ´_s = 1.2 mm^-1. The optical properties were validated by measuring large cylindrical samples of each layer. Appropriate mixtures of Intralipid and India ink were used to achieve the desired optical parameters of the tumor. Data was acquired using a clinical NIR system[35

] where the fibers were marked and photographed to extract region information (analogous to MRI images). This regional information was used to label the corresponding regions in the FEM mesh[15

]. A mesh of 1785 nodes (corresponding to 3418 linear triangular elements) was used for the diffusion model calculations and a mesh having 1360 nodes was used in the reconstruction[36

]. The meshes considered in this work have uniform nodal density across the total computation volume. NIR data was calibrated using a reference homogenous phantom to obtain initial optical properties estimates and minimize the variation between the 16 optical channels according to standard practice in human imaging studies[32

, 33

Fig. 2. (a) Reconstructed μ_a and μ´_s images from different techniques with simulated data having 5% random noise and perfect structural priors (actual images are shown in the first column of Fig. 1(a)). The first column shows the reconstruction results without the use of prior information. The middle two columns use soft prior structural information while the last row shows the result with hard prior information. In the Helmholtz case, κ = 1/8 mm^-1 (BPE) was used. (b) The mean values and standard deviations (plotted as error bars) in μ_a and μ´_s for different regions of breast geometry (labeled in actual image) with increasing noise level (1% to 10 %).

Fig. 3. Photograph for gelatin phantom (representing the idealized two-dimensional cross-sectional geometry shown as first column in Fig. 4(a)) used in the experimental studies.

4. Results

Reconstructed μ_a and μ´_s images obtained from the noisy simulated data with imperfect (7% error) glandular layer priors using the methods described in Sec. 2 are shown in Fig. 1. Using hard priors, the total number of unknowns is reduced to 6 parameters (ma and μ´_s for each of the 3 regions) and these images are presented in the last column of Fig. 1(a). The images from two different approaches of soft-priors are shown in the middle 2 columns; the first column shows the expected results. For the Helmholtz case, κ = 1/8 mm^-1 (BPE) was used, where 8 mm is the diameter of the tumor. Cross-sectional plots of the reconstructed μ_a and μ´_s distributions along the dotted line in Fig. 1(a) (first column) are provided in Fig. 1(b). Table 1(a) and (b) show the mean and standard deviation of the optical property estimations in each region of the reconstructed images. Note that the NIR reconstruction procedure without prior information (not shown) did not generate meaningful images in this complex case.

Figure 2(a) shows the reconstructed images from the data with 5% noise using all four techniques described in Sec. 2. The first column of Fig. 1(a) shows the actual distribution of optical properties. Figure 2(b) shows the mean and standard deviation values (as error bars) of reconstructed images using different techniques for different regions of the breast geometry with increasing data noise level. In the Helmholtz case, κ = 1/8 mm^-1 (BPE) was used. The actual values are also plotted for the comparison. It can be clearly seen Laplacian regularization gives lesser standard deviation for the absorption coefficient (μ_a) reconstruction compared to the Helmholtz structure. On the other hand, Helmholtz structure performs better than Laplacian in the case of scattering coefficient (μ´_s ).

A photograph of the phantom used to collect data at 785nm with 16 fibers in a single plane (giving 240 measurements) is shown in Fig. 3. Images obtained from the procedures described in Sec. 2 are presented in Fig. 4(a) along with cross-sectional plots of the reconstructed μ_a and μ´_s distributions in Fig. 4(b). Table 2(a) and (b) show the mean and standard deviation of the optical property estimations in each region of the reconstructed images. Here also, for Helmholtz case, the BPE (κ = 1/16 mm^-1) was used. Figure 5(a) gives the results for different values of κ (given on the top of the each column, true distribution is shown as first column in Fig. 4(a)) in the Helmholtz case. Corresponding cross-sections are plotted in Fig. 5(b). The mean and standard deviation values for each of this case are also given in Table 2(a) and (b).

5. Discussion

The reconstructed results (Fig. 1, 2, 4 and 5) show that the structural priors improve the reconstructed image quality dramatically. The penalized problem formulation (different type of regularizations) generates smoother images resulting in smaller standard deviations from the mean values (see Table-2) as compared to the generalized problem that does not incorporate prior information.

The hard-prior case produces significant optical property value error when the structural a priori information is imperfect in the breast geometry (Fig. 1 and Table-1). In this case, a 7% variation in the definition of the glandular layer caused false estimates of the optical properties. On the other hand, soft-priors (Laplacian and Helmholtz) yield good estimates for each layer. Hard-priors over-estimate the tumor absorption coefficient by 360% and under-estimate its scattering coefficient by 88% (Fig. 1). Soft-priors are within 6% of the expected values even with the error in the structural prior. Note that in this particular case, hard-priors failed when 7% of glandular layer made as a fatty layer, yet below this error value it gave reasonable estimates of optical properties. It is also important to realize that this is only a test case of the spatial-prior error progation. There is not enough in-vivo data for a systematic investiagion of spatial-error propagation in the reconstruction procedures, especially in the hard-priors case. With perfect structural priors, increasing the data noise level also increases the quantification error in the reconstructed images (Fig. 2(a) and(b)). Hard-Priors clearly fail in quantifying the tumor optical properties as the data noise level increases. Soft-Priors does better in the quantification than Hard-Priors. It is also evident that incorporation of structural information is key for accurate quantification of the optical properties. The experimental results (Fig. 4) from the three-layer gelatin phantom (Fig. 3) also show that incorporation of perfect structural information improves the quantification and quality of the reconstructed images. Specifically, the mean and standard deviation of the reconstructed optical property values in each region are both more accurate and more precise where the priors are included. The mean values (Table-2) show that the absorption coefficient for the Region-0 (fat) is under-estimated by 77% in the case of the Helmholtz regularization (BPE). This error can be explained by the fact that the Euclidean distance between the nodes was used rather than the distance between the nodes along the boundary of a particular region. Both Laplacian and Hard-Priors over estimate the scattering coefficient of the tumor region by a factor of 2. It is known that the photon path length is affected by the scattering coefficient. By constraining the problem based on the distance, one can expect to estimate the scattering better. In this experimental case, table-2 indicates that the Helmholtz technique produces more quantitative accuracy for the scattering coefficient estimation of the tumor and the Laplacian technique is best for the absorption coefficient estimation (as well as Fig. 2 and 4). Usage of the Helmholtz-type form for the diffusion part and the Laplacian-type form for the absoprtion part (in Eq. 8), may lead to a better estimation of both optical properties, which will be considered as a future work. Under estimation of optical properties in the regions ‘0’ and ‘1’ can be attributed to the random fluctuations in heating/cooling of the gelatin solution especially in making layered phantom, which eventually changes optical properties. Moreover, ink and TiO₂ particles settling to the bottom, making the mixture less homogenous also causes the change optical properties. These type of errors are more common for layered gelatin phantoms[37

B. Brooksby, Combining near infrared tomography and magnetic resonance imaging to improve breast tissue chromophore and scattering assessment , Ph.D. Thesis, Dartmouth College, May 2005. [PubMed]

Fig. 4. (a) Actual μ_a and μ´_s distributions (axial cross-section) of phantom (Fig. 3) case are shown in the first column. Optical properties for the region labeled ‘0’ are: μ_a = 0.0065 mm^-1 and μ´_s = 0. 65 mm^-1. Region ‘1’ values are: μ_a = 0.01 mm^-1 and μ´_s = 1.0 mm^-1. Region ‘2’ (tumor) values are: μ_a = 0.02 mm^-1 and μ´_s = 1.2 mm^-1. Reconstructed μ_a and μ´_s distribution from different techniques (discussed in Sec. 2) from the experimental phantom data. Second column of images does not use prior information. The middle rows use soft prior structural information and the last row of images were recovered with hard priors. In the Helmholtz case, κ = 1/16 mm^-1 (BPE) was used. (b) Cross-sectional plots along the dotted line in the actual image (see first column of (a)) of the true and reconstructed μ_a and μ´_s distributions.

Fig. 5. (a) Reconstructed μ_a and μ´_s images from the experimental phantom data using Helmholtz-type regularization matrix for different values of κ, which are given at the top of each column. (b) Cross-sectional plots along the dotted line of the actual images in Fig. 4(a) (first column) are shown with the data from reconstructed μ_a and μ´_s images in (a). The best prior estimate (BPE) case (κ = 1/16 mm^-1) is also presented for comparison.

Table 2. Mean and standard deviation of the reconstructed (a) μ_a and (b) μ´_s values in different regions (labeled in first column of Fig. 4(a)) recovered from the experimental phantom data. The corresponding reconstructed images are shown in Fig. 4(a) and 5(a).

Methods	Region-0	Region-1	Region-2
Actual	0.0065	0.01	0.02
No Priors	0.0025±0.0010	0.0045±0.0022	0.0120±0.0090
Laplacian	0.0031±0.0002	0.0051±0.0005	0.0174±0.0029
Helmholtz (κ = 1/16)	0.0015±0.0005	0.0058±0.0009	0.0241±0.0043
Hard Priors	0.0032	0.005	0.0213
Helmholtz (κ = 1/5)	0.0009±0.0006	0.0061±0.0008	0.0191±0.0031
Helmholtz (κ = 1/43)	0.0027±0.0003	0.0052±0.0007	0.0234±0.0043
Helmholtz (κ = 1/86)	0.0022±0.0005	0.0061±0.0032	0.0192±0.0044

Methods	Region-0	Region-1	Region-2
Actual	0.65	1.0	1.2
No Priors	0.64±0.40	0.66±0.27	0.76±0.16
Laplacian	0.38±0.03	0.63±0.07	2.37±0.41
Helmholtz (κ = 1/16)	0.46±0.02	0.59±0.02	1.08±0.12
Hard Priors	0.37	0.63	2.74
Helmholtz (κ = 1/5)	0.60±0.01	0.57±0.01	0.82±0.06
Helmholtz (κ = 1/43)	0.39±0.03	0.63±0.03	1.19±0.13
Helmholtz (κ = 1/86)	0.39±0.03	0.63±0.04	1.21±0.14

Theoretically Helmholtz and Laplacian structures are identical when κ = 0 (equivalently λ is large). Figures 5(a) and (b) (as well as Table-2) show that when κ = 1/86 mm^-1 (l = 86 mm is the diameter of the domain), Laplacian and Helmholtz structures give reasonably close reconstruction values of optical properties, which indicates the expected trend presented in this paper for the two methods. It also indicates that the BPE case (κ = 1/16 mm^-1) gives the best results, as the priors applied are close to the true structure of the feature (tumor). These results also indicate that unreasonable constraints (like κ = 1/5 mm^-1, Fig. 5(a) first column) makes the estimation problem amplify the noise resulting in physiologically invalid (scattering coefficient is greater for fatty layer compared to fibroglandular layer) estimates of optical properties. In the tumor region, as κ decreases, the estimated values of optical parameters become closer to expected values (Table-2 and Fig. 5, κ = 1/43 and κ = 1/86 cases). This is due to the correlation length becoming larger, making the covariance in the neighboring nodes larger.

In this study, it is shown that imperfect priors (commonly caused by improper image segmentation and image artifacts in MRI or X-ray) can lead to error-prone results in the hard-prior case whereas soft-priors are more immune to uncertainty in the prior data. It is also shown that the techniques used to incorporate the soft structural prior information influences the image outcome, which may lead to improvements in image accuracy if properly implemented. Srinivasan et. al[38

S. Srinivasan, B. W. Pogue, H. Dehghani, S. Jiang, X. Song, and K. D. Paulsen, “ Effect of image reconstruction bias upon spectroscopy-based quantification of chromophores in near infrared tomography,” Proc. OSA Biomedical Topical Meetings, OSA Technical Digest , WB3:1-3, Optical Society of America, Washington, DC (2004).

] have found that 5% error in optical properties introduces approximately 45% error in the chromophore concentration when a limited number of wavelengths are used for imaging. The correct “soft” inclusion of a priori information therefore can be expected to lead to a more accurate quantification of chromophore concentrations as well. It should be noted that over-weighting of the penalty term in the problem formulation may make the solution ignore the data-model misfit and emphasize smooth feature extraction. The techniques developed in this work were applied for two-dimensional test objects, and can be easily extended to threedimensional case. A more extensive study of this is left for future investigations.

6. Conclusions

This work has investigated several ways of incorporating structural information into an iterative image reconstruction. The results have been supported by gelatin phantom experiments that represent multi-layered structure which is commonly found in breast tissue with adipose (fatty) tissue on the exterior and fibroglandular tissue nearer to the interior. Soft-priors allow the tissue optical properties to vary within predefined regions, unlike hard-priors which constrain these zones to be homogeneous. Hard-priors were found to perform poorly when the prior information contained area errors as small as 7% which can easily be produced by most segmentation algorithms. True boundary extraction from MRI images introduces unavoidable segmentation and discretization errors that are better tolerated when the structural information is encoded through the soft-prior approach involving a penalty matrix.

The results reported here indicate that the optical properties of different tissue types can be quantified more accurately when their estimation is properly guided by “soft” structural a priori information. The problem formulation and results presented in this work indicate that data from other imaging modalities such as ultrasound or x-ray tomography, could also be used as the source of the structural prior. In the experimental cases investigated, the Helmholtz structure using best priors gave a better estimation of scattering coefficient of the target (tumor). However, the Laplacian type of regularization leads to more superior absorption coefficient estimate. So it is reasonable to conclude that Laplacian structure gives the best estimates of total hemoglobin concentration (Hb_T ), hemoglobin oxygen saturation (S_tO₂%) and water fraction (H₂O) (which are the main absorbers). Helmholtz structure gives the best estimates of the scattering power and scattering amplitude (scattering parameters). The framework presented here can also be extended to other regularization terms such as total variation minimization or spectral prior constraints, which may be studied in future work.

Acknowledgments

Authors are grateful to Professor Daniel R. Lynch for the useful discussions and comments on this paper. P.K.Y. acknowledges the DOD Breast Cancer predoctoral fellowship (BC050309). This work has been sponsored by the National Cancer Institute through grants RO1CA78734, PO1CA80139, and DAMD17-03-1-0405.

References and links

1.	D. A. Boas, D. H. Brooks, E. L. Miller, C. A. DiMarzio, M. Kilmer, R. J. Gaudette, and Q. Zhang, “Imaging the body with diffuse optical tomography,” IEEE Sig. Proc. Mag. 18, 57–75 (2001). [CrossRef]
2.	B. W. Pogue, M. Testorf, T. McBride, U. Osterberg, and K. D. Paulsen, “Instrumentation and design of a frequency-domain diffuse optical tomography imager for breast cancer detection,” Opt. Express 1, 391–403 (1997). [CrossRef] [PubMed]
3.	S. Srinivasan, B. W. Pogue, S. Jiang, H. Dehghani, C. Kogel, S. Soho, J. J. Gibson, T. D. Tosteson, S. P. Poplack, and K. D. Paulsen, “Interpreting Hemoglobin andWater Concentration, Oxygen Saturation and Scattering Measured in Vivo by Near-Infrared Breast Tomography,” Proc. Nat. Acad. Sci. U.S.A. 100, 12349–12354 (2003). [CrossRef]
4.	Q. Zhu, N. G. Chen, and S. C. Kurtzman, “Imaging tumor angiogenesis by use of combined near-infrared diffusive light and ultrasound,” Opt. Lett. 28, 337–339 (2003). [CrossRef] [PubMed]
5.	Q. Zhang, T. J. Brukilacchio, A. Li, J. J. Stott, T. Chaves, E. Hillman, T. Wu, M. Chorlton, E. Rafferty, R. H. Moore, D. B. Kopans, and D. A. Boas, “Coregistered tomographic x-ray and optical breast imaging: initial results,” J Biomed. Opt. 10, 024033:1-9 (2005). [CrossRef] [PubMed]
6.	B. Brooksby, S. Jiang, C. Kogel, M. Doyley, H. Dehghani, J. B. Weaver, S. P. Poplack, B. W Pogue, and K. D Paulsen, “Magnetic resonance guided near infrared tomography of the breast,” Rev. Sci. Instrum. 75, 5262–5270 (2004). [CrossRef]
7.	A. Li, E. L. Miller, M. E. Kilmer, T. J. Brukilacchio, T. Chaves, J. Stott, Q. Zhang, T. Wu, M. Chorlton, R. H. Moore, D. B. Kopans, and D. A. Boas, “Tomographic optical breast imaging guided by three-dimensional mammography,” Appl. Opt. 42, 5181–5190 (2003). [CrossRef] [PubMed]
8.	A. Li, G. Boverman, Y. Zhang, D. Brooks, E. L. Miller, M. E. Kilmer, Q. Zhang, E. M. C. Hillman, and D. A. Boas, “Optimal linear inverse solution with multiple priors in diffuse optical tomography,” Appl. Opt. 44, 1948–1956 (2005). [CrossRef] [PubMed]
9.	V. Ntziachristos, A. G. Yodh, M. Schnall, and B. Chance, “Concurrent MRI and diffuse optical tomography of breast after indocyanine green enhancement,” Proc. Nat. Acad. Sci. U.S.A. 97, 2767–2772 (2000). [CrossRef]
10.	G. Boverman, E. L. Miller, A. Li, Q. Zhang, T. Chaves, D. H. Brooks, and D. Boas, “Quantitative spectroscopic diffuse optical tomography of the breast guided by imperfect a priori structural information,” Phys. Med. Biol. 50, 3941–3956 (2005). [CrossRef] [PubMed]
11.	R. L. Barbour, H. L. Graber, J. W. Chang, S. L. S. Barbour, P. C. Koo, and R. Aronson, “MRI-guided optical tomography: Prospects and computation for a new imaging method,” IEEE Comp. Sci. Eng. 2, 63–77 (1995). [CrossRef]
12.	A. H. Hielscher and S. Bartel, “Use of penalty terms in gradient-based iterative reconstruction schemes for optical tomography,” J. Biomed. Opt. 6, 183–192 (2001). [CrossRef] [PubMed]
13.	V. Ntziachristos, X. H. Ma, and B. Chance, “Time-correlated single photon counting imager for simultaneous magnetic resonance and near-infrared mammography,” Rev. Sci. Instrum. 69, 4221–4233 (1998). [CrossRef]
14.	B. Brooksby, H. Dehghani, B. W. Pogue, and K. D. Paulsen, “Near infrared (NIR) tomography breast image reconstruction with a priori structural information from MRI: algorithm development for reconstructing heterogeneities,” IEEE J. Sel. Top. Quantum Electron. 9, 199–209 (2003). [CrossRef]
15.	B. Brooksby, S. Jiang, H. Dehghani, B. W. Pogue, K. D. Paulsen, J. Weaver, C. Kogel, and S. P. Poplack, “Combining near infrared tomography and magnetic resonance imaging to study in vivo breast tissue: implementation of a Laplacian-type regularization to incorporate magnetic resonance structure,” J. Biomed. Opt. 10, 051504:1-10 (2005). [CrossRef] [PubMed]
16.	B. Brooksby, B. W. Pogue, S. Jiang, H. Dehghani, S. Srinivasan, C. Kogel, T. D. Tosteson, J. Weaver, S. P. Poplack, and K. D. Paulsen, “Imaging breast adipose and fibroglandular tissue molecular signatures using Hybrid MRI-Guided Near-Infrared Spectral Tomography,” Proc. Nat. Acad. Sci. U.S.A. 103, 8828–8833 (2006). [CrossRef]
17.	M. F. Ernst and J. A. Roukema, “Diagnosis of non-palpable breast cancer: a review,” The Breast 11, 13 (2002). [CrossRef]
18.	H. Dehghani, B. W. Pogue, S. P. Poplack, and K. D. Paulsen, “Multiwavelength three-dimensional near-infrared tomography of the breast: initial simulation, phantom, and clinical results,” Appl. Opt. 42, 135–145 (2003). [CrossRef] [PubMed]
19.	X. Wang, B.W. Pogue, S. Jiang, X. Song, K.D. Paulsen, C. Kogel, S.P. Poplack, and W.A. Wells, , “Approximation of Mie scattering parameters in near-infrared tomography of normal breast tissue in vivo,” J. Biomed. Opt. 10, 051704:1-8 (2005). [CrossRef] [PubMed]
20.	S. R. Arridge, “Optical tomography in medical imaging,” Inv. Problems 15, R41–R93 (1999). [CrossRef]
21.	S. Srinivasan, B. W. Pogue, H. Dehghani, S. Jiang, X. Song, and K. D. Paulsen, “Improved quantification of small objects in near-infrared diffuse optical tomography,” J. Biomed. Opt. 9, 1161–1171 (2004). [CrossRef] [PubMed]
22.	H. Dehghani, B.W. Pogue, J. Shudong, B. Brooksby, and K. D. Paulsen, “Three-dimensional optical tomography: resolution in small-object imaging,” Appl. Opt. 42, 3117–3126 (2003). [CrossRef] [PubMed]
23.	B. Brooksby, S. Srinivasan, S. Jiang, H. Dehghani, B. W. Pogue, K. D. Paulsen, J. Weaver, C. Kogel, and S. P. Poplack, “Spectral priors improve near-infrared diffuse tomography more than spatial priors,” Opt. Lett. 30, 1968–1970 (2005). [CrossRef] [PubMed]
24.	H. Jiang, K. D. Paulsen, U. Osterberg, B. W. Pogue, and M. S. Patterson, “Optical image reconstruction using frequency domain data: simulations and experiments,” J. Opt. Soc. Am. A 13, 253–266 (1996). [CrossRef]
25.	M. Schweiger, S. R. Arridge, M. Hiroaka, and D. T. Delpy, “The Finite Element Model for the Propagation of Light in Scattering Media: Boundary and Source Conditions,” Med. Phys. 22, 1779–1792 (1995). [CrossRef] [PubMed]
26.	K. D. Paulsen and H. Jiang, “Spatially varying optical property reconstruction using a finite element diffusion equation approximation,” Med. Phys. 22, 691–701 (1995). [CrossRef] [PubMed]
27.	A. N. Tikhonov, “Regularization of mathematically incorrectly posed problems,” Soviet Math. Dokl. 4, 1624–1627 (1963).
28.	P. K. Yalavarthy, B. W. Pogue, H. Dehghani, and K. D. Paulsen, “Weight-matrix structured regularization provides optimal generalized least-squares estimate in diffuse optical tomography,” Med. Phys. 34(6), 2085–2098 (2007). [CrossRef] [PubMed]
29.	K. Levenberg, “A method for the solution of certain nonlinear problems in least squares,” Q. Appl. Math. 2, 164–168 (1944).
30.	D. W. Marquardt, “An algorithm for least squares estimation of nonlinear parameters,” J. Soc. Ind. Appl. Math. 11, 431–441 (1963). [CrossRef]
31.	H. Jiang, K. D. Paulsen, U. Osterberg, B.W. Pogue, and M. S. Patterson, “Simultaneous reconstruction of optical absorption and scattering maps in turbid media from near-infrared frequency-domain data,” Opt. Lett. 20, 2128–2130 (1995). [CrossRef] [PubMed]
32.	B. W. Pogue, K. D. Paulsen, H. Kaufman, and C. Abele, “Calibration of near-infrared frequency-domain tissue spectroscopy for absolute absorption coefficient quantitation in neonatal head-simulating phantoms,” J. Biomed. Opt. 5, 185–193 (2000). [CrossRef] [PubMed]
33.	S. Jiang, B. W. Pogue, T. O. McBride, M. M. Doyley, S. P. Poplack, and K. D. Paulsen, “Near-infrared breast tomography calibration with optoelastic tissue simulating phantoms,” J. Electron. Imaging 12, 613–620 (2003). [CrossRef]
34.	D. R. Lynch, Numerical partial differential equations for environmental scientists and engineers - A first practical course , (Springer, Edition 1, 2005). [PubMed]
35.	T. O. Mcbride, B. W. Pogue, S. Jiang, U. L. Osterberg, and K. D. Paulsen, “A parallel-detection frequencydomain near-infrared tomography system for hemoglobin imaging of the breast in vivo,” Rev. Sci. Instrum. 72, 1817–1824 (2001). [CrossRef]
36.	P. K. Yalavarthy, H. Dehghani, B. W. Pogue, and K. D. Paulsen, “Critical computational aspects of near infrared circular tomographic imaging: Analysis of measurement number, mesh resolution and reconstruction basis,” Opt. Express 14, 6113–6127 (2006). [CrossRef] [PubMed]
37.	B. Brooksby, Combining near infrared tomography and magnetic resonance imaging to improve breast tissue chromophore and scattering assessment , Ph.D. Thesis, Dartmouth College, May 2005. [PubMed]
38.	S. Srinivasan, B. W. Pogue, H. Dehghani, S. Jiang, X. Song, and K. D. Paulsen, “ Effect of image reconstruction bias upon spectroscopy-based quantification of chromophores in near infrared tomography,” Proc. OSA Biomedical Topical Meetings, OSA Technical Digest , WB3:1-3, Optical Society of America, Washington, DC (2004).

OCIS Codes
(100.3190) Image processing : Inverse problems
(170.3010) Medical optics and biotechnology : Image reconstruction techniques
(170.3880) Medical optics and biotechnology : Medical and biological imaging
(170.6960) Medical optics and biotechnology : Tomography

ToC Category:
Medical Optics and Biotechnology

History
Original Manuscript: March 26, 2007
Revised Manuscript: June 11, 2007
Manuscript Accepted: June 11, 2007
Published: June 13, 2007

Virtual Issues
Vol. 2, Iss. 7 Virtual Journal for Biomedical Optics

Citation
Phaneendra K. Yalavarthy, Brian W. Pogue, Hamid Dehghani, Colin M. Carpenter, Shudong Jiang, and Keith D. Paulsen, "Structural information within regularization matrices improves near infrared diffuse optical tomography," Opt. Express 15, 8043-8058 (2007)
http://www.opticsinfobase.org/oe/abstract.cfm?URI=oe-15-13-8043

Sort: Year | Journal | Reset

References

D. A. Boas, D. H. Brooks, E. L. Miller, C. A. DiMarzio, M. Kilmer, R. J. Gaudette, and Q. Zhang, "Imaging the body with diffuse optical tomography," IEEE Sig. Proc. Mag. 18, 57-75 (2001). [CrossRef]
B. W. Pogue, M. Testorf, T. McBride, U. Osterberg, and K. D. Paulsen, "Instrumentation and design of a frequency-domain diffuse optical tomography imager for breast cancer detection," Opt. Express 1, 391-403 (1997). [CrossRef] [PubMed]
S. Srinivasan, B. W. Pogue, S. Jiang, H. Dehghani, C. Kogel, S. Soho, J. J. Gibson, T. D. Tosteson, S. P. Poplack and K. D. Paulsen, "Interpreting hemoglobin and water concentration, oxygen saturation and scattering measured in Vivo by Near-Infrared Breast Tomography," Proc. Nat. Acad. Sci. U.S.A. 100, 12349-12354 (2003). [CrossRef]
Q. Zhu, N. G. Chen, and S. C. Kurtzman, "Imaging tumor angiogenesis by use of combined near-infrared diffusive light and ultrasound," Opt. Lett. 28, 337-339 (2003). [CrossRef] [PubMed]
Q. Zhang, T. J. Brukilacchio, A. Li, J. J. Stott, T. Chaves, E. Hillman, T. Wu, M. Chorlton, E. Rafferty, R. H. Moore, D. B. Kopans, and D. A. Boas, "Coregistered tomographic x-ray and optical breast imaging: initial results," J Biomed. Opt. 10, 024033:1-9 (2005). [CrossRef] [PubMed]
B. Brooksby, S. Jiang, C. Kogel, M. Doyley, H. Dehghani, J. B. Weaver, S. P. Poplack, B. W Pogue, and K. D Paulsen, "Magnetic resonance guided near infrared tomography of the breast," Rev. Sci. Instrum. 75, 5262-5270 (2004). [CrossRef]
A. Li, E. L. Miller, M. E. Kilmer, T. J. Brukilacchio, T. Chaves, J. Stott, Q. Zhang, T. Wu, M. Chorlton, R. H. Moore, D. B. Kopans, and D. A. Boas, "Tomographic optical breast imaging guided by three-dimensional mammography," Appl. Opt. 42, 5181-5190 (2003). [CrossRef] [PubMed]
A. Li, G. Boverman, Y. Zhang, D. Brooks, E. L. Miller, M. E. Kilmer, Q. Zhang, E. M. C. Hillman, and D. A. Boas, "Optimal linear inverse solution with multiple priors in diffuse optical tomography," Appl. Opt. 44, 1948-1956 (2005). [CrossRef] [PubMed]
V. Ntziachristos, A. G. Yodh, M. Schnall, and B. Chance, "Concurrent MRI and diffuse optical tomography of breast after indocyanine green enhancement," Proc. Nat. Acad. Sci. U.S.A. 97, 2767-2772 (2000). [CrossRef]
G. Boverman, E. L. Miller, A. Li, Q. Zhang, T. Chaves, D. H. Brooks, and D. Boas, "Quantitative spectroscopic diffuse optical tomography of the breast guided by imperfect a priori structural information," Phys. Med. Biol. 50, 3941-3956 (2005). [CrossRef] [PubMed]
R. L. Barbour, H. L. Graber, J. W. Chang, S. L. S. Barbour, P. C. Koo, R. Aronson, "MRI-guided optical tomography: Prospects and computation for a new imaging method," IEEE Comp. Sci. Eng. 2, 63-77 (1995). [CrossRef]
A. H. Hielscher and S. Bartel, "Use of penalty terms in gradient-based iterative reconstruction schemes for optical tomography," J. Biomed. Opt. 6, 183-192 (2001). [CrossRef] [PubMed]
V. Ntziachristos, X. H. Ma, and B. Chance, "Time-correlated single photon counting imager for simultaneous magnetic resonance and near-infrared mammography," Rev. Sci. Instrum. 69, 4221-4233 (1998). [CrossRef]
B. Brooksby, H. Dehghani, B. W. Pogue, and K. D. Paulsen, "Near infrared (NIR) tomography breast image reconstruction with a priori structural information from MRI: algorithm development for reconstructing heterogeneities," IEEE J. Sel. Top. Quantum Electron. 9, 199-209 (2003). [CrossRef]
B. Brooksby, S. Jiang, H. Dehghani, B. W. Pogue, K. D. Paulsen, J. Weaver, C. Kogel and S. P. Poplack, "Combining near infrared tomography and magnetic resonance imaging to study in vivo breast tissue: implementation of a Laplacian-type regularization to incorporate magnetic resonance structure," J. Biomed. Opt. 10, 051504:1-10 (2005). [CrossRef] [PubMed]
B. Brooksby, B. W. Pogue, S. Jiang, H. Dehghani, S. Srinivasan, C. Kogel, T. D. Tosteson, J. Weaver, S. P. Poplack and K. D. Paulsen, "Imaging breast adipose and fibroglandular tissue molecular signatures using Hybrid MRI-Guided Near-Infrared Spectral Tomography," Proc. Nat. Acad. Sci. U.S.A. 103, 8828-8833 (2006). [CrossRef]
M. F. Ernst and J. A. Roukema, "Diagnosis of non-palpable breast cancer: a review," The Breast 11, 13 (2002). [CrossRef]
H. Dehghani, B. W. Pogue, S. P. Poplack, and K. D. Paulsen, "Multiwavelength three-dimensional near-infrared tomography of the breast: initial simulation, phantom, and clinical results," Appl. Opt. 42, 135-145 (2003). [CrossRef] [PubMed]
X. Wang, B.W. Pogue, S. Jiang, X. Song, K.D. Paulsen, C. Kogel, S.P. Poplack, andW.A.Wells, "Approximation of Mie scattering parameters in near-infrared tomography of normal breast tissue in vivo," J. Biomed. Opt. 10, 051704:1-8 (2005). [CrossRef] [PubMed]
S. R. Arridge, "Optical tomography in medical imaging," Inv.Problems 15, R41-R93 (1999). [CrossRef]
S. Srinivasan, B. W. Pogue, H. Dehghani, S. Jiang, X. Song, and K. D. Paulsen, "Improved quantification of small objects in near-infrared diffuse optical tomography," J. Biomed. Opt. 9, 1161-1171 (2004). [CrossRef] [PubMed]
H. Dehghani, B.W. Pogue, J. Shudong, B. Brooksby, and K. D. Paulsen, "Three-dimensional optical tomography: resolution in small-object imaging," Appl. Opt. 42, 3117-3126 (2003). [CrossRef] [PubMed]
B. Brooksby, S. Srinivasan, S. Jiang, H. Dehghani, B. W. Pogue, K. D. Paulsen, J. Weaver, C. Kogel, and S. P. Poplack, "Spectral priors improve near-infrared diffuse tomography more than spatial priors," Opt. Lett. 30,1968-1970 (2005). [CrossRef] [PubMed]
H. Jiang, K. D. Paulsen, and U. Osterberg, B. W. Pogue, and M. S. Patterson, "Optical image reconstruction using frequency domain data: simulations and experiments," J. Opt. Soc. Am. A 13, 253-266 (1996). [CrossRef]
M. Schweiger, S. R. Arridge, M. Hiroaka, and D. T. Delpy, "The Finite Element Model for the Propagation of Light in Scattering Media: Boundary and Source Conditions," Med. Phys. 22, 1779-1792 (1995). [CrossRef] [PubMed]
K. D. Paulsen and H. Jiang, "Spatially varying optical property reconstruction using a finite element diffusion equation approximation," Med. Phys. 22, 691-701 (1995). [CrossRef] [PubMed]
A. N. Tikhonov, "Regularization of mathematically incorrectly posed problems," Soviet Math. Dokl. 4, 1624- 1627 (1963).
P. K. Yalavarthy, B. W. Pogue, H. Dehghani, and K. D. Paulsen, "Weight-matrix structured regularization provides optimal generalized least-squares estimate in diffuse optical tomography," Med. Phys. 34, 2085-2098 (2007). [CrossRef] [PubMed]
K. Levenberg, "A method for the solution of certain nonlinear problems in least squares," Q. Appl. Math. 2, 164-168 (1944).
D. W. Marquardt, "An algorithm for least squares estimation of nonlinear parameters," J. Soc. Ind. Appl. Math. 11, 431-441 (1963). [CrossRef]
H. Jiang, K. D. Paulsen, U. Osterberg, B. W. Pogue, and M. S. Patterson, "Simultaneous reconstruction of optical absorption and scattering maps in turbid media from near-infrared frequency-domain data," Opt. Lett. 20, 2128- 2130 (1995). [CrossRef] [PubMed]
B. W. Pogue, K. D. Paulsen, H. Kaufman, and C. Abele, "Calibration of near-infrared frequency-domain tissue spectroscopy for absolute absorption coefficient quantitation in neonatal head-simulating phantoms," J. Biomed. Opt. 5, 185-193 (2000). [CrossRef] [PubMed]
S. Jiang, B. W. Pogue, T. O. McBride, M. M. Doyley, S. P. Poplack, and K. D. Paulsen, "Near-infrared breast tomography calibration with optoelastic tissue simulating phantoms," J. Electron. Imaging 12, 613-620 (2003). [CrossRef]
D. R. Lynch, Numerical partial differential equations for environmental scientists and engineers - A first practical course, (Springer, Edition 1, 2005). [PubMed]
T. O. Mcbride, B. W. Pogue, S. Jiang, U. L. Osterberg, and K. D. Paulsen, "A parallel-detection frequencydomain near-infrared tomography system for hemoglobin imaging of the breast in vivo," Rev. Sci. Instrum. 72,1817-1824 (2001). [CrossRef]
P. K. Yalavarthy, H. Dehghani, B. W. Pogue, and K. D. Paulsen, "Critical computational aspects of near infrared circular tomographic imaging: Analysis of measurement number, mesh resolution and reconstruction basis," Opt. Express 14, 6113-6127 (2006). [CrossRef] [PubMed]
B. Brooksby, Combining near infrared tomography and magnetic resonance imaging to improve breast tissue chromophore and scattering assessment, Ph.D. Thesis, Dartmouth College, May 2005. [PubMed]
S. Srinivasan, B. W. Pogue, H. Dehghani, S. Jiang, X. Song, and K. D. Paulsen, " Effect of image reconstruction bias upon spectroscopy-based quantification of chromophores in near infrared tomography," Proc. OSA Biomedical Topical Meetings, OSA Technical Digest, WB3:1-3, Optical Society of America, Washington, DC (2004).

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

Figures


Fig. 1.	Fig. 2.	Fig. 3.


Fig. 4.	Fig. 5.

« Previous Article | Next Article »

OSA is a member of CrossRef.

Select a Conference
Session or Paper #
Year

OpticsInfoBase

Optics Express

Structural information within regularization matrices improves near infrared diffuse optical tomography

Article Tools

Article Navigation

Abstract

1. Introduction

2. Methods: mathematical formulation

2.1. Diffusion-based Light Transport Model

2.2. Standard image reconstruction

2.3. Inclusion of a priori information: Soft-Priors

2.3.1. a Laplacian-type structured regularization matrix

2.3.2. a Helmholtz-type structured regularization matrix

2.4. Inclusion of a priori information: hard-priors

3. Methods: simulations and experiments

3.1. Breast geometry - Effect of imperfect a priori information

3.2. Breast geometry - effect of data noise level

3.3. Phantom studies

4. Results

5. Discussion

6. Conclusions

Acknowledgments

References and links

References

Appl. Opt.

IEEE Comp. Sci. Eng.

IEEE J. Sel. Top. Quantum Electron.

IEEE Sig. Proc. Mag.

J Biomed. Opt.

J. Biomed. Opt.

J. Electron. Imaging

J. Opt. Soc. Am. A

J. Soc. Ind. Appl. Math.

Med. Phys.

Opt. Express

Opt. Lett.

Phys. Med. Biol.

Problems

Proc. Nat. Acad. Sci. U.S.A.

Q. Appl. Math.

Rev. Sci. Instrum.

Soviet Math. Dokl.

The Breast

Other

Cited By

Figures