Tutorials › Joint Modeling of Longitudinal and Time-to-Event Data › Gradient Checks At The Reported Estimate

Convergence and Identifiability · Tutorial 350 of 697

Gradient Checks At The Reported Estimate

Use scaled scores and central finite differences to check whether a reported joint-model estimate is locally stationary and whether its computed gradient is credible.

Advanced 11 min read

What You'll Learn

  • Define score components for the objective actually optimized by a joint-model fit.
  • Scale gradient components so their magnitudes reflect comparable parameter changes.
  • Compute and interpret central finite-difference checks at a reported estimate.
  • Recognize why small objective changes alone do not establish convergence.
  • Check numerical stability across step sizes and identify common gradient-check failures.
  • Distinguish interior stationarity checks from constraints and boundary solutions.

Why Check the Gradient at the Reported Estimate?

An optimizer can report success because its stopping rule was met, because successive objective values changed very little, or because its internal gradient appeared small. Those statements are not interchangeable. In a joint model, a numerically small change in the objective can arise from a tiny step, a flat likelihood direction, or limited numerical precision—even when a further feasible change in a parameter would still improve the objective. A gradient check asks a more direct question: at the reported estimate, does the objective have a non-negligible local slope in any parameter direction?

This check concerns the objective used by the fitting procedure. It might be the observed-data log likelihood, a penalized log likelihood, or an approximation to one of these. A score computed for a different objective, or under different integration settings, does not verify the reported fit. The earlier material on multiple starting values and competing optima addresses whether different runs locate different solutions; a gradient check examines local stationarity at the particular estimate being reported.

Definition: Let \(\ell(\theta)\) be the differentiable objective being maximized, with parameter vector \(\theta=(\theta_1,\ldots,\theta_p)\). Its score is \(U(\theta)=\nabla\ell(\theta)\), with component \(U_j(\theta)=\partial\ell(\theta)/\partial\theta_j\). At an interior local maximum \(\hat{\theta}\), a necessary condition is \(U_j(\hat{\theta})=0\) for every \(j\). For a minimized objective, such as \(-\ell(\theta)\), the corresponding gradient must be zero at an interior optimum.

In a shared random-effects joint model, \(\theta\) may include longitudinal fixed effects, survival-model coefficients, association parameters, baseline-hazard parameters, and random-effects covariance parameters. The score is with respect to the objective as a function of all of these parameters. A small score in one component does not excuse a large score in another: the full vector matters. The zero-gradient condition is necessary for a differentiable interior optimum, not sufficient to establish a maximum or to establish that the model is identifiable.

Scaling Score Components

Raw score components have units that depend on the parameterization. For example, a derivative with respect to a slope per year and a derivative with respect to a log standard deviation are not directly comparable by magnitude. Scale each component using a meaningful parameter increment \(a_j>0\), expressed in the units of \(\theta_j\). Define the scaled score component as \(a_j U_j(\hat{\theta})\). It represents the first-order change in the objective for a change of about \(a_j\) in parameter \(j\), with the other parameters held fixed.

$$ \theta_j^{*}=\frac{\theta_j-\hat{\theta}_j}{a_j}, \qquad \left.\frac{\partial \ell}{\partial \theta_j^{*}}\right|_{\hat{\theta}} =a_j U_j(\hat{\theta}). $$

The scale \(a_j\) should be stated and chosen to make comparisons meaningful—for example, a scientifically interpretable increment or a representative parameter scale. It is not a universal constant, and scaled scores are not automatically comparable across analyses that use different choices. They are also not standard errors or uncertainty measures. A practical assessment compares each scaled component with a stated numerical tolerance and with the precision of the objective evaluation. There is no single tolerance that applies to every model and optimizer.

Worked Example: Finding a Non-negligible Scaled Score

Suppose a hypothetical fitted model has three interior parameters: a longitudinal slope \(\beta\) measured per year, a log-scale covariance parameter \(\gamma\), and a survival association parameter \(\alpha\). At the reported estimate, suppose the analytic score is \(U=(1.8,-0.015,0.42)\), in the parameter order \((\beta,\gamma,\alpha)\). Choose illustrative increments \(a=(0.10,10,0.50)\) in the corresponding parameter units.

The scaled components are \(0.10(1.8)=0.18\), \(10(-0.015)=-0.15\), and \(0.50(0.42)=0.21\). Thus the scaled score vector is \((0.18,-0.15,0.21)\). With an illustrative componentwise tolerance of \(0.01\), none of these components is small enough to support a claim that the estimate is stationary at that tolerance. In particular, the raw covariance score \(-0.015\) might look small beside \(1.8\), but its scaled value is \(-0.15\).

This conclusion is a diagnostic, not a statement that any coefficient is statistically significant. The scales and tolerance must be reported, and a large score can indicate incomplete optimization, a coding or derivative problem, or a poorly behaved objective. This calculation alone does not identify which explanation applies.

Checking the Score by Finite Differences

An analytic or automatically differentiated score should agree with a numerical approximation to the derivative of the same objective. For parameter \(j\), the central finite-difference estimate at the reported estimate is

$$ U_{j,\mathrm{FD}}(\hat{\theta};h_j) = \frac{ \ell(\hat{\theta}+h_j e_j)-\ell(\hat{\theta}-h_j e_j) }{2h_j}, $$

where \(e_j\) is the vector with a 1 in position \(j\) and zeros elsewhere, and \(h_j>0\) is a small step in the units of that parameter. The check compares \(U_{j,\mathrm{FD}}\) with the analytic score \(U_j(\hat{\theta})\), preferably on the scaled parameter scale. One useful relative discrepancy is

$$ D_j= \frac{ \left|a_j U_j(\hat{\theta})-a_j U_{j,\mathrm{FD}}(\hat{\theta};h_j)\right| }{ \max\left\{1,\left|a_j U_j(\hat{\theta})\right|, \left|a_j U_{j,\mathrm{FD}}(\hat{\theta};h_j)\right|\right\} }. $$

A small \(D_j\) indicates agreement for that step size, but it is not enough to use one step blindly. If \(h_j\) is too large, curvature makes the finite difference a poor approximation to the derivative at the center. If it is too small, rounding error or noise in numerical integration can dominate the difference between two objective evaluations. Repeat the calculation at several sensible step sizes, often including steps that differ by factors of about ten. A credible check shows agreement over a range, not just a fortuitous match at one step.

Why the Central Difference Works

The justification follows from a Taylor expansion. If \(\ell\) has sufficiently smooth third derivatives near \(\hat{\theta}\), expanding in coordinate \(j\) gives

$$ \begin{aligned} \ell(\hat{\theta}+h_j e_j) &=\ell(\hat{\theta})+h_j U_j(\hat{\theta}) +\frac{h_j^2}{2}\frac{\partial^2\ell}{\partial\theta_j^2}(\hat{\theta}) +O(h_j^3),\\ \ell(\hat{\theta}-h_j e_j) &=\ell(\hat{\theta})-h_j U_j(\hat{\theta}) +\frac{h_j^2}{2}\frac{\partial^2\ell}{\partial\theta_j^2}(\hat{\theta}) +O(h_j^3). \end{aligned} $$

Subtracting cancels both the common objective value and the second-order term. Dividing by \(2h_j\) leaves \(U_j(\hat{\theta})+O(h_j^2)\). Thus the central difference has second-order truncation error under this smoothness condition. In computation, however, the two objective evaluations also have numerical error; as \(h_j\) shrinks, that error is magnified by division by \(2h_j\). This explains why a range of step sizes is informative and why the objective must be evaluated consistently.

Worked Example: Comparing Analytic and Numerical Scores

Consider a hypothetical reported value \(\hat{\alpha}=0.30\) for an association parameter, with scale \(a_\alpha=0.50\). The fitted objective is evaluated at \(\hat{\alpha}-0.01\) and \(\hat{\alpha}+0.01\), giving \(\ell(0.29)=812.3448\) and \(\ell(0.31)=812.34624\). Suppose the analytic score at \(0.30\) is \(0.0718\).

The central finite difference is \((812.34624-812.3448)/(2(0.01))=0.00144/0.02=0.0720\). On the scaled scale, the analytic and finite-difference values are \(0.50(0.0718)=0.0359\) and \(0.50(0.0720)=0.0360\). Their discrepancy is \(|0.0359-0.0360|/\max(1,0.0359,0.0360)=0.0001\). They agree closely at this step size. The scaled score is also small relative to the illustrative \(0.01\) tolerance used here, so this component supports local stationarity, subject to checks at other step sizes and for the other parameters.

For example, repeating with \(h_\alpha=0.001\) and \(0.05\) can reveal whether the agreement persists. The same integration approximation and convergence settings must be used for all three objective evaluations; otherwise the finite difference can reflect changes in numerical settings rather than changes in \(\alpha\).

Small Objective Changes Are Not a Stationarity Test

An optimizer may report that the objective barely changed between iterations. That fact is ambiguous unless the step size and parameter scale are also considered. For a small step \(\Delta\theta\), the first-order change is approximately \(U(\theta)^\mathsf{T}\Delta\theta\). Even a sizeable score can produce a tiny change if the step is tiny or nearly perpendicular to the score. Flat curvature can also yield little objective change over a visible distance while leaving weakly determined directions. These possibilities are reasons to inspect the gradient; they do not by themselves prove that the estimate is wrong.

Worked Example: A Tiny Change with a Large Slope

Suppose a hypothetical joint-model objective at a reported association estimate \(\hat{\alpha}\) is \(-1200.000000\). At steps of \(h=0.000001\), its values are \(\ell(\hat{\alpha}+h)=-1199.999996\) and \(\ell(\hat{\alpha}-h)=-1200.000004\). The objective changes by only a few millionths around its reported value, which can look negligible when viewed on a scale near 1200.

But the central difference is \(\{(-1199.999996)-(-1200.000004)\}/(2(0.000001))=0.000008/0.000002=4\). With an illustrative scale \(a_\alpha=0.10\), the scaled score is \(0.10(4)=0.40\), which is not small relative to \(0.01\). The objective changes are small in absolute size because the step is tiny; the derivative indicates a meaningful local slope on the chosen parameter scale. Reporting only the small objective change would therefore give an incomplete and potentially misleading convergence assessment.

Practical Checks, Failure Modes, and Limits

A useful gradient assessment checks every free parameter, not merely the parameters of greatest scientific interest. For a joint model, inspect longitudinal, survival, association, and covariance parameters, including parameters represented internally on transformed scales. If the optimizer works with a transformed parameter, such as a log variance, differentiate with respect to that optimization coordinate when assessing its stopping condition. A derivative in the original coordinate is related by the chain rule and may have a different magnitude.

  • Use the exact fitted objective. For an approximate marginal likelihood, use the same approximation and settings in the analytic and finite-difference evaluations.
  • Check step-size stability. A derivative that changes greatly across reasonable \(h_j\) values may reflect curvature, numerical noise, or a coding problem.
  • Check signs and parameter ordering. Maximizing \(\ell\) and minimizing \(-\ell\) reverse the gradient sign. A mismatch can also arise if score components are paired with the wrong parameters.
  • Respect constraints. At a boundary, a two-sided perturbation may be infeasible. The interior zero-score condition need not hold there; assess feasible one-sided directions and the relevant constrained optimality conditions instead.
  • Separate stationarity from identification. A near-zero score says the objective is locally flat to first order at the reported point. It does not show that the point is a maximum, that it is unique, or that its parameters are precisely estimated.

A finite-difference mismatch can come from an incorrect analytic derivative, a parameter transformation, an overly small or large step, noisy numerical integration, or an objective with limited smoothness. Repeatable disagreement should be investigated rather than averaged away. Conversely, agreement between the two gradients verifies the derivative calculation more than it verifies the fit: both can agree at a point that is not a maximum. The earlier material on likelihood profiles and identifiability addresses complementary questions about weakly informed directions. Hessian-based curvature and uncertainty diagnostics are a separate check and are not established by a small gradient.

Key takeaway: At an interior reported estimate, assess all score components on stated parameter scales, then verify the computed derivatives by finite differences over multiple step sizes. Small iteration-to-iteration objective changes alone do not establish stationarity.

Check Your Understanding

Answer each question using the objective, parameter scale, and constraint status specified.

  1. An interior parameter has score \(U_j=0.06\) and scale \(a_j=0.2\). What is its scaled score, and what first-order change in the maximized objective does it represent for an increase of \(0.2\) in that parameter?
  2. Why can a central finite difference be inaccurate when its step is extremely small?
  3. A minimization routine optimizes \(Q(\theta)=-\ell(\theta)\). How should its gradient relate to the score for \(\ell\)?
  4. A constrained variance estimate is at its lower boundary. Why might a two-sided finite difference and the ordinary zero-score rule be inappropriate?

Solutions. (1) The scaled score is \(0.2(0.06)=0.012\). To first order, increasing the parameter by \(0.2\) changes \(\ell\) by about \(0.012\), with other parameters fixed. (2) The difference of the two objective values can be dominated by rounding or numerical evaluation error; dividing that error by a very small \(2h_j\) can make the derivative unstable. (3) Since \(Q=-\ell\), \(\nabla Q=-\nabla\ell=-U\); both gradients are zero at the same interior stationary point, but their signs are opposite away from it. (4) A negative perturbation may violate the variance constraint, and an optimum at a boundary need not have zero derivative in the unconstrained coordinate. Use feasible one-sided directions and constrained optimality conditions instead.

Clinical Biostats