Tutorials › Joint Modeling of Longitudinal and Time-to-Event Data › Separating Parameter and Random-Effects Uncertainty

Dynamic Prediction · Tutorial 525 of 697

Separating Parameter and Random-Effects Uncertainty

Learn how random-effects uncertainty and parameter uncertainty enter a dynamic risk prediction, and how to distinguish their contributions.

Advanced 12 min read

What You'll Learn

  • Write a landmark risk conditional on the patient being event-free at the landmark.
  • Derive the random-effects distribution conditional on the patient’s observed history.
  • Integrate risk over random effects and then over uncertain model parameters.
  • Distinguish Bayesian parameter integration from a frequentist approximation.
  • Separate parameter and random-effects contributions to variation in predicted risk.
  • Recognize when plug-in predictions can conceal important uncertainty.

Two Different Sources of Uncertainty

A dynamic prediction for a patient uses both a fitted joint model and that patient’s observed biomarker history. Even if the model form is appropriate, uncertainty remains about the patient’s latent trajectory, represented by shared random effects. There is also uncertainty about the model parameters, which were estimated from finite data. These are distinct sources: learning more about a patient’s biomarker history can reduce uncertainty about that patient’s random effects, while a larger or more informative model-fitting data set can reduce parameter uncertainty.

A prediction that substitutes one estimated parameter vector and one estimated random effect into a risk formula hides both sources. The earlier material on why plug-in random effects can mislead motivates integrating over the random-effects distribution. Here we extend that calculation to also account for uncertainty in the fitted parameters. The central point is that the two integrations are nested: the distribution of a patient’s random effects depends on the parameters.

Definition: Let \(t_0\) be a landmark time, \(u>0\) a prediction horizon, and \(H_i(t_0)\) the observed history for patient \(i\) through \(t_0\), including the fact that the patient is event-free at \(t_0\). Let \(\theta\) denote the joint-model parameters and \(b_i\) the patient’s shared random effects. The conditional risk is the probability of an event in \((t_0,t_0+u]\), given \(H_i(t_0)\). Random-effects uncertainty is uncertainty about \(b_i\) given the history and model parameters; parameter uncertainty is uncertainty about \(\theta\) given the data used to fit the model.

The Landmark Risk Before Integration

Write \(S_i(t\mid b_i,\theta)\) for the patient’s survival probability through time \(t\), conditional on the random effects and parameters. For fixed values of \(b_i\) and \(\theta\), the probability of an event by \(t_0+u\), conditional on survival through \(t_0\), is

$$ r_i(b_i,\theta;t_0,u) = 1-\frac{S_i(t_0+u\mid b_i,\theta)} {S_i(t_0\mid b_i,\theta)}. $$

The ratio is essential: the prediction starts from a patient known to be event-free at the landmark, rather than from the beginning of follow-up. Assume \(S_i(t_0\mid b_i,\theta)>0\). In a shared-random-effects joint model, the longitudinal and event processes are modeled as conditionally independent given \(b_i\) and \(\theta\), with the random effects linking the marker trajectory to the event hazard. The survival probabilities in this expression therefore depend on the event submodel and its association with the latent trajectory.

At a fixed \(\theta\), uncertainty about \(b_i\) is represented by its conditional distribution given the observed history. In a standard shared-effects formulation, its density is proportional to the random-effects density, the likelihood of the observed longitudinal history, and the probability of surviving to the landmark:

$$ p(b_i\mid H_i(t_0),\theta) \ \propto\ p\{Y_i(t_0)\mid b_i,\theta\} S_i(t_0\mid b_i,\theta) p(b_i\mid\theta). $$

Here \(Y_i(t_0)\) denotes the biomarker measurements observed through \(t_0\). The survival factor matters because observing a patient event-free at the landmark provides information about the random effects when those effects are associated with the event process. This expression assumes the longitudinal and event likelihood contributions have the stated conditional structure and that the observed measurement history and event-free status are handled consistently by the joint model.

Integrating the fixed-parameter risk over this conditional random-effects distribution gives

$$ q_i(\theta;t_0,u) = \int r_i(b_i,\theta;t_0,u) p(b_i\mid H_i(t_0),\theta)\,db_i. $$

This is the patient-specific risk at parameter value \(\theta\), averaged over plausible latent effects for that patient. To account for uncertainty in \(\theta\), average this quantity over the parameter distribution conditional on the data and history:

$$ R_i(t_0,u) = \int q_i(\theta;t_0,u) p(\theta\mid D,H_i(t_0))\,d\theta. $$

Here \(D\) denotes the data used to fit the model, and the conditioning notation means the parameter distribution is based on the available model-fitting information and the patient’s observed prediction history. In a Bayesian analysis, this is a posterior predictive probability. The inner distribution of \(b_i\) is conditional on \(\theta\); treating the two distributions as unrelated would generally be incorrect.

Why the Integrations Must Be Nested

The nested expression follows by conditioning on the unknown model parameters and then on the patient’s random effects. For an event indicator \(Z_i\) that equals one if the patient experiences the event in \((t_0,t_0+u]\), conditional expectation gives

$$ \begin{aligned} P(Z_i=1\mid D,H_i(t_0)) &= \int P(Z_i=1\mid b_i,\theta,H_i(t_0))\\ &\qquad \times p(b_i\mid\theta,H_i(t_0),D) p(\theta\mid D,H_i(t_0))\,db_i\,d\theta. \end{aligned} $$

Under the landmark survival model, the first conditional probability is \(r_i(b_i,\theta;t_0,u)\). Integrating first over \(b_i\) yields \(q_i(\theta;t_0,u)\), and integrating that result over \(\theta\) yields \(R_i(t_0,u)\). This is an application of iterated expectation; it does not require the random effects and parameters to be independent. Indeed, the conditional distribution \(p(b_i\mid\theta,H_i(t_0),D)\) can change with \(\theta\).

The derivation also clarifies what must be conditioned on. The history includes the measurements available by the landmark and survival to that time. The model must describe the relevant measurement and event processes well enough for these conditional distributions to be meaningful. If the data collection process or censoring is informative beyond what the model accounts for, the stated conditional distribution may not represent the patient’s actual information.

Worked Example: Integrating Over Random-Effects Uncertainty

Suppose that, at a fixed parameter value, the conditional cumulative hazard from \(t_0\) to \(t_0+u\) is \(0.10\exp(b_i)\). For illustration, assume the patient’s random-effect distribution given the history assigns probability \(0.75\) to \(b_i=0\) and probability \(0.25\) to \(b_i=\log(4)\). Calculate the integrated risk and compare it with substituting the mean random effect.

For a cumulative hazard \(h\), the conditional event risk is \(1-\exp(-h)\). At \(b_i=0\), the cumulative hazard is \(0.10\), so the risk is \(1-\exp(-0.10)=0.09516\). At \(b_i=\log(4)\), it is \(0.10(4)=0.40\), so the risk is \(1-\exp(-0.40)=0.3297\).

The random-effects-integrated risk is the probability-weighted average:

$$ 0.75(0.09516)+0.25(0.3297) = 0.07137+0.08242 = 0.1538. $$

The mean random effect is \(0.75(0)+0.25\log(4)=0.3466\). Substitution of that mean gives a cumulative hazard \(0.10\exp(0.3466)\approx0.1414\), and hence a risk \(1-\exp(-0.1414)\approx0.1319\). This differs from \(0.1538\): averaging a nonlinear risk function is not generally the same as evaluating it at an average random effect. The example assumes the stated two-point distribution is the appropriate conditional distribution given the patient’s history and landmark survival.

Adding Parameter Uncertainty

The outer integration averages the already random-effects-integrated risk across plausible parameter values. The next example uses a deliberately small two-point parameter distribution so the arithmetic is visible. Assume the patient’s conditional random-effect weights are \(0.60\) and \(0.40\), with hazard multipliers \(1\) and \(3\), respectively. Suppose the baseline cumulative hazard over the horizon is \(0.08\) or \(0.12\), each with weight \(0.50\). For simplicity, assume these weights and multipliers are fixed in this illustration.

Worked Example: Integrating Over Both Sources

At baseline cumulative hazard \(0.08\), the risks for the two random-effect values are \(1-\exp(-0.08)=0.07688\) and \(1-\exp(-0.24)=0.2134\). Thus the risk averaged over random effects is

$$ q(0.08) = 0.60(0.07688)+0.40(0.2134) = 0.04613+0.08535 = 0.1315. $$

At baseline cumulative hazard \(0.12\), the two risks are \(1-\exp(-0.12)=0.1131\) and \(1-\exp(-0.36)=0.3023\). Therefore

$$ q(0.12) = 0.60(0.1131)+0.40(0.3023) = 0.06785+0.1209 = 0.1888. $$

Averaging over the two parameter values gives

$$ R = 0.50(0.1315)+0.50(0.1888) = 0.1601. $$

For comparison, substituting the mean baseline cumulative hazard, \(0.10\), before calculating risk gives \(0.60\{1-\exp(-0.10)\}+0.40\{1-\exp(-0.30)\}=0.1608\). The close values here are a feature of these particular numbers, not an identity. In a real model, changing \(\theta\) may also change the random-effects distribution conditional on the patient’s history, so the inner calculation generally needs to be redone for each parameter value.

Separating Contributions to Variation

A nested calculation gives a point prediction, but the nesting also shows how variation in the predicted risk can be attributed to the two sources. In the following simple example, treat the conditional risk itself as the quantity whose variation is being described. Suppose two equally weighted parameter values, \(\theta_A\) and \(\theta_B\), each have two equally weighted random-effect values. The four conditional risks are:

Parameter valueRisk at first random-effect valueRisk at second random-effect valueMean risk at this parameter value
\(\theta_A\)0.100.300.20
\(\theta_B\)0.200.500.35

Worked Example: Decomposing Risk Variation

The overall mean risk is \(0.50(0.20)+0.50(0.35)=0.275\). At \(\theta_A\), the random-effects variance of the risk is \(\{(0.10-0.20)^2+(0.30-0.20)^2\}/2=0.0100\). At \(\theta_B\), it is \(\{(0.20-0.35)^2+(0.50-0.35)^2\}/2=0.0225\). Averaging these conditional variances gives the random-effects contribution \(0.50(0.0100)+0.50(0.0225)=0.01625\).

The parameter contribution is the variance of the two parameter-specific mean risks: \(\{(0.20-0.275)^2+(0.35-0.275)^2\}/2=0.005625\). The total variance of the four equally weighted conditional risks is therefore \(0.01625+0.005625=0.02188\). Direct calculation gives the same answer: the mean of the four squared risks is \((0.10^2+0.30^2+0.20^2+0.50^2)/4=0.09750\), and \(0.09750-0.275^2=0.02188\).

This variance describes variation in the conditional risk across the assumed parameter and random-effect distributions. It is not the variance of the future event indicator, which also includes whether the event occurs. Nor does the decomposition imply independence: it is the law of total variance, with the random-effects contribution calculated conditionally on parameters.

Frequentist Parameter Uncertainty and Practical Checks

In a Bayesian analysis, the outer distribution is a posterior distribution for \(\theta\), and its integration has a direct probability interpretation conditional on the observed data and model. In a frequentist analysis, \(\theta\) is fixed. Parameter uncertainty instead refers to the sampling variability of the estimator \(\hat{\theta}\). One practical approximation is to use a distribution of parameter values centered near \(\hat{\theta}\), often informed by its estimated sampling covariance, and recompute the random-effects-integrated risk for each value. A parametric bootstrap can serve the same broad purpose by refitting the model to simulated data sets and recalculating the prediction.

Averaging over such frequentist draws is an uncertainty-aware approximation, not a statement that the fixed parameter is literally random. The average may differ from the plug-in estimate because risk is nonlinear in the parameters; it is not automatically a better-calibrated estimator in every setting. The spread of the resulting risks can help describe sensitivity to parameter estimation, but should not be confused with the probability that a fixed parameter lies in a particular interval.

Useful checks include whether the conditional random-effects distribution respects the patient’s observed history and event-free status, whether parameter draws lead to valid model parameters, and whether the risk remains between zero and one. Large changes in risk across plausible parameter values indicate meaningful parameter sensitivity. A marked difference between a plug-in random-effect result and the inner average signals that latent-effect uncertainty matters. These checks diagnose instability; they do not establish that the joint-model assumptions are correct.

The main alternatives answer different questions. A fully plug-in calculation uses \(\hat{\theta}\) and a single estimated random effect, and is simple but suppresses both uncertainties. Integrating over random effects while fixing \(\theta=\hat{\theta}\) accounts for patient-level latent uncertainty but not estimation uncertainty in the model. Integrating over both sources gives the nested prediction described here. A population-average risk instead averages over the population distribution of random effects rather than conditioning on one patient’s history; it is not a substitute for this patient-specific calculation.

Common Pitfall and Limitations

  • Averaging effects instead of risks. Substituting a posterior mean random effect, or a mean parameter vector, into a nonlinear risk formula generally does not equal integrating the risk over the corresponding distribution.
  • Ignoring landmark survival. The random-effects distribution for a patient known to be event-free at \(t_0\) should reflect that information when the event process depends on the shared effects.
  • Assuming the two uncertainties are independent. The distribution of \(b_i\) conditional on the patient’s history depends on \(\theta\); the inner integral should be conditional on each outer parameter value.
  • Calling frequentist parameter draws a posterior. Sampling variability in \(\hat{\theta}\) is not a probability distribution for a fixed parameter unless a Bayesian model is being used.

The calculation cannot correct a misspecified longitudinal trajectory, event hazard, association structure, or observation process. It can also be computationally demanding when the random-effects dimension is high or the parameter distribution is strongly non-Gaussian. The appropriate integration method is a separate computational question; regardless of method, the target remains the nested expectation, not a plug-in value.

Key takeaway: First average the landmark event risk over the patient’s random effects conditional on the observed history and parameters. Then average that result over parameter uncertainty. The order makes the conditioning explicit and preserves the dependence between the two sources.

Check Your Understanding

Answer each question, then compare with the brief solutions below.

  1. Why does a landmark risk use \(S_i(t_0+u\mid b_i,\theta)/S_i(t_0\mid b_i,\theta)\), rather than simply \(S_i(t_0+u\mid b_i,\theta)\)?
  2. At one fixed parameter value, suppose two equally likely random-effect values imply risks 0.12 and 0.28. What is the integrated risk?
  3. Suppose two equally weighted parameter values give random-effects-integrated risks of 0.10 and 0.30. What is the outer-averaged risk?
  4. In a frequentist analysis, what does a distribution of parameter draws represent, and what does it not represent?

Solutions. (1) The ratio conditions survival to \(t_0\), so the risk concerns the interval after the landmark among patients event-free at that time. (2) The risk is \(0.50(0.12)+0.50(0.28)=0.20\). (3) The outer average is \(0.50(0.10)+0.50(0.30)=0.20\). (4) The draws approximate sampling uncertainty in the estimator or arise from a resampling procedure; they do not mean that a fixed frequentist parameter is literally random.

Clinical Biostats