Tutorials › Nonlinear Mixed-Effects Models (NLME) › Stopping Rules for SAEM Estimation

The EM and SAEM Algorithms · Tutorial 223 of 797

Stopping Rules for SAEM Estimation

Learn to combine scaled change measures and stable-window checks with the SAEM gain schedule to decide when numerical stopping is defensible.

Advanced 11 min read

What You'll Learn

  • Define scaled parameter-change and sufficient-statistic-change measures for SAEM iterations.
  • Explain why a small change at one iteration is not sufficient evidence of stabilization.
  • Combine gain-schedule status with parameter and statistic checks in an operational stopping rule.
  • Choose tolerances and windows that account for parameter scales and stochastic variation.
  • Diagnose misleading stopping signals and distinguish numerical stability from convergence to a global maximum.
  • Apply a stopping rule to hypothetical SAEM traces and interpret its limitations.

Why a Stopping Rule Needs More Than One Signal

SAEM produces a sequence of parameter estimates, but stochastic variation can make that sequence appear stable for a few iterations and then move again. Conversely, a parameter can move slightly in raw numerical units while the movement is important relative to its scale. A defensible stopping rule therefore combines several pieces of evidence: whether the gain schedule is in the intended phase, whether the sufficient statistics have stabilized, and whether parameter changes remain small over repeated windows.

The rule developed here is a practical termination criterion, not a proof that the likelihood has been maximized or that the estimate is the global maximum. Even a trace that passes all numerical checks can be near an undesirable stationary point or reflect weak identification. The purpose is to avoid stopping while the algorithm is visibly unsettled, not to replace model assessment or comparisons across runs.

Definition: Let \(\theta_k\) be the vector of population and error parameters after SAEM iteration \(k\), \(S_k\) the vector of recursively updated sufficient statistics, and \(\gamma_k\) the gain used at iteration \(k\). A practical stopping rule specifies (1) a gain-schedule gate, (2) change measures for \(\theta_k\) and \(S_k\), and (3) a requirement that those measures satisfy chosen tolerances across repeated checking windows.

The schedule matters because SAEM updates combine new simulated information with accumulated statistics. In a common representation, \(S_k=(1-\gamma_k)S_{k-1}+\gamma_k s_k\), where \(s_k\) is the current simulated sufficient-statistic contribution. A large gain can make one iteration heavily dependent on the latest simulation; a smaller gain reduces that iteration's influence but does not itself establish stability. The earlier material on Robbins–Monro gain sequences and gain-schedule choice supplies the standard conditions for diminishing-gain schemes, such as \(\sum_k\gamma_k=\infty\) and \(\sum_k\gamma_k^2<\infty\), under the relevant regularity assumptions. A stopping decision should respect the chosen schedule and burn-in transition rather than treat a particular small value of \(\gamma_k\) as a universal stopping threshold.

Scaled Changes and Stable Windows

Raw differences are difficult to compare across parameters measured on different scales. For each parameter \(\theta_j\), choose a positive reference scale \(r_j\), reflecting a meaningful scale for that parameter, and a dimensionless tolerance \(\epsilon_{\theta,j}\). A denominator combining the parameter's magnitude and its reference scale avoids division by a value near zero:

$$ d_{\theta,j}(k;h)= \frac{|\theta_{j,k}-\theta_{j,k-h}|} {\max(|\theta_{j,k-h}|,r_j)}. $$

Here \(h\) is the number of iterations in a comparison window. The parameter-change check passes only if \(d_{\theta,j}(k;h)\leq\epsilon_{\theta,j}\) for every monitored parameter. A parameter near zero is judged relative to \(r_j\), not by an unstable relative percentage. The scale and tolerance must be selected for the application; there is no universal cutoff that makes a change negligible in every NLME model.

Sufficient statistics can have different units and magnitudes as well. For each statistic \(S_\ell\), choose a positive reference scale \(q_\ell\) and tolerance \(\epsilon_{S,\ell}\), then compute

$$ d_{S,\ell}(k;h)= \frac{|S_{\ell,k}-S_{\ell,k-h}|} {\max(|S_{\ell,k-h}|,q_\ell)}. $$

Use component-specific tolerances or a justified common tolerance after scaling. A single unscaled norm can conceal an important change in a small-magnitude component or be dominated by a large-magnitude one. Because a window-endpoint difference can miss movement within the window, it is also useful to inspect the range of each scaled parameter and statistic inside each window. A trace that rises and falls back to its starting value could have a small endpoint difference without being stable.

A repeated-window criterion addresses temporary quiet periods. Divide the post-burn-in portion into nonoverlapping windows of length \(h\); at each window end, compare its endpoint with its starting point and inspect within-window ranges. Require both parameter and statistic checks to pass in \(R\) consecutive windows, where \(R\) is specified in advance. Consecutive windows provide repeated evidence of stability; they do not make the chance of a misleading signal zero, particularly when adjacent iterations are dependent.

Worked Example: Scaled Parameter Changes

Suppose two hypothetical parameters are monitored at the end of a 20-iteration window. A rate parameter changes from \(0.50\ \mathrm{h}^{-1}\) to \(0.51\ \mathrm{h}^{-1}\), and its reference scale is \(r_1=0.10\ \mathrm{h}^{-1}\). A near-zero correlation parameter changes from \(0.020\) to \(0.024\), with \(r_2=0.10\). For illustration, the chosen tolerance is \(0.03\) for both scaled changes.

For the rate parameter,

$$ d_{\theta,1}= \frac{|0.51-0.50|}{\max(0.50,0.10)} =\frac{0.01}{0.50}=0.020. $$

Thus its scaled change is 0.020, below the 0.03 tolerance. For the correlation parameter,

$$ d_{\theta,2}= \frac{|0.024-0.020|}{\max(0.020,0.10)} =\frac{0.004}{0.10}=0.040. $$

This change exceeds the tolerance, so the parameter check fails for this window. Using the raw relative change \(0.004/0.020=0.20\) would give a different, scale-sensitive assessment; the specified reference scale makes the rule explicit. Neither result alone determines whether 0.004 is scientifically important: that is why the reference scale and tolerance must be justified for the model.

Constructing and Applying a Stopping Rule

A useful rule has a schedule gate and multiple stability checks. The gain-schedule gate should confirm that the planned burn-in or gain transition has finished and that the algorithm is operating in the intended estimation phase. For a diminishing-gain scheme, the schedule should be progressing as planned; a small \(\gamma_k\) is not sufficient by itself. In a constant-gain scheme, persistent stochastic fluctuation may remain, so a rule designed for diminishing-gain convergence cannot simply be transferred unchanged.

After the gate passes, monitor sufficient statistics and parameters over the same windows. A statistics check matters because the M-step maps sufficient statistics to parameter estimates: parameters may barely change because the mapping is locally insensitive even while the statistics still drift. Conversely, statistics may be stable while a numerical M-step or parameter constraint causes a parameter to keep moving. Requiring both checks reduces these blind spots.

1
Respect the schedule.
Do not evaluate stopping windows until the planned burn-in and gain transition are complete. Verify that the gain sequence is following its specified form.
2
Set scales and tolerances.
Specify \(r_j\), \(q_\ell\), and component-specific tolerances before interpreting the trace. Record the units and rationale for each reference scale.
3
Check both kinds of change.
At each window, evaluate scaled parameter and sufficient-statistic endpoint changes, and inspect within-window ranges for reversals or excursions.
4
Require repetition and review.
Stop only after both checks pass for the prespecified number \(R\) of consecutive windows, with no visible trend or other instability that contradicts the numerical summaries.

Worked Example: Statistics Pass, Parameters Do Not

Consider a hypothetical sufficient statistic \(S_1\) with reference scale \(q_1=10\) and tolerance \(0.02\). Across a 10-iteration window it changes from 50.0 to 50.6. A second statistic \(S_2\), with \(q_2=1\), changes from 2.00 to 2.05. The corresponding parameter values change from \(\theta_1=1.00\) to \(1.01\), and from \(\theta_2=0.20\) to \(0.24\). Use reference scales \(r_1=0.50\) and \(r_2=0.10\); use a parameter tolerance of 0.03.

The statistic changes are

$$ d_{S,1}=\frac{0.6}{\max(50.0,10)}=0.012, \qquad d_{S,2}=\frac{0.05}{\max(2.00,1)}=0.025. $$

The first statistic passes, but the second fails the 0.02 tolerance. For the parameters,

$$ d_{\theta,1}=\frac{0.01}{\max(1.00,0.50)}=0.010, \qquad d_{\theta,2}=\frac{0.04}{\max(0.20,0.10)}=0.20. $$

The first parameter passes and the second fails the 0.03 tolerance. Therefore neither the full statistics check nor the full parameter check passes, and this window cannot count as stable. A check based only on \(\theta_1\), or only on \(S_1\), would incorrectly suggest that the whole fit is ready to stop.

Worked Example: Repeated Stable Windows and the Gain Gate

Suppose a hypothetical run has completed its planned gain transition. The gain values at three consecutive window ends are \(0.020\), \(0.018\), and \(0.016\), consistent with the specified diminishing schedule. The planned stopping rule requires both change checks to pass in three consecutive 25-iteration windows. The maximum scaled parameter and statistic changes at the three ends are respectively \(0.018, 0.024, 0.019\) and \(0.012, 0.017, 0.015\); their respective tolerances are \(0.03\) and \(0.02\). Assume that all within-window ranges are also below their prespecified limits.

At each endpoint, the parameter maximum is at most \(0.024<0.03\), and the statistic maximum is at most \(0.017<0.02\). Both checks pass in all three consecutive windows, the transition is complete, and the gain values follow the planned schedule. Under this operational rule, the run meets the stopping criterion at the third window end, 75 iterations after the first of these windows began. The conclusion is that the monitored sequence met the chosen numerical stability requirements, not that convergence to a global maximum has been demonstrated.

Justification, Diagnostics, and Alternatives

The reason to require stable windows can be seen directly from the update structure. If \(\gamma_k\) is appreciable, one new simulated statistic \(s_k\) can move \(S_k\), and the M-step can then move \(\theta_k\). A single small endpoint difference says only that the net change over that particular interval was small. It does not rule out a subsequent movement, nor does it reveal large changes that canceled within the interval. Repeating the checks over nonoverlapping windows and examining their ranges addresses both issues: the rule asks for sustained small net changes and no large excursions during the windows.

This reasoning justifies the structure of the rule, but no finite collection of windows proves asymptotic convergence. Simulation draws can be dependent, and a trace may remain temporarily quiet by chance. Window length and repetition count are practical design choices: short windows react quickly but may be dominated by noise; longer windows can reveal slower drift but cost more iterations and can obscure localized movement if only endpoints are used.

Compare the combined rule with common alternatives. A fixed iteration cap guarantees a computational budget, but it gives no evidence that a particular run has stabilized. A parameter-only threshold is easy to compute, yet it can miss drifting sufficient statistics or be distorted by scales. A statistic-only threshold can miss movement introduced by numerical maximization in the M-step. Monitoring an observed-data likelihood or another objective can provide a useful additional diagnostic, but a noisy or approximate objective may fluctuate, and a flat objective need not imply stable estimates. These methods can supplement the combined rule; none makes its limitations disappear.

Warning signs include a parameter trend that persists across windows, large within-window excursions, one component repeatedly failing its tolerance, a gain schedule that has not completed its planned transition, or a numerical update that changes despite stable statistics. In these cases, do not relax tolerances merely to force a stop. Investigate scaling, the gain schedule, the statistic updates, and the numerical M-step. For nonmonotone traces, compare windows and retain the trace rather than relying on a single aggregate number.

Common pitfall: “The parameter changed by less than 1% in the last iteration, so SAEM has converged.” This reasoning treats one iteration as representative, ignores parameter units and near-zero values, and says nothing about sufficient-statistic drift or the gain phase. Instead, use prespecified scaled tolerances, check statistics as well as parameters, inspect within-window behavior, and require repeated passes after the planned schedule transition. Passing this rule supports a numerical stopping decision; it does not certify a global optimum.

Check Your Understanding

Answer each question using the scales and tolerances stated. Full solutions follow the questions.

  1. A parameter changes from 4.0 to 4.1 over a window. Its reference scale is 1.0 and its tolerance is 0.03. Does the scaled parameter check pass?
  2. Why can a small endpoint change fail to show that a statistic was stable throughout a window?
  3. A run passes its parameter check in three consecutive windows, but a sufficient statistic fails its tolerance in the third. Does the combined rule permit stopping?
  4. What can be concluded if all numerical checks pass after the gain transition? Name one important conclusion that cannot be drawn from those checks alone.

Solutions. (1) The scaled change is \(|4.1-4.0|/\max(4.0,1.0)=0.1/4.0=0.025\), which is below 0.03, so it passes. (2) The statistic may have moved substantially and then returned near its starting value; the endpoint difference measures net change only. Inspect within-window values or ranges as well. (3) No. The combined rule requires both checks to pass in every required consecutive window, so the failed statistic check interrupts the sequence of passes. (4) One may conclude that the monitored changes met the prespecified numerical stability criteria under the chosen schedule. These checks alone do not prove convergence to a global maximum.