Analysis of biexponential decay signals in the analog mean-delay fluorescence lifetime measurement method

Analysis of biexponential decay signals in the analog mean-delay fluorescence lifetime measurement method

Optics Communications 443 (2019) 136–143 Contents lists available at ScienceDirect Optics Communications journal homepage: www.elsevier.com/locate/o...

2MB Sizes 0 Downloads 43 Views

Optics Communications 443 (2019) 136–143

Contents lists available at ScienceDirect

Optics Communications journal homepage: www.elsevier.com/locate/optcom

Analysis of biexponential decay signals in the analog mean-delay fluorescence lifetime measurement method Wonsang Hwang a , Dongeun Kim a , Seungrag Lee b , Young Jae Won b , Sucbei Moon c , Dug Young Kim a ,∗ a

Department of Physics, Yonsei University, 50 Yonsei-ro, Seodaemun-gu, Seoul 120-749, South Korea Medical Device Development Center, Osong Medical Innovation Foundation, Cheongju, Chungbuk 361-951, South Korea c Department of Physics, Kookmin University, Seoul 136-702, South Korea b

ABSTRACT Signals obtained in fluorescence lifetime imaging microscopy (FLIM) often contain biexponential decay components, and unmixing of them is a critical issue in FLIM especially for fluorescence resonance-energy-transfer (FRET) applications. Even though the analog mean-delay (AMD) method is the fastest fluorescence lifetime measurement scheme for point-scanning confocal microscopy, its ability to analyze biexponential decay signals has not been demonstrated, yet. In this paper we have demonstrated that biexponential decay signals can be effectively analyzed with the AMD method. Effects of various measurement conditions on the accuracy of biexponential decay signal analysis are presented with numerical simulations and experimental demonstrations.

1. Introduction Fluorescence lifetime imaging microscopy (FLIM) has received considerable attention recently because of its functional imaging capability of monitoring physical, chemical, or environmental conditions near fluorescence probes. Many powerful applications of FLIM have been demonstrated in biomedical science and engineering for the measurement of various functional parameters, e.g., the inter-dot distance in a quantum-dot material [1], the existence of defects in materials [2,3], the pH distribution in a cell or tissue [4–6], the temperature of liquids in microfluidic channels [7,8], the concentration of specific ions in a cell [9–11], and the early detection of cancer [12–15]. In most practical cases, the measured fluorescence lifetime is not a mono-exponential decay curve but a weighted sum of biexponential decay curves. For example, when multiple fluorescence probes are used to monitor the distributions of many target molecules simultaneously in fluorescence microscopy, resolving the composition of mixed fluorophores is important. The unmixing of multiple fluorescence probes in the spectral domain is often limited by their broad spectra and overlap in the spectral domain. FLIM is an alternative tool for extracting the concentration ratio of multiple probes [16,17] that can be used if all the lifetimes of used fluorophores are known. Another popular application of FLIM is observing protein–protein interactions. Because contact-based protein–protein interactions control most functions in a cell, monitoring protein–protein interactions is very important in many biomedical sciences. Thus far, Förster resonance energy transfer (FRET) has been the only tool for identifying the physical contact of two proteins. If we use a proper donor–acceptor pair, the fluorescence lifetime of the donor is drastically reduced when ∗

the distance between the donor and acceptor is less than a critical value, which is on the order of 10 nm [18]. Therefore, FRET between a donor and an acceptor can be verified unambiguously by analyzing FLIM images for the donor. This FLIM-FRET scheme is known to be the most accurate method of identifying protein–protein interactions [19– 21]. A biexponential fluorescence decay model is normally used to analyze FLIM-FRET data, where a measured lifetime at a given image pixel is assumed to be composed of two different lifetime components: one with FRET and the other without FRET [21,22]. Various techniques and algorithms for deconvolving dual lifetimes and their composition from measured FLIM data have been proposed. A point-scanning confocal microscope with the time-correlated single photon counting (TC-SPC) lifetime-measurement scheme is known to be the most effective system for measuring fluorescence lifetimes of a few nanoseconds [23]. Nonlinear least-squares fitting [24,25] and maximum likelihood estimation [26,27] algorithms are widely used to fit the parameters of the biexponential decay data measured by TC-SPC-based FLIM systems. Other algorithms employing Fourier [28] and Laplace transforms [29] have been used to analyze multidimensional fluorescence lifetime data. The fluorescence lifetimes can be determined in the frequency domain by modulating an excitation light source sinusoidally and measuring the phase and the amplitude of the emitted fluorescence light [30–32]. The main advantages of frequency-domain FLIM (FD-FLIM) are the simple measurement systems and possibly fast measurement speeds. Biexponential decay components and their concentration ratio can be extracted via multi-frequency FLIM [33–35]. It was recently demonstrated that single-frequency FD-FLIM can be used to deconvolve biexponential decay lifetime components by using the

Corresponding author. E-mail address: [email protected] (D.Y. Kim).

https://doi.org/10.1016/j.optcom.2019.02.059 Received 31 October 2018; Received in revised form 27 December 2018; Accepted 24 February 2019 Available online 5 March 2019 0030-4018/© 2019 Elsevier B.V. All rights reserved.

W. Hwang, D. Kim, S. Lee et al.

Optics Communications 443 (2019) 136–143

Then, we calculate the summation of the squared errors between the measured data and the modeled data using the following equation:

polar plot representation of FD-FLIM data [32,36]. The photobleaching of fluorophores due to excessive exposure of a sample to excitation light and low photon economy are two major disadvantages of FD-FLIM. We previously reported another time-domain FLIM system based on a high-speed lifetime-measurement scheme called the analog meandelay (AMD) method [11,37–39]. In this scheme, fluorophores are excited with a picosecond laser pulse, and the time-domain waveform of the emitted fluorescence light intensity is measured using slow electronic devices whose 3-dB bandwidths are on the order of a few hundreds of megahertz. A fluorescence lifetime is obtained by calculating the expectation value of photon arrival time, treating the measured fluorescence intensity data as the probability density function of photon arrival in the time domain [38]. The definition of the mean-delay in the AMD method is same as that of the center of mass in the center of mass method (CMM), which is the expectation value of time with respect to a measured probability function of photon distribution [40– 42]. The AMD method can be considered as an analog version of the CMM algorithm, which is used to analyze time-domain FLIM data [40– 42]. Unlike the TC-SPC detection scheme, the AMD method can detect multiple fluorescence photons simultaneously for each excitation laser pulse [38,39]. We demonstrated that the measurement speed of the AMD-based FLIM can be 100 times higher than that of the conventional TC-SPC based FLIM, while its photon economy is comparable to that of the TC-SPC based FLIM [11]. Thus far, the major drawback of the AMD scheme is its inability to measure biexponential decay curves [23]. In this study, we demonstrate that biexponential decay signals can be effectively measured by the AMD method. We generate many biexponential decay data sets with various measurement conditions of the AMD method by using Monte-Carlo simulation. Two lifetime components and their mixing ratio are extracted from the computer-generated data sets by using an iterative least-squares curve-fitting algorithm. We investigate the effects of the sampling rate and the number of digitizing bits of the AMD lifetime-measurement system on the accuracy of the proposed deconvolution method. We show that the accuracy of the extracted parameters increases with the number of digitizing bits of the AMD method. Remarkably, our simulation results indicate that there exists a critical sampling frequency under which successful determination of two lifetime components is impossible. We investigate the detailed relationship between the critical sampling frequency and the experimental conditions of the AMD method.

n ( ) ∑ ( ) 𝜒 2 𝜏1 , 𝜏2 , 𝛼 = [ie tk − imodel (tk ; 𝜏1 , 𝜏2 , 𝛼)]2 ,

(4)

k=1

where n is the total number of measured data points for ie (t). We define a three-element vector X whose components are the three parameters of the model function f(t): ⎡𝜏1 ⎤ X = ⎢𝜏2 ⎥ . ⎢ ⎥ ⎣𝛼⎦

(5)

We define another vector R with n elements: ⎤ ⎡ ( ) ⎡r1 (t1 ; X)⎤ ⎢ie t1 − imodel (t1 ; X)⎥ ⎥ ⎢ ⎥ ⎢ ( ) r (t ; X)⎥ ⎢ie t2 − imodel (t2 ; X)⎥ . R (X) = ⎢ 2 2 = ⎥ ⎢ ⋮ ⎥ ⎢ ⋮ ⎢r (t ; X)⎥ ⎢ ⎥ ⎣n n ⎦ ⎢ ( ) ⎥ ⎣ie tn − imodel (tn ; X)⎦

(6)

The summation of squared errors 𝜒 2 can be written in terms of R(X) as 𝜒 2 (X) = |R(X)|2 = RT R. If we consider X and R(X) as the input and the output of a system, the Jacobian of the system becomes ⎡ 𝜕r1 ⎢ 𝜕𝜏1 ⎢ 𝜕r2 ⎢ Jr (X) = ⎢ 𝜕𝜏1 ⎢ ⋮ ⎢ 𝜕rn ⎢ ⎣ 𝜕𝜏1

𝜕r1 𝜕𝜏2 𝜕r2 𝜕𝜏2 ⋮ 𝜕rn 𝜕𝜏2

𝜕r1 ⎤ 𝜕𝛼 ⎥ 𝜕r2 ⎥ ⎥ 𝜕𝛼 ⎥ . ⋮ ⎥ 𝜕rn ⎥ ⎥ 𝜕𝛼 ⎦

𝜒 2 decreases as the modeled function imodel (t) becomes closer to the measured intensity data ie (t). We can find the best fitting parameters X that make 𝜒 2 the local minimum by using an iterative process known as the Levenberg–Marquardt algorithm [44,45]: [ ( )]−1 T ( K ) XK+1 = XK − Jr T Jr + 𝜇K diag Jr T Jr Jr R X , (7) where XK is the initial value of X, and XK+1 is the next best estimation of X in the iterative Levenberg–Marquardt algorithm. For a given measured data set ie (t) of length n, we start with an arbitrary initial vector X1 to obtain the next optimum vector X2 by using Eq. (7). The iterative index K is 1 in this case. We repeat the iterative calculation of Eq. (7) until 𝜒 2 = RT R becomes smaller than a critical value that we define. 𝜇K is a positive damping factor that controls the convergence speed of the iteration process. It can be varied in each iteration process. We apply this iterative method to find the optimum values of the biexponential decay model function for a given fluorescence intensity data set: two lifetimes 𝜏1 and 𝜏2 , along with their concentration ratio 𝛼.

2. Nonlinear least-squares curve-fitting algorithm When the impulse response function (IRF) of an AMD measurement system and the time-domain fluorescence intensity function emitted from a fluorophore are Iirf (t) and D(t), respectively, the measured electronic fluorescence intensity data ie (t) can be represented by the convolution of these two functions: ∞

ie (t) = Iirf (t) ⊗ D(t) ≡

∫−∞

Iirf (t − x) D (x) dx.

(1)

3. AMD Waveform data generation via Monte-Carlo simulation

The IRF of an AMD measurement system, Iirf (t), can be measured with a short impulse laser input. To analyze the fluorescence intensity function D(t), we use a nonlinear least-squares algorithm, which is widely used for the analysis of biexponential decay signals [24,43,44]. D(t) is modeled with a function f(t) that is composed of two exponentially decay components. When 𝜏1 and 𝜏2 are the time constants of the two exponential decay functions and 𝛼 is their mixing ratio, we can write f(t; 𝜏1 , 𝜏2, 𝛼) as follows: ( ) − t − t f t; 𝜏1 , 𝜏2 , 𝛼 = 𝛼e 𝜏1 + (1 − 𝛼)e 𝜏2 .

To demonstrate the feasibility of our proposed method for deconvolving biexponential decay lifetime components, we generated various biexponential decay intensity data sets via Monte-Carlo simulation and used them to test our proposed deconvolution method based on the least-squares curve-fitting algorithm. One of the most distinctive properties of the AMD method is its slow electronic detection system, whose bandwidth is on the order of nanoseconds, while the resolution of the lifetime is on the order of picoseconds [38,39]. We generated noisy biexponential decay fluorescence intensity data by using Monte-Carlo simulation. It is well-known that fluorescence intensities are limited by the shot noise of the photons. Because the shot noise can be reduced by increasing the number of detected photons in any photon-detection system, we neglected the shot noise of the photons. We investigated the effects of the effective number of bits (ENOB) and the sampling

(2)

For the three given parameters 𝜏1 , 𝜏2 , and 𝛼, we can generate model electronic fluorescence intensity data imodel (t) by convolving f(t; 𝜏1 , 𝜏2, 𝛼) with a measured IRF: ( ) imodel (t) = Iirf (t) ⊗ f t; 𝜏1 , 𝜏2 , 𝛼 . (3) 137

W. Hwang, D. Kim, S. Lee et al.

Optics Communications 443 (2019) 136–143

frequency of the analog-to-digital converter (ADC) used in the AMD method on the accuracy of our deconvolution method. Our proposed deconvolution method is tested with computer generated data with the following procedures. First, the lifetimes of the biexponential decay intensity curve 𝜏1 , 𝜏2 and their mixing ratio 𝛼 are set to for practical cases. Second, numerically generated data were prepared with a Monte-Carlo simulation technique by using Eq. (1) for a given sampling frequency and a fixed ENOB of an ADC. Third, the generated data are fit with the Levenberg–Marquardt algorithm to extract 𝜏1 , 𝜏2 , and 𝛼. The extracted values are compared to the original parameters for various ENOB and sampling frequencies of an ADC used to generate numerical data. Numerical data sets were generated using Eq. (1). The analytic biexponential decay function of f(t; 𝜏1 , 𝜏2 , 𝛼) given in Eq. (2) was used as the fluorescence intensity function D(t). The IRF of an AMD measurement system Iirf (t) is determined by the Gaussian low-pass filter (GLPF) used in the AMD system. The GLPF was designed to have 60-dB attenuation at the Nyquist frequency of a sampled data set [39]. If the sampling frequency of an ADC is fs , the GLPF can be written √ 2 as g(f) = e−(f∕f0 ) , where f0 = fs ∕ 24 ⋅ ln (10). For example, if fs = 1 GHz, the transmission amplitude of the GLPF becomes –60 dB at f = 0.5 GHz, which is the Nyquist frequency of the sampled data. The IRF Iirf (t) in Eq. (1) is simply the inverse Fourier transform of g(f). For a given sampling frequency of an AMD system, we can generate an equally spaced sampled data set by taking the analytic values of Eq. (1). This data set is an ideal one without any quantization error. In the aforementioned case, the number of bits of the ADC used in the measurement is infinite. To obtain a quantized digital value, we must use a do loop for each data point in programming, and repeated do loops require extensive computational time. The effects of quantization noise induced by a digitizer in experiments are simulated by adding Uniform distribution noise to ideal data with a Monte Carlo method. We used the relation between the standard deviation of a Uniform distribution noise and the ENOB of a digitizer given in Eq. (8). Denoting the full-scale range of the ADC as FS and the ENOB of the digitizer as N, the root-mean-square error 𝜎 in quantized data can be expressed as [46] 1 𝐹𝑆 . (8) 𝜎= √ 𝑁 12 2

4. Effect of number of digitizing bits

For a given ENOB N, we added random noise having a uniform distribution with a standard deviation of 𝜎 given by Eq. (8) to take into account the quantization noise of the ADC. Since errors associated with the sampling rate and the ENOB of a digitizer do not decrease with averaging, we have investigated these two effects on the unmixing of biexponential components in the AMD method. From Eq. (8) we can easily obtain the signal-to-quantization-noise ratio (SQNR) of measured √ data: SQNR ≡ FS∕𝜎 = 2N 12.’’ Fig. 1 shows typical deconvolution results and their accuracies for repeated numerical data generations and fitting procedures for a given biexponential decay model function with 𝜏1 = 8 ns, 𝜏2 = 25 ns, and 𝛼 = 0.7. These parameters were set to match experimental results for a quantum-dot sample described in Section 7. We generated 1000 data sets. The generated data correspond to an AMD setup equipped with an excitation laser with a repetition rate of 10 MHz and an 8-bit ADC with a 250-MHz sampling frequency. Therefore, each generated data set consists of 25 data points and covers a 100-ns period. each data set was fitted with the aforementioned iterative nonlinear leastsquares Levenberg–Marquardt curve-fitting algorithm. Fig. 1 shows typical statistical properties of retrieved lifetime 𝜏1 calculated from simulated data sets. Fig. 1(a) and (b) show the calculated mean lifetime 𝜏1 and its standard deviation with respect to the number of generated data sets. When the number of generated data sets increases from 2 to 1000, the mean lifetime 𝜏1 converges to 8.02 ns, whereas its real value is 8.0 ns. The standard deviation of 𝜏1 converges to 0.14 ns as the number of generated data sets increases. Fig. 1(c) shows histograms of retrieved lifetime 𝜏1 when the number of generated data sets are 50, 100, 200, and 1000, respectively. We do not show the statistical properties of 𝜏2 , because they are similar to those of 𝜏1 .

5. Effect of sampling frequency

We have investigated the relationship between the accuracy of the deconvolved parameters of a biexponential decay signal and the resolution of the ADC used in the AMD lifetime-measurement system. We assumed that the quantization error corresponding to the ENOB of the ADC is the major limiting factor in our AMD measurement method. A biexponential decay fluorescence signal with 𝜏1 = 8 ns, 𝜏2 = 25 ns, and 𝛼 = 0.7, as indicated by Eq. (1), was considered. The ENOB of the simulated data N was varied from 7 to 12. Uniform random noises with a standard deviation of 𝜎 given by Eq. (8) were added to the generated waveforms. For a given ENOB, we generated 200 waveforms, each having 25 data points covering a 100-ns period. Each generated waveform of 25 data points was interpolated to a waveform of 10,000 data points before the nonlinear least-squares curve-fitting process. Interpolation of measured data helps to increase the accuracy of numerical integration in the AMD method. We made the interval between neighboring data points to be 0.01 ns before doing numerical integration by using a cubic spline interpolation method. The Levenberg–Marquardt curvefitting algorithm described in Section 2 was used to find the best fitting parameters 𝜏1 , 𝜏2 , and 𝛼 for a given waveform. The initial values of 𝜏1 , 𝜏2 , and 𝛼 for the Levenberg–Marquardt algorithm are set to be 5 ns, 10 ns, and 0.1, respectively. The number of iterations in the Levenberg–Marquardt algorithm was fixed as 500. Fig. 2 shows the statistics of the calculated parameters obtained from the generated waveforms. The top three figures present the changes in the average values of these parameters for various ENOBs of the ADC used in the AMD method. These graphs show that the average value does not change significantly with the increase of the ENOB. The bottom three figures clearly show a decrease in the standard deviations of the calculated parameters when the ENOB is changed. As the ENOB increases from 7 to 12, the standard deviations of 𝜏1 and 𝜏2 decrease from 270 to 12 ps and from 720 to 32 ps, respectively. The standard deviation of the calculated 𝛼 is reduced from 0.017 to 0.008. Our results show that the measurement errors in deconvolved parameters for a biexponential decay model function can be effectively reduced by increasing the ENOB of the ADC used in the AMD method.

We investigated the effect of the sampling frequency of the ADC used in the AMD method on the accuracy of the retrieved parameters of a biexponential decay fluorescence signal. For a given sampling frequency fs and the parameters of the biexponential decay function 𝜏1 , 𝜏2 , and 𝛼, we generated 200 waveforms containing the Uniform quantization noises of an ADC. Each generated waveform was fitted with the Levenberg–Marquardt curve-fitting algorithm to obtain 𝜏1 , 𝜏2 , and 𝛼. The statistical properties of these calculated parameters were obtained for a given sampling frequency. We repeated this process as the sampling frequency fs was varied from 150 to 700 MHz with a step size of 50 MHz. We used the Gaussian low-pass filter to prevent aliasing in the sampled waveform data. Fig. 3(a) shows the effect of the sampling frequency on the accuracy of the deconvolved lifetimes obtained via the AMD method. A biexponential decay fluorescence signal with 𝜏1 = 2 ns, 𝜏2 = 6 ns, and 𝛼 = 0.7 in Eq. (1) was considered. The ENOB of the simulated data N was fixed to 7.2. Uniform random noises with a standard deviation of 𝜎 determined by Eq. (8) were included in the generated waveforms. For a given sampling frequency, we generated 200 digital waveforms, each covering a 100-ns period. For example, when the sampling frequency is 150 MHz, each generated waveform has only 15 data points. Each generated waveform was interpolated to have 10,000 data points before the Levenberg–Marquardt curve-fitting algorithm was applied to find the best fitting parameters 𝜏1 , 𝜏2 , and 𝛼. The initial values of 𝜏1 , 𝜏2 , and 𝛼 for the Levenberg–Marquardt algorithm are set to be 5 ns, 10 ns, and 0.1, respectively. 138

W. Hwang, D. Kim, S. Lee et al.

Optics Communications 443 (2019) 136–143

Fig. 1. Statistical properties of the deconvolved lifetime 𝜏1 for various numbers of generated data sets. The real values of 𝜏1 and 𝜏2 are 8 and 25 ns, respectively. (a) Average value of the retrieved 𝜏1 and (b) its standard deviation with respect to the number of deconvolution trials. (c) Histogram of the retrieved 𝜏1 for various numbers of generated data sets. Each time-domain data set contains 25 data points. The simulated data are assumed to be digitized by an 8-bit ADC with a 250-MHz sampling frequency.

Fig. 2. Statistical properties of the deconvolved parameters (𝜏1 , 𝜏2 , 𝛼) of the biexponential decay function calculated from generated 200 data sets. The standard deviations 𝜎 of these parameters decrease with the increasing ENOB of the ADC used in the AMD measurement method.

Fig. 3 shows that the two deconvolved lifetimes 𝜏1 and 𝜏2 well match the actual lifetimes when the sampling frequency is greater than or equal to 300 MHz. The two lifetime components are not retrieved correctly when the sampling frequency is less than 300 MHz. When the sampling frequency is 250 MHz, we have 𝜏1 = 0.05 ns, 𝜏2 = 7.10 ns, and 𝛼 = 0.55. All of these retrieved parameters are far from the real values. When the sampling frequency is further reduced to 200 MHz, the two retrieved lifetimes exhibit very large standard deviations. The average value of the shorter lifetime (𝜏1 ) is 1.18 ns, whereas that of the longer

lifetime (𝜏2 ) is 5.60 ns. When the sampling frequency is 150 MHz, we have 𝜏1 = 0.15 ns and 𝜏2 = 5.20 ns, with 𝛼 = 0.40. Even though the standard deviations of 𝜏1 and 𝜏2 are small, the retrieved lifetimes differ significantly from the actual ones. The main idea of the AMD lifetime-measurement method is that the mean value of the measured lifetime is the summation of the meandelay of the IRF Iirf (t) and the mean-delay of the exponential decay fluorescence intensity function. The mean-delay of a single exponential decay function with a lifetime 𝜏 is simply 𝜏. In previous reports, we 139

W. Hwang, D. Kim, S. Lee et al.

Optics Communications 443 (2019) 136–143

Fig. 3. Deconvolved lifetimes with respect to the sampling frequency of the ADC used in the AMD method; (a) and (b) correspond to a double exponential decay signal and a single exponential decay signal, respectively.

Fig. 4. Pseudo color plot of the sum of squared errors with respect to the three fitting parameters 𝜏1 , 𝜏2 , and 𝛼. (a) Sum of squared errors for a digital waveform generated at a sampling frequency of 200 MHz. (b) Sum of squared errors for a digital waveform generated at a sampling frequency of 1 GHz.

showed that the mean-delay of a single exponential decay function can be measured accurately even with a very slow electronic system with a broad IRF [11,37–39]. In theory, the mean-delay of a single exponential decay function can be measured accurately in the AMD method regardless of the bandwidth of the measurement system. The meandelay of the double exponential decay function f(t; 𝜏1 , 𝜏2, 𝛼) is {𝛼𝜏1 + (1-𝛼)𝜏2 }. The results shown in Fig. 3 indicate that the Levenberg– Marquardt algorithm cannot deduce the two lifetime components and their concentration ratio correctly when the sampling frequency of the measurement system is lower than the critical sampling frequency of 300 MHz. When the sampling frequency is lower than this critical value, the sampled data are broadened too much to be resolved into two different lifetime components. For a given measured mean-delay {𝛼𝜏1 + (1-𝛼)𝜏2 }, there are infinite combinations of (𝜏1 , 𝜏2 , 𝛼) that yield the same mean-delay. This is why the least-squares curve-fitting algorithm has multiple solutions for the fitting parameters 𝜏1 , 𝜏2 , and 𝛼 when the sampling frequency is lower than the critical value. We repeated similar numerical simulations for a single exponential decay signal with 𝜏 = 1 ns. Fig. 3(b) shows that the retrieved lifetimes obtained with the least-squares fitting algorithm agree well with the real lifetimes regardless of the sampling frequency of the ADC for a single exponential decay signal. Fig. 4 shows the sum of squared errors with a pseudo color map near the optimum values of the three fitting parameters. Fig. 4(a) corresponds to a sampling frequency of fs = 200 MHz, indicating that the three parameters do not exhibit a clear converging minimum point for the sum of squared errors. As the iterative nonlinear least-squares algorithm cannot find a converging point in the three-dimensional space of (𝜏1 , 𝜏2 , 𝛼), it provides multiple solutions of 𝜏1 , 𝜏2 , and 𝛼 with the repeated iterative least-squares fitting process. Fig. 4(b) shows the sum of squared errors for a generated waveform with a sampling frequency of fs = 1 GHz. A clear converging point is observed, indicating that the Levenberg–Marquardt algorithm provides a unique solution for 𝜏1 , 𝜏2 , and 𝛼 in this case. Fig. 4(b) shows that the sum of squared errors changes slowly with respect to the variations of 𝛼. We have repeated

Fig. 5. Minimum sampling frequency required for a successful least-squares curve fitting. We have 𝜏1 = 2 ns, 𝜏2 = 𝜏1 + 𝛥𝜏. The ENOB was 12, and the period of the sampled waveform was 100 ns.

these simulations for different values of 𝛼 when the two lifetimes are fixed as 𝜏1 = 2 ns and 𝜏2 , = 6 ns. The minimum sampling frequency for a successful retrieval of 𝛼 becomes lowest when 𝛼 = 0.5, and it increases as 𝛼 is either increased or decreased from 0.5. To determine the required critical sampling frequency for a given difference of the two lifetimes of a biexponential decay function, we repeated simulations similar to those shown in Fig. 3 for various values of 𝜏2 . Fig. 5 shows the relationship between the minimum sampling frequency and the lifetime difference 𝛥𝜏 between 𝜏1 and 𝜏2 . Here, the shorter lifetime 𝜏1 was fixed as 2 ns, while the longer lifetime 𝜏2 (= 𝜏1 + 𝛥𝜏) was varied from 3 to 12 ns. The ENOB was fixed as 12. For a given lifetime difference 𝛥𝜏, we repeatedly generated an AMD waveform and calculated 𝜏1 , 𝜏1 , and 𝛼 from that waveform while increasing the sampling frequency from 150 MHz to 1 GHz until the calculated lifetimes and their concentration ratio were within 10% errors of their real values. We took the sampling frequency within which the errors 140

W. Hwang, D. Kim, S. Lee et al.

Optics Communications 443 (2019) 136–143 Table 1 Comparison of parameters obtained using the TC-SPC and the AMD methods for a CdSe/ZnS core/shell type quantum dot sample.

of 𝜏1 , 𝜏1 , and 𝛼 are within 10% of their original values as the required sampling frequency to obtain the lifetime resolution 𝛥𝜏. We repeated this process 30 times to obtain the 30 required sampling frequencies for a given 𝛥𝜏. Fig. 5 shows the mean required sampling frequencies and their standard deviations for various 𝛥𝜏 values. 𝛥𝜏 can be considered as the resolution of our deconvolution scheme of the AMD method for two mixed fluorescence probes. The required sampling frequency of the AMD system exhibits a large fluctuation around the mean sampling frequency of 800 MHz when 𝛥𝜏 is 1 ns. This indicates that a sampling frequency of at least 900 MHz (= mean frequency + standard deviation) is needed to separate the 2- and 3-ns lifetime components via the AMD method. As is shown in Fig. 4(b) our least-square fitting algorithm has the lowest sensitivity with respect to the variation of 𝛼. We have repeated simulations shown in Fig. 5 with different values of 𝛼 from 0.1 to 0.9 with a step size of 0.1. Since all of these graphs overlap with each other and show similar trends, we show for the case of 𝛼 = 0.7 in this figure.

TC-PSC AMD method Difference

𝜏1 [ns]

𝜏2 [ns]

𝛼

8.9 9.1 2.2%

25.7 25.3 1.6%

0.68 0.69 1.5%

Table 2 Mean fluorescence lifetimes and their standard deviations for acridine orange and HPTS measured with the AMD method.

Lifetime ( 𝜏 ) Std. ( 𝜎 )

Acridine orange

HPTS

2.27 ns 0.18 ns

5.85 ns 0.09 ns

in measured waveforms obtained by the AMD method, the least-square fitting algorithm finds good repeatable results. In order to compare our results with those of the conventional method, we have measured the biexponential fluorescence intensity signal of the same CdSe/ZnS core/shell type QDs by using a commercially available TC-SPC-based lifetime measurement system (PicoQaunt fluotime-300). The width of the measurement window and of the temporal resolution of the TC-SPC system were set to 350 ns and 100 ps, respectively. We have used the same number of 7 × 107 photons for the measurements. Table 1 shows measured parameters obtained with our AMD lifetime measurement system and a commercially available TC-SPC system. Differences in the table are percentage mismatches of the results measured by the AMD method out of those obtained by the TC-SPC method. The second biexponential decaying fluorescence sample we have tested was a mixture of acridine orange and 8-Hydroxypyrene-1,3,6trisulfonic acid (HPTS) in water. The fluorescence lifetimes of acridine orange and HPTS are measured with our AMD setup, and the results are shown in Table 2. Mean lifetimes (𝜏) and standard deviations (𝜎) are calculated from 20 repeated measurements. Because of differences in absorption cross section, the average photon numbers collected for each measurement of acridine orange was only about 520, while that of HPTS was about 14,100. This is the reason why the standard deviation of acridine orange is much large than that of HPTS. Mean lifetimes in Table 2 are used in determining the accuracy of our deconvolution scheme for a mixed sample of acridine orange and HPTS. Acridine orange and HPTS solutions was mixed with 50 to 1 ratio to make measured fluorescence powers of the two probes become compatible. The average fluorescence lifetime of the mixed sample was 4.13 ns. Based on the measured fluorescence lifetime of acridine orange (2.27 ns) and that of HPTS (5.85 ns) shown in Table 2, the concentration ratio 𝛼 is obtained to be 0.72. As we know the two lifetimes and their concentration ratio for the mixed sample of acridine orange and HPTS, we could test the effect of the measurement bandwidth of the AMD method on the accuracy of the deconvolution scheme we have proposed. A series of measurements were made for the mixed sample of acridine orange and HPTS with the AMD method while the bandwidth of the measurement setup was changed from 300 MHz to 750 MHz with a step size of 50 MHz. The bandwidth of the AMD system was varied by changing the sampling frequency of the data acquisition (DAQ) board (Gage Eon express) and the bandwidth of electronic components. The bandwidth of electronic devices was controlled by using homemade Gaussian low-pass filters. The 70-dB cut-off frequencies of low-pass filters were set to 1/5 of the sampling frequencies of the DAQ board. The effective number of bits of the DAQ board was measured to be 9.1. Fig. 7 shows the mean values and their standard deviations of the two deconvolved lifetimes and their concentration ratio as a function of the sampling frequency of the DAQ used in the AMD method. The black solid squares and red solid circles indicate the average lifetimes

6. Minimum resolvable lifetime difference Most of the fluorescence probes used in biomedical sciences have lifetimes between 1 and 6 ns. To test the ability of the proposed method to distinguish two similar lifetime components, we performed the following two simulations. First, we calculated the lifetimes 𝜏1 and 𝜏2 from computer-generated waveforms for the AMD method when the real value of 𝜏2 was fixed as 2.0 ns and the real value of 𝜏1 was varied from 1.0 to 1.9 ns with a step size of 0.1 ns. We generated 30 data sets for each real value of 𝜏1 , and 𝜏1 and 𝜏2 were calculated from these using our proposed fitting technique for the AMD method. The sampling frequency and the ENOB were fixed as 1 GHz and 12 bits, respectively. Fig. 6(a) shows the mean retrieved lifetimes 𝜏1 and 𝜏2 with respect to the real value of 𝜏1 . The separation of the two lifetime components is possible as long as the lifetime difference is larger than 0.2 ns. The fractional concentration 𝛼 is 0.7. Second, the lifetimes 𝜏1 and 𝜏2 were retrieved from the computer-generated waveforms when the real value of 𝜏2 was fixed as 6.0 ns, and the real value of 𝜏1 was scanned from 4.0 to 5.9 ns with a step size of 0.4 ns. Fig. 6(b) shows the average values of 30 retrieved lifetimes 𝜏1 and 𝜏2 with respect to the real value of 𝜏1 . The two lifetime components could be extracted well until the lifetime difference became smaller than 0.4 ns. These two graphs show that two lifetime components can be distinguished well by the AMD lifetime-measurement method within 10% error with the nonlinear least fitting algorithm if we have an AMD system with a 1-GHz sampling frequency and a 12-bit digitizing resolution. 7. Experimental results CdSe/ZnS core/shell type quantum dot (QD) is the first sample we have used to test the feasibility of analyzing a biexponential decay signal with the AMD method. It is well known that CdSe/ZnS core/shell quantum QDs have a biexponential emission intensity profile. We used 1 mg of the QD diluted with 1 ml of toluene at room temperature. The fluorescence signal from CdSe/ZnS core/shell type QDs were measured with our AMD setup using an 8-bit digitizer (NI PCI-5114) operating at a sampling frequency of 250 MHz. In order to reduce the effect of photon shot noise, we have averaged 107 waveforms, which contain approximately 7 × 107 photons. We used the Levenberg–Marquardt algorithm to extract two lifetime components and their concentration ratio. The average values of biexponential decay parameters for 20 repeated measurements are 𝜏1 = 9.1 ns, 𝜏2 = 25.3 ns, and 𝛼 = 0.69, where 𝜏1 and 𝜏2 are the short and the long lifetime components of a biexponential decay function, and 𝛼 is the concentration ratio of the short lifetime component in a measured confocal volume. The standard deviations of measured results are 𝜎 (𝜏1 ) = 0.22 ns, 𝜎 (𝜏2 ) = 0.63 ns, and 𝜎(𝛼) = 0.02, respectively. Even though it is not easy to notice the existence of two different exponentially decay functions 141

W. Hwang, D. Kim, S. Lee et al.

Optics Communications 443 (2019) 136–143

Fig. 6. Retrieved lifetimes 𝜏1 and 𝜏2 with respect to the real value of 𝜏1 . The real value of 𝜏2 is fixed as 2.0 ns, and 𝜏1 is varied from 1.0 to 1.9 ns with a step size of 0.1 ns. (b) Retrieved lifetimes 𝜏1 and 𝜏2 with respect to the real value of 𝜏1 when 𝜏2 is fixed as 6.0 ns. The fractional concentration 𝛼 is fixed as 0.7.

Fig. 7. (a) Retrieved lifetime components of a mixed sample of acridine orange (black) and HPTS (red), and (b) their concentration ratios with respect to the sampling frequency of a digitizer used in the AMD setup.

minimum separable lifetime difference of a biexponential decay signal is determined by the bandwidth of an AMD system. Simulation results are verified with two experimental measurements. CdSe/ZnS core/shell type QD was the first sample we have tested. Two lifetime components and their concentration ratio were obtained for the QD sample with the AMD method and the Levenberg– Marquardt least square curve fitting algorithm. Measured results are compared with those obtained with the conventional TC-SPC method. Another biexponential fluorescence sample we have measured was a mixture of acridine orange and HPTS in water. Two lifetime components and their concentration ratio of the mixed sample were obtained. Repeated experiments were performed while changing the measurement bandwidth of the setup to verify the bandwidth effect on the accuracy of AMD method in analyzing biexponential decay signals.

obtained from 30 measurements each of which contains approximately 4.4 × 105 photons for a given sampling frequency. The two blue dotted lines represent 10% error boundaries of the three retrieved parameters: 𝜏1 for acridine orange (2.27 ±0.227ns), 𝜏2 for HPTS (5.85 ±0.585ns), and 𝛼 (0.72 ± 0.072). As the sampling frequency is decreased, the standard deviations of the measured lifetimes do not increase significantly while the mean lifetimes retrieved from measured data deviate off from the actual lifetimes. We define the critical sampling frequency under which the mean retrieved lifetimes are more than 10% off from the real lifetimes. Fig. 7 shows that the critical sampling frequency for this case is approximately 400 MHz. These measured results are consistent with the simulation results shown Fig. 3. 8. Conclusions

Acknowledgments

The AMD measurement method is one of the fastest lifetimemeasurement schemes that can be applied in point-scanning confocal FLIM. Even though it has good photon economy and a high measurement speed, its practical usage in FLIM or FLIM-FRET has been limited due to its lack of ability to measure biexponential lifetime components. In this study, we have demonstrated that biexponential decay waveforms measured via the AMD system can be effectively deconvolved into two lifetime components by using the least-sq"-uares curve-fitting algorithm. First, we generated a series of time-domain biexponential decay waveforms for various experimental conditions and different samples by using Monte-Carlo simulations. Then, the generated waveforms were fitted with the Levenberg–Marquardt leastsquares curve-fitting algorithm to determine two different lifetime components and their concentration ratio. With these numerically generated data sets and the least-squares curve-fitting algorithm, we showed that the ENOB of the digitizer and the bandwidth of electronics used in the AMD method are two important parameters for successful analysis of biexponential decay signals. The ENOB of the AMD system determines the precision of the deconvolved lifetime components. The

National Research Foundation of Korea (NRF-2013R1A1A2062448) through the Basic Science Research Program, Ministry of Trade, Industry and Energy of Korea (10062417) through Technology Innovation Program. References [1] K. Tai, W. Lü, I. Umezu, A. Sugimura, Inter-dot distance dependence of photoluminescence properties in CdSe quantum dot systems, Appl. Phys. Express 3 (2010) 035202. [2] G. Liaugaudas, G. Davies, K. Suhling, R. Khan, D. Evans, Luminescence lifetimes of neutral nitrogen-vacancy centres in synthetic diamond containing nitrogen, J. Phys.: Condens. Matter 24 (2012) 435503. [3] D.G. Monticone, F. Quercioli, R. Mercatelli, S. Soria, S. Borini, T. Poli, M. Vannoni, E. Vittone, P. Olivero, Systematic study of defect-related quenching of NV luminescence in diamond with time-correlated single-photon counting spectroscopy, Phys. Rev. B 88 (2013) 155201. [4] T. Nakabayashi, H.-P. Wang, M. Kinjo, N. Ohta, Application of fluorescence lifetime imaging of enhanced green fluorescent protein to intracellular pH measurements, Photochem. & Photobiol. Sci. 7 (2008) 668–670. 142

W. Hwang, D. Kim, S. Lee et al.

Optics Communications 443 (2019) 136–143 [25] D.F. Eaton, Recommended methods for fluorescence decay analysis, Pure Appl. Chem. 62 (1990) 1631–1648. [26] M. Maus, M. Cotlet, J. Hofkens, T. Gensch, F.C. De Schryver, J. Schaffer, C. Seidel, An experimental comparison of the maximum likelihood estimation and nonlinear least-squares fluorescence lifetime analysis of single molecules, Anal. Chem. 73 (2001) 2078–2086. [27] J. Enderlein, P.M. Goodwin, A. Van Orden, W.P. Ambrose, R. Erdmann, R.A. Keller, A maximum likelihood estimator to distinguish single molecules by their fluorescence decays, Chem. Phys. Lett. 270 (1997) 464–470. [28] D. O’Connor, W. Ware, J. Andre, Deconvolution of fluorescence decay curves. A critical comparison of techniques, J. Phys. Chem. 83 (1979) 1333–1343. [29] A. Gafni, R.L. Modlin, L. Brand, Analysis of fluorescence decay curves by means of the Laplace transformation, Biophys. J. 15 (1975) 263–280. [30] R.D. Spencer, G. Weber, Measurements of subnanosecond fluorescence lifetimes with a cross-correlation phase fluorometer, Ann. New York Acad. Sci. 158 (1969) 361–376. [31] E. Gratton, M. Limkeman, A continuously variable frequency cross-correlation phase fluorometer with picosecond resolution, Biophys. J. 44 (1983) 315–324. [32] J.R. Lakowicz, B.P. Maliwal, Construction and performance of a variablefrequency phase-modulation fluorometer, Biophys. chem. 21 (1985) 61–78. [33] J.-C. Brochon, A.K. Livesey, J. Pouget, B. Valeur, Data analysis in frequencydomain fluorometry by the maximum entropy method—recovery of fluorescence lifetime distributions, Chem. Phys. Lett. 174 (1990) 517–522. [34] A. Esposito, H.C. Gerritsen, F.S. Wouters, Fluorescence lifetime heterogeneity resolution in the frequency domain by lifetime moments analysis, Biophys. J. 89 (2005) 4286–4299. [35] A. Squire, P.J. Verveer, P. Bastiaens, Multiple frequency fluorescence lifetime imaging microscopy, J. Microsc. 197 (2000) 136–149. [36] G.I. Redford, R.M. Clegg, Polar plot representation for frequency-domain analysis of fluorescence lifetimes, J. Fluoresc. 15 (2005) 805. [37] Y.J. Won, S. Moon, W.-T. Han, D.Y. Kim, Referencing techniques for the analog mean-delay method in fluorescence lifetime imaging, J. Opt. Soc. Amer. A 27 (2010) 2402–2410. [38] S. Moon, Y. Won, D.Y. Kim, Analog mean-delay method for high-speed fluorescence lifetime measurement, Optics Express 17 (2009) 2834–2849. [39] Y.J. Won, W.-T. Han, D.Y. Kim, Precision and accuracy of the analog mean-delay method for high-speed fluorescence lifetime measurement, J. Opt. Soc. Amer. A 28 (2011) 2026–2032. [40] D.-U. Li, B. Rae, R. Andrews, J. Arlt, R. Henderson, Hardware implementation algorithm and error analysis of high-speed fluorescence lifetime sensing systems using center-of-mass method, J. Biomed. Opt. 15 (2010) 017006-017006-017010. [41] D.D.-U. Li, S. Ameer-Beg, J. Arlt, D. Tyndall, R. Walker, D.R. Matthews, V. Visitkul, J. Richardson, R.K. Henderson, Time-domain fluorescence lifetime imaging techniques suitable for solid-state imaging sensor arrays, Sensors 12 (2012) 5650–5669. [42] S.P. Poland, A.T. Erdogan, N. Krstajić, J. Levitt, V. Devauges, R.J. Walker, D.D.U. Li, S.M. Ameer-Beg, R.K. Henderson, New high-speed centre of mass method incorporating background subtraction for accurate determination of fluorescence lifetime, Optics Express 24 (2016) 6899–6915. [43] W.R. Ware, L.J. Doemeny, T.L. Nemzek, Deconvolution of fluorescence and phosphorescence decay curves. Least-squares method, J. Phys. Chem. 77 (1973) 2038–2048. [44] A. Grinvald, I.Z. Steinberg, On the analysis of fluorescence decay kinetics by the method of least-squares, Anal. Biochem. 59 (1974) 583–598. [45] J.J. Moré, The Levenberg–Marquardt Algorithm: implementation and theory, in: Numerical Analysis, Springer, 1978, pp. 105–116. [46] J.J. Blair, T.E. Linnenbrink, Corrected rms error and effective number of bits for sine wave ADC tests, Comput. Stand. Interfaces 26 (2004) 43–49.

[5] C. Hille, M. Berg, L. Bressel, D. Munzke, P. Primus, H.-G. Löhmannsröben, C. Dosche, Time-domain fluorescence lifetime imaging for intracellular ph sensing in living tissues, Anal. Bioanal. Chem. 391 (2008) 1871. [6] A. Orte, J.M. Alvarez-Pez, M.J. Ruedas-Rama, Fluorescence lifetime imaging microscopy for the detection of intracellular pH with quantum dot nanosensors, ACS Nano 7 (2013) 6387–6395. [7] M.A. Bennet, P.R. Richardson, J. Arlt, A. McCarthy, G.S. Buller, A.C. Jones, Optically trapped microsensors for microfluidic temperature measurement by fluorescence lifetime imaging microscopy, Lab Chip 11 (2011) 3821–3828. [8] E.M. Graham, K. Iwai, S. Uchiyama, A.P. de Silva, S.W. Magennis, A.C. Jones, Quantitative mapping of aqueous microfluidic temperature with sub-degree resolution using fluorescence lifetime imaging microscopy, Lab Chip 10 (2010) 1267–1273. [9] A.V. Agronskaia, L. Tertoolen, H.C. Gerritsen, Fast fluorescence lifetime imaging of calcium in living cells, J. Biomed. Opt. 9 (2004) 1230–1237. [10] J. Lakowicz, H. Szmacinski, K. Nowaczyk, M. Johnson, Fluorescence lifetime imaging of calcium using Quin-2, Cell Calcium 13 (1992) 131–147. [11] Y. Won, S. Moon, W. Yang, D. Kim, W.T. Han, D.Y. Kim, High-speed confocal fluorescence lifetime imaging microscopy (FLIM) with the analog mean delay (AMD) method, Optics Express 19 (2011) 3396–3405. [12] D.K. Bird, L. Yan, K.M. Vrotsos, K.W. Eliceiri, E.M. Vaughan, P.J. Keely, J.G. White, N. Ramanujam, Metabolic mapping of MCF10A human breast cells via multiphoton fluorescence lifetime imaging of the coenzyme NADH, Cancer Res. 65 (2005) 8766–8773. [13] X. Gao, Y. Cui, R.M. Levenson, L.W. Chung, S. Nie, In vivo cancer targeting and imaging with semiconductor quantum dots, Nature Biotechnol. 22 (2004) 969. [14] P.J. Tadrous, J. Siegel, P.M. French, S. Shousha, E.N. Lalani, G.W. Stamp, Fluorescence lifetime imaging of unstained tissues: early results in human breast cancer, J. Pathol. 199 (2003) 309–317. [15] Y. Sun, J. Phipps, D.S. Elson, H. Stoy, S. Tinling, J. Meier, B. Poirier, F.S. Chuang, D.G. Farwell, L. Marcu, Fluorescence lifetime imaging microscopy: in vivo application to diagnosis of oral carcinoma, Opt. Lett. 34 (2009) 2081–2083. [16] R. Cao, P. Jenkins, W. Peria, B. Sands, M. Naivar, R. Brent, J.P. Houston, Phasor plotting with frequency-domain flow cytometry, Optics Express 24 (2016) 14596–14607. [17] R. Pepperkok, A. Squire, S. Geley, P.I. Bastiaens, Simultaneous detection of multiple green fluorescent proteins in live cells by fluorescence lifetime imaging microscopy, Curr. Biol. 9 (1999) 269–274. [18] A. Grinvald, E. Haas, I. Steinberg, Evaluation of the distribution of distances between energy donors and acceptors by fluorescence decay, Proc. Natl. Acad. Sci. 69 (1972) 2273–2277. [19] M. Tramier, M. Zahid, J.C. Mevel, M.J. Masse, M. Coppey-Moisan, Sensitivity of CFP/YFP and GFP/mCherry pairs to donor photobleaching on FRET determination by fluorescence lifetime imaging microscopy in living cells, Microsc. Res. Tech. 69 (2006) 933–939. [20] H. Wallrabe, A. Periasamy, Imaging protein molecules using FRET and FLIM microscopy, Curr. Opin. Biotechnol. 16 (2005) 19–27. [21] R. Duncan, A. Bergmann, M. Cousin, D. Apps, M. Shipston, Multi-dimensional time-correlated single photon counting (TCSPC) fluorescence lifetime imaging microscopy (FLIM) to detect FRET in cells, J. Microsc. 215 (2004) 1–12. [22] P.J. Verveer, A. Squire, P.I. Bastiaens, Frequency-domain fluorescence lifetime imaging microscopy: a window on the biochemical landscape of the cell, in: Methods in cellular imaging, Springer, 2001, pp. 273–294. [23] W. Becker, Fluorescence lifetime imaging–techniques and applications, J. Microsc. 247 (2012) 119–136. [24] N. Periasamy, Analysis of fluorescence decay by the nonlinear least squares method, Biophys. J. 54 (1988) 961–967.

143