Why Numerical Diagnostics Matter
A joint model combines a longitudinal model with a time-to-event model, often linking them through subject-specific random effects. In maximum-likelihood fitting, the likelihood typically involves integrating over those random effects. The software must therefore optimize a likelihood that may be difficult to evaluate accurately, particularly when parameters are weakly identified, near a boundary, or associated with sharply varying integrands.
A message such as “converged” describes an optimizer’s stopping rule; it does not establish that the reported estimates are accurate, that the likelihood was evaluated adequately, or that the fitted model is appropriate. Conversely, a warning does not always mean the fit is unusable: some warnings reflect a conservative numerical tolerance or a parameter near a legitimate boundary. The goal is to identify the warning’s source, investigate it with a targeted check, and report the remaining uncertainty honestly.
The diagnostics below concern numerical behavior while fitting a specified model. They do not replace checking whether the software fit the intended model, which is a separate task. The basic distinction is between optimization error (the parameters may not be at a solution), integration error (the computed likelihood may be inaccurate), and weak identification or boundary behavior (the data may not support precise estimates for some parameters).
Gradients, Stopping Messages, and Curvature
Let \(\ell(\theta)\) be the log-likelihood as a function of parameter vector \(\theta\). For an interior maximum, a necessary condition is that the gradient, the vector of first derivatives, is zero. In practice, software stops when a numerical measure of the gradient or parameter change meets a tolerance. Because parameters can have different units and scales, inspect a scaled gradient when the software provides one; a raw gradient of a given size has no universal interpretation independent of parameterization.
Near a regular interior maximum, the Hessian matrix of second derivatives of the log-likelihood should be negative definite. Equivalently, the observed information, the negative Hessian, should be positive definite. A non-positive-definite or nearly singular information matrix can indicate a saddle point, a flat likelihood direction, parameters that are difficult to distinguish, or numerical error. These are warning signs to investigate, not diagnoses on their own.
This local quadratic approximation explains the checks. If \(g(\theta)\) is not close to zero, a small movement may still increase the log-likelihood. If \(H(\theta)\) has a nearly zero eigenvalue, the likelihood changes very little along some direction, so small data or numerical changes may move the estimates substantially. The approximation is local: it does not prove that the optimizer found the global maximum.
Worked Example: A “Converged” Fit With a Large Gradient
Suppose software reports convergence for a hypothetical joint model, but its scaled gradient components for three parameters are \(0.0002\), \(0.018\), and \(-0.006\). The fitting procedure used a stopping tolerance of \(0.001\) for the largest absolute scaled component.
Calculate the largest component: \( \max(0.0002, |0.018|, |-0.006|)=0.018 \). Since \(0.018>0.001\), the stated gradient criterion is not met, even though the software printed “converged.” The second parameter is the main contributor. The result is evidence that the stopping message alone is insufficient; it does not establish that the estimate is wrong by a particular amount.
A focused response is to refit from the reported estimates with a larger iteration limit and, if available, a more reliable gradient calculation or tighter convergence control. Then compare the log-likelihood, estimates, gradient, and curvature. If the gradient falls below tolerance and the estimates stabilize, the original run likely stopped prematurely. If the gradient remains large or the estimates move materially, the fit remains numerically unresolved and should not be treated as a settled maximum.
A useful check is whether independent starts lead to similar solutions. Start values should be plausible and should vary enough to reveal whether the optimizer is sensitive to its starting point. Similar estimates and log-likelihoods, together with acceptable gradients and curvature, support numerical stability. Different solutions may reflect multiple local maxima, a flat likelihood, or optimization failure; matching output from two starts is reassuring but not a proof of global optimality.
Boundary Estimates and Nearly Singular Covariance Matrices
Variance parameters are constrained to be nonnegative. An estimate at or very near zero can be a legitimate boundary estimate, but it may also signal that the data provide little evidence for that variance component or that the model is difficult to estimate. Correlations in a random-effects covariance matrix are constrained to the interval from \(-1\) to \(1\). A correlation extremely close to either limit can make the covariance matrix nearly singular and can make standard errors or optimization unstable.
For a two-dimensional random-effects covariance matrix with variances \(\sigma_1^2,\sigma_2^2\) and correlation \(\rho\), its determinant is \(\sigma_1^2\sigma_2^2(1-\rho^2)\). If this determinant is close to zero, the matrix has a direction with very little variation, so calculations involving its inverse can be sensitive to rounding and small changes in estimates.
Worked Example: A Correlation Near the Boundary
In a hypothetical fit, the estimated random-intercept variance is \(0.40\), the estimated random-slope variance is \(0.25\), and their estimated correlation is \(0.999\). The determinant of the estimated covariance matrix is \(0.40(0.25)(1-0.999^2)\).
First, \(0.40(0.25)=0.10\), and \(1-0.999^2=1-0.998001=0.001999\). Therefore the determinant is \(0.10(0.001999)=0.0001999\). This very small determinant indicates near-singularity; it is not, by itself, a formal cutoff for rejecting a fit.
A targeted check is to inspect the reported covariance matrix and its eigenvalues, then refit from multiple reasonable starts and compare the correlation, likelihood, and standard errors. If a simpler covariance structure is scientifically defensible, it can be fitted as a separate model, but the change should be justified rather than used merely to force a warning to disappear. A zero variance or boundary correlation should be reported as such, with appropriate caution about ordinary interior-based standard errors and tests.
Do not interpret a boundary estimate as automatic proof that the corresponding random effect is absent. The estimate can be affected by limited information, a small sample, highly correlated random effects, or numerical instability. Equally, do not repeatedly adjust optimizer controls to push a variance away from zero without a substantive reason. The objective is to determine whether the fitted solution is stable and interpretable under the stated model.
Integration Warnings and Accuracy Checks
In many maximum-likelihood joint models, the likelihood contribution for a subject is obtained by integrating the conditional contribution over the random-effects distribution. Numerical integration approximates this integral, for example using adaptive Gaussian quadrature. An integration warning may indicate that the approximation did not meet its requested accuracy, that the integrand was difficult to evaluate, or that an internal calculation encountered an extreme value. The exact meaning depends on the software; record the message and consult its documentation rather than assuming all warnings have the same cause.
A practical accuracy check is to repeat the fit with increased integration accuracy, such as a higher quadrature order when the software supports it. Compare the log-likelihood and key estimates across settings. A change in the last displayed decimal place may be inconsequential; a material shift in an association estimate or standard error calls the original approximation into question. The comparison is informative only if other settings are held fixed and the fits have otherwise converged.
Worked Example: Checking Quadrature Stability
Consider a hypothetical maximum-likelihood fit using the same data, model, optimizer controls, and starting values. At quadrature order 7, the log-likelihood is \(-812.436\) and an association estimate is \(0.58\). At order 15, the corresponding values are \(-812.401\) and \(0.59\).
The log-likelihood difference is \(|-812.401-(-812.436)|=0.035\), and the association estimate changes by \(0.59-0.58=0.01\). Those changes are small in the displayed values, but their practical importance depends on the inferential scale and the precision required; the numbers alone do not establish an acceptable universal threshold.
A focused follow-up is to increase the order again, for example to 25, while keeping the other settings fixed. If the log-likelihood and key estimates remain stable to the precision relevant for the analysis, that supports the adequacy of the approximation. If they continue to shift materially, investigate the software’s integration guidance, parameter scaling, and whether the optimizer is converging at each order. Do not attribute a difference to integration accuracy if the fits also stopped at different parameter values.
A Targeted Troubleshooting Sequence
A reliable response to warnings changes one relevant feature at a time. This makes it possible to tell whether a refit addressed the suspected issue. Keep the original fit and its full output so that the comparison is auditable.
Record whether it concerns optimizer termination, gradients or Hessian, a boundary estimate, or numerical integration. Save the exact message, software version, and fitting controls.
Inspect convergence status, scaled gradients, curvature or information-matrix warnings, parameter estimates, and any reported boundary constraints. A single reassuring item does not override a conflicting diagnostic.
For a large gradient, try more iterations or an alternative optimizer; for start sensitivity, use several plausible starts; for integration warnings, increase integration accuracy while holding other settings fixed.
Compare log-likelihoods, estimates, standard errors, gradients, and relevant warnings. State what changed and whether the evidence supports a stable numerical solution.
Changing the optimizer, starting values, quadrature order, parameterization, and model structure all at once may yield a different answer, but it does not reveal why. If a targeted refit still fails, report the unresolved issue and consider whether the available data support the complexity being estimated. A successful computation is not a reason to suppress a warning that remains relevant.
Alternatives, Limitations, and Common Pitfalls
For maximum-likelihood fitting, an alternative optimizer can help determine whether a particular algorithm is struggling; it is a cross-check, not a guarantee of correctness. Multiple starts check sensitivity to initialization. Increased quadrature accuracy probes integration error. These checks address different failure modes and are not interchangeable. Bayesian MCMC fitting uses different diagnostics, such as chain mixing and convergence assessments; an optimizer gradient is not an MCMC convergence diagnostic. A reassuring result from one computational approach does not erase concerns in another.
Numerical diagnostics have limitations. Thresholds depend on parameter scaling, software implementation, and the required precision. A Hessian warning may result from weakly identified parameters, boundary behavior, or inaccurate numerical derivatives. Stable estimates across several runs do not prove that the model is scientifically appropriate or that all uncertainty has been captured.
Common mistake: “The software said converged, so the fit is valid.” This reasoning confuses an algorithm’s stopping condition with evidence about the quality of the solution. Instead, check the gradient and curvature, examine boundary estimates, investigate any integration warning, and compare appropriately targeted refits. The reverse mistake is to treat every warning as proof that all estimates are unusable; determine what the message concerns and whether the relevant diagnostic changes the result.
Check Your Understanding
Answer each question, then compare your reasoning with the solutions.
- A fit reports convergence, but the largest absolute scaled gradient is \(0.012\) and the stated tolerance is \(0.001\). What should you conclude, and what is a focused next step?
- Why can a random-effects correlation near \(1\) be a numerical concern even when both estimated variances are positive?
- A higher quadrature order changes both the log-likelihood and estimates. What should be held fixed, and what should be compared before drawing a conclusion?
- Does a variance estimate of zero prove that the corresponding random effect is scientifically absent? Explain.
- Why is “use another optimizer” not a complete response to an integration warning?
Solutions. (1) The stated gradient criterion is not met because \(0.012>0.001\); convergence status alone is insufficient. Refit from the reported estimates with more iterations or a more reliable gradient calculation, then inspect whether the gradient and estimates stabilize. (2) The covariance determinant includes the factor \(1-\rho^2\), which approaches zero as \(\rho\) approaches \(1\); the matrix can therefore be nearly singular and sensitive to small changes. (3) Hold the model, data, optimizer controls, and starting values fixed as far as possible. Compare log-likelihoods, estimates, standard errors, and convergence diagnostics at each order; otherwise changes may have more than one cause. (4) No. Zero is a boundary estimate that may reflect limited information or numerical and identification difficulties; it is not by itself proof of scientific absence. (5) An optimizer change targets the optimization algorithm, whereas an integration warning concerns likelihood approximation. Increase integration accuracy and compare results while maintaining adequate optimization diagnostics.