Spline Trajectories for Binary Biomarkers

Why Link a Binary Trajectory to Event Risk?

A binary biomarker records whether a condition or threshold is present at each assessment. A logistic mixed model can describe how the probability of a positive result changes over time, but a simple linear time trend may be too restrictive. A spline lets the trajectory bend. In a joint model, subject-level random effects describe differences between individuals in these trajectories and also enter the time-to-event model. This creates a statistical link between the repeated binary measurements and event risk, while accounting for the fact that measurements from the same person are related.

The model has two scales to keep distinct. The longitudinal spline describes the biomarker’s log odds, which are converted to probabilities by the logistic function. The event submodel describes a hazard, conditional on the subject’s random effects. A spline coefficient is not itself a probability, and a hazard ratio is not a difference in probabilities.

Model Specification and Assumptions

For subject \(i\), let \(Y_{ij}\) be the binary biomarker measured at time \(t_{ij}\), where \(Y_{ij}=1\) denotes a positive result and \(Y_{ij}=0\) a negative result. Let \(b_i\) be a \(q\)-dimensional vector of subject-level random effects, with \(b_i\sim N_q(0,D)\), where \(D\) is a positive-definite covariance matrix. A spline-based logistic mixed model is

Here \(p_{ij}=P(Y_{ij}=1\mid b_i)\), \(x_i\) contains baseline covariates, and \(B(t)\) is a fixed spline basis specified by its basis functions and knots. The coefficient vector \(\gamma\) describes the population-level spline trajectory on the log-odds scale. The vector \(Z(t)\) selects or constructs the spline components that vary by subject, and \(b_i\) gives that subject’s deviations. The logistic link is \(\operatorname{logit}(p)=\log\{p/(1-p)\}\), so \(p=1/\{1+\exp(-\operatorname{logit}(p))\}\). The basis and random-effects structure must be specified before interpreting fitted coefficients.

A proportional-hazards event submodel that shares these random effects is

In this expression, \(h_0(t)\) is the baseline hazard, \(w_i\) contains baseline event covariates, and \(\alpha\) describes how the shared random effects are associated with event risk. Conditional on covariates and \(b_i\), the hazard ratio for two subjects with random effects \(b_i\) and \(b_k\) and identical \(w\) is \(\exp\{\alpha^\mathsf{T}(b_i-b_k)\}\). This is an association conditional on the model; it does not establish that changing a biomarker trajectory would cause an event-risk change.

The usual assumptions include a correctly specified spline basis and link, an appropriate distribution for \(b_i\), and conditional independence of the repeated binary outcomes given \(b_i\) and covariates. The event submodel must also describe the hazard adequately. The visit schedule is treated as fixed or noninformative under the model; if assessment times depend on unobserved health status, that visit process may need its own model. Censoring is assumed noninformative conditional on the modeled information, so the censoring mechanism does not add an unmodeled dependence on the event time. These assumptions are substantive, not consequences of using splines.

Subject Likelihood and the Shared-Effects Link

Let \(T_i\) be the observed event or censoring time, and let \(\delta_i=1\) indicate an observed event and \(\delta_i=0\) indicate right censoring. Conditional on \(b_i\), the event-time contribution is the hazard at the observed time when an event occurs, multiplied by the probability of surviving event-free through the observed time. Define \(J_i\) as the set of observed biomarker visits. Under the stated conditional-independence assumptions, the subject’s marginal likelihood contribution is

Here \(p_{ij}(b)\) and \(h_i(t\mid b)\) are the longitudinal probability and hazard evaluated at the same candidate random-effects value \(b\), and \(\phi_q(b;0,D)\) is the multivariate normal density. For a censored subject, \(\delta_i=0\), so the event contribution is survival through \(T_i\), with no hazard multiplier. The population likelihood is the product of these subject contributions, or equivalently its log is the sum of their logarithms.

This construction also explains why numerical integration is needed. The binary outcomes contribute Bernoulli probabilities, while the event term depends nonlinearly on \(b_i\) through the hazard and cumulative hazard. The integral is generally not available in closed form. Numerical integration approximates it; the accuracy of that approximation must be adequate for the number and dimension of the random effects.

Derivation, Interpretation, and Practical Checks

The likelihood follows by conditioning first. Given \(b_i\), the Bernoulli observations have product likelihood \(\prod_{j\in J_i}p_{ij}^{Y_{ij}}(1-p_{ij})^{1-Y_{ij}}\). The event-time density contribution is \(h_i(T_i\mid b_i)^{\delta_i}S_i(T_i\mid b_i)\), where

Conditional independence makes the two contributions multiply. Since \(b_i\) is unobserved, average their product over its density. This gives the integral in the likelihood above. If instead one multiplied the separately averaged marker and event contributions, one would use the product of two integrals rather than the integral of the product. Those are generally different because both conditional contributions vary with \(b_i\).

The longitudinal fixed effects and spline coefficients describe the population pattern in log odds, subject to the chosen covariates and basis. The random effects describe deviations around that pattern. The event association parameter \(\alpha\) quantifies how those deviations relate to the hazard. For instance, if a component of \(b_i\) represents a higher subject-specific log-odds level and its corresponding association is positive, higher values of that component are associated with a higher conditional hazard, holding the other modeled terms fixed. Because the random effects can be correlated, interpretation of one component should acknowledge the model’s joint distribution.

Useful checks include comparing observed positive-result proportions over time with fitted probabilities, examining whether the fitted spline behaves plausibly over the observed assessment range, and checking that event-time residual diagnostics do not reveal substantial disagreement with the hazard model. Inspect convergence warnings, estimated random-effects variances and correlations, and sensitivity to reasonable integration accuracy settings. A near-zero estimated variance may indicate that a random component is unsupported; a boundary estimate or unstable covariance estimate can also make the association parameters difficult to estimate.

The spline is flexible but does not guarantee a good trajectory. Sparse measurements in parts of follow-up, very few positive or negative results, or near-separation can make spline coefficients unstable. Predictions beyond the observed time range deserve particular caution. If numerical integration changes materially when its accuracy is increased, the fitted likelihood and resulting inferences may not be reliable.

Alternatives and Limitations

A linear-time logistic mixed model is simpler and may be adequate when the log-odds trajectory is close to linear, but it cannot represent meaningful curvature without additional terms. A spline-based logistic model without an event submodel can describe the repeated binary outcomes but does not use the shared-effects structure to connect them with event risk. Conversely, fitting separate longitudinal and survival models does not account for dependence through latent subject characteristics. These alternatives answer different questions and rely on different assumptions.

Sharing random effects is also a modeling choice, not a universal description of dependence. The particular components placed in \(b_i\), and the form \(\alpha^\mathsf{T}b_i\), determine how the processes are linked. This model does not automatically make the hazard depend on the current fitted biomarker probability; that would be a different association structure. Nor does shared dependence alone correct bias from informative visit timing, misspecified censoring assumptions, or an inadequate baseline hazard. Model interpretation should therefore stay tied to the specified conditional structure.