Profile Likelihood for an Association Parameter

Why Profile the Association Parameter?

In a joint model, the association parameter describes how the underlying longitudinal trajectory relates to the event hazard. Its estimate is rarely independent of the other model parameters: changing the association value may be offset, in part, by changes in the longitudinal mean, random-effects distribution, or survival parameters. A profile likelihood makes this dependence visible by asking: for each fixed association value, how well can the rest of the model fit?

The focus here is the construction of that profile, not a test or interval based on it. The association parameter is denoted by \(\alpha\), and all other unknown parameters are grouped into \(\eta\). Depending on the model, \(\eta\) includes longitudinal fixed effects and residual variance, random-effects covariance parameters, survival regression coefficients, and baseline-hazard parameters.

A profile is therefore not a likelihood evaluated at one fixed set of nuisance estimates. It is a sequence of constrained maximizations. The overall maximum-likelihood estimate of \(\alpha\) maximizes \(\ell_p(\alpha)\), provided the global maximum is attained within the permitted parameter space.

The Joint Likelihood Being Profiled

Let \(y_i\) denote subject \(i\)'s longitudinal measurements, \(t_i\) the observed event or censoring time, and \(\delta_i\) the event indicator, equal to 1 for an observed event and 0 for right censoring. Let \(b_i\) be the subject-specific random effects. A common shared-random-effects specification has longitudinal density \(f(y_i\mid b_i;\eta)\), random-effects density \(f(b_i;\eta)\), and conditional hazard

where \(w_i\) contains survival covariates, \(\gamma\) is their coefficient vector, and \(m_i(t;b_i)\) is the longitudinal model's latent mean at time \(t\). In this parameterization, \(\alpha\) is the log-hazard change per one-unit increase in the latent mean, conditional on the other terms in the hazard model. Its numerical interpretation depends on the units used for that mean.

Under conditional independence of the longitudinal measurements and event time given \(b_i\) and covariates, the subject's likelihood contribution is

Here \(S_i(t\mid b_i)=\exp\{-\int_0^t h_i(u\mid b_i)\,du\}\) is the conditional event-free survival function. This likelihood relies on a correctly specified longitudinal model, random-effects distribution, hazard and survival models, and an appropriate censoring assumption. In particular, censoring must be noninformative for the event process under the conditioning used by the model. The integration over \(b_i\) is essential: the profile is formed from the marginal likelihood, not by treating estimated random effects as observed data.

Constructing the Profile

Choose a value of \(\alpha\), hold it fixed, and maximize the integrated log-likelihood over every component of \(\eta\). Repeat this operation across a grid of association values. The nuisance estimates will generally differ from one grid point to another. A reliable implementation uses the solution at a nearby grid point as a starting value, but should also use additional starts when local maxima are plausible.

Why Reoptimizing Nuisance Parameters Matters

The reason for reoptimization can be expressed through derivatives. Suppose that at each interior value \(a\), the nuisance maximum \(\hat\eta(a)\) is unique and smooth, and that the nuisance score is zero there. Then \(\ell_p(a)=\ell(a,\hat\eta(a))\). Differentiating by the chain rule gives

because the nuisance score in braces is zero at the constrained maximum. Thus an interior maximum of the profile satisfies the association score equation evaluated at the nuisance estimates for that value of \(\alpha\). If \(\alpha\) is on a boundary, a one-sided condition replaces the interior zero-score condition.

The profile's local curvature also accounts for nuisance-parameter adjustment. If the full log-likelihood Hessian is partitioned into blocks for \(\alpha\) and \(\eta\), then, at a regular interior maximum,

This follows by differentiating the nuisance score equation, which gives \(d\hat\eta/d\alpha=-\ell_{\eta\eta}^{-1}\ell_{\eta\alpha}\), and substituting into the derivative of the profile score. The second term represents the adjustment from nuisance parameters being allowed to move. Holding \(\eta\) fixed removes that adjustment and produces a different, conditional slice of the likelihood.

A Profile Is Not a Fixed-Nuisance Slice

At the joint maximum, the nuisance estimates are optimized for one particular association value. Evaluating the likelihood at those same nuisance estimates for all other values of \(\alpha\) may be useful as a diagnostic, but it is not the profile likelihood. The distinction can materially affect the apparent peak and curvature.

Refining the Peak and Checking Computation

A coarse grid helps reveal the shape of the profile but does not, by itself, identify the maximum precisely. A local quadratic approximation can suggest where to refine the search. Consider another hypothetical analysis in which nuisance parameters are reoptimized at \(\alpha=0.50, 0.60,\) and \(0.70\) per unit, yielding \(-287.41,-287.28,\) and \(-287.39\). Set \(x=\alpha-0.60\), and fit \(q(x)=c+bx+ax^2\) through these three values. Since \(q(0)=-287.28\), the other two points imply

The fitted parabola has its vertex at \(x=-b/(2a)=0.00417\), suggesting a peak near \(\alpha=0.6042\). Its predicted log-likelihood gain over \(q(0)\) is \(b^2/(-4a)=0.000208\). This is only a local interpolation: the appropriate check is to reoptimize the nuisance parameters at the refined value and confirm that the actual profile value and nearby shape agree. A small apparent improvement should not be reported as meaningful if it is below the numerical accuracy of the likelihood calculation.

For a joint model, numerical integration over random effects is part of every likelihood evaluation. If integration error changes as nuisance parameters or \(\alpha\) change, it can create a jagged profile or a false maximum. Check that integration settings are sufficiently accurate across the profile, not just at the final estimate. Also inspect convergence codes, gradients, and nuisance estimates at grid points near the peak. Abrupt jumps in nuisance estimates or irregular profile values can indicate local optimization failures, weak identification, or inadequate integration accuracy.

Alternatives, Limitations, and Common Mistakes

Profiling does not remedy model misspecification, weak identification, or a poorly computed integrated likelihood. It can be costly because nuisance optimization must be repeated many times. A profile can also be flat, asymmetric, or have multiple local peaks; in such cases, a single local quadratic summary may conceal important features. Constraints on nuisance parameters must be enforced at every fixed value of \(\alpha\), and a boundary in \(\alpha\)'s permitted range should be reported rather than silently treated as an interior estimate.

Common mistake: “I estimated all parameters once, so I can vary \(\alpha\) while holding the other estimates fixed; that curve is the profile.” The correction is that a profile is defined by a new maximization over all nuisance parameters for every fixed \(\alpha\). A fixed-nuisance curve is just a slice. Confusing the two can give the wrong profile peak and an incorrect picture of how well nearby association values fit.