MemExp: Fits by Distributed and Discrete Kinetic Models

The MEM

The MemExp program makes synergistic use of both the maximum entropy method (MEM) and maximum likelihood (ML) fitting to interpret kinetics data from the distributed and discrete perspectives. Its hybrid algorithm has been described in the literature [1]. MemExp version 2.0 extended the user's ability to manipulate the default model for kinetics involving overlapping phases [3]. Versions 3.0 and beyond support the rigorous ML analysis of Poisson data in addition to the treatment of Gaussian noise via nonlinear least-squares (NLS) fitting. Deconvolution of an instrument response function is supported for data that are equally spaced in time. A correction for scattered light was also first implemented in MemExp 3.0. New to version 4.0 was the easier-to-use 'auto' mode and a simple way to filter artifacts from the lifetime distribution HREF="node15.html#AB">2] by setting IBIGF = 3. Additional algorithmic details and applications of the MEM can be found elsewhere [4,5,6].

MemExp fits kinetics using one or two liftetime distributions:

$${\cal F}_i = D_0 \int_{-\infty}^{\infty} dlog\tau [g(log \ta... ...e^{ -t_i/\tau } + \sum_{k=0}^3 (b_k - c_k) (t_i/t_{max})^k,$$ (1)

where $g(log \tau)$ and $h(log \tau)$ are the distributions describing decaying and rising kinetics, respectively, and the polynomial describes the baseline. Baseline coefficients are scaled using the constant $t_{max}$ so they remain comparable in magnitude.

When deconvolving an instrument response function R(t), corresponding to experimentally measured values $R_i$, the data are described by:

$${\cal F}_i = b + \xi R_i + D_0 \int_{-\infty}^{\infty} dlog\... ..._0^t dt^\prime R(t^\prime + \delta) e^{ -(t-t^\prime)/\tau }.$$ (2)

Here, the data are assumed equally spaced in time. The convolution integral is calculated in version 6.0 using a cubic approximation of the count function for the instrument response [7]. Scattered excitation light can be accounted for with positive values of the parameter $\xi$. The zero-time shift parameter $\delta$ is optimized during a preliminary series of MEM runs using a reduced resolution in log $\tau$, before being fixed for the final MEM run.

The maximum entropy method (MEM) is used to obtain this fit ${\cal F}$ by distributions. A MEM calculation is a 'tug of war' between two or three terms, depending on whether or not the normalization is to be constrained. These terms are summed in the function Q,

$$Q \equiv S - \lambda C - \alpha I,$$ (3)

that is to be maximized. S is the entropy of the 'image' f,

$$S(\vec{f},\vec{F}) = \sum_{j=1}^M [f_j - F_j - f_j ln(f_j/F_j)].$$ (4)

Here, the image f includes both the g and h lifetime distributions that describe kinetic processes of opposite sign as well as any baseline parameters used. In this documentation, f is alternatively referred to as the image or the distribution obtained by the MEM. F is the MEM 'prior' distribution used to incorporate prior knowledge into the solution. Unconstrained maximization of S without regard for the fit or normalization ( $\lambda = \alpha = 0$) yields $f_j = F_j$.

The quality of the fit ${\cal F}$ to the data $D$ produced by f is measured by C, which is a minimum (zero) when the data and fit are everywhere equal. When analyzing data with normally distrbuted noise, C takes a familiar form:

$$C \equiv \chi^2 = \frac {1} {N} \sum_{i=1}^{N} \left( \frac { {\cal F}_i-D_i} {\sigma_i} \right) ^{2}.$$ (5)

When analyzing Poisson-distributed data, the appropriate function for C is the Poisson deviance:

$$C \equiv P = \frac {2} {N} \sum_{i=1}^{N} \left[ D_i ln(D_i/{\cal F}_i) + {\cal F}_i - D_i \right].$$ (6)

For normally distributed noise, all residuals plotted are given by:

$$r_i = (D_i - {\cal F}_i)/\sigma_i,$$ (7)

whereas for Poisson noise, all residuals are defined as:

$$r_i = (D_i - {\cal F}_i)/\sqrt{{\cal F}_i}.$$ (8)

The image normalization is given by

$$I = \delta \left( \sum_{j=1}^{M_1} g_j - \sum_{j=1}^{M_2}h_j \right),$$ (9)

where $\delta$ is the spacing in $log \tau$ and $M_1$ and $M_2$ points are used to discretize the g and h distributions, respectively. The competing goals of high entropy, good fit, and appropriate image normalization are controled by $\lambda$ and $\alpha$, the Lagrange multipliers used to constrain the goodness-of-fit C and normalization. A Newton-Raphson optimization of Q is performed in MemExp. The Lagrange multipliers are adjusted automatically.

NLS and ML

Upon convergence of the MEM calculation, MemExp performs a series of fits to the kinetics by $n_e$ discrete exponentials using either nonlinear least-squares (NLS) or maximum likelihood (ML) fitting:

$${\cal F}_i = D_0 \sum_{j=1}^{n_e} A_j e^{ -t_i/\tau_j } + \sum_{k=0}^3 B_k (t_i/t_{max})^k,$$ (10)

where $A_j$ and $\tau_j$ are the amplitude and lifetime of the $j^{th}$ exponential, respectively, and the baseline is again approximated as a polynomial. In MemExp 5.0, initial values of the fit parameters are assigned based on the recommended MEM distribution, adding exponentials, one at a time, to represent MEM peaks in decreasing order of peak area.

When deconvolving a known instrument response function, R, these discrete fits take the form:

$${\cal F}_i = D_0 \sum_{j=1}^{n_e} A_j \int_{t_0}^{min(t_i,t_... ...(t_i-t^\prime)/\tau_j } + \sum_{k=0}^3 B_k (t_i/t_{max})^k.$$ (11)

The value of C given by either equation 5 or 6 is calculated for both the continuous and the discrete fits, facilitating a comparison of the distributed and discrete fits. The discrete fits are also characterized by either

$$\chi^2_r = \frac {1}{N-N_p} \sum_{i=1}^{N} \left( \frac { {\cal F}_i-D_i} {\sigma_i} \right) ^{2},$$ (12)

or

$$P_r = \frac {2}{N-N_p} \sum_{i=1}^{N} \left[ D_i ln(D_i/{\cal F}_i) + {\cal F}_i - D_i \right].$$ (13)

where $N_p$ is the number of fit parameters (amplitudes, lifetimes, scattering parameter, baseline coefficients) in the discrete fit. Note: The value of $\chi^2_r$ or $P_r$ (not $\chi^2$ or P) is used in determining which fit is to be recommended by MemExp as the optimal discrete kinetic description.

Some Definitions

The graphical summaries output by MemExp (.out.ps and .exp.ps files) are labeled by several quantities, including the values of C given by equations 5 (or 6) and 12. The correlation length of the residuals $\tau_c$ is also reported for both the MEM and NLS/ML fits. It is estimated by summing the absolute value of the autocorrelation function of the residuals over a span that depends on how slow the slowest process is. The normalization of the fits are also given: I (equation 9) for the MEM fits and the sum of exponential amplitudes, $\Sigma A$, for the discrete fits.

For the MEM fits, the value of Test [8] is also given (denoted as T). Test is zero when the gradients of S and C are parallel. Test is therefore a good measure of convergence when only the Lagrange multiplier $\lambda$ is used (equation 3); small values of Test imply that the maximum-entropy solution has indeed been found for the given value of C. Ratio is another measure of convergence that is recorded in the .out files. Ratio is the norm of the gradient of Q divided by the norm of a unit vector [9]. Unlike Test, Ratio accounts for the effects of constraining the normalization I. Ratio is therefore a good measure of convergence whether or not the Lagrange multiplier $\alpha$ is used (equation 3). If Test and Ratio both approach 1.0, the convergence should be slowed down by reducing EPSL and EPSU.