Tomographic fluorescence lifetime multiplexing in the spatial frequency domain

The ability to simultaneously recover multiple fluorophores within biological tissue (multiplexing) can have important applications for tracking parallel disease processes in vivo. Here we present a novel method for rapid and quantitative multiplexing within a scattering medium, such as biological tissue, based on fluorescence lifetime contrast. This method employs a tomographic inversion of the asymptotic (late) portion of time-resolved spatial frequency (SF) domain measurements. Using Monte Carlo simulations and phantom experiments, we show that the SF-asymptotic time domain (SF-ATD) approach provides a several-fold improvement in relative quantitation and localization accuracy over conventional SF-time domain inversion. We also show that the SF-ATD approach can exploit selective filtering of high spatial frequencies to dramatically improve reconstruction accuracy for fluorophores with subnanosecond lifetimes, which is typical of most near-infrared fluorophores. These results suggest that the SF-ATD approach will serve as a powerful new tool for whole-body lifetime multiplexing. measurements, such as rapid whole-body imaging and high resolution, with quantitative multiplexing capabilities using TD measurements. We validate SF-ATD using Monte Carlo (MC) simulations with MCX, a GPU-accelerated high-speed MC computing platform [7] that can accept arbitrary spatial have presented a novel approach for tomographic multiplexing using time-resolved SF domain measurements that provide high accuracy for quantitative tomography of multiple fluorophores simultaneously present in a turbid medium. We have also shown that SF filtering allows quantitative tomographic imaging of shorter lifetimes than possible with point measurements. This result has high importance for whole-body molecular imaging that the lifetimes of most near-infrared fluorophores are in the subnanosecond range Further studies will focus on in vivo applications and on optimizing experimental parameters, including time gates and spatial frequencies, to further improve the imaging performance of SF-ATD.

patterns as input. We show that SF-ATD provides a dramatic improvement in localization and relative quantitation accuracy compared to conventional direct inversion of SF-TD fluorescence data. We also show that by allowing selective filtering of high spatial frequency data, SF-ATD can quantitatively recover shorter lifetimes than possible using point-scanning methods. We validate the feasibility and advantages of SF-ATD over SF-DTD using experimental measurements with tissue-mimicking phantoms.
Consider a turbid medium obeying the radiative transport equation with absorption μ a x, m and scattering μ s x, m at the excitation (λ x ) and emission (λ m ) wavelengths, and with N embedded fluorophores, each having a lifetime τ n = 1/Γ n and yield distribution η n (r) (product of the quantum yield and extinction coefficient). It is straightforward to show (Supplement 1) that when τ n > τ a , where τ a = (νμ a x, m (r)) −1 , the TD fluorescence signal in the asymptotic region, defined as t > τ a , is factorized into functions of spatial frequency and time as where k s and k d are the source and detector spatial frequencies (SF) on the boundary of the medium, respectively, and a n are time-independent decay amplitudes, which are related to corresponding yield distributions a n (k s , k d ) = ∫ d 3 rW n (k s , k d , r)η n (r), (2) where W n = G n x (k s , r)G n m (r, k d ) is a continuous wave (CW) fluorescence weight matrix evaluated as the product of transport Green's functions, G n x, m , computed at a reduced absorption μ a x, m (r) − Γ n ν. The reduced absorption is crucial for the accurate relative quantitation of multiple lifetime components. Equations (1)-(2) compactly and rigorously describe asymptotic TD fluorescence and are valid for arbitrary complex-shaped heterogeneous transport media, provided τ n > τ a (≈0.45 ns for μ a x, m = 0.1/cm). We will show that the use of high spatial frequencies extends the applicability of Eq. (1) to even shorter lifetimes by effectively lowering τ a .
Equation (2) separates the TD fluorescence from a turbid medium containing multiple fluorophores into independent forward problems for each lifetime. To see the advantage of this separation, consider the inverse problem of recovering η n from TD measurements U F . Discretizing the problem into V medium voxels, L time points, and M pairs of sources and detector SFs at the medium boundaries, the full TD forward problem and its asymptotic limit [Eq. (1)] can be expressed as a matrix equation: Here, y is a (ML × 1) data vector, W is the full (ML × NV) TD weight matrix, η = [η 1 , η 2 , …, η N ] T is a (NV × 1) vector of unknown fluorescence yields for both lifetimes, A = [exp(−t/τ 1 ) ⊗ I, exp(−t/τ 2 ) ⊗ I, …, exp(−t/τ N ) ⊗ I] is a (ML × MN), over-determined "basis" matrix containing Kronecker products of exponential decay terms and the (M × M) identity matrix, I, and W = diag(W 1 , W 2 , …, W N ) is a (NM × NV) block diagonal matrix containing reduced absorption CW weight matrices for each lifetime. The inversion of Eq.
(3) to recover η can now proceed along two alternate ways. In the first way, called the SFasymptotic TD (SF-ATD) method, we first invert the well-conditioned matrix A using its Moore-Penrose pseudo inverse A † , followed by a Tikhonov inversion of the matrix W [6].
The inverse problem takes the form where a = A † y = [a 1 , …, a N ] T is a (MN × 1) vector of decay amplitudes for the N lifetimes and M pairs of source-detector spatial frequencies. The conventional approach for spatial frequency domain time-resolved fluorescence tomography [3], which we call the spatial frequency-direct TD (SF-DTD) approach, does not exploit the two-step nature of the inverse problem but rather performs a direct inversion of the full TD weight matrix W. It can be shown that in the asymptotic limit, we can still exploit the factorization, W = AW, implicit in Eq. (3), to arrive at the following expression for the SF-DTD inverse problem [6]: where C a = (A T A σ y 2 ) −1 , C y = σ y 2 I is the data covariance matrix, and λ = (σ y /σ η ) 2 , with σ y and σ η as the variances in the measurement and yield. Note that the size of the matrix to be inverted in Eq. (5) is NM × NM (for WW T ), which is much smaller than the matrix for full TD inversion (LM × LM for WW T ), since typically N ≪ L. Thus, Eq. (5) offers a dramatic reduction in the computation of SF-DTD, which can otherwise be prohibitive for more than a few gates. The key distinction between SF-ATD [Eq. (4)] and SF-DTD [Eq. (5)] is that the inverse problem in Eq. (4) is block diagonal, implying zero cross talk between multiple lifetimes, whereas Eq. (5) is not block diagonal. We will see that this difference results in significant performance difference between the two methods.
We first demonstrate the advantages of SF-ATD [Eq. (4)] over SF-DTD [Eq. (5)] using Monte Carlo (MC) simulations, performed using "MCX," a GPU-based accelerated MC computing platform [7]. We consider a slab medium (6 cm × 4 cm × 1.9 cm, 1 mm 3 voxels) with absorption μ a x = μ a m = 0.1 cm −1 , scattering μ s x = μ s m = 10 cm −1 , anisotropy of 0.01, and refractive index of 1.37. Two 1 mm 3 fluorescent inclusions of equal concentration are placed at a height of z = 13 mm, with varying separations of 10 mm, 2 mm, and 0 mm. We first consider a case where the inclusions have lifetimes of 0.8 ns and 1.2 ns, which are longer than the intrinsic absorption timescale [(νμ a ) −1 ≈ 0.45 ns]. The sources were 2D sinusoidal patterns of the form S n = 0.5[1 + cos(2π(k x x + k y y) + ϕ n ))], where the three phases ϕ n = 2nπ/3 (n = 0, 1, 2) are used to eliminate the D.C. background [1,8]. k x and k y ranged from 0 to 0.5/cm in six steps, resulting in 36 sources at z = 0. The detectors consisted of the same 36 patterns at z = 1.9 cm. The source and detector MC TD Green's functions G n x, m were computed for 2 × 10 9 photons, and subsequently, the CW weight functions W 1 and W 2 were calculated for the 36 source and detector patterns and 64 time points between 0 and 6.3 ns (resulting in 1296 × 64 measurements). Separately, the decay amplitudes for the two lifetimes, a 1 and a 2 were recovered from the decay portion of the TD signal for all M measurements, using a linear bi-exponential fit with the known lifetimes. The linear fit employed 12 gates ranging from 1.9 ns to 6.3 ns (use of more gates did not significantly improve the quality of the fits). The a n s were finally used in Eqs. (4) and (5) to recover the yield distributions for the two lifetime components using the SF-ATD and SF-DTD methods.
The regularization λ was chosen to provide the least reconstruction error. The optimal regularization (least error) for SF-ATD was higher than that for SF-DTD, resulting in broader distributions (low precision) with highly accurate centroids and yield values (high accuracy) compared to SF-DTD. We note that 2 mm separation is shorter than previously used separations for fluorescence tomography. Also note that the simulations do not consider possible model errors, such as incorrect estimation of optical properties and medium geometry, which affect SF-DTD more than SF-ATD, as will be seen with experimental data.
While SF-ATD thus provides superior performance for lifetimes longer than the absorption timescale, τ a , a more interesting scenario is when lifetimes are comparable to, or even shorter than, τ a . In this case, the late time evolution of the TD fluorescence is dominated by intrinsic photon scattering timescales and no longer reflects the true lifetime of the fluorophore. Thus, the benefits of the asymptotic factorization [Eq. (3)] can no longer be exploited for quantitative multiplexing when using all SFs (or equivalently, point illumination and detection). However, the use of high SFs can alleviate the influence of scattering and allow detection of short lifetimes in turbid media. This can be understood based on an increase in the "effective" absorption of the medium with increasing SF as μ a + k 2 D, where k is the spatial frequency and D is the diffusion coefficient [8]. The increased μ a implies that the highly scattered, late-arriving photons are attenuated more, thereby decreasing τ a and resulting in a faster approach toward the asymptotic limit and also lowering the limit for the shortest lifetimes that can be detected in the asymptotic potion. This phenomenon is illustrated in Fig. 2, which shows the error in lifetime recovered from single exponential fits to the asymptotic TD signal for a single inclusion with a lifetime of either 0.3 ns [ Fig. 2(a)] or 0.5 ns [ Fig. 2(b)] for the same geometry as in Fig. 1. It is clear that for low frequencies, the error in lifetime is significant for the 0.3 ns case but decreases toward higher SF. At very high SFs, the lifetime error increases again due to a decreased signal-to-noise ratio (SNR) with increasing SF, which is in turn due to the medium acting as a low pass filter. high SF data improves quantitation, even for lifetimes comparable to or slightly longer than τ a (Supplement 1). SF-ATD therefore enables quantitative multiplexing using a wider range of near-infrared fluorophores, which typically have short (subnanosecond) lifetimes.
We evaluated the feasibility of SF-ATD using experiments with a tissue phantom. Two parallel tubes (0.965 mm outer diameter, 0.58 mm inner diameter) separated by 4 mm were embedded at a height of 12 mm in a circular petri dish (8.8 cm diameter) filled with an Intralipid+nigrosin mixture (μ a = 0.3 cm −1 and μ s′ = 10 cm −1 ) to a height of 1.7 cm [ Fig.   4(a)]. The tubes contained IRdye800CW (Licor biosciences, τ 1 = 0.45 ns) and 3, 3' diethylthiatricarbocyanine iodide (τ 2 = 0.65 ns). The sample was excited at 745 nm (Ti:sapphire, Spectra Physics MaiTai, 150 fs, 80 MHz) and the transmitted TD fluorescence was detected with a λ > 800 nm filter (Chroma) attached to a intensified CCD camera (LaVision Picostar HRI, 400 ps gatewidth, 500 ms exposure, 600 V gain). To generate spatial patterns, we modified a Pico projector (Aaxa P300) with a digital micromirror device (DMD) to accept collimated fiber input from the Ti:sapphire laser. Sinusoidal patterns at 19 frequencies ranging from k x = 0 to 0.51 cm −1 along the axis perpendicular to the tubes (x) were projected under the dish. Figure 4(a) shows the lifetime maps of the tubes without Intralipid for k x = 0. Figure 4(b) shows examples of the input patterns (integrated over all time gates) for zero (top) and the highest SF (0.51 cm −1 ) used, measured directly with a white paper placed on the imaging plate. These patterns were used directly as the input source patterns in MCX (rather than relying on perfect sinusoids), thereby implicitly accounting for experimental factors, such as variation in the response across the DMD, laser illumination, and the ICCD sensitivity. For simplicity of analysis, we used point detectors, assigned as 26 pixels (2 × 2 binning) on the CCD images with a 2 mm spacing. The decay amplitudes, a 1 and a 2 , for the 0.45 ns and 0.65 ns components were recovered for each measurement from linear bi-exponential fits [ Fig. 4(c)] to the asymptotic portion of the TD data using 9 time gates.
The SF-DTD [ Fig. 4 The true ratio of the fluorescence yields was η 1 :η 2 = 1.36, as estimated directly from the CW fluorescence of the dye-filled tubes without Intralipid [as in Fig. 4(a)]. While SF-ATD provided 3.5% error in relative quantitation and perfect lateral localization, the depth localization of η 2 had an error of 1 mm. SF-DTD resulted in 44% quantitation error and was unable to delineate the lifetimes as spatially separate. The poorer localization of SF-DTD with experimental data compared to the simulations could possibly be attributed to model errors not accounted for in our reconstructions, including errors in optical properties and positioning, incorrect noise model, or poorer SNR of the experimental data. The SF-ATD is, however, more robust to these unknown parameters and provides accurate relative quantitation and spatial localization under the same conditions.
We have presented a novel approach for tomographic multiplexing using time-resolved SF domain measurements that provide high accuracy for quantitative tomography of multiple fluorophores simultaneously present in a turbid medium. We have also shown that SF filtering allows quantitative tomographic imaging of shorter lifetimes than possible with point measurements. This result has high importance for whole-body molecular imaging given that the lifetimes of most near-infrared fluorophores are in the subnanosecond range [9]. Further studies will focus on in vivo applications and on optimizing experimental parameters, including time gates and spatial frequencies, to further improve the imaging performance of SF-ATD. Dependence of fluorescence lifetimes on source and detector spatial frequencies, simulated for a 1.9 cm thick slab with an inclusion of lifetime 0.3 ns (left) or 0.5 ns (right) placed at a height of 1.3 cm. The images show the error, |τ fit -τ actual |/τ actual × 100, of the lifetimes recovered from a single exponential fit to the asymptotic data.