Skip to main content
Advertisement
Browse Subject Areas
?

Click through the PLOS taxonomy to find articles in your field.

For more information about PLOS Subject Areas, click here.

  • Loading metrics

Blind deconvolution estimation by multi-exponential models and alternated least squares approximations: Free-form and sparse approach

  • Daniel U. Campos-Delgado ,

    Roles Conceptualization, Funding acquisition, Methodology, Supervision, Validation, Visualization, Writing – original draft, Writing – review & editing

    ducd@fciencias.uaslp.mx

    Affiliations Facultad de Ciencias, Universidad Autonoma de San Luis Potosi, San Luis Potosi, Mexico, Instituto de Investigacion en Comunicacion Optica, Universidad Autonoma de San Luis Potosi, San Luis Potosi, Mexico

  • Omar Gutierrez-Navarro,

    Roles Software, Validation, Writing – review & editing

    Affiliation Department of Biomedical Engineering, Universidad Autonoma de Aguascalientes, Aguascalientes, Mexico

  • Ricardo Salinas-Martinez,

    Roles Conceptualization, Investigation, Software, Writing – original draft, Writing – review & editing

    Affiliation Facultad de Ciencias, Universidad Autonoma de San Luis Potosi, San Luis Potosi, Mexico

  • Elvis Duran,

    Roles Data curation, Software

    Affiliation Department of Biomedical Engineering, Texas A&M University, College Station, Texas, United States of America

  • Aldo R. Mejia-Rodriguez,

    Roles Validation, Writing – review & editing

    Affiliation Facultad de Ciencias, Universidad Autonoma de San Luis Potosi, San Luis Potosi, Mexico

  • Miguel J. Velazquez-Duran,

    Roles Software, Visualization, Writing – review & editing

    Affiliation Facultad de Ciencias, Universidad Autonoma de San Luis Potosi, San Luis Potosi, Mexico

  • Javier A. Jo

    Roles Conceptualization, Funding acquisition, Supervision, Validation, Writing – review & editing

    Affiliation School of Electrical and Computer Engineering, University of Oklahoma, Norman, Oklahoma, United States of America

Abstract

The deconvolution process is a key step for quantitative evaluation of fluorescence lifetime imaging microscopy (FLIM) samples. By this process, the fluorescence impulse responses (FluoIRs) of the sample are decoupled from the instrument response (InstR). In blind deconvolution estimation (BDE), the FluoIRs and InstR are jointly extracted from a dataset with minimal a priori information. In this work, two BDE algorithms are introduced based on linear combinations of multi-exponential functions to model each FluoIR in the sample. For both schemes, the InstR is assumed with a free-form and a sparse structure. The local perspective of the BDE methodology assumes that the characteristic parameters of the exponential functions (time constants and scaling coefficients) are estimated based on a single spatial point of the dataset. On the other hand, the same exponential functions are used in the whole dataset in the global perspective, and just the scaling coefficients are updated for each spatial point. A least squares formulation is considered for both BDE algorithms. To overcome the nonlinear interaction in the decision variables, an alternating least squares (ALS) methodology iteratively solves both estimation problems based on non-negative and constrained optimizations. The validation stage considered first synthetic datasets at different noise types and levels, and a comparison with the standard deconvolution techniques with a multi-exponential model for FLIM measurements, as well as, with two BDE methodologies in the state of the art: Laguerre basis, and exponentials library. For the experimental evaluation, fluorescent dyes and oral tissue samples were considered. Our results show that local and global perspectives are consistent with the standard deconvolution techniques, and they reached the fastest convergence responses among the BDE algorithms with the best compromise in FluoIRs and InstR estimation errors.

Introduction

Fluorescence microscopy has become a powerful tool to characterize the chemical properties of tissue samples [14]. In fluorescence lifetime imaging microscopy (FLIM), a sample is excited by an electromagnetic source, usually an ultraviolet (UV) narrow laser-pulse. Given a sample with synthetic dyes or endogenous fluorophores, the sample will emit light at a longer wavelength than the excitation source [3, 4]. This fluorescent response contains information of the chemical environment in the sample, its fluorescent components and their concentrations. In multi-spectral FLIM (mFLIM) several spectral channels are recorded simultaneously, and the resulting time-responses are concatenated to construct a characteristic time-frequency signature [5, 6]. The mFLIM methodology has shown to be a powerful medical resource for early and non-invasive diagnosis for different pathologies such as cardiovascular and dermatological diseases [5, 79], oral pre-cancer and cancer conditions [10, 11], colonic dysplasia [12], skin cancer [13, 14], and to assess therapeutic responses of anticancer drugs [15].

The fluorescence response measured by FLIM can be modeled as the convolution between the instrument response (InstR) and the particular fluorescence impulse response (FluoIR) of the tissue sample. In order to identify the FluoIR of the sample and provide quantitative information of the FLIM data, a deconvolution stage needs to isolate the InstR from the fluorescence decay (FluoD) [1620]. There are different strategies to solve this inverse problem, usually the InstR is assumed known or measured a priori, and then carefully aligned with the FluoIRs to avoid bias in the estimations. Other strategies quantify FLIM data by analyzing the FluoDs with a linear unmixing approach [2125], or in a lower-dimensional domain using the phasor approach [2628].

In the literature, there are two principal tendencies for the deconvolution process for FLIM datasets: 1) multi-exponential models [1618], and 2) Laguerre-basis approximation [19, 20, 29]. The former assumes that the FluoIR of each fluorophore in the sample is a linear combination of multi-exponential functions, characterized by the time constants for the exponential functions and the scaling coefficients for each sample position (local approach). On the other hand, a global approach on the time constants is defined when these parameters are assumed constant throughout the sample, i.e. just the scaling coefficients of the exponential functions change for every sample position. Given the non-linear dependence on the time constants, a non-linear least squares approach is assumed to estimated these parameters [30], while the scaling coefficients can be solved using standard least squares routines for each spatial point. In other approaches, the FluoIR is modeled as a linear combination of discrete-time Laguerre basis functions. In this way, the objective is to estimate the scaling coefficients of the Laguerre functions, whose values can be computed through a linear least squares estimation. The disadvantage of this approach is that for some cases, the resulting FluoIRs might not have a monotonic decay, which is a characteristic of the fluorescence decays, or the fitting could be poor for extremely short or long lifetimes [4]. Therefore, a constrained optimization method is applied with restrictions on the second or third order time derivatives of the resulting FluoIRs to overcome this drawback. Finally, fluorescence lifetimes can be computed from the estimated FluoIRs to provide information about the chemical composition of the sample [31].

An alternative approach is blind deconvolution estimation (BDE), which aims to simultaneously estimate the InstR and the FluoIRs in the sample [32]. Under this perspective, the resulting InstR will be automatically aligned with the FluoDs to minimize the bias estimation of the FluoIRs. To the best of the author’s knowledge, there is only one BDE algorithm reported for FLIM data, which considers a Laguerre-basis approach [32]. However, a key hyper-parameter of this proposal has to be carefully tuned for the given input data, and it requires a constrained optimization scheme to achieve a monotonically decaying FluoIR. An equivalent formulation to BDE was presented in [33] for spike train inference from fluorescence measurements to evaluate neurons activity. In this formulation, the generative model is linear auto-regressive, and the spike inference is obtained by an approximated maximum a posterior formulation with a sparse constraint. This same formulation is further extended in [34] by an active set method to solve the sparse non-negative deconvolution problem (non-negative least absolute shrinkage and selection operator, LASSO, formulation) motivated by isotonic regression. Recently, in [35], the authors introduced the short-and-sparse deconvolution methodology which is based on a bilinear LASSO problem, sphere constraints and data-driven initialization.

In this context, this paper introduces two new BDE algorithms based on modeling the FluoIR by a linear combination of multi-exponential functions. The first BDE algorithm seeks for the characteristic parameters of the exponential functions (time constants and scaling coefficients) with a local perspective in each spatial point of the sample, i.e. pixel-by-pixel. On the other hand, in the global approach, the exponential functions in the FluoIR model are assumed common for all spatial points of the dataset, while their contributions change across the sample, i.e. the time constants of the exponentials are common to all pixels, but the scaling coefficients in each pixel are different. Given the monotonically decaying nature of the exponential functions, there is no need to include a shape-constraint during the estimation process in contrast to [32]. In addition, the only hyperparameters are the lower and upper bounds for the time constants of the exponential functions, which can be easily defined based on prior information of the fluorophores expected in the sample. To overcome the nonlinear interaction between FluoIRs and InstR variables, an alternating least squares (ALS) methodology, iteratively solves both estimation problems. In fact, due to a linear dependence on the InstR parameters in the observation model, a global search is made to estimate and align the shape of the UV laser-pulse through a linear non-negative least squares approximation. The InstR estimation is assumed with a free-form structure and a sparcity condition. In this way, our BDE methods jointly provides an estimation of the InstR and FluoIRs in the sample. The synthetic and experimental results, including synthetic datasets, fluorescence dyes, and oral tissue samples, show that the proposals are robust to different noise types and levels, and achieve a low computational time compared to other strategies in the state of the art. An initial version of the BDE algorithm based on the global approach was presented in [36]. Contrary to [33, 34], in our approach, the InstR does not have a spike train shape, and the observation model relies on a convolution with a multi-exponential structure. Furthermore, our formulation does not assume a short kernel nor a sparse activation map, as in comparison to [35]; and also our multi-exponential kernel is not restricted to a sphere constraint.

The notation used in this paper is described next. Scalars are denoted by italic letters, while vectors and matrices by lower and upper-case bold letters, respectively. and represent the real and integer numbers, respectively, and N-dimensional real vectors. For a real vector x, the transpose operation is denoted by x, the l-th element by x[l], the Euclidean norm by , and x ≥ 0 (x > 0) represents that each element in the vector is non-negative (positive). For a square matrix , Xi,j represents the element in the i-th row and j-th column (i, j ∈ [1, N]), Tr(X) = ∑i Xi,i denotes the trace operation, and the Frobenius norm. X ≥ 0 represents that each element in the matrix is non-negative. For all vectors and , denotes a Toeplitz matrix with x and y as its first column and row, respectively, where the first element in x and y must be equal. An N-dimensional vector for which all elements are ones (zeros) is represented by 1N (0N), and IN denotes the identity matrix of order N. For a random variable x, represents that x is uniformly distributed in the interval [a, b] (b > a), and that x is normally distributed with zero mean and variance σ2.

Blind deconvolution estimation

We assume that the FLIM dataset is sampled regularly with a period T over a spatial domain of K points in the field of view (FOV) of the instrument [13]. Hence, the FluoDs at any spatial location are L-th dimensional vectors that satisfy zk ≥ 0 ∀k ∈ [0, K − 1]. The set of FLIM measurements is denoted as , where their spatial location is not relevant so the ordering in is indistinct. Considering an initial pre-processing stage, all the measured FluoDs are normalized to sum-to-one, therefore: (1) such that the dataset of scaled FluoDs is given by and . The FLIM observation model for the l-th time sample in the k-th spatial position is given by (2) where yk[l], u[l] and hk[l] denote the measured FluoD, InstR, and FluoIR, respectively, ⋆ represents the convolution operator, and vk[l] stands for random noise related to measurement uncertainty. As a result, the scaled FluoDs represent the causal convolution between the InstR and the FluoIR at each spatial location [37]. The FLIM observation model is depicted in Fig 1. In our formulation, the FluoIR time sample at k-th position is characterized by the linear combination of N exponential functions. Thus, in the FluoIR model, the free parameters are the scaling coefficients of the exponential functions, and their time constants. Nonetheless, according to the nature of the exponential functions, we adopt two strategies to model the l-th FluoIR time sample:

  • Local Approach: the time constants of the exponential functions are different at each k-th spatial location: (3)
  • Global Approach: common exponential functions with time constants are scaled at each k-th spatial location: (4)
thumbnail
Fig 1. FLIM observation model.

The instrument response u[l] is common to all K spatial points, and the FluoDs result from the convolution between the FluoIRs and the InstR u[l].

https://doi.org/10.1371/journal.pone.0248301.g001

The scaling coefficients n ∈ [0, N] are selected different for each spatial location k, such that the estimated FluoD matches the actual measurement [19, 20, 29]. In this way, the local perspective in the FluoIR construction allows a major diversity in the fitting of the measured FluoDs, since at k-th spatial location there are 2N + 1 parameters and to be estimated. On the other hand, the global perspective assumes the same exponential functions in the whole dataset, so at k-th spatial location just the N + 1 scaling parameters are obtained, resulting in a faster fitting procedure compared to the local approach at the expense of limited diversity. Nonetheless, the fitting accuracy of the measured FluoDs by the local and global approaches will depend on the studied FLIM dataset and the order selection of the multi-exponential models in Eqs (3) and (4).

The observation model in Eq (2) can be expressed in vector notation for the local approach as (5) and similarly, for the global approach as (6) where (7) (8) (9) (10) (11) (12) (13)

By following a similar perspective to [32], the InstR is assumed common to all K spatial points in the FLIM dataset, and with a free-form. The InstR samples are also considered non-negative and normalized to sum-to-one to avoid scaling ambiguity, such that (14) In FLIM applications, the InstR is sparse, so there is no need to estimate every sample. Hence, we consider only the first terms to represent the InstR. As a result, the input matrix U in Eq (7) can be parametrized as a linear combination of Toeplitz matrices: (15) where the parameter θl = u[l] represents l-th sample in the InstR, and

With this mathematical description, the blind deconvolution estimation (BDE) problem is formulated as jointly obtaining the InstR components , the scaling coefficients and the time constants, either or just τ for a given FLIM dataset . Hence assuming Gaussian noise and independent measurements in Eq (5), we can formulate the local BDE as a maximum likelihood estimation, which is equivalent to the following nonlinear least-squares (NLS) approximation problem [30]: (16) such that (17) (18) A similar approximation problem is defined for the global BDE: (19) with restrictions in Eq (17), ck ≥ 0 ∀k, and τ > 0.

The inverse problems in Eqs (16) and (19) involve nonlinear interactions among the free parameters of the InstR , and those involved in the FluoIR or τ. To tackle these problems, we applied ALS approaches similar to [21] and [38, 39], where we estimate a solution for the FluoIR components while fixing the InstR, and vice versa until convergence.

Alternated least-squares methodology for blind estimation

The proposed ALS methodology solves iteratively the local and global BDE problems in Eqs (16) and (19), and these two mathematical formulations have distinctive features that are discussed next. We observe that the InstR parameters have a linear dependence on the approximation error, as well as the scaling coefficients in the FluoIRs. However, the exponential time constants or τ show a nonlinear interaction. Furthermore, the scaling coefficients for each spatial point are always non-negative, while the time constants are restricted to be positive values. Moreover, in our formulation for FLIM datasets, the InstR is a narrow pulse without repetitions, so each measured FluoD will exhibit a sharp increase to its peak value, followed by a monotonic decrease. Moreover, all the measured FDs are scaled to sum-to-one, and also the InstR is limited to sum-to-one. As a result, the scaled shift symmetry described in [35] will not hold in our BDE formulation. In addition, at each iteration of the ALS scheme in the local and global approaches, a quadratic approximation problem is solved by either a non-linear least squares or linear least squares. So, at each iteration, the estimation error is reduced or at least maintained. Consequently, convergence is guaranteed in the iterative scheme, but only to a local minimum. For this reason, the initialization based on processing the FLIM dataset is a crucial step in our formulation to obtain meaningful results.

We set initial conditions for the InstR parameters, and for the time constants of the exponential functions by considering a raw estimation of the overall lifetime present in the dataset, this process is fully described in the following sections. To speed up the estimation of the InstR parameters, and considering that the InstR is common for all spatial points in the dataset, only a random subset of FluoDs () was used for the initial estimation. The overall process for ALS assumes fixing one unknown parameter while optimizing for the other, i.e. fixing the InstR and estimating the FluoIRs or vice versa. In this sense, we first fixed the InstR to its initial condition and estimated the FluoIRs, for each spatial position of the reduced subset, by NLS following either the local or the global perspective [30]. Next, the estimated FluoIRs are fixed and the InstR is computed by constrained linear least squares. This alternated process is repeated until convergence for the InstR is reached. Finally, the FluoIRs are estimated taking into account the estimated InstR from the last iteration. The mathematical derivations for each optimization methodology are presented in the appendix, and the block diagrams of the local and global implementations are shown in Figs 2 and 3, respectively.

thumbnail
Fig 2. Blind deconvolution estimation: Local approach.

Block diagram of the general methodology with local approach to estimate FluoIRs.

https://doi.org/10.1371/journal.pone.0248301.g002

thumbnail
Fig 3. Blind deconvolution estimation: Global approach.

Block diagram of the general methodology with global approach to estimate FluoIRs.

https://doi.org/10.1371/journal.pone.0248301.g003

Estimation of FluoIR parameters

At this stage, the components of the InstR are assumed fixed, i.e. matrix U is known in Eqs (5) and (6). The local FluoIR estimation considers the estimation of the parameters {ck, τk} at k-th spatial point by NLS with a Levenberg-Marquardt (LM) approach [30]. The damping term in the LM approach is adapted according to a similarity metric ρk [40]. The mathematical derivations of the local approach are described in S2 Appendix in S1 File (see Fig 2). For this purpose, we implement two iterative blocks for the k-th measurement to estimate the FluoIR parameters (ck, τk):

  1. A convergence loop that computes the ALS for the optimal parameters () indexed by h until the percentage approximation error converges with respect to the previous iteration h − 1: (20) where denotes the k-th estimated FluoD at h-iteration, i.e., the convergence criterion is given by (21) where ϵ1 > 0 is the stopping threshold.
  2. An internal loop that adjusts the damping factor in the LM approach until , where ϵ2 > 0. After this condition is achieved, the damping factor is updated for subsequent iterations.

The global FluoIR estimation takes into consideration a non-negative least-squares (LS) approximation to estimate the scaling parameters {ck} for each k-th pixel in the dataset, which is followed by a NLS optimization for the time constants τ over the reduced dataset until convergence [30]. Similar to the local approach, the NLS optimization is implemented by a LM approach [40]. The mathematical derivations of the global approach are described in S3 Appendix in S1 File (see Fig 3). As previously outlined, we implement two iterative blocks to estimate the global FluoIR parameters:

  1. An outer-loop that computes the ALS for the optimal scalings and time constants τh indexed both by h until the percentage approximation error eh over the reduced dataset converge with respect to the previous iteration: (22) where i.e. (23) where ϵ1 > 0 is the stopping threshold.
  2. An inner-loop that adjusts the damping factor λh in the LM approach to estimate the time constants τ until ρh+1 > ϵ2, where ϵ2 > 0. After this condition is achieved, the damping factor λh is updated for subsequent iterations.

Estimation of InstR parameters

At this step, the scaling coefficients {ck}, and time constants {τk} or τ are assumed known and fixed in the cost functions in Eqs (16) or (19), and the optimization is implemented with respect to the samples of the InstR . In this case, a closed-form solution can be calculated, as will be shown next by using the reduced dataset . In this formulation, we are not assuming a pre-defined form or time-pattern for the InstR samples, just that are always non-negative values and sum-to-one. The optimization problem is re-formulated in a compact structure as (24) such that where (25) i.e. where . The mathematical derivations of the InstR parameter estimation are described in S4 Appendix in S1 File (see Figs 2 and 3). An iterative scheme following the ALS philosophy is implemented between the FluoIR estimation at each spatial pixel, and the InstR parameters over the reduced dataset. The iterative scheme is stopped after convergence of the overall parameters in the FluoIRs and InstR. If this iterative structure is indexed by t, i.e. the convergence performance is evaluated by a normalized metric: (26) and the next stoppage criterion is considered as (27) where ϵ3 > 0 is the stopping threshold in this outer-loop.

Initial conditions

An initial condition is required in the ALS scheme for InstR parameters θ0, and time constants in the FluoIRs estimation for the local or global τ0 approaches. We start by estimating an initial condition for θ0, which is defined as follows: (28) The proposed initial InstR represented by θ0 is a narrow pulse, like the one used in fluorescent measurements applications [32, 41], which will be refined in each subsequent iteration. Another advantage is that the resulting estimated InstR is aligned with the FluoDs, resulting in a better approximation of the FluoIRs.

An initial estimation for the time constants τ0 or is computed from the average lifetime in the whole dataset, which is defined as: (29) where the peak value of each FluoD is defined as: (30) By defining the mean and standard deviation of the values in : (31) (32) the elements in τ0 and are selected in the range . By setting the initial condition on for the local approach, the matrices U and HL(τ0) can be constructed from Eqs (7) and (8), respectively. To set the initial scaling coefficients in the NLS scheme for k-th spatial pixel, a LS problem is formulated from Eq (16) as: (33) which can be solved using standard numerical methods efficiently [30].

Synthetic and experimental validation

To validate the proposed BDE algorithms, we consider synthetic and experimental FLIM datasets. In both cases, the performance of the proposals will be tested by measuring the estimation errors on the InstR and FluoDs; in addition, for the synthetic evaluation, we can also quantify the errors on the estimated FluoIRs. We also evaluate different scenarios of Gaussian and Poisson noise in the FluoDs for the synthetic datasets [23]. Furthermore, we quantify the shape of the InstR by the full-width at half-maximum (FWHM) parameter Δufwhm: (34) where the time indexes 0 < l1 < Imax and Imax < l2 < L − 1 satisfy (35) where (36) By assuming that , and denote the FluoD, FluoIR and InstR estimations by a BDE algorithm, respectively, and yk, hk and u the corresponding ground-truths, we employ four estimation performance metrics: (37) (38)

The metrics (Eu, Eh, Ey) provide information of the estimation performance with respect to the estimated InstR, FluoIRs and FluoDs in a percentage fashion. For the synthetic FluoIR hk at k-th spatial point, the average lifetime (ALT) τk is a key parameter related to the fluorescent property of the sample: (39) where represents a vector of the sampling times. These ALTs are used to generate a representative image with spatial information of the chemical composition of the sample. Hence, the error between the ALTs generated through the synthetic FluoIR and the estimated one can be quantified by using the next normalized metric: (40)

In our evaluations, we compare the proposed BDE local and global algorithms based on multi-exponential models (BDELME and BDEGME) to the BDE by a Laguerre-basis (BDELB) [32], and a blind extension of the exponentials library deconvolution (BDEEL) approach [42]. In addition, we also compare the standard local and global deconvolution methodologies with a multi-exponential model that assume available the InstR, and they are denoted as: DELME and DEGME, respectively [17, 18]. As suggested in [32], BDELB was implemented with a 8th order approximation and shape parameter 0.85. Meanwhile, for BDEEL, the exponential library contains 25 elements and a weight factor of 0.25, as suggested in [41] and [42]. All the MATLAB implementations of the methodologies: DELME, DEGME, BDELME, BDEGME, BDELB and BDEEL are freely available in the website http://galia.fc.uaslp.mx/~bde.

Synthetic evaluation

The proposals were first validated by using synthetic datasets under different types and noise levels. The synthetic datasets were generated considering a measured InstR [5], with a sampling interval T = 0.25 ns and a length of 186 samples (L = 186). This measured InstR is a positive time-signal with a sharp rising time and exponential decay with Δufwhm = 1.53 ns. The k-th synthetic FluoIR is modeled as a sum of Nsynth ∈ {2, 3, 4} exponential functions: (41) where the magnitudes ak,i were selected to have uniform regions of high concentration of all components in the dataset [41], ak,0 = 0.001 ∀k, and the characteristic times are selected randomly, but in a limited interval ns, ns, ns, and ns ∀k. With these definitions of the synthetic FluoIRs, the BDE local and global perspectives are both viable.

Next, the synthetic noise-free FluoDs are obtained by applying the convolution operator in (2), i.e. . In our evaluation, we included Gaussian and Poisson noise to the FluoDs to take into account uncertainty in the equipment according to the following model [23]: (42) where and represent normal random variables, and the variances , are selected with respect to a desired signal-to-noise ratio (SNR) and peak-to-noise signal ratio (PSNR): (43) (44)

Since the construction of the synthetic datasets involved uniformly distributed random samples, we carried out a Monte Carlo evaluation with 10 repetitions by implementing the BDE algorithms according to the parameters listed in Table 1, whose selection is explained next. In the synthetic evaluation, a spatial domain of 100 × 100 = 10, 000 samples were generated by Eqs (41) and (42) at different values of SNR ∈ {40, 45, 50, 55} dB, and PSNR ∈ {10, 15, 20, 25} dB, and the resulting datasets were analyzed by the BDE algorithms, i.e. K = 10, 000. The random spatial sampling retained only to estimate the InstR, i.e. 20% of the complete dataset. As previously mentioned, we considered three different values for the number of exponential functions in the synthetic FluoIRs Nsynth ∈ {2, 3, 4}. Since for DELME, DEGME, BDELME and BDEGME, we could not know a priori the number of exponentials in the FluoIRs model, we assume a fourth order model in the synthetic datasets, i.e. N = 4.

thumbnail
Table 1. Parameters of synthetic dataset and BDE implementation during the synthetic evaluation.

https://doi.org/10.1371/journal.pone.0248301.t001

First, we have included a numerical convergence evaluation of the BDE schemes (BDELME and BDEGME) with the synthetic datasets for different noise levels, and orders in the multi-exponential model of the synthetic impulse response. The stopping threshold is set to ϵ3 = 1 × 10−3. The resulting normalized metric Υt in the alternated least squares iterations is described below in Fig 4, where it is observed that in either BDE scheme, the convergence is always monotonic for any noise combination and model order, and just its final steady-state value depends on the noise level. From Fig 4, we conclude that the local approach (BDELME) converges slower than the global scheme (BDEGME), which could be intuitively expected, since BDELME has roughly twice the free parameters to optimize. For BDELME, the most drastic reductions are achieved in the first five iterations, and for BDEGME, in the first two.

thumbnail
Fig 4. Normalized metric (Υt) vs iteration (t) in the ALS scheme for (top row) BDELME and (bottom row) BDEGME by considering synthetic FLIM datasets at different noise levels (SNR,PSNR) ∈ {(10, 40), (15, 45), (20, 50), (25, 55)} dB, and impulse response orders (2nd, 3rd and 4th).

https://doi.org/10.1371/journal.pone.0248301.g004

In the following evaluations, to evaluate the convergence in the iterative process of the BDE algorithms (see Figs 2 and 3), we use a 5% tolerance in all iterations, ϵ1 = ϵ2 = ϵ3 = 0.05, as a good balance between precision and complexity in our evaluations. All the data processing was carried out in MATLAB, and the computational time was measured using the function “tic-toc”. For the experiments, we used a MacBook Pro with an Intel Dual Core i5 CPU at 2.3 GHz, and 16 GB of RAM. In addition, our evaluation also considered DELME and DEGME with the measured InstR by applying the procedures in S1 and S3 Appendices in S1 File. Hence, the performance indexes Ey, Eh and Ealt can examine the estimation accuracy in the FluoDs, FluoIRs and ALTs by a direct deconvolution process, as well as the computational time. Meanwhile, the metrics Eu and Efwhm can also quantify the accuracy in the InstR estimations by the BDE strategies.

Fig 5 illustrates the computational time and performance metrics (Ey, Eh, Ealt) as a function of ascending SNR/PSNR pairs for 2nd, 3rd and 4th order synthetic FluoIR models. As expected, the lowest computational time is achieved by the standard deconvolution techniques, DELME and DEGME. In fact, the results in the computational time are independent of the order in the synthetic model. Now, from the BDE strategies, BDEGME obtained the lowest computational time consistently, and BDELME was just surpassed by BDEEL in all scenarios. For the Ey metric, the results had the same tendency for all algorithms, as the pair SNR/PSNR increased, the error metric decreased. In the FluoIR estimations, BDELB reached the worst performance in all scenarios. For all BDE algorithms, the results for Eh show a small variability with respect to the SNR/PSNR pairs. The proposed BDELME and BDEGME had always Eh errors lower than 20%, in contrast with state of the art approaches. Finally, the errors in ALT Ealt are small (<10%) for all methodologies, despite all different noise types and levels. Fig 6 shows the estimation performance in the InstR by the metrics Eu and Efwhm. It is remarkable that the proposed blind techniques BDELME and BDEGME reached the best performance for all noise pairs SNR/PNSR in 2nd and 4th order synthetic models. Meanwhile, for the 3rd order model, their performance was in between the worst by BDELB and the minimum by BDEEL.

thumbnail
Fig 5. Performance metrics in synthetic evaluation with 2nd, 3rd and 4th order FluoIR models: (top row) Computational time, (middle-top row) Ey, (middle-bottom row) Eh, and (bottom row) Ealt.

https://doi.org/10.1371/journal.pone.0248301.g005

thumbnail
Fig 6. Performance metrics in synthetic evaluation with 2nd, 3rd and 4th order FluoIR models: (top-row) Eu, and (bottom row) Efwhm.

https://doi.org/10.1371/journal.pone.0248301.g006

To illustrate the estimation performance, Fig 7 shows the ALT map for one realization of the synthetic dataset with the most severe noise conditions, i.e. SNR = 40 dB and PSNR = 10 dB, and 3rd order synthetic model. As shown in Fig 6, all the ALT maps are consistent with the ground-truth, since the errors are less than 10% for all deconvolution techniques. Finally, Fig 8 presents the resulting time responses for k = 1, 000 spatial points in the same testing scenario. Thus, the top plot highlights the heavy noise condition in the synthetic FluoD measurement, and the accurate estimation by all techniques, since the errors are lower than 7%, as shown in the middle column for Ey. The middle plot illustrates that the closest response to the ground-truth is achieved by DELME, BDELME, DEGME and BDEGME. Finally, the bottom plot shows that the InstR is accurately estimated by all BDE methodologies, since the errors are always lower than 20% (see Fig 6).

thumbnail
Fig 7. Average lifetime map for one realization in synthetic evaluation with 3rd order FluoIR model, SNR = 40 dB and PSNR = 10 dB: (a) Ground-truth, (b) DELME, (c) BDELME, (d) DEGME, (e) BDEGME, (f) BDELB, and (g) BDEEL.

https://doi.org/10.1371/journal.pone.0248301.g007

thumbnail
Fig 8. Estimated time responses for k = 1, 000 spatial position for one realization with 3rd order FluoIR model, SNR = 40 dB and PSNR = 10 dB: (top) FluoDs, (middle) FluoIRs, and (bottom) InstRs.

https://doi.org/10.1371/journal.pone.0248301.g008

As conclusions of the synthetic evaluation, BDEGME reached the fastest convergence in the blind techniques, independently to the noise scenario. BDELME required more computational time, but not significantly more to BDEGME. In the FluoIR and InstR estimations, BDEGME and BDELME reached the best performance in the 2nd and 4th synthetic order models, and for the 3rd order model, their results were between the best and worst of all the approaches considered.

Experimental evaluation with fluorescence dyes

The proposals were first validated experimentally with fluorescence dyes: POPOP, FAD and NADH. The lifetimes of these dyes are reported in the literature: ∼1.3 ns, 2.0-2.5 ns and 0.3-0.6 ns, respectively [1, 4]. The FLIM datasets were collected at the wavelength channels: 390±40 nm (channel 1), and 494±41 nm (channel 2), with a sampling time of 0.25 ns (T = 0.25 ns). The FAD dye has an emitting response only in channel 1, the POPOP only in channel 2, while the NADH response to both channels. For each channel measurement, there are defined 80 time samples (L = 80) over a spatial resolution of 1,000 pixels (K = 1, 000). Since the number of spatial samples is relatively small, there is no random subsampling to estimate the InstR, i.e. . The InstR was measured at both channels, where the FWHM of the UV laser-pulse (355 nm) was Δufwhm = 1.78 ns and 1.97 ns for channel 1 and 2, respectively. Since we expect two fluorophores per channel, we set N = 2 for BDELME and BDEGME. The estimated ALT by using BDELME, BDEGME, BDELB and BDEEL are presented in Table 2 with parameters τmin = 0.25 ns and τmax = 4 ns (based on the literature results [1, 4]). The computational times for the BDE algorithms in channels 1/channel 2 were 1.78/2.76 s, 1.75/2.07 s, 2.95/2.80 s, and 2.45/2.92 s for BDELME, BDEGME, BDELB and BDEEL, respectively. Thus, the lowest computational times were obtained by the proposed algorithms: BDELME and BDEGME.

thumbnail
Table 2. Estimated average lifetimes (ns) for synthetic dyes.

https://doi.org/10.1371/journal.pone.0248301.t002

These results show a good agreement with the literature by considering that there is an uncertainty factor in the estimation due to the sampling interval of 0.25 ns. In addition, the performance of the proposed BDE algorithms is visualized by using the Bland-Altman (B&A) methodology with the estimated ALTs [43]. The B&A plot measures the similarity between two types of data, by setting limits of agreement in terms of the mean and the standard deviation (SD) of the differences between the two sources. In this case, the B&A analysis was computed with the ALTs generated by pairs of the studied BDE algorithms: BDELME, BDEGME, BDELB and BDEEL, and Figs 9 and 10 illustrate the B&A plots for all pairs. In addition at the top of each subplot in Figs 9 and 10, the correlation coefficient ρ was computed between both estimated ALTs. As a result, we observe that in all scenarios, two distinctive ALTs are obtained per channel by the two characteristic clusters in the B&A plots. Also, the plots illustrate that the estimated ALTs are mostly contained in the 95% confidence interval defined by the red dashed-lines (mean ± 1.96 SD), and the correlation coefficients are always greater than 0.99 for both channels.

thumbnail
Fig 9. Bland-Altman plots for estimated average lifetimes with fluorescence dyes in channel 1 (390±20 nm) and correlation coefficients: (a) BDELME vs BDEGME, (b) BDELME vs BDELB, (c) BDELME vs BDEEL, (d) BDEGME vs BDELB, (e) BDEGME vs BDEEL, and (f) BDELB vs BDEEL.

https://doi.org/10.1371/journal.pone.0248301.g009

thumbnail
Fig 10. Bland-Altman plots for estimated average lifetimes with fluorescence dyes in channel 2 (494±41 nm) and correlation coefficients: (a) BDELME vs BDEGME, (b) BDELME vs BDELB, (c) BDELME vs BDEEL, (d) BDEGME vs BDELB, (e) BDEGME vs BDEEL, and (f) BDELB vs BDEEL.

https://doi.org/10.1371/journal.pone.0248301.g010

Experimental evaluation with oral tissue samples

From the study in [44], oral tissue samples were used for validation in this section that belong to different regions in the oral cavity. Dysplastic and cancerous oral lesions were analyzed by in vivo clinical endogenous mFLIM images. The imaging protocol was approved by the Institutional Review Board at Texas A&M University. The clinical diagnosis and more detailed information about the samples are presented in Table 3. The temporal resolution of the measurements is 0.16 ns (T = 0.16 ns). All the measurements included the fluorescent responses to three wavelength bands: 390 ± 20, 452 ± 22.5, and 550 ± 20 nm, that correspond to channel 1, channel 2 and channel 3, respectively. Only 186 time samples were considered for each channel (L = 186). The spatial dimensions of the tissues are approximately 10 mm × 10 mm, divided equiespatially in 160 × 160 pixels (K = 25, 600). Due to low SNR, some pixels in the images were manually masked, so their information will not misguide the BDE strategies. During the experimental evaluations, the following parameters were set for BDELME and BDEGME as τmin = 0.64 ns, τmax = 16 ns, and ϵ1 = ϵ2 = ϵ3 = 0.05. The experimental InstRs used to generate the fluorescent measurements in channels 1, 2 and 3 have a FWHM of Δufwhm = 1.70 ns, 1.88 ns and 2.21 ns, respectively. For each oral tissue the dataset contains two samples, lesion and reference. The reference sample was taken from the symmetrical location of the sagittal plane.

thumbnail
Table 3. Detailed information of the oral tissue lesion samples.

https://doi.org/10.1371/journal.pone.0248301.t003

Fig 11 shows the resulting computational time and performance metrics (Ey, Eu, Efwhm) for all the datasets describe in Table 3. We can observe that the computational time is always the lowest for BDELME and BDEGME for all three channels. With respect to the estimation error in the FluoDs, Ey exhibits the same tendency for all datasets and channels, where the differences are lower than 1%. Meanwhile, for the InstR estimations, no particular BDE methodology achieved consistently the lowest errors Eu and Efwhm. Nonetheless, BDELME and BDEGME presented more regular responses in Eu and Efwhm, especially in the second and third channels. In general, the largest errors were obtained by BDELB in the three channels.

thumbnail
Fig 11. Performance metrics in experimental evaluation with lesion and normal samples for the three spectral channels: (top row) Computational time, (middle-top row) Ey, (middle-bottom row) Eu, and (bottom row) Efwhm.

https://doi.org/10.1371/journal.pone.0248301.g011

To illustrate a specific response for the estimated FluoIRs, Fig 12 presents the estimated ALT maps for lesion sample No. 82 in the 3rd channel with the four BDE algorithms. The subplots illustrate the same morphological patterns with just small differences in the minimum and maximum ALTs. Finally, Fig 13 presents the B&A plots for all pairs of ALT estimations, and the corresponding correlation coefficients. In all cases, the correlation coefficients are larger than 0.94, which highlights high consistency among all BDE techniques.

thumbnail
Fig 12. Average lifetime maps for lesion sample No. 82 in 3rd spectral channel: (a) BDELME, (b) BDEGME, (c) BDELB, and (d) BDEEL.

https://doi.org/10.1371/journal.pone.0248301.g012

thumbnail
Fig 13. Bland-Altman plots for estimated average lifetimes for lesion sample No. 82 in 3rd spectral channel and correlation coefficients: (a) BDELME vs BDEGME, (b) BDELME vs BDELB, (c) BDELME vs BDEEL, (d) BDEGME vs BDELB, (e) BDEGME vs BDEEL, and (f) BDELB vs BDEEL.

https://doi.org/10.1371/journal.pone.0248301.g013

Conclusions

In this work, we introduced two new BDE algorithms based on a linear combination of multi-exponential functions for the FluoIRs modeling. Hence, the proposed algorithms estimate the FluoIRs in all the spatial points of the dataset and, simultaneously, the InstR for the fluorescence excitation. The InstR is assumed with a free-form and a sparcity condition. The local perspective of the BDE methodology (BDELME) assumes that the characteristic parameters of the exponential functions (time constants and scaling coefficients) are estimated for each spatial point of the dataset, i.e. pixel-by-pixel. On the other hand, for the global perspective, the exponential functions are assumed common to all the points in the dataset, and just their scaling coefficients are updated for each spatial point. By using a convolution modeling between the FluoIRs and InstR for the measured FluoDs, the time samples of the InstR and the scaling coefficients of the exponential functions exhibit a linear dependence in the observation model, but for the exponentials time constants, the dependence is nonlinear. To overcome the nonlinear interaction on the decision variables, an alternating least squares (ALS) methodology iteratively solves both estimation problems based on non-negative and constrained optimizations. The validation stage considered synthetic datasets at different noise types and levels, and a comparison with the standard deconvolution techniques: DELME and DEGME, as well as, two more BDE methodologies in the state of the art: BDELB and BDEEL. In the validation with experimental datasets, fluorescent dyes and oral tissue samples were considered. Our results show that BDELME and BDEGME reached the fastest convergence with the best compromise in FluoIRs and InstR estimation errors compared to BDELB and BDEEL. Also, the estimation performance of BDELME and BDEGME was consistent with the standard deconvolution techniques that assume the InstR is available: DELME and DEGME. For future work, we will develop parallel implementations of all the BDE methodologies to reduce computational time, as well as design a more comprehensive comparison with diverse tissue samples.

References

  1. 1. Marcu L, French PM, Elson DS. Fluorescence lifetime spectroscopy and imaging: principles and applications in biomedical diagnostics. CRC Press; 2014.
  2. 2. Anthony N, Guo P, Berland K. Principles of fluorescence for quantitative fluorescence microscopy. In: FLIM microscopy in biology and medicine. Chapman and Hall/CRC; 2009. p. 67–96.
  3. 3. Periasamy A, Clegg RM. FLIM microscopy in biology and medicine. Chapman and Hall/CRC; 2009.
  4. 4. Datta R, Heaster TM, Sharick JT, Gillette AA, Skala MC. Fluorescence lifetime imaging microscopy: fundamentals and advances in instrumentation, analysis, and applications. Journal of Biomedical Optics. 2020;25(7):1–43.
  5. 5. Park J, Pande P, Shrestha S, Clubb F, Applegate BE, Jo JA. Biochemical characterization of atherosclerotic plaques by endogenous multispectral fluorescence lifetime imaging microscopy. Atherosclerosis. 2012;220(2):394–401.
  6. 6. Shrestha S, Applegate BE, Park J, Xiao X, Pande P, Jo JA. High-speed multispectral fluorescence lifetime imaging implementation for in vivo applications. Optics Letters. 2010;35(15):2558–2560.
  7. 7. O’Donnell M, McVeigh ER, Strauss HW, Tanaka A, Bouma BE, Tearney GJ, et al. Multimodality cardiovascular molecular imaging technology. Journal of Nuclear Medicine: Official Publication, Society of Nuclear Medicine. 2010;51(Suppl 1):38S. pmid:20457794
  8. 8. De Giorgi V, Massi D, Sestini S, Cicchi R, Pavone F, Lotti T. Combined non-linear laser imaging (two-photon excitation fluorescence microscopy, fluorescence lifetime imaging microscopy, multispectral multiphoton microscopy) in cutaneous tumours: first experiences. Journal of the European Academy of Dermatology and Venereology. 2009;23(3):314–316.
  9. 9. Galletly N, McGinty J, Dunsby C, Teixeira F, Requejo-Isidro J, Munro I, et al. Fluorescence lifetime imaging distinguishes basal cell carcinoma from surrounding uninvolved skin. British Journal of Dermatology. 2008;159(1):152–161. pmid:18460029
  10. 10. Jabbour JM, Cheng S, Malik BH, Cuenca R, Jo JA, Wright J, et al. Fluorescence lifetime imaging and reflectance confocal microscopy for multiscale imaging of oral precancer. Journal of Biomedical Optics. 2013;18(4):046012. pmid:23595826
  11. 11. Duran-Sierra E, Cheng S, Cuenca-Martinez R, Malik B, Maitland KC, Lisa Cheng YS, et al. Clinical label-free biochemical and metabolic fluorescence lifetime endoscopic imaging of precancerous and cancerous oral lesions. Oral Oncology. 2020;105:104635. pmid:32247986
  12. 12. Mycek MA, Schomacker KT, Nishioka NS. Colonic polyp differentiation using time-resolved autofluorescence spectroscopy. Gastrointestinal Endoscopy. 1998;48(4):390–394.
  13. 13. Liu L, Yang Q, Zhang M, Wu Z, Xue P. Fluorescence lifetime imaging microscopy and its applications in skin cancer diagnosis. Journal of Innovative Optical Health Sciences. 2019;12(5):1930004.
  14. 14. Romano RA, Teixeira Rosa RG, Salvio AG, Jo JA, Kurachi C. Multispectral autofluorescence dermoscope for skin lesion assessment. Photodiagnosis and Photodynamic Therapy. 2020;30:101704.
  15. 15. Walsh AJ, Cook RS, Sanders ME, Aurisicchio L, Ciliberto G, Arteaga CL, et al. Quantitative optical imaging of primary tumor organoid metabolism predicts drug response in breast cancer. Cancer Research. 2014;74(18):5184–5194. pmid:25100563
  16. 16. Lee KB, Siegel J, Webb S, Leveque-Fort S, Cole M, Jones R, et al. Application of the stretched exponential function to fluorescence lifetime imaging. Biophysical Journal. 2001;81(3):1265–1274. pmid:11509343
  17. 17. Pelet S, Previte M, Laiho L, So P. A fast global fitting algorithm for fluorescence lifetime imaging microscopy based on image segmentation. Biophysical Journal. 2004;87(4):2807–2817.
  18. 18. Warren SC, Margineanu A, Alibhai D, Kelly DJ, Talbot C, Alexandrov Y, et al. Rapid global fitting of large fluorescence lifetime imaging microscopy datasets. PLoS One. 2013;8(8):e70687. pmid:23940626
  19. 19. Jo JA, Fang Q, Papaioannou T, Marcu L. Fast model-free deconvolution of fluorescence decay for analysis of biological systems. Journal of Biomedical Optics. 2004;9(4):743–753.
  20. 20. Dabir AS, Trivedi C, Ryu Y, Pande P, Jo JA. Fully automated deconvolution method for on-line analysis of time-resolved fluorescence spectroscopy data based on an iterative Laguerre expansion technique. Journal of Biomedical Optics. 2009;14(2):024030.
  21. 21. Campos-Delgado DU, Gutierrez-Navarro O, Rico-Jimenez JJ, Duran-Sierra E, Fabelo H, Ortega S, et al. Extended Blind End-Member and Abundance Extraction for Biomedical Imaging Applications. IEEE Access. 2019;7:178539–178552. pmid:31942279
  22. 22. Qin B, Hu C, Huang S. Target/Background Classification Regularized Nonnegative Matrix Factorization for Fluorescence Unmixing. IEEE Transactions on Instrumentation and Measurement. 2016;65(4):874–889.
  23. 23. Gutierrez-Navarro O, Campos-Delgado DU, Arce-Santana ER, Maitland KC, Cheng S, Jabbour J, et al. Estimation of the number of fluorescent end-members for quantitative analysis of multispectral FLIM data. Optics Express. 2014;22(10):12255–12272. pmid:24921344
  24. 24. Pengo T, Muñoz-Barrutia A, Zudaire I, Ortiz-de Solorzano C. Efficient Blind Spectral Unmixing of Fluorescently Labeled Samples Using Multi-Layer Non-Negative Matrix Factorization. PLoS ONE. 2013;8(11):e78504.
  25. 25. Xu H, Rice BW. In-vivo fluorescence imaging with a multivariate curve resolution spectral unmixing technique. Journal of Biomedical Optics. 2009;14(6):064011.
  26. 26. Fereidouni F, Bader AN, Gerritsen HC. Spectral phasor analysis allows rapid and reliable unmixing of fluorescence microscopy spectral images. Optics Express. 2012;20(12):12729–12741.
  27. 27. Fereidouni F, Blab GA, Gerritsen HC. Blind unmixing of spectrally resolved lifetime images. Journal of Biomedical Optics. 2013;18(8):086006.
  28. 28. Campos-Delgado DU, Navarro OG, Arce-Santana E, Jo JA. Extended output phasor representation of multi-spectral fluorescence lifetime imaging microscopy. Biomedical optics express. 2015;6(6):2088–2105.
  29. 29. Liu J, Sun Y, Qi J, Marcu L. A novel method for fast and robust estimation of fluorescence decay dynamics using constrained least-squares deconvolution with Laguerre expansion. Physics in Medicine & Biology. 2012;57(4):843.
  30. 30. Nocedal J, Wright SJ. Numerical optimization. 2nd ed. Springer; 2006.
  31. 31. Kaye B, Foster PJ, Yoo TY, Needleman DJ. Developing and Testing a Bayesian Analysis of Fluorescence Lifetime Measurements. PLoS ONE. 2017;12(1):e0169337.
  32. 32. Campos-Delgado DU, Gutierrez-Navarro O, Arce-Santana ER, Skala MC, Walsh AJ, Jo JA. Blind deconvolution estimation of fluorescence measurements through quadratic programming. Journal of Biomedical Optics. 2015;20(7):075010.
  33. 33. Vogelstein JT, Packer AM, Machado TA, Sippy T, Babadi B, Yuste R, et al. Fast Nonnegative Deconvolution for Spike Train Inference From Population Calcium Imaging. Journal of Neurophysiology. 2010;104(6):3691–3704. pmid:20554834
  34. 34. Friedrich J, Zhou P, Paninski L. Fast online deconvolution of calcium imaging data. PLOS Computational Biology. 2017;13(3):1–26.
  35. 35. Lau Y, Qu Q, Kuo HW, Zhou P, Zhang Y, Wright J. Short and Sparse Deconvolution—A Geometric Approach. In: 8th International Conference on Learning Representations; 2020. p. 1–45.
  36. 36. Salinas-Martinez R, Campos-Delgado DU, Jo JA. Global Blind Deconvolution of Fluorescence Lifetime Imaging Microscopy. In: Latin America Optics and Photonics Conference. Optical Society of America; 2018. p. W2C.4.
  37. 37. Proakis JG, Manolakis DG. Digital Signal Processing. 4th ed. Prentice Hall; 2006.
  38. 38. Young FW, De Leeuw J, Takane Y. Regression with qualitative and quantitative variables: An alternating least squares method with optimal scaling features. Psychometrika. 1976;41(4):505–529.
  39. 39. Jaumot J, Gargallo R, de Juan A, Tauler R. A graphical user-friendly interface for MCR-ALS: a new tool for multivariate curve resolution in MATLAB. Chemometrics and Intelligent Laboratory Systems. 2005;76(1):101–110.
  40. 40. Madsen K, Nielsen H, Tingleff O. Methods for non-linear least squares problems, Informatics and Mathematical Modeling, Technical University of Denmark. URL Available online: http://www.imm.dtu.dk. 2004;.
  41. 41. Campos-Delgado DU, Navarro OG, Arce-Santana E, Walsh AJ, Skala MC, Jo JA. Deconvolution of fluorescence lifetime imaging microscopy by a library of exponentials. Optics Express. 2015;23(18):23748–23767.
  42. 42. Campos-Delgado DU, Gutierrez-Navarro O, Mejia-Rodriguez A, Jo JA. Blind Deconvolution Estimation by an Exponentials Library. In: 2020 IEEE International Autumn Meeting on Power, Electronics and Computing (ROPEC 2020). Ixtapa, Mexico: IEEE; 2020. p. 1–6.
  43. 43. Giavarina D. Understanding bland altman analysis. Biochemia medica: Biochemia medica. 2015;25(2):141–151.
  44. 44. Duran-Sierra E, Cheng S, Cuenca-Martinez R, Malik B, Maitland KC, Lisa Cheng YS, et al. Clinical label-free biochemical and metabolic fluorescence lifetime endoscopic imaging of precancerous and cancerous oral lesions. Oral Oncology. 2020;105:104635. pmid:32247986