Fluorescence Lifetime Imaging Microscopy ( FLIM ) Data Analysis with TIMP

Fluorescence Lifetime Imaging Microscopy (FLIM) allows fluorescence lifetime images of biological objects to be collected at 250 nm spatial resolution and at (sub-)nanosecond temporal resolution. Often ncomp kinetic processes underlie the observed fluorescence at all locations, but the intensity of the fluorescence associated with each process varies per-location, i.e., per-pixel imaged. Then the statistical challenge is global analysis of the image: use of the fluorescence decay in time at all locations to estimate the ncomp lifetimes associated with the kinetic processes, as well as the amplitude of each kinetic process at each location. Given that typical FLIM images represent on the order of 10 timepoints and 10 locations, meeting this challenge is computationally intensive. Here the utility of the TIMP package for R to solve parameter estimation problems arising in FLIM image analysis is demonstrated. Case studies on simulated and real data evidence the applicability of the partitioned variable projection algorithm implemented in TIMP to the problem domain, and showcase options included in the package for the visual validation of models for FLIM data.


Introduction
This paper describes the utility of the TIMP package for the R language and environment for statistical computing (R Development Core Team 2006) for the analysis of images collected by Fluorescence Lifetime Imaging Microscopy (FLIM) experiments.FLIM experiments allow the collection of fluorescence lifetime images of biological objects at 250 nm lateral resolution and at (sub-) nanosecond temporal resolution.FLIM has been widely applied in cell biology to detect interactions between fluorescently labeled biological molecules such as proteins, lipids, DNA and RNA.One experimental technique, which is sensitive at these small length scales of 1-10 nm, is the detection of Förster Resonance Energy Transfer (FRET).FRET is a bimolecular process in which the excited-state energy of a donor fluorophore is non-radiatively transferred to a ground-state acceptor molecule by dipole-dipole coupling.The FRET efficiency varies with the inverse 6th power of the distance between donor and acceptor and is usually negligible when the distance is larger than 10 nm.FRET is a fluorescence quenching process and can be identified by a shorter fluorescence lifetime of the donor.Intracellular proteins of interest can be genetically tagged with variants of the green fluorescent protein (GFP) (Tsien 1998).Spectral variants cyan fluorescent protein (CFP, donor) and yellow fluorescent protein (YFP, acceptor) have been proven a valuable FRET-pair.FRET as measured by FLIM can therefore be used as a "spectroscopic ruler" to map protein-protein interactions inside cells (for recent applications see Barber et al. (2005) and Grailhe et al. (2006)).Suhling et al. (2005) have comprehensively reviewed different FLIM methods, FLIM and FRET examples and other FLIM applications.
Well-designed data analysis techniques are required to process the measured lifetime images to unravel photophysical phenomena in complex molecular and cellular systems.It is often the case that the dynamics of the system are well-described by a model in which a small number of exponential decays with equal decay rates across all pixels underly the measured fluorescence, with amplitude parameters for the exponentials varying per-pixel.Then estimates for the amplitude parameters are conditionally linear on estimates for the decay rate parameters, allowing application of the variable projection algorithm (Golub and LeVeque 1979), which has been shown to have many desirable properties (Golub and Pereyra 2003;Mullen et al. 2007) for problems of this form.A disadvantage of the variable projection method that has prevented its application in this problem domain (Verveer et al. 2000) is that large memory resources are required.The TIMP package contains an implementation of a partitioned variable projection algorithm that returns the same results as the standard variable projection algorithm but requires much less memory (Mullen and van Stokkum 2007).The ability to apply the variable projection functional to estimation problems in the absence of large memory resources is a primary advantage of the application of TIMP in the FLIM image analysis problem domain.A further primary advantage is the support the package provides for visual interpretation and validation of the results of model fit.
The organization of the paper is as follows.Section 2 describes the sum-of-exponentials model that is often used to describe FLIM images and the optimization problem associated with fitting the parameters of this model.Section 3 discusses approaches to the parameter estimation task associated with the analysis of FLIM images, including the partitioned variable projection approach employed by TIMP.Section 4 describes extensions to TIMP implemented to support FLIM image analysis.Section 5 contains a simulation study in the application of the package to the analysis of datasets inspired by measured FLIM data.Section 6 describes the use of TIMP to fit a measured FLIM image.Section 7 contains conclusions.

Statistical model for FLIM data
FLIM images represent the decay of fluorescence in time at many different locations in the underlying system.Each location is represented by a pixel x, so that the image may be represented as a matrix Each column of Ψ represents a fluorescence decay in time at a given pixel x.
The decay of fluorescence data in time ψ(t) can often be satisfactorily modeled as a sum of n comp first-order kinetic processes convolved with an instrument response function (IRF) g(t), so that where c l represents the contribution to the data from process l in time t, a l represents the amplitude of decay l, and is the convolution operator.The model parameters to be fit are then the lifetimes τ l and their associated linear coefficients a l representative of intensity.
When the same kinetic processes underly the fluorescence at all n locations, Equation 2 can be applied globally to the image Ψ, so that where C is a matrix in which column l represents the time-profile of the lth kinetic process, and E is a matrix in which column l represents the intensity of kinetic process l across pixels.Then the parameter estimation task is global analysis: estimation of the n comp lifetimes τ associated with the image as a whole and the n comp amplitude parameters a l associated with each pixel (so that n * n comp amplitude parameters are estimated in total).Under least-squares criteria this is min This is an instance of the multi-exponential analysis problem, which is common in physics applications.Its difficulty is well-known, as Istratov and Vyvenko (1999) review.
Ψ represents the number of photons fluorescing from the location represented by pixel x at time t, and is therefore count data, the noise associated with which is assumed to be Poisson distributed.

Methods for FLIM data collection and analysis
FLIM data is collected by exciting a sample (such as a cell) to fluoresce (i.e., emit photons) using a laser pulses.After a laser pulse the time t until the first arrived photon is detected at pixel x is measured; The process of applying a laser pulse and recording the time of arrival of the first photon is repeated many times.The resulting histograms of arrival times represent fluorescence decays per pixel x.
Since samples of interest are in vivo, the power of the laser light used must be low if the sample is to remain alive throughout the measurement.The use of laser light of low power results in the arrival of few photons at location x.The requirement to keep the sample alive and in the same condition also means that the acquisition time cannot be long, (i.e., the sample cannot be subjected to too many laser pulses).
To improve signal-to-noise ratio (SNR), the time resolution may be decreased, since under a wider histogram channel more photons will be collected.The FLIM experiment is thus always compromising between time resolution and SNR.FLIM experiments that measure process with sub-nanosecond time resolution often have low SNRs (6-15 is typical).Further discussion of the methodology of FLIM experiments may be found in, e.g., Becker and Bergmann (2003).
Methods for the analysis of FLIM data commonly applied (Becker et al. 2001(Becker et al. , 2002) ) fit the model given by Equation 2 independently to each of the n pixels in an image.This yields estimates for the decay rates of the n comp kinetic process as well as n comp estimates of the amplitude of each process, for each pixel in the image.This is seldom desired, for the assumption that the same n comp kinetic processes underly measurements at all pixels is usually valid.
Data analysis methods that acknowledge the desirability of global analysis, which assumes that the underlying kinetic processes have the same lifetimes but different amplitudes across all pixels measured, as in Equation 3, often restrict themselves to the bi-exponential instance of the model (Barber et al. 2005;Pelet et al. 2000).Verveer et al. (2000) acknowledge that the global analysis associated with Problem 4 is a separable nonlinear least-squares problem that may be solved using the variable projection approach, though the authors state that for the large number of variables involved in fitting typical FLIM data, the memory requirements prohibit the approach.
TIMP allows an arbitrary number of exponentials may be fit to the data, though under experimentally realistic SNRs it is most often possible to resolve only one or two components.TIMP applies a partitioned variable algorithm to the global analysis problem associated with the analysis of FLIM images.This algorithm forms the residual vector prescribed by the variable projection functional without the need to store and operate on prohibitively large matrices, as is described in detail in Mullen and van Stokkum (2007).The present study of the application of TIMP to FLIM data is to the best of the authors' knowledge the first application of variable projection to FLIM data in the literature.

Extension of TIMP for FLIM data analysis
Several new capabilities were added to the TIMP package to facilitate the analysis of FLIM data.As described in Section 4.1, a file format was defined for the input of FLIM images into TIMP.A method for numerical convolution of an exponential decay with a measured IRF was added to the kinetic model options of the package, as Section 4.2 elaborates.New options to visually validate the results of fitting were also added, and are discussed in Section 4.3.

Data format
Given a FLIM image, it is often desirable to select those pixels associated with the subject of interest for modeling.For example, given a FLIM image of a cell, only those pixels interior to the cell wall are typically representative of the fluorescence decay of interest, and accordingly only these pixels are usually selected for modeling.Pre-processing dedicated to pixel selection is currently performed outside of TIMP.The indices of selected pixels are then included in the ASCII input file.
The format of the input file is as follows.
line 1: reserved for comments, not read line 2: reserved for comments, not read line 3: the character string "FLIM Image" line 4: dimension of image as x y (space-delimited) line 5: number of timepoints t in image line 6: number of pixels x selected for analysis line 7: vector of the times at which measurements were made line 8: index of selected pixel and decay trace from this pixel Ψ(, p) . . . . . .Line (8+number of selected pixels): the character string "Intensity map" Remaining lines: FLIM intensity image as matrix of dimension x × y Examples of this file format are included with this contribution.Files in this format may be read into R using the TIMP function readData.Section 4.3 defined the meaning of the intensity image; Section 6 contains an example of the use of the readData function.

Options for numerical convolution with a measured IRF
The analysis of FLIM data typically employs a measured IRF g(t) in fitting the exponential decay model contained in Equation 3. Evaluation of Equation 3 requires the numerical convolution of g(t) with an exponential decay.Methods to perform this convolution have been addressed in the literature at least since the seminal paper of Grinvald and Steinberg (1974), as Bajzer et al. (1995) discuss.For FLIM data (in which g(t) and the exponential decay are very often represented by 256 or less time points), we have found that methods based on a Fourier transformation are problematic, and that iterative methods give better results.Since an iterative method for the convolution of a vector and an exponential decay was not found by the authors implemented in R or in openly available scientific programming libraries, a method based on an iterative technique suggested in Grinvald and Steinberg (1974) was implemented in the shared C library used by TIMP, as the function Conv1.
To validate that the implementation of this iterative convolution technique returns an unbiased result, we considered its operation on the convolution of an IRF g(t) simulated as a Gaussian with location µ and full width half maximum ∆ parameters inspired by values occurring in FLIM experiments.The convolution of a Gaussian with an exponential decay is determined by the analytical expression where ∆ = ∆/(2 2 log(2)) and erf is the error function.Note that Equation 4.2 uses the decay rate k l (which is also the parameter estimated), whereas its reciprocal τ l = 1/k l is commonly reported.For times and decay rates k l inspired by values in measured data of interest, the results determined by the implementation of the iterative technique are unbiased as compared to results obtained using the analytical expression.

Model validation
Model validation in the FLIM image analysis application domain is ideally largely visual.The magnitude of residuals and fitted parameter estimated are possible to map per-pixel onto the modeled image as colors, allowing the results of fitting to be quickly evaluated.Several options for this display are implemented in TIMP.The analysis of a FLIM image with the TIMP function fitModel results in a multipanel summary plot as shown in Figure 1, whose components will be explained in turn.
First histograms of the estimated amplitudes associated with each component, with the corresponding global lifetime estimate on the bottom are displayed.In Figure 1 these are the two plots contained in row 1, columns 1 and 2. These plots allow for an impression of the absolute contributions of the components across all pixels of the image.
The following plots in the summary figure are n comp −1 histograms of the relative contribution p l of component l, where In Figure 1 this is the plot contained in row 1, column 3.These plots allow for an impression of the relative contribution of the component l across all pixels of the image.
A plot of the intensity image is then given.This intensity image includes those pixels not selected for modeling, and represents the number of photons per pixel measured over the course of all times t represented by the dataset.From this intensity image only some pixels are typically selected for modeling.The selected pixels are shown in the next plot in blue.
The intensity image and the intensity image with selected pixels in blue are contained in row 1, column 4 and row 2, column 1, respectively, in Figure 1.
The plot titled "< τ >" presents the average lifetime for each pixel from the selected region, where the average lifetime is given as (row 2, column 2 of Figure 1).The average lifetime may allow insight into the rate of energy transfer in processes on a per-pixel basis.The next n comp plots show normalized amplitudes in a color code mapped to the associated image, for each component l.In Figure 1 these are the plots contained in row 2, columns 3 and 4. The normalized amplitude plots allow insight into spatial patterns in the contribution of components.For example, these plots may allow identification of specific structures in a cell where the contribution of a given component is large.
Next the residuals associated with each pixel from the selected region are given as a color image, providing information on the quality of the fit both spatially and temporally, (row 3, column 1 of Figure 1).The first left singular vector of the residuals as results from a singular value decomposition (SVD) is plotted next (row 3, column 2 of Figure 1).This plot allows insight into structure in the residuals in time.For typical FLIM experiments, this structure is large around time 0, where the exponentially decaying components and the IRF contribute most.Structure in the left singular vectors after time 0 may be indicative of an inadequacy in the applied model.The next plot shows the first right singular vector associated with the SVD of the residuals mapped to the pixels selected for analysis, which provides information on the quality of the fit per pixel, and allows determination of whether the lack of fit is spatially structured.The last plot shows the singular values associated with the SVD of the data.The number of singular values that stand out in this plot indicate how many spatially and temporally independent components are present in the data.Further discussion of the use of the rank of the data in the estimation of the number of components can be found in e.g., Henry (1997).

A simulation study
A study of the application of TIMP to the analysis of simulated FLIM images was made in order to investigate the capabilities of the package in the problem domain.The study was designed in two parts.
The first part, described in Section 5.1, examines the ability of the software to estimate the lifetimes associated with bi-exponential decays in which the decay of fluorescence in time was measured over 64 and 256 timepoints (which we refer to as channels throughout).64 and 256 channel data is commonly collected in FLIM experiments, and thus was of particular interest.Simulation of bi-exponential decays was performed because Gratton et al. (2003) have shown that resolution of more than two components is not possible over this number of channels for experimentally realistic lifetime values and signal-to-noise ratios.
The second part of the simulation study, described in Section 5.2, smoothly varies the two amplitude parameters associated with bi-exponential decays across columns of the image for the purpose of examining whether the software is able to accurately estimate the relative contribution of the components.
Images Ψ(t, x) were simulated using Equation 2, as shown in Figure 2.Each pixel is associated with a decay in the time window 12.5 ns, over either 64 channels or 256 channels (equidistant in the interval 0-12.5 ns).The IRF g(t) was simulated as a Gaussian with mean 9 and 34 and standard deviation .4 and 1 for the 64 channel and 256 channel cases, respectively, in units of channels.Note that non-zero contribution of the IRF in both the 64 channel and the 256 channel case is represented by very few channels (3-8), as is commonly the case in FLIM experiments.Poisson noise was added to each decay trace ψ x (t) to obtain data of the desired signal-to-noise ratio (SNR) (using the R function rpois).The result may be considered as count data where Ψ(t, x) represents the number of photons collected at a given pixel x and time t, as in measured time-correlated single photon counting data (Maus et al. 2001).The SNRs of simulated images were chosen to reflect those commonly obtained in FLIM experiments.

A simulation study in the resolution of bi-exponential decays
This part of the simulation study examines the ability of TIMP to recover satisfactory estimates for the lifetimes underlying simulated images representative of two components.Images simulated with three pairs of lifetimes (in nanoseconds) collated in Table 1 were studied.For each pair of lifetimes studied, the relative contribution of the two components was varied between .1 and .9,so that 9 different images were simulated using each pair of lifetime values.The lifetime values are experimentally motivated (Borst et al. 2005).The images were simulated for both the SNR 8 and the SNR 15 case; the SNR 8 case is average for typical FLIM experiments, while the SNR 15 is higher than average.
A bi-exponential model was fit to the images, with the relative contribution of the two components being estimated as conditionally linear on values for the nonlinear lifetime estimates.The results are shown for images simulated with the pairs of lifetimes on row 1 and 2 of Table 1 in Figure 3.Note that each boxplot describes the variance in lifetime estimates as the relative amplitude of the components is varied.Our criteria for a satisfactory lifetime estimate is that the estimate is ±5% of the lifetime value used in simulation for data containing 256 channels, and within ±10% of the lifetime value used in simulation for data containing 64 channels.Under this criteria, the lifetime estimates obtained and shown graphically in Figure 3 are satisfactory.The small bias is attributed to using the number of photons at the leftmost point of each bin of times comprising a time-channel as representative of the average lifetime within the bin; because the data is exponentially decaying, there are always more photons to the left of the bin than to the right, and the average lifetime is thereby underestimated.The bias disappears when the number of channels is increased (for example, for data containing 1024 channels and the same SNR and lifetime values, it is insignificant).For the third pair of lifetimes studied, with τ 2 = .2ns and τ 2 = .5ns, it is impossible to determine satisfactory estimates even for data with SNR 15.The very short lifetimes are represented by only a few channels, so that there is not sufficient information.

5
Table 1: Parameter values in nanoseconds used in simulation of bi-exponential images.Instances of each group were simulated with contributions from the component with the longer lifetime τ 2 as 10%, 20%, 30%, 40%, 50%, 60%, 70%, 80% and 90% of the total intensity.Figure 3: Boxplots of lifetime estimates for each of two components given datasets simulated with lifetimes described in the first two rows of Table 1.Each boxplot is of lifetime estimates determined from fitting 9 different images, simulated with different relative amplitudes between the lifetimes.Lifetimes used in simulation are marked as red lines.Green represents results obtained on images in which the decay was represented by 256 channels, whereas blue represents results obtained on images in which the decay was represented by 64 channels.(a) shows results on images simulated with the lifetimes given in row 1 of Table 1 and (b) shows results on images simulated with the lifetimes given in row 2 of Table 1.
We found that for the lifetime values examined, for cases in which the contribution of one component was lower than 20% and the SNR was 8, lifetime estimates were not satisfactorily estimable.For SNR 15, lifetimes were not satisfactorily estimable for cases in which the contribution of one component was less than 10%.
This part of the simulation study was also repeated using an IRF measured on a FLIM set-up (as opposed to using a simulated IRF with a Gaussian distribution) to check that noise present in the IRF does not significantly decrease the accuracy of lifetime and amplitude estimates.
The obtained lifetime and amplitude estimates were very similar to those reported for the Gaussian IRF case, validating that the parameter estimation methodology is robust to an experimentally realistic amount of noise in the IRF.
We consider this part of the simulation study to demonstrate some limits of the resolvability of bi-exponential lifetimes on images inspired by measured data, and that for cases of practical interest TIMP lifetime estimates returned by TIMP are satisfactory.

A simulation study in the estimation of relative amplitudes of biexponential decays
A simulation study was made on instances of the image shown in Figure 4.The decay curve associated with each pixel is bi-exponential, with the two components having lifetimes of 0.6 and 2.5 ns respectively.The amplitude of the contribution a 1 from the first component varies from 0 to 1 across each column of the image, while the contribution a 2 from the second component varies from 1 to 0. Fitting a bi-exponential model to such images allows examination of whether the software is able to accurately estimate the two amplitude parameters a 1 and a 2 associated with each pixel.This part of the simulation study is inspired by a similar study by Pelet et al. (2000).The size of each analyzed image was 64 × 64 pixels (4096 pixels).the interval 0-12.5 ns (this is the 48 ps/channel case described in Section 5.1).Images were simulated with both SNR 8 and SNR 25.
Table 2 shows that the lifetime estimates well-approximate the values used in simulation of the images in both the SNR 25 and the SNR 8 case.The deviations from the values of the amplitudes a l used in simulation are small and unbiased, as shown in Figure 6 (c) graphically.Furthermore the lifetime estimates collated in Table 2 are also satisfactory.We conclude that this part of the simulation study demonstrates the ability of TIMP to return satisfactory estimates of the amplitude parameters a l determining the relative contribution of components.

Case study on measured CFP data
We were interested in investigating the capabilities of TIMP for FLIM image analysis of measured data.In cell biology studies FRET-FLIM is often used to demonstrate molecular interactions in vivo.For this purpose the fluorescent proteins cyan fluorescent protein (CFP) and yellow fluorescent protein (YFP) are the most widely used as donor-acceptor FRET pairs (Grailhe et al. 2006).However, the fluorescent decay of CFP is bi-exponential, making quantitative analysis of an interacting FRET population challenging (Russinova et al. 2004;Peter et al. 2005).
Time-correlated single photon counting experiments with a very high SNR (unattainable in FLIM experiments) described in Borst et al. (2005) established the lifetimes of CFP in a solution.We performed an experiment to collect FLIM images of the same sample in a microcapillary, using the experimental set-up described in Borst et al. (2003).Note that FLIM images of proteins in solution are not usually measured (the study of protein conformational dynamics in situ being the goal of most FLIM experiments), but that this experiment offers an opportunity to validate the ability of the software to estimate the lifetimes associated with the fluorescence decay of this important donor.
The SNR of the FLIM experiment was approximately 9. The time resolution was 48 ps/channel (over 256 channels).The fluorescence intensity image and region selected for analysis are shown in Figure 7.
To convey how the package is used to analyze a FLIM image, we describe all commands used to perform this part of the study.

Reading FLIM data into TIMP and preprocessing
The package is loaded.

R> library("TIMP")
Data is read into R using the readData function of TIMP.
R> cfp_data_sel <-preProcess(serT, sel_time=c(33,230)) A measured IRF is then read in and the same time points as selected in the data are chosen.
R> mea_IRF <-scan ("xtetoh_256_060822-bg_int.txt")[33:230]6.2.Initial model for CFP in solution: Mono-exponential decay The first model applied is based on a mono-exponential decay.The starting value for the decay rate given as 0.3, and is constrained positivity.The model is specified using the TIMP function initModel.

Fitting and validation of initial mono-exponential model
The TIMP function fitModel is used to fit the mono-exponential model.R> mono_result <-fitModel(list(cfp_data_sel), + mono_cfp_model, + opt=list(iter=20, linrange = 20, + makeps ="cfp_mono", + notraces = TRUE, xlabel = "time (ns)", + ylabel = "pixel number", FLIM=list())) The plot of the residuals returned is shown in Figure 8.The image plot of the residuals in the upper left hand corner show that there is a pattern of misfit around time 3.5 ns.This pattern of misfit is also indicated in the large upward trend of the left singular vector of the residuals shown in the upper right plot of Figure 8, which peaks at 3.5 ns.The root mean square (RMS) error associated with the fit is 5.2.We conclude that a mono-exponential decay model for CFP is not sufficient.

Refined model for CFP in solution: Bi-exponential decay
Based on the inadequacy of the fit of the mono-exponential model as evidenced by analysis of the residuals, the initModel function was used to specify a bi-exponential model for the measured CFP image.The first right singular vector of an SVD of the residuals mapped to the associated pixels on the intensity image shows that the residuals remain homogeneous in space.The RMS error associated with this fit is 4.9, less than for the mono-exponential model fit.

Fitting and validation of initial bi-exponential model
R> bi_result <-fitModel(list(cfp_data), bi_cfp_model, + opt=list(iter=20, linrange = 20, + makeps ="cfp_bi", + notraces = TRUE, + xlabel = "time (ns)", + ylabel = "pixel number", FLIM=list())) An image plot of the residuals under the bi-exponential model fit shows less structure around time 3.5 ns as compared to the same plot for the mono-exponential fit in Figure 8. Table 3: Parameters estimates obtained using TIMP on a measured CFP dataset analyzed with a bi-exponential model, and values in the literature for a dataset collected under similar experimental conditions analyzed using the same bi-exponential model.Note that variability in the experimental set-up, laser power and sample preparation limit the degree to which the results are directly comparable.left singular vector resulting from an SVD of the residuals also shows less structure around 3.5 ns.Note that we are not concerned about the structure in the SVD around 0 ns because misfit at this time results from the large contribution of the IRF and the peak in the amplitude of components at the start of their decay at this timepoint.The first right singular vector of an SVD of the residuals mapped to the associated pixels on the intensity image shows that the residuals remain homogeneous in space.Furthermore, the RMS square error has decreased to 4.9 from the RMS error of 5.2 under the fit of the monoexponential model.
The lifetime estimates under the bi-exponential model agree well with values published in Borst et al. (2005) for analysis of a dataset collected under similar experimental conditions, as tabulated in Table 3. Figure 10 (a) shows that the estimate for the average lifetime per pixel over the course of the decay (as determined with Equation 7) has no spatial structure, as is expected since the measured image represents a homogenous solution.Figure 10 (b) and (c) show that the normalized amplitudes of the components are also spatially homogenous, also as expected from the homogeneity of the solution.

Conclusions
A feasibility study has been made to investigate the use of the TIMP package of R for the analysis of FLIM data.In the course of the study new options for the fitting and validation of FLIM images with the package were developed.
In a simulation study the package was shown to return satisfactory estimates of both lifetime and amplitude parameters, the latter of which are estimated as conditionally linear parameters.On a real dataset it was possible to resolve the contributions of two components known to exist in terms of lifetime and amplitude estimates known from the literature, which further confirms the applicability of the partitioned variable projection fitting algorithm that TIMP implements to modeling FLIM images.
Future work will apply TIMP to the analysis of further experimentally collected FLIM data.Energy transfer between components will be modeled using the compartmental modeling options for TIMP described in Mullen and van Stokkum (2007).Implementation of a graphical user interface (GUI) to facilitate interactive model validation is also planned, along with a study to benchmark and optimize the package for speed on problems in FLIM analysis.

Figure 1 :
Figure 1: An example multipanel summary plot of residuals and fit of a bi-exponential model for measured FLIM data.Individual plots are explained in the text of Section 4.3.The image is taken from a fixed BHK (baby hamster kidney) cell with CFP expressed.

Figure 2 :
Figure 2: (Left) Intensity image of simulated data comprised of 1600 pixels in a 40 × 40 pixel arrangement, where intensity means the total photons summed over all channels.(Right) A fluorescence decay trace over 256 channels in the interval 0-12.5 ns is associated with each of the 1600 pixels comprising an image.

Figure 4 :Figure 5 :
Figure 4: A 64 × 64 pixel simulated image at one timepoint, in which the relative contribution of the first component increases linearly from 0 to 1 and the contribution of the second component decreases linearly from 1 to 0 along each column of the image.Each simulated dataset is comprised of 256 such images, representing the 256 timepoints (channels) simulated.
Figure 7: (a) Intensity image of a measured image of CFP in solution, where color represents the number of photons detected in a given pixel (b) Intensity image with pixels selected for analysis in blue.

Figure 8 :
Figure 8: An image plot of the residuals under the mono-exponential model fit shows structure in time before 3.5 ns.The first left singular vector resulting from an SVD of the residuals also shows this structure.The first right singular vector of an SVD of the residuals mapped to the associated pixels on the intensity image shows the residuals are relatively homogeneous in space.The RMS error associated with this fit is 5.2.

Figure 9 :
Figure9: An image plot of the residuals under the bi-exponential model fit shows less structure in time as compared to the same plot for the mono-exponential fit in Figure8.The first left singular vector resulting from an SVD of the residuals also shows less structure.The first right singular vector of an SVD of the residuals mapped to the associated pixels on the intensity image shows that the residuals remain homogeneous in space.The RMS error associated with this fit is 4.9, less than for the mono-exponential model fit.
Figure 10: (a) Distributions of the average lifetimes per location determined with Equation 7. Normalized amplitudes for component 1 (b) and component 2 (c) as a color on the associated image.

Figure 11 :
Figure 11: Histograms of amplitudes of components (A,B) and normalized amplitude of component 2 calculated with Equation 6 (C)

Table 2 :
The decay of the intensity at each pixel was represented by 256 timepoints equidistant in Lifetime values in nanoseconds used in simulation and estimated lifetimes for simulated images with smoothly varying contributions from two components.