A method to unmix multiple fluorophores in microscopy images with minimal a priori information

The ability to quantify the fluorescence signals from multiply labeled biological samples is highly desirable in the life sciences but often difficult, because of spectral overlap between fluorescent species and the presence of autofluorescence. Several so called unmixing algorithms have been developed to address this problem. Here, we present a novel algorithm that combines measurements of lifetime and spectrum to achieve unmixing without a priori information on the spectral properties of the fluorophore labels. The only assumption made is that the lifetimes of the fluorophores differ. Our method combines global analysis for a measurement of lifetime distributions with singular value decomposition to recover individual fluorescence spectra. We demonstrate the technique on simulated datasets and subsequently by an experiment on a biological sample. The method is computationally efficient and straightforward to implement. Applications range from histopathology of complex and multiply labelled samples to functional imaging in live cells. © 2009 Optical Society of America OCIS codes: (170.3650) Lifetime-based sensing; (180.2520) Fluorescence microscopy. References and links 1. J. H. Frank, A. D. Elder, J. Swartling, A. R. Venkitaraman, A. D. Jeyasekharan, and C. F. Kaminski, “A white light confocal microscope for spectrally resolved multidimensional imaging,” J Microsc. 227, 203–215 (2007). 2. C. F. Kaminski, R. S. Watt, A. D. Elder, J. H. Frank, and J. Hult, “Supercontinuum radiation for applications in chemical sensing and microscopy,” Appl. Phys. B: Lasers and Optics 92(3), 367–378 (2008). 3. S. Schlachter, A. D. Elder, A. Esposito, G. S. Kaminski, J. H. Frank, L. K. van Geest, and C. F. Kaminski, “mhFLIM: Resolution of heterogeneous fluorescence decays in widefield lifetime microscopy,” Opt. Express 17(3), 1557–1570 (2009). 4. A. Esposito, T. Tiffert, J. M. A. Mauritz, S. Schlachter, L. H. Bannister, C. F. Kaminski, and V. L. Lew, “FRET Imaging of Hemoglobin Concentration in Plasmodium falciparum-Infected Red Cells,” PLoS ONE 3(11), e3780 (2008). 5. T. Zimmermann, J. Rietdorf, and R. Pepperkok, “Spectral imaging and its applications in live cell microscopy,” FEBS Lett. 546(1), 87–92 (2003). #116166 $15.00 USD Received 2 Sep 2009; revised 9 Nov 2009; accepted 14 Nov 2009; published 30 Nov 2009 (C) 2009 OSA 7 December 2009 / Vol. 17, No. 25 / OPTICS EXPRESS 22747 6. A. D. Elder, A. Domin, G. S. Kaminski Schierle, C. Lindon, J. Pines, A. Esposito, and C. F. Kaminski, “A quantitative protocol for dynamic measurements of protein interactions by Frster resonance energy transfersensitized fluorescence emission,” Journal of The Royal Society Interface 6(Suppl 1), S59–S81 (2009). 7. H. Shirakawa and S. Miyazaki, “Blind Spectral Decomposition of Single-Cell Fluorescence by Parallel Factor Analysis,” Biophys. J. 86(3), 1739–1752 (2004). 8. D. Chorvat, J. Kirchnerova, M. Cagalinec, J. Smolka, A. Mateasik, and A. Chorvatova, “Spectral Unmixing of Flavin Autofluorescence Components in Cardiac Myocytes,” Biophys. J. 89(6), L55–L57 (2005). 9. R. A. Neher, M. Mitkovski, F. Kirchhoff, E. Neher, F. J. Theis, and A. Zeug, “Blind Source Separation Techniques for the Decomposition of Multiply Labeled Fluorescence Images,” Biophys. J. 96(9), 3791–3800 (2009). 10. J. R. Lakowicz, Principles of Fluorescence Spectroscopy (Springer, 2006). 11. P. Verveer, A. Squire, and P. Bastiaens, “Global Analysis of fluorescence lifetime imaging microscopy data,” Biophys. J. 78, 2127–2137 (2000). 12. G.-J. Kremers, E. B. van Munster, J. Goedhart, and T. W. J. Gadella, “Quantitative Lifetime Unmixing of Multiexponentially Decaying Fluorophores Using Single-Frequency Fluorescence Lifetime Imaging Microscopy,” Biophys. J. 95(1), 378–389 (2008). 13. G. Strang, Linear algebra and its applications (Thomson, 2006). 14. A. H. A. Clayton, Q. S. Hanley, and P. J. Verveer, “Graphical representation and multicomponent analysis of single-frequency fluorescence lifetime imaging microscopy data,” J. Microscopy 213(1), 1–5 (2004). 15. F. S. Wouters and A. Esposito, “Quantitative analysis of fluorescence lifetime imaging made easy,” HFSP Journal 2(1), 7–11 (2008). 16. M. A. Digman, V. R. Caiolfa, M. Zamai, and E. Gratton, “The phasor approach to fluorescence lifetime imaging analysis,” Biophys. J. 94, L14–L16 (2007). 17. A. D. Elder, S. C. Schlachter, and C. F. Kaminski, “Theoretical investigation of the photon efficiency in frequency-domain FLIM,” J. Opt. Soc. Am. A 25, 452–462 (2008). 18. D. M. Owen, E. Auksorius, H. B. Manning, C. B. Talbot, P. A. A. de Beule, C. Dunsby, M. A. A. Neil, and P. M. W. French, “Excitation-resolved hyperspectral fluorescence lifetime imaging using a UV-extended supercontinuum source,” Opt. Lett. 32(23), 3408–3410 (2007). 19. P. Serge, J. R. P. Michael, K. Daekeun, K. Ki Hean, J. S. Tsu-Te, and T. C. S. Peter, “Frequency domain lifetime and spectral imaging microscopy,” Microsc. Res. Tech. 69(11), 861–874 (2006). 20. Q. S. Hanley, D. J. Arndt-Jovin, and T. M. Jovin, “Spectrally Resolved Fluorescence Lifetime Imaging Microscopy,” Appl. Spectrosc. 56, 155–166 (2002). 21. J. A. Levitt, D. R. Matthews, S. M. Ameer-Beg, and K. Suhling, “Fluorescence lifetime and polarization-resolved imaging in cell biology,” Curr. Opin. Biotechnol. 20(1), 28–36 (2009). 22. N. P. Galletly, J. McGinty, C. Dunsby, F. Teixeira, J. Requejo-Isidro, I. Munro, D. S. Elson, M. A. A. Neil, A. C. Chu, P. M. French, and G. W. Stamp, “Fluorescence lifetime imaging distinguishes basal cell carcinoma from surrounding uninvolved skin.” British Journal of Dermatology 159, 152–161 (2008). 23. A. D. Elder, S. M. Matthews, J. Swartling, K. Yunus, J. H. Frank, C. M. Brennan, A. C. Fisher, and C. F. Kaminski, “The application of frequency-domain fluorescence lifetime imaging microscopy as a quantitative analytical tool for microfluidic devices,” Opt. Express 14(12), 5456–5467 (2006). 24. S. M. Matthews, A. D. Elder, K. Yunus, C. F. Kaminski, C. M. Brennan, and A. C. Fisher, “Quantitative kinetic analysis in a microfluidic device using frequency-domain fluorescence lifetime imaging,” Anal. Chem. 79(11), 4101–4109 (2007). 25. X. Dai, Z. Yue, M. E. Eccleston, J. Swartling, N. K. H. Slater, and C. F. Kaminski, “Fluorescence intensity and lifetime imaging of free and micellar-encapsulated doxorubicin in living cells,” Nanomedicine: Nanotechnology, Biology and Medicine 4(1), 49–56 (2008). 26. X. W. Dai, M. E. Eccleston, Z. L. Yue, N. K. H. Slater, and C. F. Kaminski, “A spectroscopic study of the self-association and inter-molecular aggregation behaviour of pH-responsive poly(L-lysine iso-phthalamide),” Polymer 47(8), 2689–2698 (2006).


Introduction
It is becoming increasingly common to use multiple fluorescent species to label biological samples.This permits one to unravel complex biological pathways within living systems, because multiple, interconnected events can be tracked simultaneously and their spatial and temporal dynamics can be followed.However, the interpretation of measurements from multiply labeled samples presents a complex challenge because it is difficult to unmix the contributions of individual fluorophores from the intensity of a given image pixel.Most biological fluorophores have excitation and emission spectra with bandwidths of approximately 30 to 60 nm, and one always faces the problem of spectral and spatial overlap between them.The likelihood of inadvertently exciting multiple fluorescent species at once is thus large, and this makes the interpretation of the resulting images prone to error.
Ultimately, the goal is to recover the molecular fraction (i.e., the relative concentrations) of each fluorescent species in every pixel location of the image.This is however not possible in multiply labelled samples without additional independent measurements.Here we show how one can combine several attributes of the fluorescence process, including lifetime, excitation and emission spectra to permit the quantitative separation, or unmixing of a multiply labeled sample with a minimum of a priori information.Our method makes no assumption on the shape of the fluorophore species' excitation or emission spectra, or their degree of overlap.The only assumption is that their lifetimes differ.
We have recently shown how supercontinuum radiation can be used for spectrally resolved confocal microscopy by providing excitation light that is continuously tunable across the visible spectrum and beyond into the near infra-red region [1] [2].We shall show here that such excitation-wavelength resolved microscopy is ideally suited to the present unmixing problem.
The pulse width of most supercontinuum sources is in the fs to ps range.Combined with their high repetition rates (on the order of tens of MHz) this means that they are ideal sources for lifetime resolved imaging [3] [4].Here we use supercontinuum excited, spectrally resolved Time Correlated Single Photon Counting (TCSPC) to provide multi-dimensional data as input to our unmixing algorithm.We acquire 4-dimensional datasets (two spatial dimensions, excitation spectrum, and lifetime) to greatly enhance the contrast between multiple labels.
Traditionally, unmixing has relied on measuring fluorescence emission using several channels, each covering a segment of the emission spectrum.Prior to imaging a multiply labeled sample, a 'calibration' is performed from which the contribution of each fluorophore species to each emission channel can be calculated.A sample labelled with several fluorophore species can then be quantitatively unmixed using the calibration information.This method is known as linear unmixing [5] and requires at least N spectral emission channels to resolve N fluorophores.Linear unmixing is also used in FRET quantification [6].
In contrast to linear unmixing, so called 'blind unmixing' algorithms seek to quantitatively unmix a sample labeled with unknown fluorophores.Several methods have been proposed for this purpose, such as PARAFAC [7], principle component analysis [8], and non negative matrix factorization (NMF) [9].These methods use various algorithms to unmix multiple fluorophores in spectrally resolved microscopy images.
Here we present a novel method for blind unmixing that uses both spectral and lifetime differences between fluorophores to differentiate them.Though demonstrated here for confocal TCSPC image data, the method is equally valid for widefield, frequency domain FLIM methods, for emission resolved microscopy (in contrast to the excitation resolved microscopy we use here), or combinations between the different approaches.
The method offers a powerful tool for biomedical research.We envisage applications ranging from histopathology, where one seeks to differentiate between classes of tissue via multiple stains, to physiology, where multiple fluorophores could be used to unravel complex signalling cascades.The paper begins with an outline of the theory behind our method followed by its validation using simulated data.The method is also validated experimentally before conclusions are drawn and an outlook is given.

The fractional intensity, α
We start with the only assumption that the fluorescent species have different lifetimes with monoexponential decays.As an example, consider a sample containing just two fluorescent species.Each pixel of a FLIM image can thus be fitted with a biexponential decay in the form: where p is the pixel index, τ i are the lifetimes, and a and b are the preexponential coefficients.It is straightforward to determine the fractional contribution α i of each fluorescent species to the overall fluorescence intensity.For α 1 we thus obtain: It must be noted here, that α(p), referred to as the 'fractional intensity' [10], is a function of excitation and emission wavelengths (as, in general, the fluorophores have different spectra) and thus α(p) is not necessarily representative of the molecular fraction, which is wavelength invariant and represents the quantity we seek to estimate.In Section 2.2 we show that by taking FLIM images at two or more excitation (or emission) wavelengths, it is possible to produce a wavelength-independent estimate for the molecular fraction.
A robust method for determining the fractional intensity is global analysis [11].In a global analysis one assumes that every pixel in a FLIM image may contain contributions from different lifetime components, but, importantly, the component lifetimes are not themselves a function of spatial location.Only the proportion of fluorescence belonging to each of the pixel invariant lifetimes (i.e. the fractional intensity) is allowed to vary.This assumption allows α(p) to be calculated with a higher signal-to-noise ratio (SNR) than available to pixel-by-pixel fitting methods, because a large number of data points contribute to the determination of the component lifetimes [11].

Algorithm for unmixing F fluorophores
Consider now the most general case where a sample is labeled with F different fluorescent species, each having a different lifetime.Every pixel p in this image contains a mixture of fluorescent species (with excitation spectra that may, but need not, differ from one another).In this case, the intensity detected in a given pixel p in response to a very short pulse of light varies according to the lifetimes, number of fluorophores and excitation efficiencies of each fluorescent species in that pixel: where B i represents the brightness factor (a combination of the Einstein absorption and spontaneous emission coefficients, excitation laser intensity, spectral overlap of the absorption and emission spectra with the excitation bandwidth and detection bandwidth respectively, and other microscope efficiency factors) of fluorophore species i; N i represents the number of fluorophores present in pixel p from fluorophore species i; τ i is the lifetime of fluorophore species i and λ is the excitation wavelength.Importantly for what follows, one should note that B i contains the spectral properties of the fluorescent species and is thus a function of wavelength which may differ between different molecular species.
The time averaged intensity, as would be detected in response to continuous excitation, is calculated by integrating over all time: The fractional intensity of each fluorophore species (α i ) can then be written as: α i varies with excitation wavelength, because each fluorescent species has a different excitation spectrum.The fractional intensities for a multiply labeled sample can easily be measured from a global analysis.However, we are not interested in α i , which changes with excitation wavelength, emission window, laser power, etc. but instead seek a more fundamental quantity, the molecular fraction ( f i ).The molecular fraction denotes the number of fluorophores belonging to species i in a particular pixel divided by the total number of fluorophores.It represents the ultimate goal of all unmixing algorithms and is defined as follows: How can one determine f i (p) for each fluorescent species?If the brightness spectra of the fluorophore species are identical, then f i (p) = α i (p) [12].Alternatively, if the brightness spectra of the fluorophore species (B i (λ )) were known from a priori measurements, one need only measure < I(λ , p) > at F wavelengths, and then f i (p) could easily be determined from Eq. ( 5).This is the approach taken by the traditional linear unmixing method [5].In linear unmixing, each pixel contains F unknowns which need to be determined to obtain the molecular fractions f i (p).It is thus clear that one needs to perform at least Λ = F measurements, where Λ is the number of spectral positions at which each pixel is sampled, to recover the molecular fraction in each pixel.
Let us now assume, however, that we do not have any a priori knowledge of B i (λ ).It is now clearly not sufficient to perform Λ = F spectral measurements to unmix the signal contributions to each pixel.Additional information is required, which, as we shall show, can be obtained from lifetime measurements and global analysis.The latter yields the intensity fractions α i (λ , p).Let us assume then, that we have available data for α i (λ , p) and < I(λ , p) >, each measured at Λ excitation wavelengths for each image consisting of P pixels.
To develop the method by which this unmixing is achieved, it is useful to recast Eq. ( 6) in vector form: where I is a matrix containing the total fluorescence intensity and α i α i α i is the fractional intensity of species i.At each pixel, α i and I are determined for Λ wavelength positions and the resulting matrices (α i α i α i , I) are P × Λ in size.N i is a P × 1 vector, B i is a Λ × 1 vector, '•' is the entrywise (Hadamard) product and '⊗' is the outer product.Eq. ( 8) contains the essence of the unmixing algorithm demonstrated in this paper: Because the orthogonal vectors B i and N i define the intensity matrix α i α i α i • I in a unique fashion, information on molecular fraction (contained in N i ) and the individual spectra of the fluorophores (contained in B i ) can in principle be recovered from a decomposition of This goal can be reached using Singular Value Decomposition (SVD).To understand how this is accomplished, recall that SVD can factor any matrix into a product of three matrices UΣ Σ ΣV T [13], hence: where U (P × P) and V (Λ × Λ) are orthogonal matrices and Σ Σ Σ (P × Λ) is a diagonal matrix containing the r singular values (σ j ) (in descending order) of α i α i α i • I, where r is the matrix rank where u i,j and v i,j are the j th column vectors of U i and V i , respectively (also called the left and right singular vectors).Referring to Eq. ( 8), one can see that u i,1 and v i,1 are estimators ( Bi = u i , Ni = v i ) for N i and B i respectively.
There is, however, the problem of how to divide σ i,1 between the two singular vectors.This is not specified, and thus there exists an unresolved proportionality constant between the estimators and the physical values (i.e., n i Ni = N i and m i Bi = B i where n i m i = σ i,1 ).
Thus the molecular fraction estimator, fi (p): is related to the true physical value f i (p) by where c i is a correction factor based on either the ratio Thus if the true molecular fraction is known at any one pixel in the image, or if the ratio of the brightness spectra is known at one wavelength, then the estimated molecular fraction can be corrected to yield true molecular fractions.
The algorithm hinges on the orthogonality between B i and N i .This is only the case when the brightness of a given fluorophore species changes with wavelength (not pixel number) and the number of fluorophores of a given species changes only with pixel number (and not wavelength).In reality, of course, there will be variability in the spectral properties with pixel location as a result of molecular interactions of the fluorophores with their environment.We take this into account by allowing some fluctuations in lifetime and spectral properties for each fluorophore from pixel to pixel and the algorithm is robust enough to deal with such abberations.However, the method is no longer valid for cases where the fluorophores are subject to very large property changes (such as may be caused by FRET).

Simulations
In order to demonstrate the unmixing algorithm described in Section 2 and to explore its limitations we have constructed a model to simulate experimental data for various signal-to-noise ratios (SNR), fluorophore excitation and lifetime properties.
The simulated datasets, shown in Fig. 1, mimic data acquired by confocal microscopy with TCSPC hardware and a wavelength flexible excitation source.
We show that the algorithm is capable of unmixing the fluorophores in the simulated sample, and then continue to further characterize the robustness of the algorithm for faintly emitting samples with highly overlapping spectra and/or lifetimes.

Generation of test data
Two-dimensional patterns were generated as test sets with each pixel containing spectral and lifetime information.We assume that every spatial pixel in the simulated image contains multiple fluorophores, drawn from two separate fluorophore species.Each of the two species was assigned its own (monoexponential) lifetime and associated brightness spectrum B(λ ).B(λ ) was assumed to be a Gaussian (although the shape of the assumed spectrum is immaterial), with peak excitation wavelength λ μ1 and bandwidth λ σ 1 for species 1 (and likewise for species 2): In order to simulate the variability of real fluorophores, each fluorophore from a given species is modeled as having a slightly different peak excitation wavelength and bandwidth, drawn from a Gaussian distribution with mean Λ μ 1 and standard deviation Λ σ 1 for species 1.These concepts are exemplified by the simulated excitation spectra shown in Fig. 1(d).
The fluorescence decays were furthermore simulated assuming that the lifetimes of fluorescent species 1 and 2 were distributed around 2 ns and 1 ns as shown in Fig. 1(e) (similar to lifetime differences found in fluorophores of biological interest, see Section 5.For a two component system, the instantaneous intensity at each pixel p and wavelength λ can be written as follows (see section 2): Equations ( 13) and ( 14) were used to create a 4 dimensional (x, y,t, λ ) dataset with in-house developed software, written in Matlab (Mathworks, Inc.).The impulse response of the TCSPC system is modeled by a 300 ps FWHM Gaussian, close to the value measured for our hardware.Poisson noise is then added to the dataset to model photon noise (See Fig. 1(f)).We assumed mean signal levels of around 800 photons per pixel (see Fig. 1(g)), mimicking a typical signal measurement.Intensity images were obtained by integrating the time resolved signals in each pixel for each excitation wavelength.

Global analysis of lifetime data
The purpose of the simulated datasets is to assess the capability of our algorithm to retrieve the original fluorophore distributions in the phantom images in the presence of noise from various sources.For this we assume that we do not know anything about the fluorophores' spectral, lifetime, and concentration distributions across the image field.The only assumption we make is that the mean lifetimes τ 1 and τ 2 for fluorophores 1 and 2 are measurably different.For a two component mixture the global analysis method is a particularly powerful and simple way to recover τ 1 and τ 2 and the fractional intensities as defined in Section 2.1 [14].We start the global analysis by calculating A and B parameters from the simulated decay data according to [11]: The AB parameters are then plotted on a phasor diagram (AB plot) for each excitation wavelength and pixel (see Fig. 2) [15] [16].From this, we extract two global lifetime estimates, τ1 and τ2 , from which, in turn we are now able to calculate fractional intensities α1 (λ , p) [11] using the following relations (see Fig. 2(c)): Section 2 outlined how SVD is used to convert the estimated fractional intensities into estimated molecular fractions.
The results are shown in Fig. 3. Clearly, the spatial distribution of fluorophores 1 and 2 are recovered from the data and the estimated molecular concentrations N1 and N2 (Fig. 3(b)) show a similar distribution to the original (Fig. 3(a)).Likewise, the shape of the estimated fluorophore spectra, B1 and B2 bear a close resemblance to the true spectra (Fig. 3(c): normalised estimated spectra shown in green and yellow, original spectra shown in red and blue).As discussed in section 2 the estimators represent molecular concentration and spectra, subject to the constraint n i m i = σ i,1 .
If the true molecular fraction is known in one pixel (e.g. through a calibration experiment or knowledge of expected fluorophore distributions at a certain point in the sample) then Ni Ni Ni can be corrected to yield true molecular fractions.The results of such a correction are shown in Fig. 4.

(d).
There is no spatial structure discernible in Fig. 4(d) proving that a single proportionality factor accounts for the difference between the original and estimated molecular fraction and no spatial bias is evident.

Figure of merit calculations
Only two excitation wavelengths are required to unmix two fluorophore species.Excitation at additional wavelength positions simply increases the SNR.
A Figure-of-Merit (FoM) was devised to characterize the fidelity of the unmixing algorithm.In an ideal situation, the estimator for Ni differs from the true distribution N i by a factor n i which is constant across the entire image.Deviations from this are either caused by photon noise, or because the unmixing algorithm is ineffective in separating the fluorophore distributions.We define the FoM as the coefficient of variation of n i compared to the shot noise limit (where N ph is the number of photons collected) [17]: Figure 5 shows the FoM as a function of both spectral separation (from 1 nm to 40 nm) and lifetime separation (from 0.2 ns to 1.1 ns).The lower left hand corner of the FoM plot corresponds to fluorophores with small spectral and lifetime separations.In this situation the fluorophores are difficult to unmix and the FoM is correspondingly large.At the other extreme, in the upper right hand corner, spectra and lifetimes are well separated, and unmixing is correspondingly easier.This is reflected by a FoM close to one, the theoretical limit.The insets show simulations corresponding to the FoMs at these two extremes: In the top right panel, the estimated distributions of the two fluorophores are shown (left hand column).The middle column shows the distribution of n i across the image and the right hand column the spectral and lifetime separations of the fluorophore species.There is little spatial bias evident in n i , and the fluorophores and the channels are easily separated by the method.In contrast, the lower inset shows data with the highest represented FoM value.This situation corresponds to spectra which are completely overlapping for the two fluorophores and a lifetime separation of only 200 ps.The images show that there is significant crosstalk and the two fluorophore channels are poorly unmixed.The FoM plot shows that even at complete spectral overlap, a good degree of unmix-ing can be obtained for lifetime separations exceeding around 1 ns.Likewise, at low lifetime separations, good FoMs can be obtained when the spectral separation is large.The lifetime of Fluorophore 1 is held at 2 ns and peak excitation at 500 nm.The lifetime of Fluorophore 2 is then varied from 1.2 ns to 2.2 ns and spectral peak is varied from 500 nm to 540 nm.The insets show the estimators for the unmixed images ( N1 and N2 ), ratios to the true images (n 1 and n 2 ) and the associated excitation spectra and lifetime distributions.See text for details.

Instrumentation
A commercial white light supercontinuum source (SC450, Fianium Inc.) was used for excitation.The source features a repetition rate of approximately 40 MHz with pulse widths on the order of tens of picoseconds.An acousto optic tunable filter, AOTF (AA Optoelectronique AOTFnc-vis, France), was used to select the excitation wavelengths, with a bandwidth of approximately 1 nm.The supercontinuum source is part of a custom-built confocal scanning microscope, details of which are given in [1].

Image acquisition
All data presented in this paper were acquired using a 60x objective (NA=1.20,Olympus Inc.).FLIM images were taken at eight excitation wavelengths: 500 nm, 520 nm, 540 nm, 550 nm, 555 nm, 560 nm, 565 nm, and 570 nm using the supercontinuum source described above.Fluorescence emission at each excitation wavelength was collected through a longpass filter, with a 590 nm cutoff.Lifetimes were resolved using a TCSPC system (SPC-830 Becker and Hickl, Germany).Each FLIM image took 60 s to acquire, during which time the photon count rate was kept well below 1% of the laser repetition rate to prevent pulse pile-up.Photobleaching was negligible.All TCSPC images were binned into 2 × 2 pixels and subsequently processed using Matlab software (all software is available upon request to the authors).

Experimental results
Figure 6 shows the lifetime-excitation dataset from our biological sample on an AB plot.The AB cloud is distributed between two lifetimes, which are extracted via graphical analysis (see Section 3.2).Fluorophore 1 is estimated to have a longer lifetime (2.94 ns) compared to Fluorophore 2 (1.20 ns).Combining these lifetimes with the intensity images at each wavelength (Fig. 6), we can calculate α images that are then unmixed using SVD to extract the contributions of the two fluorophores.
The results of the SVD unmixing are shown in Fig. 7.As seen in Fig. 7(a), the extracted spectra compare well to literature values for Alexa 546 and 555, especially considering the highly overlapping nature of these excitation spectra.The spectra have been corrected for changes in the supercontinuum laser power at each wavelength and excitation-emission efficiencies of the microscope as detailed in [1].From the spectra, we can deduce that Fluorophore 1 is the estimator for Alexa 546 dye distribution, which is further confirmed by its longer lifetime.As expected from the behaviour of the laminin IgG antibody (to which the Alexa 546 dye is conjugated), the largest signal contributions originate from the nucleus and nuclear envelope of the cells.Likewise the estimated distribution for Fluorophore 2 is mostly confined to the cytoplasm, which is also expected as Fluorophore 2's antibody (Map 2 IgG, to which the Alexa 555 dye is conjugated) targets the cytoplasm.These distributions can be combined to form the estimated molec-ular fraction, shown in Fig. 7(d).As expected, the estimated molecular fraction shows a high proportion of Fluorophore 1 (Alexa 546) in the nuclear region, and a low proportion in the cytoplasmic region.This image has not been calibrated and is thus not a true molecular fraction, but related to the true fraction via a constant factor.The estimated molecular fraction -intensity overlay shown in Fig. 7 (e) shows the nuclear localization of Fluorophore 1.

Conclusion
In conclusion, we have developed an unmixing algorithm that takes advantage of two aspects of fluorescence emission to unmix doubly labeled samples.A summary of the method consists of the following steps: 1.A wavelength-resolved lifetime image dataset is acquired using at least 2 wavelength positions for a two fluorophore sample.Wavelength resolution can be either excitation or emission resolved; here we measure at different excitation wavelengths with the advantage that signal can be collected over the entire emission bandwidth.
2. The obtained datasets are represented on AB plots and a graphical global analysis is performed to extract the lifetimes and fractional intensities in each pixel.
3. SVD is used to perform a quantitative unmixing of the the two fluorophore channels in each image.
There are several advantages to our method.It is computationally very efficient, requires no initial guesses, and does not require minimisation.Associated computation speed is compatible with the FLIM image acquisition rates.The method is applicable to both widefield and confocal FLIM techniques and works independently of whether the data are collected using wavelength resolved excitation or wavelength resolved emission imaging.Here we report on the use of excitation resolved FLIM, which has the great advantage that fluorescence can be detected over the entire emission bandwidth, using long pass filters.We show that fluorophores can be quantitatively unmixed even if no a priori information of spectra and lifetimes of the fluorophores are available.The only constraint is that the fluorophores lifetimes be discernible from one another so that a global analysis can be performed.Although we have presented the method for samples containing two fluorophores, it can be adapted for the unmixing of samples containing three or more fluorophores.The multidimensional datasets one needs to acquire in such situations are very large, containing information in multiple spatial, temporal and spectral dimensions.The technology to perform multidimensional microscopy is already widely available [18] [19] [20] [21].The opportunities this generates are great, e.g. for the study of multiplexed pathways, to provide contrast between healthy and diseased tissue [22], for the analysis of multicomponent fluidic assays [23] [24] and for the study of functional aspects at the molecular scale [25], [26].Speedy and accurate algorithms are needed for the intelligent reduction and representation of multidimensional datasets if these applications are to be realized.

Fig. 2 .
Fig. 2. (A) AB plot showing simulated dataset.AB clouds of different colour represent data obtained at different excitation wavelengths.Cloud centroids are indicated by yellow dots.(B) Graphical global analysis showing linear fit to signals from all pixels and excitation wavelengths.The global lifetimes are described by the intersection of the linear fit with the red monoexponential circle.(C) Fractional intensity images calculated using Eq.(9) and global lifetime estimates.

Figure 4 (
Figure 4(a) shows the original molecular fraction across the image.The estimated, uncorrected molecular fraction f1 (p) is shown in Fig. 4(b) and is calculated using Eq.(11).The corrected fraction fcorrected (p) is shown in Fig.4(c) and the mean squared difference error between true molecular fraction and estimated corrected molecular fraction is shown in Fig.4(d).There is no spatial structure discernible in Fig.4(d)proving that a single proportionality factor accounts for the difference between the original and estimated molecular fraction and no spatial bias is evident.

Fig. 5 .
Fig.5.Figure-of-Merit (FoM) for unmixing algorithm at a range of spectral -lifetime separations.The lifetime of Fluorophore 1 is held at 2 ns and peak excitation at 500 nm.The lifetime of Fluorophore 2 is then varied from 1.2 ns to 2.2 ns and spectral peak is varied from 500 nm to 540 nm.The insets show the estimators for the unmixed images ( N1 and N2 ), ratios to the true images (n 1 and n 2 ) and the associated excitation spectra and lifetime distributions.See text for details.

Fig. 6 .
Fig. 6.AB plot showing biological sample FLIM images taken at several excitation wavelengths (mean of each AB cloud is shown by a yellow dot) and (B) Intensity images at each excitation wavelength.Images are approximately 150 μm × 150 μm Figure 7(b) shows N1 , the estimator for the distribution of Fluorophore 1.As expected from the behaviour of the laminin IgG antibody (to which the Alexa 546 dye is conjugated), the largest signal contributions originate from the nucleus and nuclear envelope of the cells.Likewise the estimated distribution for Fluorophore 2 is mostly confined to the cytoplasm, which is also expected as Fluorophore 2's antibody (Map 2 IgG, to which the Alexa 555 dye is conjugated) targets the cytoplasm.These distributions can be combined to form the estimated molec-