From a Study Description to a Joint Model
A joint model is not fully specified by saying that a biomarker is “related” to an event time. The specification must identify what is measured, the scale of the latent trajectory, the response distribution conditional on that trajectory, the random effects shared across outcomes, and the way the trajectory enters the hazard. These choices determine both the scientific interpretation and the likelihood being fitted.
Consider participants followed over time, with repeated measurements \(Y_{ij}\) at times \(t_{ij}\), and an event or right-censoring time \(T_i\). Let \(\delta_i=1\) indicate that the event was observed and \(\delta_i=0\) indicate right censoring. Let \(X_i\) denote baseline covariates. The model below uses a flexible spline trajectory and an association between the latent biomarker and the event hazard. It specifies the conditional model; it does not by itself imply a causal effect of the biomarker.
Specify the Trajectory, Response, and Hazard
One convenient trajectory uses a fixed cubic regression-spline basis and participant-specific random intercept and slope:
For example, a truncated-power basis with fixed internal knots \(\kappa_1,\kappa_2\) is \(B(t)=(1,t,t^2,t^3,(t-\kappa_1)_+^3,(t-\kappa_2)_+^3)^\mathsf{T}\), where \(x_+=\max(x,0)\). This gives a piecewise cubic function with continuity through its second derivative at the knots. The knots, time origin, and time units must be stated. The vector \(b_i\) allows individual trajectories to depart from the population trajectory in level and linear trend; it does not make every spline coefficient participant-specific.
The linear predictor \(\eta_i(t)\) is linked to a conditional response mean \(\mu_i(t)\) through a chosen link \(g\): \(g\{\mu_i(t)\}=\eta_i(t)\). The likelihood must then suit the observed response. For counts, a negative-binomial model may use a log link; for a binary response, a Bernoulli model may use a logit link. The hazard specification must say which trajectory scale is associated with event risk. In the examples below, it uses the current link-scale value \(\eta_i(t)\):
Here \(h_0(t)\) is the baseline hazard, \(\gamma\) is a vector of covariate coefficients, and \(\alpha\) is the association parameter. With this convention, a one-unit increase in \(\eta_i(t)\) at time \(t\), holding the other terms fixed, multiplies the hazard by \(\exp(\alpha)\). Because \(\eta_i(t)\) is on the link scale, this is not automatically a one-unit increase in the biomarker’s response-scale mean.
Assumptions and Likelihood Structure
A standard shared-random-effects specification assumes that, conditional on \(b_i\) and the modeled covariates, repeated responses are independent with their stated response distributions, and the event process depends on the longitudinal history through the specified latent trajectory and hazard. The random effects are commonly modeled as \(b_i\sim N(0,D)\), with a positive-definite covariance matrix \(D\). This permits correlation between the random intercept and slope. Conditional independence is an assumption about the model after accounting for the latent effects; it does not say the observed measurements or event times are marginally independent.
The event contribution is conditional on \(b_i\). For a hazard \(h_i(t\mid b_i,X_i)\), the conditional survival function is \(S_i(T_i\mid b_i,X_i)=\exp\{-\int_0^{T_i}h_i(u\mid b_i,X_i)\,du\}\). Under conditionally independent censoring, the event contribution is \(h_i(T_i\mid b_i,X_i)^{\delta_i}S_i(T_i\mid b_i,X_i)\). Combining it with the repeated-response likelihood and integrating over the random-effects distribution gives the participant’s marginal likelihood contribution:
This expression follows by multiplying the conditional response densities, the event density or censoring survival contribution, and the random-effects density, then integrating out the unobserved \(b_i\). It makes clear why fitting separate biomarker and survival models generally does not reproduce a shared-random-effects joint model: those separate fits do not integrate the two sources of information through the same latent effects.
The specification also requires appropriate handling of entry and censoring. The displayed likelihood assumes follow-up begins at time zero and that censoring is independent of the event process conditional on modeled information. Delayed entry requires conditioning on survival to entry; informative censoring may require a censoring model. Neither feature is repaired merely by adding a spline to the biomarker submodel.
Worked Example: Repeated Overdispersed Counts
Suppose a hypothetical study records a count biomarker at irregular visits, in years since enrollment, and follows participants to a clinical event. Counts are nonnegative and show more variation than a Poisson model would allow. Specify a model with knots at 3 and 6 years, random intercept and slope, and a link-scale current-value hazard association.
Use the basis \(B(t)=(1,t,t^2,t^3,(t-3)_+^3,(t-6)_+^3)^\mathsf{T}\). Let the conditional count mean be \(\mu_i(t)=\exp\{\eta_i(t)\}\), with negative-binomial variance \(\mu_i(t)+\mu_i(t)^2/\phi\), where \(\phi>0\) is a dispersion parameter. Let \(b_i\sim N(0,D)\). These choices suit overdispersed counts; the mean must be positive, as ensured by the log link.
For a numerical illustration, take \(\beta=(0.2,0.03,0,0,-0.002,0.001)^\mathsf{T}\), \(b_i=(0.10,-0.01)^\mathsf{T}\), and evaluate at \(t=4\) years. The basis is \(B(4)=(1,4,16,64,1,0)^\mathsf{T}\), so
Thus the conditional expected count at that time is about 1.459 counts for the stated random effects. It is not the observed count and does not assert that a count of 1.459 can be observed.
For the event model, suppose \(h_0(t)=\lambda\rho t^{\rho-1}\), with \(\lambda=0.02\), \(\rho=1.2\), a binary baseline covariate \(X_i=1\), \(\gamma=-0.3\), and \(\alpha=0.4\). At 4 years, \(h_0(4)=0.02(1.2)4^{0.2}\approx0.03167\) per year. The hazard is
This is the model’s conditional hazard at year 4 for the given covariate and random effects. Holding the other terms fixed, doubling the conditional mean count changes \(\eta\) by \(\log 2\), so the hazard ratio is \(\exp(0.4\log 2)\approx1.320\). That is an association under the specified model, not a causal effect.
Worked Example: Binary Repeated Biomarker
Suppose a hypothetical biomarker is whether a test is positive at each visit. Use the same spline basis and random intercept and slope, but set \(Y_{ij}\mid b_i\sim\operatorname{Bernoulli}(p_i(t_{ij}))\) and \(\operatorname{logit}\{p_i(t)\}=\eta_i(t)\). The response likelihood is then a product of Bernoulli probabilities, not a Gaussian density. This is appropriate if each measurement is a single binary outcome and the conditional success probability is the target mean.
At \(t=2\), the two knot terms vanish. For illustration, let the nonzero fixed coefficients be \(\beta_0=-1\) and \(\beta_1=0.2\), and let \(b_i=(0.30,-0.05)^\mathsf{T}\). Then
The fitted conditional probability is about 0.4013. If the hazard uses \(\alpha\eta_i(t)\) with \(\alpha=0.5\), an increase of one unit in log odds multiplies the hazard by \(\exp(0.5)\approx1.649\). This is a link-scale contrast; it is not the hazard ratio for a one-unit increase in probability, which is bounded between zero and one.
Worked Example: Continuous Biomarker Bounded Between Zero and One
Suppose a hypothetical biomarker is continuous and strictly between zero and one at every recorded visit. A beta regression specification can use \(Y_{ij}\mid b_i\sim\operatorname{Beta}(\mu_i(t_{ij})\phi,\{1-\mu_i(t_{ij})\}\phi)\), with \(\operatorname{logit}\{\mu_i(t)\}=\eta_i(t)\) and precision \(\phi>0\). The resulting conditional mean is between zero and one, and the conditional variance is \(\mu_i(t)\{1-\mu_i(t)\}/(1+\phi)\).
At a visit where \(\eta_i(t)=\log(3)\), \(\mu_i(t)=3/(1+3)=0.75\). If \(\phi=20\), the conditional variance is \(0.75(0.25)/21\approx0.008929\). This likelihood is unsuitable for exact zeros or ones without an extension that models those boundary values. For a current link-scale association with \(\alpha=0.4\), doubling the odds of the conditional mean adds \(\log 2\) to \(\eta_i(t)\), giving hazard ratio \(\exp(0.4\log 2)\approx1.320\), conditional on other terms.
Checks, Alternatives, and Limitations
Before fitting, check that the response likelihood matches the data support and that its link yields valid means. Check that spline knots lie within the observed time range and are scientifically plausible; a flexible basis cannot create information where measurements are sparse. The chosen random-effects structure should be estimable from the number and timing of repeated observations. In particular, random slopes are weakly informed when almost all participants have only one measurement or visits concentrated at nearly identical times.
Useful diagnostics include comparing observed and fitted response distributions over time, checking residual patterns appropriate to the response family, assessing whether the spline leaves systematic temporal structure, and checking convergence, integration accuracy, and sensitivity to knot placement or random-effects structure. For the event submodel, examine calibration and residual diagnostics appropriate to the chosen baseline hazard, and check whether the fitted association is unduly sensitive to influential participants. A flexible trajectory can still be misspecified, and a fitted association does not prove that the latent trajectory captures all dependence between biomarker and event.
| Specification choice | When it may be appropriate | Key consequence |
|---|---|---|
| Shared random effects with a link-scale current value | The hazard is hypothesized to depend on the current latent level on the response-link scale | Association is interpreted per unit of the link-scale predictor |
| Association with the response-scale mean | A unit change in expected biomarker response is the intended contrast | The hazard association has a different scale and may be nonlinear in the link predictor |
| Separate longitudinal and event models | A joint latent-dependence assumption is not intended or needed | Does not account for shared latent effects in a joint likelihood |
These are alternatives, not interchangeable labels for the same model. For example, a hazard term \(\alpha\mu_i(t)\) and a term \(\alpha\eta_i(t)\) impose different relationships between biomarker level and hazard. The scale must be chosen and stated before interpreting \(\alpha\). A model may also be limited by omitted time-varying predictors, misspecified random-effects distributions, informative visit timing, or a baseline hazard that cannot represent the event-time pattern.
A Common Specification Mistake
A frequent error is to write \(h_i(t)=h_0(t)\exp\{\alpha\mu_i(t)\}\), then interpret \(\exp(\alpha)\) as the hazard ratio for a one-unit increase in \(\eta_i(t)\). That interpretation is wrong: the model includes the response-scale mean, not the link-scale predictor. With a log link, \(\mu=\exp(\eta)\), so a one-unit increase in \(\eta\) changes the hazard multiplier from \(\exp\{\alpha\exp(\eta)\}\) to \(\exp\{\alpha\exp(\eta+1)\}\), which depends on the starting value of \(\eta\). State the hazard association exactly as intended, then derive its contrast on that scale.
A second practical warning is to avoid treating every observed response as if it were the latent trajectory. In a non-Gaussian model, \(\eta_i(t)\) is generally a link-scale predictor, while \(\mu_i(t)\) is a conditional response mean and \(Y_{ij}\) is a realized observation. Keeping these three quantities distinct prevents scale errors in both the likelihood and the hazard interpretation.
Check Your Understanding
Answer each question, then compare with the solutions.
- A count model has a log link and hazard term \(\alpha\eta_i(t)\), with \(\alpha=0.4\). What is the hazard ratio for doubling the conditional mean count?
- For a binary response with \(\operatorname{logit}(p)=\eta\), what probability corresponds to \(\eta=0\)?
- Why does a beta response model with support strictly between zero and one need modification when exact zero values occur?
- What does conditional independence of repeated responses mean in a shared-random-effects model?
Solutions. (1) Doubling the mean under a log link changes \(\eta\) by \(\log 2\), so the hazard ratio is \(\exp(0.4\log 2)=\exp(0.2773)\approx1.320\). (2) \(p=\exp(0)/\{1+\exp(0)\}=0.5\). (3) A beta distribution assigns probability only to values in the open interval \((0,1)\), so it cannot give positive probability to an exact zero; a model for boundary values is needed. (4) Conditional on the random effects and modeled covariates, the repeated responses have the specified independent response distributions; marginally, they can be dependent because they share the random effects.