Monitoring Sufficient-Statistic Stabilization

Why Monitor Sufficient Statistics?

SAEM estimates population parameters using simulated latent effects and a stochastic approximation to conditional sufficient statistics. The M-step uses those statistics to update parameters. Consequently, parameter traces are important, but they are not the only useful record of the algorithm’s progress: a parameter trace can look nearly flat while a statistic feeding the M-step is still changing.

This distinction matters especially when statistics have different scales or when a parameter depends only weakly on a changing statistic. A careful check examines the statistic components themselves, their changes over iterations, and whether the changes are small for meaningful reasons. A flat curve alone does not establish stabilization.

Definition: Let $Z$ denote the latent effects and let $T(y,Z)$ be the vector of complete-data sufficient statistics for observed data $y$. At SAEM iteration $k$, let $S_k$ be the stochastic-approximation estimate of the conditional expectation of $T(y,Z)$ given $y$. The sequence $S_0,S_1,\ldots$ is the sufficient-statistic trajectory. Its components are monitored because the M-step uses $S_k$ to update the population parameters.

For a complete-data model in exponential-family form, its log likelihood can be written as a parameter-dependent term involving the statistic vector $T(y,Z)$, plus a term not depending on the parameters. The conditional expectation of the complete-data log likelihood therefore depends on the latent variables through the conditional expectation of that statistic vector. In this setting, tracking $S_k$ directly checks the quantities that carry the latent-data information into the M-step.

SAEM statistic update: If $T_k$ is the statistic vector computed from the latent-effect draw or draws at iteration $k$, the stochastic-approximation update is $$ S_k=S_{k-1}+\gamma_k(T_k-S_{k-1}) =(1-\gamma_k)S_{k-1}+\gamma_kT_k, $$ where $0<\gamma_k\leq 1$ is the gain for that iteration. With multiple draws, $T_k$ can be their average statistic vector. The update is applied componentwise.

The gain controls how much the current simulated statistic changes the running estimate. A larger gain gives the current draw more influence; a smaller gain gives more weight to the preceding trajectory. This is why the gain schedule must be considered when interpreting a trajectory: a slowly changing $S_k$ can result from a small $\gamma_k$, even if the newly simulated $T_k$ values have not settled near $S_{k-1}$.

What Iteration-to-Iteration Changes Show

Define the change in the statistic vector as $\Delta S_k=S_k-S_{k-1}$. Subtracting $S_{k-1}$ from both sides of the update gives the exact identity

\Delta S_k=\gamma_k(T_k-S_{k-1}).

Thus, the observed change is the gain multiplied by the difference between the current simulated statistic and the previous running estimate. A small $\Delta S_k$ does not by itself show that $T_k$ is close to $S_{k-1}$: the gain may simply be small. Inspecting the statistic trajectory alongside the gain, and where useful the innovation $T_k-S_{k-1}$, helps distinguish these explanations. The innovation itself is based on simulation and should not be treated as a noise-free measure of the conditional expectation.

For components measured on different scales, compare changes using meaningful component-specific scales rather than raw absolute differences alone. One useful descriptive measure is

r_{k,j}=\frac{|S_{k,j}-S_{k-1,j}|}{\max(|S_{k,j}|,a_j)},

where $a_j>0$ is a prespecified reference scale for component $j$, chosen to avoid an unstable ratio when that component is near zero. This is a way to summarize relative movement, not a universal convergence threshold. The scale should be recorded and scientifically interpretable; changing it can change the reported relative change.

Worked Example: Calculating a Running Statistic

Suppose a hypothetical SAEM run has a scalar statistic estimate $S_0=8.00$, uses gain $\gamma_k=0.20$ for three iterations, and produces simulated statistics $T_1=11.0$, $T_2=9.00$, and $T_3=10.5$. Calculate the running values and changes.

At iteration 1, $$ S_1=8.00+0.20(11.0-8.00)=8.60, $$ so $\Delta S_1=8.60-8.00=0.600$.

At iteration 2, $$ S_2=8.60+0.20(9.00-8.60)=8.68, $$ so $\Delta S_2=8.68-8.60=0.0800$.

At iteration 3, $$ S_3=8.68+0.20(10.5-8.68)=9.044, $$ so \(\Delta S_3=9.044-8.68=0.364. $$

The change fell sharply at iteration 2, but increased at iteration 3. The calculated trajectory is not steadily flattening over these three updates. In particular, $T_3=10.5$ is well above $S_2=8.68$, and the gain turns that difference into a change of $0.364$. A single small change would not justify describing the statistic as stabilized.

Monitor Every Informative Component

A statistic vector may contain components for different features of the latent effects: for example, sums or averages related to their means, second moments related to their variability, and cross-products related to covariances. Which components are required depends on the complete-data model and its M-step. Monitoring only one component—or only the population parameters that are easiest to plot—can miss movement in another component that matters to the update.

Worked Example: A Steady Mean Does Not Establish Stable Variability

Consider hypothetical running summaries for 100 latent effects. The mean statistic remains $S_{\text{mean}}=0.10$, while the second-moment statistic moves from 1.20 to 1.24 to 1.28 to 1.32. For this illustration, suppose the corresponding variance update is $S_{\text{second}}-S_{\text{mean}}^2$.

At the first displayed point, the variance estimate is $$ 1.20-(0.10)^2=1.19. $$ At the last point, it is $$ 1.32-(0.10)^2=1.31. $$ The implied standard deviations are approximately $\sqrt{1.19}=1.091$ and $\sqrt{1.31}=1.145$, respectively.

The mean statistic, and therefore a mean parameter directly updated from it, is unchanged in these summaries. But the second moment and the implied variability have moved. A plot showing only the mean parameter would miss that movement; even a broad vertical axis or low display precision could make the variability trace look deceptively flat. The check is to inspect the sufficient statistics used by the update, at adequate numerical precision, rather than infer their stability from one parameter curve.

This example uses a simple moment update to show the diagnostic issue; an actual NLME model can have a more complicated or numerical M-step. The general point remains: establish which statistic components the model’s update uses, and check each one. Earlier material on closed-form and numerical M-steps provides the relevant update forms; the monitoring question here is whether their input statistics continue to move.

Scaled Changes and a Joint Check

Raw changes can be misleading when a statistic vector mixes units or magnitudes. For instance, a change of $0.01$ may be substantial for a correlation-like component near zero but negligible for a sum near 250. Componentwise scaled changes make such differences visible, while retaining the individual components prevents a small average change from hiding one problematic coordinate.

Worked Example: Finding the Component With the Largest Relative Change

Suppose three hypothetical statistic components change from $(0.500,4.00,-0.200)$ to $(0.502,4.04,-0.190)$. Choose reference scales $a=(0.100,1.00,0.100)$, respectively. Calculate $r_{k,j}$ for each component.

For component 1, $$ r_{k,1}=\frac{|0.502-0.500|}{\max(0.502,0.100)} =\frac{0.002}{0.502}\approx 0.003984, $$ or about $0.3984\%$.

For component 2, $$ r_{k,2}=\frac{|4.04-4.00|}{\max(4.04,1.00)} =\frac{0.04}{4.04}\approx 0.009901, $$ or about $0.9901\%$.

For component 3, $$ r_{k,3}=\frac{|-0.190-(-0.200)|}{\max(0.190,0.100)} =\frac{0.010}{0.190}\approx 0.05263, $$ or about $5.263\%$. Component 3 has the largest relative change under these chosen scales, despite its small absolute change. This comparison points to the component needing closer inspection; it does not, by itself, establish that the whole algorithm is or is not stabilized.

Diagnostics, Comparisons, and Limitations

A practical review combines trajectories with their changes. Inspect plots of all statistic components against iteration, paying attention to persistent direction, changes in apparent variability, abrupt shifts, and whether the behavior changes when the gain schedule changes. Also examine componentwise differences or scaled differences over successive iterations or a clearly specified window. A single endpoint comparison can hide a sustained trend or a reversal.

The comparison with parameter traces is complementary rather than competitive. Parameter traces summarize the quantities ultimately being estimated; statistic traces show the changing inputs to the M-step. If parameters look steady but a statistic continues to drift, check the statistic’s role in the update, the display precision and axis range, and whether the M-step is weakly sensitive to that component. Conversely, a statistic may fluctuate while parameter updates appear steadier because the M-step maps its inputs nonlinearly or combines several components.

Raw simulated latent effects answer a different question. Individual draws can vary considerably even when their running sufficient-statistic estimate changes little. Looking at a single draw is therefore not a substitute for examining the relevant statistic trajectory. At the same time, a smooth running trajectory does not prove that the latent-effect simulation is exploring adequately; chain-mixing diagnostics remain necessary.

Trajectory evidence has limits. A statistic’s apparent stability depends on the gain schedule, the simulated draws, the length of the displayed run, and the plotting scale. Small updates late in a run can be mechanically induced by small gains. Nor does a stable statistic trajectory establish that the fitted model is appropriate, that the estimate is a global optimum, or that the result would be reproducible from another run. These are separate questions.

Common pitfall: “The parameter plot is flat and the latest statistic change is tiny, so the latent statistics have stabilized.” This reasoning can fail because the plot may omit an informative statistic, its scale may hide meaningful movement, or a small gain may suppress the update even when the simulated statistic differs from the running estimate. Check every informative component over a sequence of iterations, relate changes to the gain, and interpret scaled differences using stated reference scales. Treat the result as diagnostic evidence, not proof from a single visual feature.

Identify the statistic vector.
List the components used in the model’s M-step, including moments or cross-products relevant to random-effect variability and covariance.

Plot the trajectories.
Use adequate numerical precision and inspect each component’s direction and variation over iterations.

Calculate changes in context.
Compare successive changes and, for components on different scales, use justified reference scales. Check the gain when a change is small.

Compare with parameter traces.
Investigate mismatches rather than assuming either display tells the whole story; verify that latent-effect simulation diagnostics are also acceptable.

Check Your Understanding

Answer each question using the update and diagnostic ideas above. Short solutions follow.

If $S_{k-1}=5.00$, $T_k=7.00$, and $\gamma_k=0.10$, calculate $S_k$ and $\Delta S_k$.
Why can a small $\Delta S_k$ occur even when $T_k$ is not close to $S_{k-1}$?
A plotted mean statistic is nearly constant, but the second-moment statistic has a sustained upward trajectory. What should be checked before describing the run as stabilized?
Why might raw absolute changes be inadequate when comparing components with different scales?
Does a smooth sufficient-statistic trajectory prove that the latent-effect simulation mixes well? Explain.

Solutions. (1) $S_k=5.00+0.10(7.00-5.00)=5.20$, so $\Delta S_k=0.20$. (2) Since $\Delta S_k=\gamma_k(T_k-S_{k-1})$, a small gain can produce a small update despite a larger difference between the simulated statistic and the running estimate. (3) Check whether the second moment feeds the M-step, inspect it at adequate precision over multiple iterations, and consider the gain and scaled changes; a flat mean alone is insufficient. (4) The same absolute change can represent very different relative movement for components with different magnitudes or units; use justified component-specific scales. (5) No. A smooth running estimate does not establish adequate exploration by the latent-effect chain, so chain-mixing diagnostics are still needed.