2.1.

Instrument Design

A diagram and photograph of our instrument are shown in Figs. 1(a) and 1(b), respectively. We used a pulsed super-continuum laser (NKT Photonics) with 80-MHz repetition rate, fitted with an infrared beam-splitting unit, so that the emitted light was broadband in the visible and NIR range from $\sim 550$ to 850 nm. Depending on the measurement, one of two optical paths was used to illuminate the sample. In the first path, fluorescence excitation light was generated by filtering the laser output with a 680-nm shortpass and $655/50\text{-}\mathrm{nm}$ bandpass filter (Chroma Technology Corp), as indicated by the red line in Fig. 1(a). In the second path, the unfiltered broadband NIR light was used to measure wavelength-dependent optical properties of the samples. The path could be selected with a set of broadband laser mirrors (Edmund Optics) on flip-mounts.

Fig. 1

(a) Diagram and (b) photograph of the time-resolved multiplexed FMT instrument. See text for component details.

We used dual motorized rotation stages similar to the configuration in our previously reported time-resolved scanner.²⁵ Samples were mounted on an inner rotation stage, and two detector fibers were mounted on an outer rotation stage, so that an arbitrary number of scanning angles could be generated. Light escaping the surface of the sample was first filtered with 700-nm longpass filters (Chroma) and collected by two 1-mm multimode optical fibers that were coupled into Spec-PML-16C detector units (Becker & Hickl, Berlin, Germany). Each Spec-PML-16C is composed of a multiwavelength detection assembly combining a PML-16C (16-channel multianode PMT) and a $1/8\mathrm{m}$ polychromator. In combination, the PMT-array and grating were configured to cover a 200-nm wavelength range from 695 to 895 nm. The detectors were controlled with a DCC card (Becker & Hickl, Berlin, Germany) and connected to two SPC130 time-correlated single photon counting (TCSPC) modules (Becker & Hickl), with 1024 time bins per measurement, covering 12.5 ns following the laser pulse. Data acquisition by the TCSPC cards was synchronized with the laser pulse by a trigger signal from the laser. Overall, this detection configuration allowed us to measure the emitted fluorescent signal with 12.7-ps temporal and 12.5-nm spectral resolution simultaneously.

2.2.

Optical Phantoms

We characterized the instrument and algorithm with custom made optical phantoms.²⁵ Briefly, the phantom is a hollowed out 25-mm diameter, 1-mm solid-walled cylinder, with 12, 1-mm diameter holes drilled in the base as shown in Figs. 2(a) and 2(b). This design allows us to flexibly insert tubes filled with fluorophore into the holes to generate an arbitrary number of spatial combinations. The hollow part of the phantom was filled with a liquid phantom to mimic the approximate optical properties of tissues. Specifically, we used 1% intralipid (Sigma-Aldrich, Michigan) with 50-ppm India ink added. Overall, the optical properties are $\sim {\mu }_{\mathrm{a}}=0.1{\mathrm{cm}}^{-1}$, ${\mu }_{\mathrm{s}}^{\prime }=10{\mathrm{cm}}^{-1}$ in the NIR region.^26,27

Fig. 2

(a) Diagram and (b) photograph of the optical phantom used in these experiments. Normalized (c) absorption and (d) emission spectra of Alexa Fluor 680, 700, 750, and 790. The vertical line in (c) and shaded rectangle in (d) represent the instrument excitation and collection wavelengths, respectively.

We used our system to image several combinations of four different commonly used fluorophores: Alexa Fluor (AF) 680, AF 700, AF 750, and AF 790 (ThermoFisher Scientific). The normalized (to maximum) absorption and emission spectra for each are shown in Figs. 2(c) and 2(d). The vertical bar in Fig. 2(c) indicates the center emission wavelength of the laser in the first configuration (see Sec. 2.1); as indicated it allows us to excite all four of the fluorophores with varying efficiency due to the relative extinction coefficient for each. The shaded rectangle in Fig. 2(d) indicates the detection wavelength range of the spectrometers, indicating that our system is also capable of detecting the emission of all four fluorophores simultaneously. The absorption and emission maxima, as well as the fluorescence lifetime for each fluorophore are summarized in Table 1.²⁸

Table 1

Photophysical properties of fluorophores used in these experiments.

Three illustrative cases were studied for this paper:

During data acquisition, phantoms were placed on the inner rotation stage. The power of the laser light at the sample was 37 mW, and the gain of PMT detectors was set to 98% (adjustable in Becker & Hickl software). Phantoms rotated through 360 deg with 5 deg increments (72 angles). After each full rotation of the sample stage, the outer stage was rotated by 10 deg, and the acquisition was repeated. In total then, we positioned the detectors at six angular positions ($-45\mathrm{deg}$, $-25\mathrm{deg}$, $-10\mathrm{deg}$, 0 deg, $+35\mathrm{deg}$, and $+45\mathrm{deg}$ relative to the optical axis) and obtained $72\times 6=432$ projective measurements. However, more or fewer measurements could easily be obtained by adjusting the rotation angles of the two stages. Each “measurement” here was 1024 temporal measurements (with 12.7-ps resolution) and 16-spectral measurements (with 12.5-nm resolution), for a total of $7.1\times {10}^6$ data points per scan.

2.3.

Signal Demixing and Image Reconstruction

We previously proposed and validated (in silico) a two-stage algorithm to efficiently utilize the large data sets generated by our system.²³ This algorithm is described briefly here for completeness. The workflow is shown in Fig. 3. In the first stage, the measured temporal and spectral signal from each projection was demixed independently and accumulated tomographic measurements were extracted for each fluorophore using a signal library. In the second stage, an FMT reconstruction was performed to retrieve the location and concentration for each fluorophore.

Fig. 3

Workflow of the two-stage demixing and reconstruction algorithm used in this study.

2.4.

Generation of Extended Libraries

To implement this demixing algorithm, it was necessary to generate an extended library of fluorophore signatures ‘$\mathbf{S}$’ [Eqs. (1) and (2), above] to capture the unknown source-target-detector pathlength variability in the tomographic problem. We used two methods to do this.

2.4.1.

Method 1: experimental measurement

As shown in Fig. 4(a), a straw filled with each fluorophore at a concentration of $1\mu \mathrm{M}$ was inserted into a fixed hole on the phantom base. By rotating the detector to three angles with labels 1, 2 and 3, we can measure the temporal and spectral signatures corresponding to three different path lengths: $L_1+L_2=23\mathrm{mm}$, $L_1+L_3=25\mathrm{mm}$, and $L_1+L_4=27\mathrm{mm}$. This process was repeated for all four fluorophores (AF680, AF700, AF750, and AF790), so that a library containing 12 signatures in total was generated. Example signatures are shown in Figs. 4(b)–4(e) for AF 680 [Fig. 4(b)], AF700 [Fig. 4(c)], AF750 [Fig. 4(d)], and AF790 [Fig. 4(e)], corresponding to $L=23\mathrm{mm}$. Here, each signature was normalized by the summation of all temporal-spectral photons, and would be transformed to a column vector in $\mathbf{S}$ defined in Eq. (2).

Fig. 4

(a) Geometry used for experimental measurement of the signature library with three source-detector pair combinations. Example measured signatures for (b) AF680, (c) AF00, (d) AF750, and (e) AF790, respectively.

2.4.2.

Method 2: calculated signature library

More generally it is desirable to be able to directly calculate (estimate) the library based on the photophysical properties of the fluorophores and mathematical models of light propagation in diffusive media. The latter requires knowledge of the optical properties of the diffusive media, in particular the absorption coefficient ${\mu }_{\mathrm{a}}$ and the reduced scattering coefficient ${\mu }_{\mathrm{s}}^{\prime }$. In general, these are unknown, so we used our instrument to make time-resolved transmission measurements to first determine the optical properties of the medium, as has been reported previously.^30,31

For this, we used the broadband emission of our laser (path # 2, Sec. 2.1) to illuminate the sample at the range of fluorescence emission wavelengths and measured the time-resolved spectral signal through the medium. We then fit the following time-resolved diffusion equation to data at each wavelength channel:

Eq. (8)

$$T({\mathbf{r}}_s,{\mathbf{r}}_d,t)=A{\int }_t\mathrm{TIRF}(t-t^{\prime })G(|{\mathbf{r}}_d-{\mathbf{r}}_s|,t-t^{\prime })\mathrm{d}{t}^{\prime },$$

where $\mathrm{TIRF}(t)$ is the temporal impulse response function of the system which could be measured experimentally. $G$ denotes the Green’s function from the source to the detector, and has the form

Eq. (9)

$$G=\frac{c}{{(4\pi Dct)}^{3/2}}\exp (-\frac{{|{\mathbf{r}}_d-{\mathbf{r}}_s|}^2}{4Dct}-{\mu }_{\mathrm{a}}ct),$$

where $D=1/[3({\mu }_{\mathrm{a}}+{\mu }_{\mathrm{s}}^{\prime })]$ and $c$ is the light velocity in the medium. Here, the parameter $A$ is a scalar constant that accounts for the system calibration, which depends on, e.g., the specific source power, detector efficiency. For curve fitting, we used the “lsqcurvefit” function in MATLAB^®.

This procedure retrieved ${\mu }_{\mathrm{a}}$ and ${\mu }_{\mathrm{s}}^{\prime }$ for each of the detection wavelength channels of our system. Example measured wavelength-dependent optical properties for our phantom are shown in Figs. 5(a) and 5(b), and these are generally in agreement with the literature values for intralipid in the NIR range.²⁶ Using these data, we then computed the signature library for each fluorophore and pathlength combination using

Eq. (10)

$$S({\mathbf{r}}_s,{\mathbf{r}}_d,t)={\int }_{t}Fl(\mathbf{r},t-t^{\prime })W({\mathbf{r}}_s,{\mathbf{r}}_d,\mathbf{r},t^{\prime })\mathrm{d}{t}^{\prime },$$

where $Fl(\mathbf{r},t)$ denotes the time-resolved emission light profile of the fluorophore, which incorporates quantum yield, lifetime information, as well as relative absorption and emission efficiency for the laser and wavelength combination. $W$ is the time-dependent weight function and has the form

Eq. (11)

$$W({\mathbf{r}}_s,{\mathbf{r}}_d,\mathbf{r},t)={\int }_t{G}^x({\mathbf{r}}_s,\mathbf{r},t^{\prime })G^m({\mathbf{r}}_d,\mathbf{r},t-t^{\prime })\mathrm{d}{t}^{\prime },$$

where $G^x$ and $G^m$ denote the excitation and emission Green’s function, respectively. To find $Fl(\mathbf{r},t)$, we first measured the time-resolved fluorescence profile of each fluorophore in a nonscattering tube, which incorporated the relative emission at each wavelength and the fluorophore lifetime. We then convolved this with computed Green’s functions for each pathlength combination (23, 25, and 27 mm as above). Example simulated signatures for AF680, AF700, AF750, and AF790 (red dotted line) are shown in Figs. 5(c)–5(f), along with the corresponding experimental signatures (blue line). We plotted six dominant wavelength channels for each fluorophore type. By inspection, the signatures generated by the two methods were very similar, and the computed cosine similarity between the simulated and experimental signatures was greater than 99% for all fluorophore types.

Fig. 5

(a) Measured wavelength-dependent absorption coefficient and (b) reduced scattering coefficient of the phantom. Experimentally measured (blue line) and computed fluorescence signatures (dotted red line) for (c) AF680, (d) AF700, (e) AF750, and (f) AF790, respectively.

3.1.

Imaging Two Fluorophores Concurrently

We first evaluated the performance of our instrument and algorithm in simultaneously imaging Alexa Fluor 680 (AF680) and AF700. AF680 and AF700 have emission peaks separated by only 20 nm and lifetimes separated by only 0.2 ns (1.2 and 1.0 ns, respectively), and are therefore in general difficult to image simultaneously with either spectral or temporal data alone. The positions of the two fluorescent targets in this example are shown in Fig. 6(a). A representative measurement from a single source-detector pair is shown in Fig. 6(b). In total we collected 432 such measurements from a single axial slice. There are significant contributions from both fluorophores in the first three detector channels (centered at 702, 714.5, and 727 nm, respectively). We plotted the data from the first and third channels as a function of rotation angle for one source-detector pair, along with the corresponding tomographic reconstructions [Figs. 6(d) and 6(e)] for these data. By inspection, significant cross-talk was observed between the channels, illustrating the challenges associated with imaging two closely overlapping fluorophores within the NIR wavelength range.

Fig. 6

(a) Placement of two fluorophores in the phantom. (b) Example measurement from a single source-detector projection. (c) Example measured CW traces corresponding to the first (green) and third (red) wavelength channels (centered at 702 and 727 nm). Image reconstructions for the same data corresponding to (d) 702 and (e) 727 nm, illustrating the significant bleed associated with pure spectral data.

We next implemented our two-stage imaging algorithm²³ and generated the extended library for demixing of the fluorophores (step 1) using either the experimental or computational methods. The demixed signals from the two fluorophores and the associated reconstruction for the experimental library are shown in Figs. 7(a)–7(c). Likewise, the demixed signals and reconstructions obtained using the computational method are shown in Figs. 7(d)–7(f).

Fig. 7

(a) Demixed signal as a function of rotation angle for AF680 (green) and A700 (red) using the experimental signature library. The corresponding reconstructions for (b) AF680 and (c) AF700 are shown. (d) Example demixed fluorescent signals for AF680 (green) and AF700 (red) using the computed signature library, as well as the (e) AF680 and (f) AF700 reconstructions.

Both methods for computing the library were able to demix the signal contributions from the two fluorophores. The image cross-talk between the two fluorophores was zero with the two methods, compared to 61.23% for the “naïve” approach shown in Fig. 6(e). In the final images, the targets were reconstructed with only 0.71, 0.5, 0.5, and 1.12 mm center-to-center error from the known positions of the targets for the cases in Figs. 7(b), 7(c), 7(e), and 7(f) respectively. In combination, this illustrates the ability of the algorithm to separate and image two highly overlapping fluorophores with experimental data.

3.2.

Computed Error in the Demixing Stage

We next studied the quantitative error introduced at the demixing stage as follows. We scanned the phantom with a tube of either AF680 or AF700 (but not both) present in the same position as above, and then summed them numerically in MATLAB^®. We then applied the demixing algorithm (stage 1) and computed the point-by-point quantitative error from the original signals using

In addition, we tested of our system in imaging multiple targets placed at positions deeper in the medium. AF680 and AF700 were positioned near the center of the phantom, and separated by 7-mm center-to-center as shown in Fig. 8(a). Analogous to the above, reconstructed images obtained using the experimental library are shown in Figs. 8(b) and 8(c), and the final relative reconstructed concentration ratio was $1.3:1$. Likewise, the reconstructions for AF680 and AF700 obtained using the computed library are shown in Figs. 8(d) and 8(e), respectively. The final relative concentration ratio here was $1.06:1$. The error in reconstruction of the position of the targets (center-to-center) was 0.5, 0, 0.71, and 0.5 mm in Figs. 8(b)–8(e), respectively.

Fig. 8

(a) Placement of fluorophores at deep-seated position in the optical phantom. Reconstruction of (b) AF680 and (c) AF700 with the experimental library. Reconstruction of (d) AF680 and (e) AF700 with the computed library.

3.3.

Four-Fluorophore Phantom

We then performed tomographic imaging of four fluorophores in our phantom model. The positions and concentrations of AF680, AF700, AF750, and AF790 are shown in Fig. 9(a). An example measured signal for a single source-detector pair is shown in Fig. 9(b). Red curves in Figs. 9(c)–9(f) show the demixed signal (after step 1) for each fluorophore where the experimentally measured library was used. The corresponding reconstructed images for each fluorophore are shown in Figs. 9(g)–9(j). As above, zero cross-talk was observed between fluorophores for each of the reconstructed images, indicating the ability of the algorithm to efficiently demix four fluorophores. The center-to-center error in the position of the reconstructed fluorophores was 0.5, 0.71, 1.58, and 1 mm for each fluorophore, respectively. The relative fluorophore concentrations for AF680: AF700: AF750: AF790 was $1.4:1.3:1:4.6$ which is in relatively good agreement with the true ratios of $1:1:1:5$. We attribute most of this quantitative error to the impact of optical properties (absorption) of the fluorescent targets on the measured signals, which were not modeled in our approach. Our extended library could be expanded to incorporate this in the future.

Fig. 9

(a) Placement of four fluorophores in the optical phantom. (b) Example measurement from a single source-detector projection. (c)–(f) Example CW traces after demixing with experimental (red solid) and computed (green x) fluorescent signature libraries corresponding to AF680, AF700, AF750, and AF790, respectively. (g)–(j) Reconstructions of AF680, 700, 750, and 790 obtained using data demixed with experimentally measured fluorescence library, respectively. (k)–(n) Reconstructions of AF680, 700, 750, and 790 obtained using data demixed with computed fluorescence signature library, respectively.

Likewise, the green curves in Figs. 9(c)–9(f) show the demixed signal for each fluorophore with the computed library, and Figs. 9(k)–9(n) show the corresponding image reconstructions. Zero cross-talk was observed among fluorophores. The center-to-center position errors of the reconstructed targets were 0.5, 0.71, 1.58, and 1 mm, respectively. The final concentration ratio of the image reconstructions was about $1.2:1.6:1:7$, which is in reasonable agreement with the true values, considering the complexity of the phantom model. We note that the error associated with the computed library was higher than that obtained with the experimentally measured library. This was expected, but nevertheless the computed library is practically preferable in many real experimental situations.

3.4.

Colocalized Fluorophores

Last, we considered the potential impact of spatially overlapping (colocalized) fluorophores on the instrument and algorithm. In this configuration, we placed AF680 and AF750 at one location and AF700 and AF790 at a second location. The rationale here was that these combinations were most likely to produce unintended bleed among the different fluorophore channels. To avoid chemical interactions between different fluorophores, we designed a structure made with seven capillary tubes, six of which were filled with one of the two types of fluorophores. The configuration of this phantom is shown in Fig. 10(a). Figures 10(b)–10(e) show the final reconstructions for each fluorophore, where the demixing was performed with the experimentally measured library. Again, zero-cross talk was observed, illustrating the robustness of the algorithm with spatially overlapping targets. The center-to-center error in the position of the reconstructed flurophores was 0.71, 0.5, 0.71, and 0.5 mm, respectively. Last, the computed concentration ratios were $1.04:1:1.52:6.07$, which is in relatively good agreement with the true ratios of $1:1:1:5$.

Fig. 10

(a) Placement of four fluorophores in the optical phantom, where AF680 and AF750, and AF700 and AF790 were arranged in pairs as shown. Image reconstructions for (b) AF680, (c) AF700, (d) AF750, and (e) AF790.