OSA's Digital Library

Biomedical Optics Express

Biomedical Optics Express

  • Editor: Joseph A. Izatt
  • Vol. 1, Iss. 5 — Dec. 1, 2010
  • pp: 1514–1531
« Show journal navigation

Singular value decomposition metrics show limitations of detector design in diffuse fluorescence tomography

Frederic Leblond, Kenneth M. Tichauer, and Brian W. Pogue  »View Author Affiliations


Biomedical Optics Express, Vol. 1, Issue 5, pp. 1514-1531 (2010)
http://dx.doi.org/10.1364/BOE.1.001514


View Full Text Article

Acrobat PDF (1854 KB)





Browse Journals / Lookup Meetings

Browse by Journal and Year


   


Lookup Conference Papers

Close Browse Journals / Lookup Meetings

Article Tools

Share
Citations

Abstract

The spatial resolution and recovered contrast of images reconstructed from diffuse fluorescence tomography data are limited by the high scattering properties of light propagation in biological tissue. As a result, the image reconstruction process can be exceedingly vulnerable to inaccurate prior knowledge of tissue optical properties and stochastic noise. In light of these limitations, the optimal source-detector geometry for a fluorescence tomography system is non-trivial, requiring analytical methods to guide design. Analysis of the singular value decomposition of the matrix to be inverted for image reconstruction is one potential approach, providing key quantitative metrics, such as singular image mode spatial resolution and singular data mode frequency as a function of singular mode. In the present study, these metrics are used to analyze the effects of different sources of noise and model errors as related to image quality in the form of spatial resolution and contrast recovery. The image quality is demonstrated to be inherently noise-limited even when detection geometries were increased in complexity to allow maximal tissue sampling, suggesting that detection noise characteristics outweigh detection geometry for achieving optimal reconstructions.

© 2010 OSA

1. Introduction

Whole-body small animal fluorescence tomography (FT) is a pre-clinical imaging technique relying on the acquisition of diffuse near-infrared (NIR) light signals transmitted through the body which are used for reconstruction within the framework of an inverse problem [1

1. V. Ntziachristos and R. Weissleder, “Charge-coupled-device based scanner for tomography of fluorescent near-infrared probes in turbid media,” Med. Phys. 29(5), 803–809 (2002). [CrossRef] [PubMed]

]. Three-dimensional (3D) images are reconstructed by minimizing the difference between acquired measurements and theoretical measurements defined by analytical or numerical solutions derived from a light transport forward model using boundary conditions consistent with the geometry of the interrogated specimen. Five important factors affect image resolution and recovered contrast of FT images: noise, forward model approximations, detection geometry, anatomical priors, and data types. Understanding the interplay and tradeoffs between these factors within the confines of currently available photodetection instrumentation is central to system optimization. This study looks at an approach to compare the effects of these factors upon the tomography image quality.

2. Background material

2.1 Imaging system and geometry

The mathematical analysis and simulations presented in this work are based on the detection geometry of a non-contact diffuse fluorescence tomography instrument, allows specimens to be interrogated with multiple laser projections [3

3. D. Kepshire, N. Mincu, M. Hutchins, J. Gruber, H. Dehghani, J. Hypnarowski, F. Leblond, M. Khayat, and B. W. Pogue, “A microcomputed tomography guided fluorescence tomography system for small animal molecular imaging,” Rev. Sci. Instrum. 80(4), 043701 (2009). [CrossRef] [PubMed]

]. In its current implementation, the system is co-registered with a microcomputed tomography (microCT) system allowing delineation of the specimen surface to be retrieved and used for light transport modeling using finite elements methods (FEM). Figure 1
Fig. 1 Optical detection geometries evaluated for diffuse fluorescence tomography simulations. In the limit that the number of detectors, Nd, is large and the angle between adjacent detection points, θ, is small, wide-field detection is achieved.
presents a schematic of the generalized detection geometry considered in this work. The detection of filtered fluorescence light is achieved in transmission using a number Nd of photodetectors, separated by an angle θ, and opposing the point of illumination. The laser and photodetectors are fixed on a gantry that can rotate 360° around the specimen with an angular separation providing a number Np of laser projections. In this theoretical study, fluorescence as well as transmission (excitation) measurements were simulated, allowing reconstructions to be performed based on fluorescence-to-transmission ratio data (Born normalization) [4

4. V. Ntziachristos and R. Weissleder, “Experimental three-dimensional fluorescence reconstruction of diffuse media by use of a normalized Born approximation,” Opt. Lett. 26(12), 893–895 (2001). [CrossRef] [PubMed]

] for a range of focused-detection data collection schemes (Table 1

Table 1. Technical specifications associated with the five representative imaging geometries for which the performance is compared. The last column in the table shows the numerical rank of the forward model fluorescence matrix associated with the geometries. The number of projections is in reference to the system described in Section 2.1, and corresponds to the number of laser positions used (manipulated by rotating the gantry). In other words, the number of measurements per 2D slice is equal to the number of projections times the number of detectors.

table-icon
View This Table
| View All Tables
). A wide-field geometry (i.e., CCD detection) was emulated by minimizing θ and maximizing Nd.

2.2 Forward modeling and image reconstruction

Diffuse fluorescence tomography relies on the resolution of an inverse problem of the form (see, e.g., Ref [5

5. F. Leblond, H. Dehghani, D. Kepshire, and B. W. Pogue, “Early-photon fluorescence tomography: spatial resolution improvements and noise stability considerations,” J. Opt. Soc. Am. A 26(6), 1444–1457 (2009). [CrossRef] [PubMed]

].),

Dρ=a=1NvAρaCa  ρ=1,2,,Nm,
(1)

where, for continuous-wave (CW) imaging, Dρ is a column vector with N m = N d × N p entries corresponding to the fluorescence-to-transmission ratio measurements obtained with the representative geometry shown in Fig. 1. The column vector C a (a = 1,2,…,N v) is composed of entries representing the concentration (e.g., in units of Molars or μg/ml) of fluorophores in each of the N v voxels contained in the discretized volume (e.g., an FEM mesh) associated with the interrogated specimen. Finally, the N m × N v matrix A is derived from the light transport model of choice and relates the CW measurements vector to the fluorophore concentration vector as shown in Eq. (1). If the diffusion approximation to the radiative-transport equation (RTE) is used, the forward model matrix takes the form,

Aρa(rs,rd;ra)=QFεFτΦx(rs,ra)Φe(ra,rd),
(2)

where Q F is the quantum efficiency of the fluorophore, ε F the extinction coefficient and τ the lifetime. The vectors r s, r d and r a (a = 1,2,…,N v) represent the location of a laser source point (labeled s) projected on the specimen, a detection point (labeled d) also projected on the specimen and the position of a voxel inside the interrogated medium, respectively. The functions Φx and Φe are fluence fields at the excitation and the fluorescence emission wavelengths, respectively. These fields are computed by solving the diffusion equation on a discretized mesh either numerically, or analytically in cases where simple imaging geometries are considered. Fluorescence tomography involves the retrieval of the vector Ca in Eq. (1) − a volumetric fluorescence image − through the resolution of an inverse problem. An analysis of the inverse problem for guiding source-detector geometry optimization is discussed in the following section.

2.3 Singular value decomposition analysis for imaging geometry optimization

Singular value analysis (SVA) is a tool that has been used in the past to derive mathematical criteria allowing relative a priori assessment of the performance of various imaging geometry setups in diffuse optical tomography (DOT) as well as in diffuse FT [6

6. J. P. Culver, V. Ntziachristos, M. J. Holboke, and A. G. Yodh, “Optimization of optode arrangements for diffuse optical tomography: a singular-value analysis,” Opt. Lett. 26(10), 701–703 (2001). [CrossRef] [PubMed]

,7

7. E. E. Graves, J. P. Culver, J. Ripoll, R. Weissleder, and V. Ntziachristos, “Singular-value analysis and optimization of experimental parameters in fluorescence molecular tomography,” J. Opt. Soc. Am. A 21(2), 231–241 (2004). [CrossRef] [PubMed]

]. Thus far, the tool has predominantly been limited to the study of the singular values. Here, the evaluation is extended to a more detailed study of the singular mode vectors and how their behavior can translate into tomography features such as spatial resolution.

In singular value analysis, the inverse problem in Eq. (1) is solved by performing singular value decomposition (SVD) of the matrix A into two sets of singular vectors and a finite number of singular values (Eq. (2).1) in Ref [8

8. P. C. Hansen, Rank-Deficient and Discrete Ill-Posed Problems (SIAM, 1998).

].),

A=UΣVT=i=1NvdisiciT,
(3)

where, for Nv > Nm, U = (c1, …,cNm) ∈ RNm x Nv and V = (d1, …,dNv) ∈ RNm x Nv are orthonormal matrices, and where the diagonal matrix Σ = diag(s1, …,sNm) has non-negative elements appearing in non-increasing order (s1s2≥…≥sm≥0). A least-squares solution to the problem then takes the form,

Ca=i=1NNrankFidiTDsi(ca)i,
(4)

where si are the singular values, di the singular data (SD) vectors and ci the singular image (SI) vectors. Each singular value is associated with one SD vector and one SI vector, all of which are labeled with the index i running from 1 to N, where N is the maximal number of modes included in a reconstruction. Equation (4) represents a truncated singular value decomposition (TSVD), for which NNrank. Nrank is the numerical rank of the matrix A, which is demarcated by the mode at which a significant and abrupt decrease in singular values is observed [8

8. P. C. Hansen, Rank-Deficient and Discrete Ill-Posed Problems (SIAM, 1998).

]. The singular value solution in Eq. (4) is roughly similar to the projection of a signal onto different spectral modes in Fourier analysis. According to this analogy, the vectors ci are the Fourier modes, the coefficients Fi are regularization filter factors allowing for the relative contribution of the spectral components to be weighted, and the coefficients,

Ki=diTDsi,
(5)

are factors weighting the contribution of each of the different spectral components according to how much overlap there is between the SD vectors and the tomography data vector D.

3. Methods

Evaluations of focused-detection and wide-field imaging geometries were performed on five potential imaging geometries described in Fig. 1 and Table 1. Geometries I to IV correspond to potential focused-detection applications, while geometry V corresponds to a wide-field application. As mentioned in Section 2.1, the focused-detection geometries were chosen as potential variants of an in-house developed non-contact diffuse FT instrument that allows specimens to be interrogated with multiple laser projections [2

2. D. S. Kepshire, S. L. Gibbs-Strauss, J. A. O’Hara, M. Hutchins, N. Mincu, F. Leblond, M. Khayat, H. Dehghani, S. Srinivasan, and B. W. Pogue, “Imaging of glioma tumor with endogenous fluorescence tomography,” J. Biomed. Opt. 14(3), 030501 (2009). [CrossRef] [PubMed]

,3

3. D. Kepshire, N. Mincu, M. Hutchins, J. Gruber, H. Dehghani, J. Hypnarowski, F. Leblond, M. Khayat, and B. W. Pogue, “A microcomputed tomography guided fluorescence tomography system for small animal molecular imaging,” Rev. Sci. Instrum. 80(4), 043701 (2009). [CrossRef] [PubMed]

]. A rigorous analysis of all given imaging geometries was preceded by and completed in the context of an investigation into the relative contributions of forward model mismatch and stochastic noise to perturbations in the raw data set (Section 4.1). This was followed by an investigation into the characteristics of both singular value (Section 4.2) and novel singular vector analyses (Section 4.3) for evaluation of imaging geometries. Finally, both the model mismatch/noise analysis and the singular value/vector analyses were evaluated in the context of image quality, specifically, spatial resolution and contrast recovery (Section 4.4). All analyses were carried out assuming known anatomical priors for FEM boundary conditions. The benefits of Born normalization, a commonly used technique for improving fluorescence data reconstructions by normalizing to transmittance [4

4. V. Ntziachristos and R. Weissleder, “Experimental three-dimensional fluorescence reconstruction of diffuse media by use of a normalized Born approximation,” Opt. Lett. 26(12), 893–895 (2001). [CrossRef] [PubMed]

], were also investigated for completeness. In the past, fluorescence-to-transmission data sets have predominantly been shown to improve reconstructions by partially negating the impact of experimental factors, such as out-of-focus data acquisition for non-contact instruments. However, the emphasis here is on the further benefits this normalization scheme provides in terms of reducing the burden on precise light transport modeling in FT, a phenomenon referred to as model mismatch throughout this work.

3.1 The difference between model mismatch errors and stochastic noise

Analysis of the effects of model mismatch and stochastic noise on raw data vectors was carried out for geometry IV applied to a mouse head mesh created from a microCT data set acquired in vivo (Fig. 3c
Fig. 3 (a) Simulated data vectors associated with the distribution of inclusions shown in (c) (right-most image) for which a 10:1 fluorescent contrast with respect to the background has been assigned. Two of the curves in (a) represent raw fluorescence data sets, with one curve being associated with an homogenous medium (10% stochastic noise added) and the other one being associated with an heterogeneous medium (no stochastic noise added) where different optical properties (absorption and reduced scattering coefficients) were associated with the different regions shown in (c) (left-most image). The third curved in (a) (normalized by multiplying with the homogenous transmittance data set) represents the fluorescence-to-transmission ratio data set. (b) Graphs showing the percentage difference between a fluorescence data set generated for a homogenous medium (0% noise) and the data sets shown in (a).
below). For fluorescence tomography simulations, two 4 mm-diameter inclusions were inserted with a 5 mm separation between their center-of-masses, as shown in Fig. 3c (right-most image). The fluorescence contrast between the inclusions and the homogeneously fluorescent background was set to 10:1. For all geometries, data vectors were simulated using two-dimensional diffusion theory using the simulation software NIRFAST [9

9. S. C. Davis, B. W. Pogue, H. Dehghani, and K. D. Paulsen, “Contrast-detail analysis characterizing diffuse optical fluorescence tomography image reconstruction,” J. Biomed. Opt. 10(5), 050501 (2005). [CrossRef] [PubMed]

], with an FEM mesh containing Nv = 1979 nodes. In this first study, five data vectors were simulated: three fluorescence data vectors, two of which were simulated assuming an optically homogeneous background with and without 10% stochastic (Gaussian) noise added and one assuming an optically heterogeneous background and 0% noise. The last two simulated data vectors were transmittance vectors, one assuming an optically heterogeneous medium and one assuming a homogeneous medium (noise-free). The left-most image in Fig. 3c illustrates the different anatomical regions imposed to provide heterogeneity (1: brain, 2: bone, 3: other soft tissue). Table 2

Table 2. Optical properties (absorption and reduced scattering coefficients) used for the light transport simulations presented in this work

table-icon
View This Table
| View All Tables
summarizes the optical properties (absorption and reduced scattering coefficients) that were assigned to each region in the scope of two different types of simulations: heterogeneous and homogeneous media.

3.2 Singular value decomposition for optimizing imaging system geometries

The theoretical formulation of singular value and singular vector analysis is introduced in Section 2.3. Typically, only the singular values are used as figures of merit for assessing imaging geometries [7

7. E. E. Graves, J. P. Culver, J. Ripoll, R. Weissleder, and V. Ntziachristos, “Singular-value analysis and optimization of experimental parameters in fluorescence molecular tomography,” J. Opt. Soc. Am. A 21(2), 231–241 (2004). [CrossRef] [PubMed]

,10

10. T. Lasser and V. Ntziachristos, “Optimization of 360° projection fluorescence molecular tomography,” Med. Image Anal. 11(4), 389–399 (2007). [CrossRef] [PubMed]

]. In this study, characteristics of singular values and both SI and SD vectors were investigated using simulated data sets from the homogeneous mouse head mesh described in Section 3.1 and for the imaging geometries summarized in Table 1. For comparison with other studies investigating singular values, the magnitude of the singular values was plotted against mode order and the number of available modes in the presence of stochastic noise and model-mismatch was evaluated.

To provide a more rigorous analysis of the imaging geometries, the shapes of both the SI and the SD vectors were also investigated as a function of mode order for all geometries. The spatial frequency of the SI vector modes typically increases with the order, i, of the singular order (Fig. 2
Fig. 2 (a) Singular image (SI) vectors associated with increasing singular order. Mode images are plotted on a log-scale to convey the idea that increasing the order of the mode typically leads to spatial frequency increases. (b) Cross-section of the modes shown in (a) along the dotted line.
, Fig. 5
Fig. 5 Figure showing the metric labeled mode resolution as a function of singular mode order. The resolution of each SI mode is defined here as the average distance (in mm) between two consecutive maxima along the dotted line shown in Fig. 2. Higher spatial frequency is therefore equivalent to smaller mode resolution. The red dots represent (for each geometry) the maximum mode order (minimum mode resolution) that can be used to reconstruct an image when 1% stochastic noise is added to the simulated data set (see Table 3). The dotted line should be used as a visual guide indicating the approximate location of the noise floor threshold for all geometries.
); whereas, the shape evolution of the SD vectors follows a more complex pattern with increasing mode order (Fig. 6
Fig. 6 (a) Representative singular data (SD) vectors associated with imaging geometry IV (see Table 1). Inspection of the modes of order i = 1, i = 10 and i = 20 show that the frequency (as a function of the source-detector pair index label measurement number) increases as a function of the singular order i. (b) Frequency of SD modes as a function of singular mode order for the five imaging geometries presented in Table 1. The curves demonstrate that the frequency, Eq. (6), tends to increase monotonically within a certain range, but that phases of sudden decrease always occur for large values of i. Red dots on the curves represent the order from which modes cannot be used for reconstruction because of noise propagation (see Table 3).
). The resolution (inverse of spatial frequency) of each individual SI vector mode was defined as the average distance, in millimeters, between two consecutive extrema (maxima or minima) along one-dimensional cross-sections of SI vector images. An example for this is the dotted line shown in Fig. 2a. More precisely, a numerical value has been assigned to each two-dimensional SI mode (Fig. 5) corresponding to the number of extrema, NexSI, along a one-dimensional cross-section divided by the length of the cross-section, L (in millimeters),

Mode resolution (mm) =(NexSIL)1.
(6)

The number of extrema was computed using an algorithm that automatically localizes and counts the total number of minima and maxima along cross-sections of the SI modes. As an example, Fig. 2b shows three such curves where the number of extrema are NexSI = 3, NexSI = 5 and NexSI = 9 for SI modes of singular order i = 1, i = 5 and i = 30, respectively. Although Fig. 5 has been plotted using the one-dimensional cross-section (dotted line) shown in Fig. 2, numerical analysis (not presented) was performed demonstrating that the final conclusions derived from the mode resolution analysis in Section 4 are unchanged for different cross-sections. For this metric, favorable detection geometries were determined as those for which the modes available for reconstruction had the smallest possible spatial resolution.

Each SD vector (e.g., see Fig. 6a) was also characterized by a numerical criterion, in this case labeled mode frequency (Fig. 6b). The frequency for each SD mode was computed by evaluating the number of extrema NexSD for each SD mode and dividing this number by the total number of measurements Nm,,

Mode frequency =NexSDNm.
(7)

The number of extrema, NexSD, was automatically computed using an algorithm similar to that described above for the SI vectors. While the SI mode resolution is an effective measure for the level of spatial resolution each image mode can contribute to a reconstructed FT image, the SD mode frequency does not lend itself to such a direct interpretation. However, as explained in Section 4, the numerical value of the SD mode frequency can be used as an indirect assessment of the propagation of stochastic noise and model mismatch error into FT images.

For the FT image reconstructions presented in this work, the number of usable modes per imaging geometry was determined by using the singular value decomposition to reconstruct images by sequentially including higher order modes in the presence of different levels of stochastic noise. This method, referred earlier, is defined as the TSVD [8

8. P. C. Hansen, Rank-Deficient and Discrete Ill-Posed Problems (SIAM, 1998).

]. The usable mode cutoff was empirically defined as the point when further inclusion of higher order modes did nothing to improve the image quality. The use of filter functions F corresponding to a Tikhonov regularization was also considered but did not provide significant image quality improvements.

Other aspects of the SVD analysis used in the a priori assessment of detector geometries included the steepness of the singular value vs. mode order and the number of singular modes above the noise floor (i.e., the modes available for image reconstruction).

3.3 Image quality: spatial resolution and recovered contrast

4. Results and discussion

4.1 The difference between model mismatch errors and stochastic noise

Figure 3a displays the fluorescence data vectors described in Section 3.1 as a function of measurement number for the homogeneous medium with 10% noise and for the heterogeneous medium with no noise added. The third curve displayed is the fluorescence-to-transmittance (Born ratio) for the heterogeneous medium. The Born ratio curve was normalized by multiplication with the simulated transmittance data set associated with the homogeneous optical properties in Table 2. Based on the similarities between the normalized data set and the data set of the homogeneous medium data set, the results suggest that if unknown heterogeneous optical properties are suspected, a homogenous forward model may still be applicable as long as Born normalization is applied. This is further demonstrated in Fig. 3b.

Figure 3b displays the percentage difference between the three vectors plotted in Fig. 2a and the simulated fluorescence data vector for a homogeneous background with no added noise (not shown in Fig. 3a) as a function of measurement number. The differences between the heterogeneous and homogeneous data sets are seen to be up to 15-20% and the measurement-to-measurement variation is significant. In addition to corroborating the observations made about Fig. 3a, this result more precisely suggests that reconstructing a heterogeneous data set using a homogeneous forward model will lead to imprecise reconstructions in terms of center-of-mass localization, as well as potentially decreasing spatial resolution and contrast recovery. Again, befitting the preliminary conclusion drawn from Fig. 3a, the minimal difference between the Born curve and the homogeneous curve suggests that reconstructions based on normalized data sets are less reliant on the veracity of the forward model, which often inaccurately represents the imaging object: in vivo values for different organs are typically not available and in cases where they would be available, significant inter-specimen variations are expected. Additionally, different systems use different excitation and emission wavelengths, preventing use of generic values in the scope of routine experiments, and accurate co-registration with anatomical imaging modalities is difficult but required.

Figure 3b also displays the differences between the noise and noiseless homogeneous vectors (labeled homogenous medium - 10% noise), illustrating the characteristics of noise propagation into the data sets. As expected, the difference between the two vectors consists of stochastic noise normally distributed around zero, with a maximum amplitude of 10% of the fluorescence measurement amplitude. This curve demonstrates that sources of noise such as photon shot noise will always contribute high-frequency components to a data vector. This is in opposition to the predominantly low frequency contributions associated with model mismatch. Therefore, the effect of both stochastic noise and model mismatch on fluorescence tomography can be empirically written as,

D=Dhomo+Dhigh-f+Dlow-f,
(8)

where Dhomo is the vector that would result for an homogenous medium while Dhigh-f and Dlow-f are vectors corresponding high-frequency (stochastic noise) and low-frequency (predominantly model mismatch) contributions. This parametric classification of data vector contributions is clearly not exact; however, it will be conceptually useful in the upcoming discussion relating to singular value decomposition of the matrix A. It should also be noted here that model mismatch contributions associated with discretization in Eq. (2) will mainly contribute high-frequency components to the data vector in situations where the size of the individual mesh components is comparable or smaller than typical light mean free paths. On the other hand, the fact that most of the low frequency contributions is associated with model mismatch errors is a direct result of the diffusive nature of the light transport sensitivity functions, which leads to averaging of signals emerging from large regions within the interrogated medium.

4.2 Singular value analysis of image geometries

Figure 4
Fig. 4 Singular values as a function of the singular mode order for the source-detection geometries presented in Table 1. The red dots represent (for each geometry) the minimum singular value (maximum mode order) that can be used to reconstruct an image when 1% stochastic noise is added to the simulated data set (see Table 3 for numerical values). The dotted line should be used as a visual guide indicating the approximate location of the noise floor threshold for all geometries.
displays the singular value vs. mode order curves associated with the five imaging geometries described in Table 1 and the red dots indicate the maximum usable mode order as described in Section 3.2. The maximum number of modes used for different reconstruction scenarios are tabulated in Table 3

Table 3. Tabulation of the number of singular modes that can be used to reconstruct images in situations where different levels of noise and model mismatch are added to simulated tomography data sets. In cases where no noise is added, all singular modes (100% of rank value) can be used to reconstruct an image. In situations where noise is added there is a dramatic decrease in the number of useful singular modes.

table-icon
View This Table
| View All Tables
. The truncation value corresponds to the largest singular mode order where stochastic noise does not corrupt the reconstructed images. The values used in Fig. 4 (red dots) correspond to the three columns associated with a 1% stochastic noise level. Detection geometries are usually deemed favorable when the slope of the singular value curve is as flat as possible and when the number of singular values or modes above a certain noise threshold is maximized [7

7. E. E. Graves, J. P. Culver, J. Ripoll, R. Weissleder, and V. Ntziachristos, “Singular-value analysis and optimization of experimental parameters in fluorescence molecular tomography,” J. Opt. Soc. Am. A 21(2), 231–241 (2004). [CrossRef] [PubMed]

]. In other words, conventional wisdom in the fluorescence tomography community is that geometries providing a large value of N in Eq. (4) are superior because more modes can be used to reconstruct images. Based on those criteria, inspection of Fig. 4 indicates that the wide-field detection geometry (V) is superior compared to the other geometries. It can also be deduced that the ranking of the most favorable geometries goes like: V, IV, II, III, I.

The maximum number of singular modes that can in principle be used to reconstruct a fluorescence image is set by the rank of the matrix A, i.e., N = N rank. However, in practice the number of modes that can actually be used to reconstruct an image is limited by, for example, the presence of high-frequency noise in the data vector D. The projection of stochastic noise contribution onto the data vector, D high-fdiT, typically causes the coefficients Ki to contribute disproportionally to reconstructed images for some value i < N rank, systematically leading to noisy reconstructed images. This effect might be perceived as a divergence in the singular coefficients caused by two competing effects: namely, monotonic decrease in singular values with increasing order and projection of high-frequency components of the data vector, e.g. stochastic noise or volume discretization noise, onto the SD vectors, which is disproportionately amplified for some values of i. For a better conditioned problem, sequentially smaller singular values (for i<Nrank) will not cause the Ki’s to diverge since the monotonic increase of 1/si is opposed by a continuously decreasing overlap of the SD and raw data vectors as i increases. However, in the presence of high-frequency noise, the large i behavior of the term D high-fdiT introduces divergences, which hide the contribution of the lower frequency components, i.e., the modes containing the physically relevant biological information. In a nutshell, regularization in optical tomography amounts to minimizing the contribution of those singular vectors contributing disproportionally to reconstructed images due to the propagation of high-frequency noise. When using singular value decomposition to solve the inverse problem, regularization can be achieved by appropriately choosing the filter function F i, or, by introducing a smoothing norm by using a generalized singular value decomposition (GSVD) [8

8. P. C. Hansen, Rank-Deficient and Discrete Ill-Posed Problems (SIAM, 1998).

].

Whereas singular value analysis can provide insight into the number of singular modes achievable for a given imaging geometry, it gives little insight into the shape of those modes, and therefore does not provide a full and fair analysis of the imaging geometry in question. In an attempt to utilize singular value decomposition to its full potential, the spatial frequency of the SI and SD vectors, novel metrics for analyzing imaging geometries, was investigated in the context of quantifying mode resolution characteristics for each of the described imaging geometries.

4.3 Singular vector analyses of image geometries

Figure 5 displays the mode resolution of the SI vectors (Eq. (6)) as a function of mode order for the five representative geometries described in Table 1. According to this figure of merit, favorable geometries will be those for which the modes used to reconstruct an image have the highest possible spatial frequency (smallest possible mode resolution). In other words, geometries with the steepest mode resolution versus mode order curves will be deemed favorable according to the SI mode resolution metric. Based on these criteria, the order of the most favorable metrics is (I, III, V, II, IV). Another criterion is the detection geometry-specific, minimum spatial resolution attainable (red dots in Fig. 5). The corresponding maximum number of modes used for different reconstruction scenarios are tabulated in Table 3. The values used in Fig. 5 correspond to the three columns associated with a 1% stochastic noise level. Using this last figure of merit and assuming equivalent levels of noise, geometry V (wide-field) seems to be favorable, although only by a small margin compared to the focused-based detection geometries. Perhaps a surprising observation derived from Fig. 5 is that the geometry with the smallest measurement density has resolution features comparable with geometries associated with a far larger number of measurements.

Figure 6b displays the frequency of SD vectors (Eq. (7)) as a function of mode order for individual imaging geometries. This metric was observed to generally increase monotonically with mode order, but was subject to interspersed dramatic decreases in resolution for some critical singular mode order values. This implies that low-frequency model mismatch errors, in addition to high-frequency stochastic noise errors, have to be taken into account at higher order modes. For a monotonically increasing SD mode frequency, it is expected that as the mode order increases, projection onto the high-frequency components in Eq. (8) will gradually make way to the projections onto the lower-frequency components. However, for geometries associated with interspersed, highly variable frequency-mode dependencies, this simple relationship between the contributions of high- and lower-frequency components is not expected to hold. An interesting observation here is that the truncation of useful modes for image reconstructions at 1% noise (red dots in Fig. 6b) always seems to occur shortly after (in terms of singular mode order) the critical points where low frequency noise starts to contribute to high order projection coefficients.

4.4 Image quality

Finally, Table 4

Table 4. Figures of merit used to assess a priori and a posteriori performance of the imaging geometries presented in Table 1. For each set of figures of merit, a ranking is provided for geometries I-V starting with the one for which a figure of merit predicts the optimal performance.

table-icon
View This Table
| View All Tables
shows a summary of the image geometry performance predictions associated with both the singular value curve and the SI mode resolution curve. The ranking that is provided in the last column of the table orders the geometries starting with the most favorable geometry as assessed with either singular values or SI vectors. The last line of the table provides a ranking where image quality was assessed by evaluating the contrast recovery and the spatial resolution of the images presented in Fig. 7. The important thing to note here is that the ranking associated with either of the singular value analysis metrics (singular values or SI vectors) alone does not entirely agree with the a posteriori evaluation of the images based on image quality. For example, based on the singular value curve, geometry II is preferable to geometry III. However, based on the spatial resolution of the modes, geometry III is expected to perform better than geometry II. Similarly, based on singular values, the high density of measurements geometry IV is expected to be superior to geometry III, which is in disagreement with the prediction derived from the spatial resolution of the SI modes. Again, the prediction from the later criteria is more in tune with the ranking based on image quality.

4.5 Summary and conclusions

In this work, singular value decomposition analysis of image reconstruction matrices for fluorescence tomography was investigated in terms of its potential to predict the performance of various imaging geometries (corresponding to focused- and wide-field detection), quantified by reconstructed image quality. In the past, analysis of the number of usable singular values has been used to pre-evaluate the performance of proposed imaging geometries, generally predicting that, for larger fields of view, more source-detector projections are always preferable [7

7. E. E. Graves, J. P. Culver, J. Ripoll, R. Weissleder, and V. Ntziachristos, “Singular-value analysis and optimization of experimental parameters in fluorescence molecular tomography,” J. Opt. Soc. Am. A 21(2), 231–241 (2004). [CrossRef] [PubMed]

]. A singular value analysis was repeated in the present study for completeness, corroborating this finding; however, an analysis of the image quality improvements corresponding to increased source-detector projections demonstrated only very marginal improvements in the presence of realistic levels of noise and model-data mismatch. In fact, inspection of Fig. 7 provides support to the conservative claim that focused-based detection can achieve image quality close to what can be achieved with wide-field detection. To further support this claim, Fig. 9
Fig. 9 Different simulation scenarios where FT is used to monitor the progress of a tumor with center-of-mass the same as top inclusion shown on the fluorescence image in Fig. 3c. Reconstruction results are presented as a cross-section of the corresponding tomography image through a horizontal line across the center-of-mass of the tumor. Scenario A shows a tumor of constant contrast-to-background ratio (10:1) but with decreasing size, Scenario B shows a tumors of decreasing contrast but with constant size, Scenario C shows tumor stages where both size and contrast vary. For each scenario, four tumor stages are shown and each image is reconstructed with three different imaging configurations, namely geometries I, II, V (Table 1).
presents simulated imaging scenarios illustrating how different imaging geometries perform when it comes to monitoring the evolution of a hypothetical fluorescent tumor. Three different evolution scenarios are presented where different tumor parameters are varied at four different stages: in Scenario A the contrast-to-background ratio remains constant (10:1) but the size of the tumor varies, in Scenario B the size of the inclusion is constant but the contrast varies, and, in Scenario C both the size and the contrast of the tumor vary. The tomography results for each scenario are presented in Fig. 9 where images were reconstructed using an homogenous matrix A but with fluorescence-to-transmission data generated with the heterogeneous optical properties found in Table 2. For all the reconstructions, 1% stochastic noise has been added to the data vector. The results presented in Fig. 9 show that under the same noise conditions, a wide-field detection geometry provides only minor imaging gains when it comes to spatial resolution and recovered contrast.

It should be mentioned here that the noise characteristics of wide-field detection, which are usually based on cooled CCD chips [11

11. F. Leblond, S. C. Davis, P. A. Valdés, and B. W. Pogue, “Pre-clinical whole-body fluorescence imaging: Review of instruments, methods and applications,” J. Photochem. Photobiol. B 98(1), 77–94 (2010). [CrossRef] [PubMed]

,12

12. B. W. Pogue, F. Leblond, V. Krishnaswamy, and K. D. Paulsen, “Radiologic and near-infrared/optical spectroscopic imaging: where is the synergy?” AJR Am. J. Roentgenol. 195(2), 321–332 (2010). [CrossRef] [PubMed]

] are expected to be different than focused-detection geometries, typically based on PMTs, when it comes to the detection of low photon fluxes. For high fluorescence signals, the noise is expected in both cases to be shot-noise limited and to some extent the levels of noise could be equivalent. However, when it comes to the detection of very low photon fluxes, non-contact CCD-based systems will be limited by the read-out noise of the chip. This is particularly important given the limited dynamical range associated with non-contact CCD-based small animal imaging. In terms of noise characteristics for low photon fluence detection, the optimal systems are time-domain PMT systems using time correlated single photon counting [3

3. D. Kepshire, N. Mincu, M. Hutchins, J. Gruber, H. Dehghani, J. Hypnarowski, F. Leblond, M. Khayat, and B. W. Pogue, “A microcomputed tomography guided fluorescence tomography system for small animal molecular imaging,” Rev. Sci. Instrum. 80(4), 043701 (2009). [CrossRef] [PubMed]

]. In this case, detection is achieved in such a way that the signal-to-noise levels achievable for low photon fluxes encountered in fluorescence tomography applications can be orders of magnitude larger than what can be achieved with cooled CCD’s. In conclusion, when comparing wide-field detection with focused-based detection it might be more realistic to compare wide-field reconstructions at 1% noise with focused-based detection with 0.1% stochastic noise added.

Even assuming equal noise considerations for all imaging geometries, the geometry evaluation ranking of the singular values did not match the image quality rankings (Table 4). These findings suggest that the number of usable singular values alone cannot fully characterize the performance of a given imaging geometry for realistic data sets. In response, we developed new metrics, mode resolution and mode frequency, to quantify the characteristics of the singular image and singular data vectors, respectively, corresponding to each singular value. Like the singular value analysis, alone, these singular vector metrics did not adequately reflect the image quality ranking; however, their ranking could be used to account for the discrepancies between the singular value ranking and the image quality ranking, suggesting that together, singular value and singular vector metrics are complementary and could conceivably be used to fully characterize the efficacy of proposed imaging geometries. Future research efforts will be devoted to findings methodologies to better characterize the relative importance and interplay between the different singular value decomposition metrics.

4.6 Limitations and future work

Singular value analysis has traditionally been used in optical tomography as an analysis tool for imaging geometry optimization. This method is rarely used to actually reconstruct diffuse fluorescence images in the scope of biological experiments. One of the reasons for this is that singular value decomposition is a direct matrix inversion method that requires significant computer resources, which become particularly time intensive for reconstructing 3D volumes in whole body, small animal imaging or for larger volumes such as breast imaging. Though the current study focused on 2D reconstructions of a mouse head, the conclusion that increasing source-detector densities has diminishing gains in terms of quality of image reconstruction is expected to hold for 3D and larger volume reconstructions although further studies are warranted. Another reason why SVD is not commonly used for image reconstruction is that the approach is not flexible in terms of allowing for the inclusion of a generic penalty function when resolving the inverse problem. Although a conventional Tikhonov penalty term can be introduced using a specific filter function F i in Eq. (4), an analysis of the problem based on a generalized SVD (GSVD) must be performed in order to allow penalty functions that would be consistent with, e.g., soft as well hard spatial priors. While the results presented in this paper provide useful information on the FT problem, they should be interpreted within the confines of those limitations. In fact, it is more commonplace in whole-body fluorescence tomography to use inverse problem resolution methods with regularization functions allowing certain aspects of the images to be favored as determined by prior knowledge (e.g., spectral and/or spatial) or pre-determined features that are expected for the images (e.g., edge preservation, absence of high spatial frequencies).

Finally, the analysis presented in this work was concerned with how different imaging configurations can affect FT results assuming that the data were from continuous-wave measurement systems. A conclusion that is reached here is that the relative gains, in terms of spatial resolution and recovered contrast, associated with using a high density of measurements, is not as impressive as one might expect based upon linear system performance. In fact, the results presented in Fig. 7 show only modest image quality gains associated with a system where the number of measurements was increased up to 10 - 50 times. Therefore, future developments in FT should be concerned with investigating imaging schemes that would allow more significant improvements to be achieved within the confine of the intrinsic limitations posed by light scattering. One way of achieving this would be to use the singular value analysis methods presented in this work to investigate how different non-redundant data types could decrease the impact of the ill-posedness of the problem. Approaches using tomography data sets derived from multi-spectral or time-domain systems constitute promising research directions.

Acknowledgments

This work has been funded by the National Institutes of Health (NIH) grants RO1CA120368 and K25 CA138578 through the National Cancer Institute (NCI).

References and links

1.

V. Ntziachristos and R. Weissleder, “Charge-coupled-device based scanner for tomography of fluorescent near-infrared probes in turbid media,” Med. Phys. 29(5), 803–809 (2002). [CrossRef] [PubMed]

2.

D. S. Kepshire, S. L. Gibbs-Strauss, J. A. O’Hara, M. Hutchins, N. Mincu, F. Leblond, M. Khayat, H. Dehghani, S. Srinivasan, and B. W. Pogue, “Imaging of glioma tumor with endogenous fluorescence tomography,” J. Biomed. Opt. 14(3), 030501 (2009). [CrossRef] [PubMed]

3.

D. Kepshire, N. Mincu, M. Hutchins, J. Gruber, H. Dehghani, J. Hypnarowski, F. Leblond, M. Khayat, and B. W. Pogue, “A microcomputed tomography guided fluorescence tomography system for small animal molecular imaging,” Rev. Sci. Instrum. 80(4), 043701 (2009). [CrossRef] [PubMed]

4.

V. Ntziachristos and R. Weissleder, “Experimental three-dimensional fluorescence reconstruction of diffuse media by use of a normalized Born approximation,” Opt. Lett. 26(12), 893–895 (2001). [CrossRef] [PubMed]

5.

F. Leblond, H. Dehghani, D. Kepshire, and B. W. Pogue, “Early-photon fluorescence tomography: spatial resolution improvements and noise stability considerations,” J. Opt. Soc. Am. A 26(6), 1444–1457 (2009). [CrossRef] [PubMed]

6.

J. P. Culver, V. Ntziachristos, M. J. Holboke, and A. G. Yodh, “Optimization of optode arrangements for diffuse optical tomography: a singular-value analysis,” Opt. Lett. 26(10), 701–703 (2001). [CrossRef] [PubMed]

7.

E. E. Graves, J. P. Culver, J. Ripoll, R. Weissleder, and V. Ntziachristos, “Singular-value analysis and optimization of experimental parameters in fluorescence molecular tomography,” J. Opt. Soc. Am. A 21(2), 231–241 (2004). [CrossRef] [PubMed]

8.

P. C. Hansen, Rank-Deficient and Discrete Ill-Posed Problems (SIAM, 1998).

9.

S. C. Davis, B. W. Pogue, H. Dehghani, and K. D. Paulsen, “Contrast-detail analysis characterizing diffuse optical fluorescence tomography image reconstruction,” J. Biomed. Opt. 10(5), 050501 (2005). [CrossRef] [PubMed]

10.

T. Lasser and V. Ntziachristos, “Optimization of 360° projection fluorescence molecular tomography,” Med. Image Anal. 11(4), 389–399 (2007). [CrossRef] [PubMed]

11.

F. Leblond, S. C. Davis, P. A. Valdés, and B. W. Pogue, “Pre-clinical whole-body fluorescence imaging: Review of instruments, methods and applications,” J. Photochem. Photobiol. B 98(1), 77–94 (2010). [CrossRef] [PubMed]

12.

B. W. Pogue, F. Leblond, V. Krishnaswamy, and K. D. Paulsen, “Radiologic and near-infrared/optical spectroscopic imaging: where is the synergy?” AJR Am. J. Roentgenol. 195(2), 321–332 (2010). [CrossRef] [PubMed]

OCIS Codes
(170.0110) Medical optics and biotechnology : Imaging systems
(170.3010) Medical optics and biotechnology : Image reconstruction techniques
(170.3660) Medical optics and biotechnology : Light propagation in tissues

ToC Category:
Image Reconstruction and Inverse Problems

History
Original Manuscript: September 7, 2010
Revised Manuscript: November 9, 2010
Manuscript Accepted: November 20, 2010
Published: November 29, 2010

Citation
Frederic Leblond, Kenneth M. Tichauer, and Brian W. Pogue, "Singular value decomposition metrics show limitations of detector design in diffuse fluorescence tomography," Biomed. Opt. Express 1, 1514-1531 (2010)
http://www.opticsinfobase.org/boe/abstract.cfm?URI=boe-1-5-1514


Sort:  Author  |  Year  |  Journal  |  Reset  

References

  1. V. Ntziachristos and R. Weissleder, “Charge-coupled-device based scanner for tomography of fluorescent near-infrared probes in turbid media,” Med. Phys. 29(5), 803–809 (2002). [CrossRef] [PubMed]
  2. D. S. Kepshire, S. L. Gibbs-Strauss, J. A. O’Hara, M. Hutchins, N. Mincu, F. Leblond, M. Khayat, H. Dehghani, S. Srinivasan, and B. W. Pogue, “Imaging of glioma tumor with endogenous fluorescence tomography,” J. Biomed. Opt. 14(3), 030501 (2009). [CrossRef] [PubMed]
  3. D. Kepshire, N. Mincu, M. Hutchins, J. Gruber, H. Dehghani, J. Hypnarowski, F. Leblond, M. Khayat, and B. W. Pogue, “A microcomputed tomography guided fluorescence tomography system for small animal molecular imaging,” Rev. Sci. Instrum. 80(4), 043701 (2009). [CrossRef] [PubMed]
  4. V. Ntziachristos and R. Weissleder, “Experimental three-dimensional fluorescence reconstruction of diffuse media by use of a normalized Born approximation,” Opt. Lett. 26(12), 893–895 (2001). [CrossRef] [PubMed]
  5. F. Leblond, H. Dehghani, D. Kepshire, and B. W. Pogue, “Early-photon fluorescence tomography: spatial resolution improvements and noise stability considerations,” J. Opt. Soc. Am. A 26(6), 1444–1457 (2009). [CrossRef] [PubMed]
  6. J. P. Culver, V. Ntziachristos, M. J. Holboke, and A. G. Yodh, “Optimization of optode arrangements for diffuse optical tomography: a singular-value analysis,” Opt. Lett. 26(10), 701–703 (2001). [CrossRef] [PubMed]
  7. E. E. Graves, J. P. Culver, J. Ripoll, R. Weissleder, and V. Ntziachristos, “Singular-value analysis and optimization of experimental parameters in fluorescence molecular tomography,” J. Opt. Soc. Am. A 21(2), 231–241 (2004). [CrossRef] [PubMed]
  8. P. C. Hansen, Rank-Deficient and Discrete Ill-Posed Problems (SIAM, 1998).
  9. S. C. Davis, B. W. Pogue, H. Dehghani, and K. D. Paulsen, “Contrast-detail analysis characterizing diffuse optical fluorescence tomography image reconstruction,” J. Biomed. Opt. 10(5), 050501 (2005). [CrossRef] [PubMed]
  10. T. Lasser and V. Ntziachristos, “Optimization of 360° projection fluorescence molecular tomography,” Med. Image Anal. 11(4), 389–399 (2007). [CrossRef] [PubMed]
  11. F. Leblond, S. C. Davis, P. A. Valdés, and B. W. Pogue, “Pre-clinical whole-body fluorescence imaging: Review of instruments, methods and applications,” J. Photochem. Photobiol. B 98(1), 77–94 (2010). [CrossRef] [PubMed]
  12. B. W. Pogue, F. Leblond, V. Krishnaswamy, and K. D. Paulsen, “Radiologic and near-infrared/optical spectroscopic imaging: where is the synergy?” AJR Am. J. Roentgenol. 195(2), 321–332 (2010). [CrossRef] [PubMed]

Cited By

Alert me when this paper is cited

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


« Previous Article

OSA is a member of CrossRef.

CrossCheck Deposited