## Abstract

The development of effective multi-modality imaging methods typically requires an efficient information fusion model, particularly when combining structural images with a complementary imaging modality that provides functional information. We propose a composition-based image segmentation method for X-ray digital breast tomosynthesis (DBT) and a structural-prior-guided image reconstruction for a combined DBT and diffuse optical tomography (DOT) breast imaging system. Using the 3D DBT images from 31 clinically measured healthy breasts, we create an empirical relationship between the X-ray intensities for adipose and fibroglandular tissue. We use this relationship to then segment another 58 healthy breast DBT images from 29 subjects into compositional maps of different tissue types. For each breast, we build a weighted-graph in the compositional space and construct a regularization matrix to incorporate the structural priors into a finite-element-based DOT image reconstruction. Use of the compositional priors enables us to fuse tissue anatomy into optical images with less restriction than when using a binary segmentation. This allows us to recover the image contrast captured by DOT but not by DBT. We show that it is possible to fine-tune the strength of the structural priors by changing a single regularization parameter. By estimating the optical properties for adipose and fibroglandular tissue using the proposed algorithm, we found the results are comparable or superior to those estimated with expert-segmentations, but does not involve the time-consuming manual selection of regions-of-interest.

©2010 Optical Society of America

## 1. Introduction

Accurate diagnosis of diseases, particularly cancers, requires more detailed, disease-specific and preferably non- or minimally-invasive imaging information. The morphological images, traditionally obtained from X-ray based imaging modalities, are becoming insufficient in many cases [1]. In the past several decades, increasing research efforts in combining multiple imaging modalities, particularly those providing complementary information, have formed the grounds and active frontiers of multi-modality imaging [2]. Dual- or triple-modality systems are becoming abundant in publications and gradually become part of the standard clinical practices. Among them, a successful example is the combined Positron Emission Tomography and Computed Tomography (PET/CT) [3]. In PET/CT, the low-resolution albeit disease-specific PET image is co-registered to a high-resolution 3D CT image providing the background tissue anatomy. By interpreting the functional PET and structural CT images simultaneously, one is able to show a diagnostic accuracy that is not achievable by each individual modality [4]. PET/MRI (magnetic resonance imaging) [5], SPECT(single photon emission computed tomography)/CT [6] and MicroCT-Bioluminescence [7,8] are among other examples of multi-modality imaging that are gaining popularity for human or animal disease model studies.

Diffuse optical tomography (DOT) is another emerging technique that can provide functional assessment about tissue status [9,10]. By measuring the scattered light at multiple wavelengths within the near-infrared frequency range, one can recover hemoglobin (oxygenated, HbO, or deoxygenated, HbR), water and lipid concentrations as well as the tissue scattering properties. These measured quantities were found to be clinically relevant and are particularly promising as biomarkers for breast cancer diagnosis [11–15]. However, similar to PET, standalone DOT images suffer from low resolution due to the diffusive nature of the photons inside the tissue and the subsequent ill-posedness in solving the inverse model. Combining DOT with a high-resolution structural-oriented imaging modality can be highly beneficial [2]. Various dual-modality breast imaging systems involving DOT have been developed in the past decade and yielded different levels of success. Zhu et al. [16] built a combined DOT/Ultrasound hand-held probe for breast tumor diagnosis. Ntziachristos et al. [17], Brooksby et al. [18] and Carpenter et al. [19] have reported efforts combining DOT with MRI by either simultaneous measurement or post-registration. At Massachusetts General Hospital, we have developed a combined DOT and digital breast tomosynthesis (DBT) system for breast cancer screening [20,21]. In the past 4 years, we have been focusing on the clinical evaluation of the system by measuring healthy or tumor-bearing breasts from nearly 200 volunteers; we have shown promising diagnostic efficiency for malignant tumors using the reconstructed optical parameters from this pilot clinical study [22].

In most multi-modality DOT studies as mentioned above, the structural imaging modalities were only used to provide the exterior boundary of the breast, or overlaid on the reconstructed DOT images and used for image interpretation. This appears to be an under-use of the available information as the internal tissue structures are not fully exploited. In fact, the under-used internal structure is the key ingredient that is missing to produce DOT images of higher spatial resolution. To incorporate this ingredient into the optical image construction, novel algorithms were proposed and tested by combining a spatial prior to “inform” the optical image recovery. The simplest approach to bind a spatial structure to the DOT reconstruction is to segment the structural image and represent nodes within the same tissue segment by a single set of optical parameters, i.e. the so-called “hard-prior” approach [18,23]. Several limitations were found for this approach. The implied piecewise-constant optical properties within each tissue type is likely inaccurate in realistic biological tissue. The highly over-determined reconstruction equation effectively cancels the information from most optical measurements, which are usually precious to get. The resulting images are also intolerable to segmentation error, as it is essentially identical to the segmentation itself. To balance the contributions from the structural prior and optical measurements, “softer” constraints, in the form of regularizations, have been investigated. Li et al. has applied different regularization coefficients for different tissue types based on a segmented structural image [24]. Brooksby et al. [23] and Yalavarthy et al. [25] used a regularization matrix, defined by a Laplacian operator for the nodes sharing the same tissue type. Although these algorithms allow certain variations from the structural prior, the underlying piecewise-constant assumption remains to be arguable in realistic tissues. Only recently, researchers begin to exploit more general prior-guided DOT reconstruction by treating a target domain as statistical mixtures of different tissue types [26,27].

Here, we present a structural-prior guided DOT reconstruction algorithm using a compositional segmentation. Different from the binary approach, in a compositional segmentation we assume each pixel/node in the structural image is a combination of two or more tissue components. Using a weighted-graph representation in the compositional space, we are able to construct a regularization matrix, being the modified Laplacian operator of the weighted-graph, and use it in the DOT image reconstruction. We investigate the control of the “softness” of the structural constraints by changing the weight of the regularization. To test this algorithm, we use a clinical data set collected in the past few years with both DOT and DBT measurements. With these clinical data, we are able to identify the potential improvement in image resolution and quality using the proposed algorithm. We also calculate the estimates for the optical properties for each component of the breast tissue, i.e. adipose and fibroglandular tissue, using the compositional priors and compare them to our previous results using prior-free reconstruction and expert segmentation.

In the Methods section of the paper, we first report the theory for compositional segmentation of a structural image and then the formulations of the regularization matrix. In the Results section, we use a subset of our healthy breasts measurements (*N* = 31) to establish an empirical correlation between the X-ray image intensities for adipose and fibroglandular tissue. Using this relationship, we apply the compositional segmentation to another 58 healthy breasts from bilateral measurements of 29 subjects and study the improvement in image resolution and accuracy. Our findings are summarized in the Discussion section.

## 2. Methods

The proposed algorithm includes two essential steps. In the first step, we decompose a structural image into a compositional map using a given set of components and their pre-defined characteristics. After the decomposition, each pixel or node in the structural image is represented by a vector recording the concentration of each component. In the second step, we build a regularization matrix based on the previously generated compositional maps and run a DOT image reconstruction to produce a solution informed by the structural constraints. These steps are explained in detail in the following sub-sections.

#### 2.1 Compositional segmentation and applications to DBT images

We consider that the imaging target is composed *of* a finite number of components, where each component is independent of the others and has differentiable and additive contrast when measured in the structural imaging modality. We assume the image intensity of the target in the structural image is a linear combination of the separate contrast from each of its components, which is proportional to its concentration. We define a compositional vector at a given location *r* in the target as

where C* _{i}*(

*r*) denotes the concentration (0≤C

*(*

_{i}*r*)≤1) of the

*i*-th component at location

*r*, Ω denotes the target domain and N

_{c}is the total number of components. By definition, Eq. (1) also implies ${\sum}_{i}{C}_{i}(r)}=1$ for

*$\forall r\in \Omega $. If we assume the imaging intensity produced by the components can be represented by a transfer function f(·) governed by the biophysics of the tissue, then we have*

*where I _{s}(r) is the image intensity in the measured structural image and f(P_{i}) is the intensity anticipated if the target at r is occupied by the pure i-th component (P
_{i}). The image segmentation problem is now identified as the task to uncover C = {C_{i}(r)} from I_{s}(r) for $r\in \Omega $. If there are only two components, i.e. N_{c} = 2, from Eq. (2) and ${\sum}_{i}{C}_{i}(r)}=1$, we are able to uniquely solve for C_{i}(r) if all f(P_{i}) are known. When there are more than 2 components, additional knowledge needs to be acquired in order to uniquely solve for C_{i}(r).*

*For a healthy breast, two major tissue types are concerned: fibroglandular (denoted with subscript “ f”) and adipose (denoted with subscript “a”) tissue. From a 3D DBT image, it is possible to identify a small region where considered as pure adipose (or fibroglandular) tissue, from which, we can estimate f(P
_{a}) (or f(P
_{f})), and the compositional vector C(r) = {C_{f},(r), C_{a}(r)} can be calculated as *

*$$\begin{array}{l}{C}_{f}={I}_{s}-f({P}_{a})/[f({P}_{f})-f({P}_{a})]\\ {C}_{a}=1-{C}_{f}\end{array}$$*

*The above approach can be extended to other types of structural images, for example, anatomical MRI. One may use the Gaussian-Mixture-Model (GMM) [28] or other data-clustering algorithms to identify the compositional vectors inside the target domain [27]. For algorithm and simulation studies, one can also derive such a compositional map from a published segmented breast data archive [29].*

*When such a compositional map is available for a given target, for any reconstructed functional image within the same domain, we can compute the functional values for each “pure” component by a least-square estimation as*

*where matrix H is a horizontal concatenation of the vertical compositional vectors for all parameter nodes and H
^{+} = (H
^{T}
H)^{−1} is the Moore-Penrose pseudo-inversion of H; μ represents a vertical vector recording the reconstructed functional values at all parameter nodes.*

*2.2 Regularized image reconstruction*

*In this subsection, we illustrate the compositional-map-guided regularization using our dual-modality data processing pipeline. This approach is general and can be easily extended to other pipelines that share a similar reconstruction scheme.*

*In a prior-guided DOT image reconstruction, we first discretize the target domain by generating a pair of meshes, a fine one for the forward solution (referred to as the forward mesh), and a coarse one for the inversion (referred to as the parameter mesh) [30]. For a co-registered dual-modality experiment, we register the meshes to the structural image space. Subsequently, we map the nodes on the forward mesh to the structural image volume, and assign a compositional vector, {C_{f}(r),C_{a}(r)} for each forward node. Using the mapping between the forward and parameter meshes, we calculate the compositional vectors at each node on the parameter mesh by weighted averages.*

*In our particular pipeline, we compute the forward solutions using an FE solver and perform an iterative Gauss-Newton reconstruction to recover the 3D profiles of the optical properties on the parameter mesh. For each iteration, the general form of the update equation looks like [31]*

*where ${J}_{k}\in {C}^{M\times N}$ is the complex Jacobian matrix at the k-th iteration, M is the number of measurements and N is the number of parameters; λ is the regularization parameter and L is the regularization matrix; $\Delta {\mu}_{k}$is the update to the parameters; A is the forward matrix and $A{\mu}_{k-1}$ is the forward solution using the estimated parameters from the k-1 iteration. Note that matrix $({J}_{k}^{T}{J}_{k}+\lambda {L}^{T}L)$ has dimension of NxN. In practice, N is typically substantially bigger than the number of measurements M. In this situation, an under-determined form [32] of the update equation is more computationally efficient:*

*$$\Delta {\mu}_{k}={\Theta}^{-1}{J}^{T}{\left[{J}_{k}{\Theta}^{-1}{J}_{k}^{T}+\lambda I\right]}^{-1}(y-A{\mu}_{k-1})$$*

*where $\Theta ={L}^{T}L\in {R}^{NxN}$ and can be pre-computed for each parameter mesh, and I is an identify matrix. Very often, when there are multiple (S>1) optical properties associated with each node, for example, HbO, HbR, scattering amplitude and scattering power, the calculations in Eq. (6) in the full-matrix form can be slow as it involves the multiplications of large matrices. This can be accelerated by the block-diagonal nature of the system. We simply expand Eq. (6) by*

*$$\begin{array}{l}{\Theta}^{-1}=\left[\begin{array}{ccc}{\left({L}^{T}L\right)}^{-1}& & \\ & \ddots & \\ & & {\left({L}^{T}L\right)}^{-1}\end{array}\right]\\ {J}_{k}{\Theta}^{-1}{J}_{k}^{T}={\displaystyle \sum _{i=1}^{S}{J}_{i,k}{({L}^{T}L)}^{-1}{J}_{i,k}^{T}}\end{array}$$*

*where the sub-matrix ${J}_{i,k}$ represents the Jacobian of the i-th optical property at iteration k. Using the above block-diagonal formula leads to a significant saving in memory and computational time.*

*2.3. Compositional-map-guided regularization*

*A non-trivial regularization matrix L ($L\ne I$
) encodes the prior knowledge about the solution. When L is defined as a spatial Laplacian or Helmholtz operator [23] based on the mesh connectivity, it essentially penalizes large gradients and ensures smoothness of the solution. In the binary “soft-prior” approach [25], it penalizes the differences between the nodes falling into the same piecewise-constant segment. Here, we further extend this notion by taking advantage of the refined compositional segmentation scheme. To illustrate the method, we first construct a weighted graph in the composition space. In Fig. 1
, we plot all nodes in the parameter mesh using its compositional vector. In fact, it is a hyper-plane in the ${R}^{{N}_{c}-1}$
space defined by
${\sum}_{i}{C}_{i}(r)}=1$. For each node, we connect it to the nodes sharing similar compositions within a predefined limit, i.e $\left|\right|{C}_{i}-{C}_{j}\left|\right|<\alpha {N}_{c}$ where $\left|\right|\cdot \left|\right|$ is the L_{2} norm. Scalar $\alpha \in (0,1)$ is a regularization parameter and was heuristically set to 0.2 in our simulations. Smaller α values appear to produce sharper images Then we add an edge between the nodes satisfying the above condition and assign a weight to this edge proportional as ${u}_{i,j}=\alpha -\left|\right|{C}_{i}-{C}_{j}\left|\right|/{N}_{c}$ where N_{c} is the number of total components. To this end, we have constructed a weighted graph in the compositional space.*

*Based on [33], the (normalized) Laplacian operator, L = [l_{i,j}], for a weighted graph is defined by*

*$${l}_{i,j}=\{\begin{array}{cc}1& \text{if}i=j\\ -\frac{{u}_{i,j}}{\sqrt{{d}_{i}{d}_{j}}}& \text{if}i\text{and}j\text{areconnected}\\ 0& \text{otherwise}\end{array}$$*

*where d_{i} is the “weighted” degree of the i-th node, defined by ${d}_{i}={\displaystyle {\sum}_{n,\text{\hspace{0.17em}}\left|\right|{C}_{i}-{C}_{n}\left|\right|<\alpha {N}_{c}}{u}_{i,n}}$. To avoid degeneracy and ensure the existence of ${\left({L}^{T}L\right)}^{-1}$, we slightly modified L by*

*$${l}_{i,j}=\{\begin{array}{cc}1& \text{if}i=j\\ -\frac{{u}_{i,j}}{\beta \sqrt{{d}_{i}{d}_{j}}}& \text{if}i\text{and}j\text{areconnected}\\ 0& \text{otherwise}\end{array}$$*

*where β>1 is a constant. We choose β = 1.2 although other values may also be used. Using the L defined in Eq. (9) effectively penalizes the differences between nodes sharing similar compositions. One can combine it with the spatial Laplacian to impose additional smoothness constraints.*

*2.4. Clinical measurements for algorithm validation*

*We test the proposed regularization algorithm with a set of healthy subjects measured with our combined DOT/DBT imaging system. The details of the experiments and protocols are described in [21] and [22]. Both radio-frequency modulated (RF) and continuous-wave (CW) lasers at 685 nm and 830nm were used for the data acquisition. From these data, we recover the 3D profiles of total hemoglobin concentration (HbT = HbO + HbR), oxygen saturation (SO _{2} = HbO/HbT), scattering amplitude and scattering power of the target breast. Water and lipids volume fractions are assumed to be constants at 23% and 58%, respectively, based on our previous studies [21].*

*3. Results*

*In this section, we show the recovered images for healthy breasts using the above described algorithm and compare the output with the results from unconstrained image reconstructions. We also estimated the optical properties for adipose and fibroglandular tissue and compare the findings with the previous results computed based on expert segmentation.*

*First, using a set of 31 DBT images of healthy breast as a learning set, we study the X-ray contrasts of different tissue types and seek possibilities for an automatic compositional segmentation approach. The intensity values in all DBT images were first rescaled to the range of 0-255. For each DBT image, we manually selected two regions-of-interest (ROI) corresponding to “pure” adipose and fibroglandular tissue. Each ROI is a 2D rectangular region about 1cm x 1cm in size from a selected image slice. From these ROIs, we estimated f(P_{f}) and f(P
_{a}) values by averaging the DBT image intensities inside the ROIs. We plot the estimated f(P_{f}) and f(P
_{a}) pairs for all selected breasts in Fig. 2
.*

*From Fig. 2, we assume a linear relationship between f(P_{f}) and f(P
_{a}), i.e.*

*and the coefficients a = 1.226 and b = 13.76 are determined by linear regression from the plot in Fig. 2. The R
^{2} of the linear fitting is 0.92. Thus, to compute the compositional map from Eq. (3), we only need to select a single ROI for adipose (or fibroglandular) tissue and compute the other from Eq. (10).*

*Using the above correlation, we apply the compositional segmentation for another 29 pairs of bilateral DOT/DBT measurements. In this case, we only manually select a single ROI corresponding to the “pure” adipose tissue, estimate f(P
_{a}) from the ROI and calculate the value for f(P_{f}) from Eq. (10). The compositional maps of the 58 breasts are then obtained by Eq. (3).*

*To use the compositional maps in our FE-based image reconstructions, we first generated forward/parameter meshes using the “iso2mesh” mesh generator [34]. In all reconstructions, we first estimated a set of bulk optical properties which served as the homogeneous initial guess for the full image reconstruction. We used both RF and CW measurements at 685 nm and 830 nm in all the reconstructions. A multi-spectral algorithm [35,21] was used to estimate HbO, HbR and scattering properties from multi-wavelength measurements simultaneously. All image reconstructions ran for 5 Gauss-Newton iterations. All calculations were performed on a single core of an Intel Xeon E5530 (2.4 GHz) CPU.*

*In Fig. 3
, we show sample reconstructed HbT images (in μM) from the right breast of a 45-year-old healthy volunteer. The node numbers for forward and parameter meshes are 14824 and 1934, respectively. We solved for the forward solutions at 13/6 RF source/detector and 26/19 CW source/detector locations at two wavelengths and built the Jacobian matrix using the adjoint method [36]. The DBT image and the compositional map slice for fibroglandular tissue are shown in Figs. 3 (a) and (b). From (g) to (j), we show the image slices reconstructed with the compositional-prior-guided algorithm with different λ values, (λ = 2.5 to 0.039 with a multiplicative step factor of 0.25). In (c), we show the reconstruction output from the “binary soft-prior” algorithm (the binary segmentation of fibroglandular tissue was achieved by thresholding at C_{f} >0.5) and (d) the result without structural prior (i.e. L = I). All images were extracted at a horizontal plane z = 2.6 cm from the bottom of the compressed breast. The prior-guided image reconstructions took about 102 seconds per iteration, while that for the prior-free reconstruction is around 73 seconds. The relative residuals (normalized by the residual from a common initial guess) for the prior-guided reconstruction range from 0.27 (λ = 0.039) to 0.34 (λ = 2.5); that for the prior-free reconstruction is 0.32.*

*Notice that the compositional map in Fig. 3(b) contains an edge artifact due to the application of an edge-enhancing algorithm in the DBT image processing. We deliberately choose this case as we want to see how the algorithm performs when the compositional map contains errors.*

*Based on Fig. 3, we chose λ around 1 as the default regularization parameter. We show SO_{2} and μ’_{s} (at 830nm) for the same breast using compositional priors (λ = 1); in comparison, those without priors are shown in Figs. 3(e) and 3(f), respectively. Then we ran reconstructions for 58 healthy breasts. Similar to the previous cases, we ran 5 Gauss-Newton iterations to recover 3D volumetric images for HbO, HbR, scattering amplitude and power, respectively. Images from 4 additional healthy subjects are shown in Fig. 4
.*

*From the recovered images, we estimate the optical properties by Eq. (4) for adipose and fibroglandular tissue. In Figs. 5 (a)-(c)
, we plot the estimated adipose HbT, SO _{2} and reduced scattering coefficient (μ
_{s}’) at 830nm for left vs. right breasts (circles). We also overlapped the results estimated by the prior-free algorithm using manually created ROIs [22] (dots). We also used these same manually generated ROIs to extract the optical properties for adipose and fibroglandular tissue from the prior-guided images. These were plotted as “pluses” in Fig. 5. In Table 1
, we summarize the correlation coefficients between the optical properties of the adipose tissue between left and right breasts, as well as the parameter values recovered from the 3 approaches. Those for the fibroglandular tissue are similar (not shown).*

*4. Discussion*

*The reconstructed images using the compositional priors in Figs. 3 (g)-(j) are very encouraging and exhibit several merits over the conventional image reconstructions. First of all, for most of the images, we can clearly identify the presence of the compositional map patterns in the recovered optical images, indicating that it has indeed informed the image reconstruction. Another finding from Figs. 3 is the robustness of the algorithm with respect to segmentation error: the edge artifacts in Fig. 3(b) are effectively removed in (g)-(j) with the incorporation of the optical data. In comparison, the result from the “binary soft-prior” approach seems more sensitive to segmentation error. Comparing the images from (g)-(j) and (d), we find that by adjusting the regularization parameter λ we are able to tune the strength of the priors: with an increasing λ, the reconstructed images blend toward the compositional map as we are giving more weight to the priors; with a small λ, on the other hand, the result approaches that of the prior-free reconstruction. By comparing (g)-(j) to (c), we find that the compositional prior approach gives more spatial details than the “binary” soft-prior approach, which is essentially a smoothed version of the binary segmentation. The images for SO_{2} and μ
_{s}’ show similar improvement in resolution. In (k), the SO_{2} for the adipose tissue is generally greater than that of the fibroglandular tissue. This is consistent with our previous findings [21]. Additional images from 4 other subjects show similar improvement as we observed in Fig. 3, except for a few localized SO_{2} contrasts in (i) and (s), which may be a result of coupling differences from the nearby optodes.*

*As we discussed in [22], the bilateral correlations serve as a metric to test the robustness of the system. From Fig. 5 and Table 1, we found that the optical properties estimated by the compositional-map-guided reconstructions generally match those estimated from prior-free reconstruction with expert segmentations. Particularly, the least-square estimates using the compositional map (circles) show a great similarity to those using the same image but with ROIs (pluses). These results render the prior-guided approach a very attractive choice as it is almost entirely automated, while the conventional method requires a manual segmentation to create the adipose and fibroglandular ROIs. Moreover, the compositional-prior-guided approach offers more objectivity than expert segmentation as the latter approach depends on the experiences of the operator. We also noticed that the standard deviations for HbT and SO2 values became larger when using priors in the reconstructions. This may be potentially a result of different regularization mechanisms in these methods and deserves further investigation. The optical properties estimated using priors are generally consistent with literature findings [11,18,21–23].*

*In the next stage of the study, we will further investigate the use of this approach for accurate quantification of breast cancers. To add cancer as the third component in the segmentation, we can simply segment the tumor and superimpose it on the adipose/fibroglandular tissue compositional maps. However, in cases where the size/location of the tumor is unknown, the application of this scheme can be more challenging. We will consider a two-step reconstruction where in the first step, we treat the breast as a healthy breast and reconstruct optical images as described here; in the second step, we use the previous result as the initial guess and reduce the strengths of the prior in order to recover the contrast from the tumor. Alternatively, we can use the HbT contrast from a prior-free reconstruction to estimate the compositional map for tumor tissue as we know tumors usually have high HbT values [11,12,22]. Combining the optically derived tumor tissue composition with the DBT-derived adipose/fibroglandular compositions, we can run the prior-guided reconstruction. In either case, we will use over 50 available tumor breast measurements in our current data set to assess the improvement using the statistical test similar to our previous work [22]. In addition, we will also explore strategies for optimal selection of the regularization parameters, λ, α and β, which are currently determined in a heuristic fashion. A deterministic selection scheme is expected to further improve the robustness of the algorithm.*

*Acknowledgments*

*The authors would like to thank Jayne Cormier and Donna Burgess for the assistance in the clinical experiments for this study. The authors also acknowledge the NIH funding support under grants R01-CA97305 and U54-CA105480. QF is partly supported by Harvard Catalyst Pilot Grant.*

*References and links*

**1. **D. W. Townsend and S. R. Cherry, “Combining anatomy and function: the path to true image fusion,” Eur. Radiol. **11**(10), 1968–1974 (2001). [CrossRef] [PubMed]

**2. **F. S. Azar, and X. Intes, eds., *Translational multimodality Optical Imaging*, Artech House, Norwood (2008)

**3. **J. Czernin, and H. R. Schelbert, eds., *PET/CT in cancer patient management.* The Journal of Nuclear Medicine, 48 (2007)

**4. **M. Charron, T. Beyer, N. N. Bohnen, P. E. Kinahan, M. Dachille, J. Jerin, R. Nutt, C. C. Meltzer, V. Villemagne, and D. W. Townsend, “Image analysis in patients with cancer studied with a combined PET and CT scanner,” Clin. Nucl. Med. **25**(11), 905–910 (2000). [CrossRef] [PubMed]

**5. **M. S. Judenhofer, H. F. Wehrl, D. F. Newport, C. Catana, S. B. Siegel, M. Becker, A. Thielscher, M. Kneilling, M. P. Lichy, M. Eichner, K. Klingel, G. Reischl, S. Widmaier, M. Röcken, R. E. Nutt, H.-J. Machulla, K. Uludag, S. R. Cherry, C. D. Claussen, and B. J. Pichler, “Simultaneous PET-MRI: a new approach for functional and morphological imaging,” Nat. Med. **14**(4), 459–465 (2008). [CrossRef] [PubMed]

**6. **Z. Keidar, O. Israel, and Y. Krausz, “SPECT/CT in tumor imaging: technical aspects and clinical applications,” Semin. Nucl. Med. **33**(33 Issue 3), 205–218 (2003). [CrossRef] [PubMed]

**7. **M. Doubrovin, I. Serganova, P. Mayer-Kuckuk, V. Ponomarev, and R. G. Blasberg, “Multimodality in vivo molecular-genetic imaging,” Bioconjug. Chem. **15**(6), 1376–1388 (2004). [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. **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 Signal Process. Mag. **18**(6), 57–75 (2001). [CrossRef]

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

**11. **A. E. Cerussi, A. J. Berger, F. Bevilacqua, N. Shah, D. Jakubowski, J. Butler, R. F. Holcombe, and B. J. Tromberg, “Sources of absorption and scattering contrast for near-infrared optical mammography,” Acad. Radiol. **8**(3), 211–218 (2001). [CrossRef] [PubMed]

**12. **D. Grosenick, K. T. Moesta, M. Möller, J. Mucke, H. Wabnitz, B. Gebauer, C. Stroszczynski, B. Wassermann, P. M. Schlag, and H. Rinneberg, “Time-domain scanning optical mammography: I. Recording and assessment of mammograms of 154 patients,” Phys. Med. Biol. **50**(11), 2429–2449 (2005). [CrossRef] [PubMed]

**13. **P. Taroni, A. Torricelli, L. Spinelli, A. Pifferi, F. Arpaia, G. Danesini, and R. Cubeddu, “Time-resolved optical mammography between 637 and 985 nm: clinical study on the detection and identification of breast lesions,” Phys. Med. Biol. **50**(11), 2469–2488 (2005). [CrossRef] [PubMed]

**14. **L. C. Enfield, A. P. Gibson, N. L. Everdell, D. T. Delpy, M. Schweiger, S. R. Arridge, C. Richardson, M. Keshtgar, M. Douek, and J. C. Hebden, “Three-dimensional time-resolved optical mammography of the uncompressed breast,” Appl. Opt. **46**(17), 3628–3638 (2007). [CrossRef] [PubMed]

**15. **S. P. Poplack, T. D. Tosteson, W. A. Wells, B. W. Pogue, P. M. Meaney, A. Hartov, C. A. Kogel, S. K. Soho, J. J. Gibson, and K. D. Paulsen, “Electromagnetic breast imaging: results of a pilot study in women with abnormal mammograms,” Radiology **243**(2), 350–359 (2007). [CrossRef] [PubMed]

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

**17. **V. Ntziachristos, A. G. Yodh, M. D. Schnall, and B. Chance, “MRI-guided diffuse optical spectroscopy of malignant and benign breast lesions,” Neoplasia **4**(4), 347–354 (2002). [CrossRef] [PubMed]

**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 by using hybrid MRI-guided near-infrared spectral tomography,” Proc. Natl. Acad. Sci. U.S.A. **103**(23), 8828–8833 (2006). [CrossRef] [PubMed]

**19. **C. M. Carpenter, S. Srinivasan, B. W. Pogue, and K. D. Paulsen, “Methodology development for three-dimensional MR-guided near infrared spectroscopy of breast tumors,” Opt. Express **16**(22), 17903–17914 (2008). [CrossRef] [PubMed]

**20. **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**(2), 024033 (2005). [CrossRef] [PubMed]

**21. **Q. Fang, S. A. Carp, J. Selb, G. Boverman, Q. Zhang, D. B. Kopans, R. H. Moore, E. L. Miller, D. H. Brooks, and D. A. Boas, “Combined optical imaging and mammography of the healthy breast: optical contrast derived from breast structure and compression,” IEEE Trans. Med. Imaging **28**(11 issue 1), 30–42 (2009). [CrossRef] [PubMed]

**22. **Q. Fang, J. Selb, and S. A. Carp, “G. Boverman G, E. L. Miller, D. H. Brooks, R. H. Moore, D. B. Kopans, D. A. Boas, “Combined optical and x-ray tomosynthesis breast imaging,” Radiology . in press.

**23. **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**(5), 051504 (2005). [CrossRef] [PubMed]

**24. **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**(25), 5181–5190 (2003). [CrossRef] [PubMed]

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

**26. **P. Hiltunen, S. J. D. Prince, and S. Arridge, “A combined reconstruction-classification method for diffuse optical tomography,” Phys. Med. Biol. **54**(21), 6457–6476 (2009). [CrossRef] [PubMed]

**27. **M. Guven, B. Yazici, X. Intes, and B. Chance, “Diffuse optical tomography with *a priori* anatomical information,” Phys. Med. Biol. **50**(12), 2837–2858 (2005). [CrossRef] [PubMed]

**28. **A. Gelman, J. B. Carlin, H. S. Stern, and D. B. Rubin, *Bayesian Data Analysis*, Boca Raton, FL, Chapman and Hall/CRC (2004)

**29. **E. Zastrow, S. K. Davis, M. Lazebnik, F. Kelcz, B. D. Van Veen, and S. C. Hagness, “Development of anatomically realistic numerical breast phantoms with accurate dielectric properties for modeling microwave interactions with the human breast,” IEEE Trans. Biomed. Eng. **55**(12), 2792–2800 (2008). [CrossRef] [PubMed]

**30. **K. D. Paulsen, P. M. Meaney, M. J. Moskowitz, and J. R. Sullivan Jr., “A dual mesh scheme for finite element based reconstruction algorithms,” IEEE Trans. Med. Imaging **14**(3), 504–514 (1995). [CrossRef] [PubMed]

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

**32. **P. K. Yalavarthy, D. R. Lynch, B. W. Pogue, H. Dehghani, and K. D. Paulsen, “Implementation of a computationally efficient least-squares algorithm for highly under-determined three-dimensional diffuse optical tomography problems,” Med. Phys. **35**(5), 1682–1697 (2008). [CrossRef] [PubMed]

**33. **F. Chung and K. Oden, “Weighted graph Laplacians and isoperimetric inequalities,” Pac. J. Math. **192**(2), 257–273 (2000). [CrossRef]

**34. **Q. Fang and D. Boas, “Tetrahedral mesh generation from volumetric binary and gray-scale images,” Proceedings of IEEE International Symposium on Biomedical Imaging2009, 1142–1145 (2009).

**35. **S. Srinivasan, B. W. Pogue, S. Jiang, H. Dehghani, and K. D. Paulsen, “Spectrally constrained chromophore and scattering near-infrared tomography provides quantitative and robust reconstruction,” Appl. Opt. **44**(10), 1858–1869 (2005). [CrossRef] [PubMed]

**36. **Q. Fang, P. M. Meaney, S. D. Geimer, A. V. Streltsov, and K. D. Paulsen, “Microwave image reconstruction from 3-D fields coupled to 2-D parameter estimation,” IEEE Trans. Med. Imaging **23**(4), 475–484 (2004). [CrossRef] [PubMed]

**37. **P. A. Yushkevich, J. Piven, H. C. Hazlett, R. G. Smith, S. Ho, J. C. Gee, and G. Gerig, “User-guided 3D active contour segmentation of anatomical structures: significantly improved efficiency and reliability,” Neuroimage **31**(3), 1116–1128 (2006). [CrossRef] [PubMed]