Marginal Likelihood With Non-Gaussian Markers

Why the Random Effects Must Be Integrated After Combining the Data

In a joint model, repeated marker values and an event time can both carry information about a subject’s latent trajectory. Conditional on that subject’s random effects, the marker observations and event outcome may be modeled as independent. But once the random effects are unobserved, the two sources of data are generally dependent: an observed marker can change how plausible different random-effect values are, which in turn changes the event contribution.

The observed-data likelihood therefore has a specific order. For each subject, first multiply the marker contributions and the event-time contribution conditional on the same random effects. Then integrate that product over the random-effects distribution. This tutorial derives that likelihood and works through discrete-marker examples, including a calculation in which integrating separate marker and event likelihoods gives a different answer.

Definition: Let \(Y_{ij}\) denote subject \(i\)’s marker value at visit \(j\), observed at time \(t_{ij}\). Let \(T_i\) be the event or censoring time and \(\delta_i\) the event indicator, equal to 1 when the event is observed and 0 when the subject is right-censored. Let \(b_i\) be a subject-specific random-effects vector with density \(f_b(b_i;\theta)\). The observed data for subject \(i\) are the available marker values and times, \(T_i\), and \(\delta_i\). The marginal, or observed-data, likelihood contribution is obtained by integrating the joint conditional contribution over \(b_i\).

Conditional Contributions and the Marginal Likelihood

Write \(f_Y(y_{ij}\mid b_i;\theta)\) for the conditional probability mass or density of a marker observation under its non-Gaussian longitudinal model. For example, it can be a Bernoulli probability, an ordinal-category probability, or a Poisson probability. The random effects can enter the model’s linear predictor, and the response distribution and link determine how that predictor maps to the marker scale.

For the event process, let \(h_i(t\mid b_i;\theta)\) be the hazard conditional on the random effects and \(S_i(t\mid b_i;\theta)\) the corresponding survival function. For an observed event, the event-time contribution is the conditional density \(h_i(T_i\mid b_i;\theta)S_i(T_i\mid b_i;\theta)\). For right-censoring, it is \(S_i(T_i\mid b_i;\theta)\): the subject is known to have remained event-free through the censoring time.

L_i(\theta) = \int \left\{\prod_{j\in\mathcal{O}_i} f_Y(Y_{ij}\mid b_i;\theta)\right\} h_i(T_i\mid b_i;\theta)^{\delta_i} S_i(T_i\mid b_i;\theta) f_b(b_i;\theta)\,db_i .

Here \(\mathcal{O}_i\) is the set of visits with observed marker values. If the event is observed, \(\delta_i=1\), so the factor involving the event process is \(h_i(T_i\mid b_i;\theta)S_i(T_i\mid b_i;\theta)\). If it is censored, \(\delta_i=0\), and that factor reduces to \(S_i(T_i\mid b_i;\theta)\). The likelihood for independent subjects is the product of these subject-level integrals, \(L(\theta)=\prod_i L_i(\theta)\).

Assumptions for this factorization: Conditional on \(b_i\) and included covariates, the marker observations and event process are independent, and the marker observations are conditionally independent across visits. The random effects follow the specified distribution \(f_b\). Censoring is non-informative for the event time under the model, conditional on the variables included, so its distribution contributes no parameter-dependent factor to this likelihood. The visit schedule is treated as fixed or conditionally non-informative, and missing marker values are handled as ignorable when their missingness mechanism can be omitted under the model’s assumptions.

These assumptions are about conditional distributions, not unconditional independence. Integrating out \(b_i\) can induce dependence between a subject’s marker values and event outcome even when the conditional independence statement is appropriate. This shared latent structure is precisely why the marker and event factors must be multiplied before integration.

Justification of the Integral

The likelihood follows from the law of total probability. Given \(b_i\), the joint density or mass of the observed marker values and event data factors into the marker term and event term under the stated conditional independence assumption. Since \(b_i\) is not observed, average that conditional contribution over its density. For a right-censored observation, the probability contribution for the event data is the probability of surviving beyond \(T_i\); for an observed event, the density at \(T_i\) is hazard times survival.

\begin{aligned} f(\mathbf{Y}_i,T_i,\delta_i;\theta) &=\int f(\mathbf{Y}_i,T_i,\delta_i\mid b_i;\theta) f_b(b_i;\theta)\,db_i\\ &=\int \left\{\prod_{j\in\mathcal{O}_i}f_Y(Y_{ij}\mid b_i;\theta)\right\} \left\{h_i(T_i\mid b_i;\theta)^{\delta_i} S_i(T_i\mid b_i;\theta)\right\} f_b(b_i;\theta)\,db_i . \end{aligned}

The second equality uses the conditional factorization and the event contribution for the observed censoring status. If the random-effects distribution is discrete, the integral is replaced by a weighted sum over its support. In either case, the integrand is the full conditional contribution for that subject, including both data sources.

Worked Example: A Poisson Marker and a Censored Event

Consider one hypothetical subject with a count marker \(Y=1\) and right-censoring at \(T=2\), so \(\delta=0\). For illustration, suppose \(b\) takes values \(-0.5\) and \(0.5\), each with probability \(0.5\). Conditional on \(b\), let \(Y\) be Poisson with mean \(\mu(b)=\exp(-0.2+b)\), and let the event hazard be constant in time at \(\lambda(b)=0.1\exp(0.3b)\). These specifications imply conditional independence of the count and event time given \(b\), with non-informative censoring assumed.

For a Poisson observation \(Y=1\), the marker contribution is \(\mu(b)\exp\{-\mu(b)\}\). Because the subject is censored, the event contribution is survival through time 2, \(\exp\{-2\lambda(b)\}\). Their product is calculated separately at each random-effect value before averaging.

At \(b=-0.5\), \(\mu=\exp(-0.7)=0.4965853\), so the marker probability is \(0.4965853\exp(-0.4965853)\approx0.3022\). The hazard is \(\lambda=0.1\exp(-0.15)=0.0860708\), and the survival contribution is \(\exp\{-2(0.0860708)\}\approx0.8419\). Thus the conditional product is approximately \(0.3022(0.8419)=0.2544\); using the unrounded model values gives approximately \(0.2544308\).

At \(b=0.5\), \(\mu=\exp(0.3)=1.3498588\), so the marker probability is \(1.3498588\exp(-1.3498588)\approx0.3500\). The hazard is \(\lambda=0.1\exp(0.15)=0.1161834\), and the survival contribution is \(\exp\{-2(0.1161834)\}\approx0.7927\). Their product is approximately \(0.2774\).

Weighting the two conditional products by their random-effect probabilities gives the marginal likelihood contribution:

L_i = 0.5(0.2544308\ldots)+0.5(0.2774\ldots) \approx 0.2659.

This is the probability contribution for the observed count and censoring information under this hypothetical model. It is a likelihood value, not a posterior probability that the model or any parameter value is true.

Repeated Discrete Markers and an Observed Event

With multiple non-Gaussian marker observations, each contributes a conditional probability or density factor. Those factors are multiplied together for a given \(b_i\); they are not first averaged separately over the random-effects distribution. The event factor joins that same conditional product.

Worked Example: Two Binary Markers and an Observed Event

Suppose, hypothetically, \(b\) takes values \(-1,0,1\) with probabilities \(0.25,0.50,0.25\). A subject has two binary marker observations, both equal to 1. At each visit, let the probability of a 1 be \(p(b)=\operatorname{expit}(-1+b)\), where \(\operatorname{expit}(x)=1/(1+\exp(-x))\). Suppose an event is observed at time \(T=1\), with conditional exponential hazard \(\lambda(b)=0.2\exp(b)\). Assume the marker visits are conditionally independent given \(b\), and the event is conditionally independent of the markers given \(b\).

The conditional marker product is \(p(b)^2\), and the event density contribution is \(\lambda(b)\exp\{-\lambda(b)\}\). At \(b=-1\), \(p=\operatorname{expit}(-2)\approx0.1192\), so \(p^2\approx0.01421\). Here \(\lambda=0.2\exp(-1)\approx0.07358\), and the event contribution is approximately \(0.07358\exp(-0.07358)=0.06836\). Their product is about \(0.000971\).

At \(b=0\), \(p=\operatorname{expit}(-1)\approx0.2689\), so \(p^2\approx0.07233\). The event contribution is \(0.2\exp(-0.2)\approx0.1637\), giving a conditional product of about \(0.01184\).

At \(b=1\), \(p=\operatorname{expit}(0)=0.5\), so \(p^2=0.25\). The event contribution is approximately \(0.2\exp(1)\exp\{-0.2\exp(1)\}=0.3156\), giving a product of about \(0.07890\). Integrating over the three-point random-effects distribution gives

L_i \approx 0.25(0.000971)+0.50(0.01184)+0.25(0.07890) \approx 0.0259.

The small contribution at \(b=-1\) reflects the low probability of two marker successes at that random-effect value; the larger contribution at \(b=1\) reflects both the higher marker probability and the event density. The likelihood averages these joint conditional contributions, not the marker and event pieces separately.

Why Separate Marginalization Changes the Model

Worked Example: A Censored Event Shows the Factorization Error

Suppose a hypothetical subject has one binary marker \(Y=1\) and is censored at time 1. Let \(b=-1\) or \(b=1\), each with probability \(0.5\). Set \(P(Y=1\mid b)=\operatorname{expit}(b)\), and let the event hazard be constant at \(\lambda(b)=0.2\exp(b)\). Given \(b\), the marker and event are independent.

At \(b=-1\), the marker probability is \(0.26894\), and survival to time 1 is \(\exp\{-0.2\exp(-1)\}\approx0.9291\). The joint conditional contribution is approximately \(0.26894(0.9291)=0.2499\). At \(b=1\), the marker probability is \(0.73106\), and survival is \(\exp\{-0.2\exp(1)\}\approx0.5806\), giving a product of about \(0.4244\). Therefore the correct marginal joint contribution is

0.5(0.2499)+0.5(0.4244)\approx0.3371.

If the marker and event contributions are marginalized separately and then multiplied, the result is instead \((0.5(0.26894)+0.5(0.73106))\{0.5(0.9291)+0.5(0.5806)\}\approx0.5(0.7548)=0.3774\). This is not the joint likelihood: it treats the two data sources as though they averaged over unrelated draws of \(b\). In the model, the same subject-specific \(b\) links both contributions.

Assumptions, Diagnostics, and Limitations

The likelihood is only as appropriate as its conditional models and sampling assumptions. Check that the chosen response distribution and link are plausible for the marker scale, that the random-effects structure represents the intended subject-level dependence, and that the event hazard and its association with the random effects are specified on the intended scale. For example, a Poisson marker model assumes its conditional mean and variance are equal; a binary model assigns probabilities through its specified link. Violations can distort the marker factors and therefore also affect the inferred random-effects distribution and event association.

Useful diagnostics include comparing observed marker frequencies, category proportions, or count distributions with model-based predictions, and checking event calibration across relevant follow-up times. Examine whether the fitted model systematically misses marker patterns or event outcomes in groups with different observed histories. These checks do not prove the conditional-independence or censoring assumptions, which require subject-matter justification as well as data review.

The integral is generally not available in closed form for non-Gaussian markers, even when the random effects have a normal distribution. Its dimension is the number of random effects. The likelihood definition remains the same, but numerical evaluation may be needed; the specific integration methods are a separate topic. Computational approximations must be accurate enough that likelihood comparisons and parameter estimates are not artifacts of integration error.

A conditional likelihood that treats \(b_i\) as known is not the observed-data likelihood when those effects are unobserved. Plugging in estimated subject-specific effects also does not generally replace integrating over their distribution. In contrast, a model with no shared random effects or with conditionally independent sources that truly have separate latent structures may permit a different factorization. The factorization must follow the model’s dependence assumptions, rather than convenience.

Common pitfall: Multiplying the marker likelihood averaged over random effects by the event likelihood averaged over random effects. This incorrectly allows each source of data to use its own latent effect. For a shared-random-effects joint model, multiply the marker and event contributions conditional on the same \(b_i\), then integrate that product.

Check Your Understanding

Answer each question, then compare your reasoning with the solutions.

For a subject with one observed event at \(T_i\), what is the conditional event contribution? What is it if the subject is censored at \(T_i\)?
A subject has two conditionally independent Bernoulli marker observations given \(b_i\), both equal to 1. If their conditional success probability is \(p(b_i)\), what marker factor enters the integrand?
Suppose the random effect has two support points \(b_1,b_2\) with probabilities \(0.3,0.7\). Write the marginal contribution when the conditional marker-and-event product is \(g(b)\).
Why is the product of separately marginalized marker and event contributions generally not the joint observed-data likelihood?

Solutions. (1) For an observed event, the contribution is \(h_i(T_i\mid b_i)S_i(T_i\mid b_i)\); for right-censoring it is \(S_i(T_i\mid b_i)\). (2) The marker factor is \(p(b_i)^2\), which is then multiplied by the event factor before integration. (3) The marginal contribution is \(0.3g(b_1)+0.7g(b_2)\). (4) Separate marginalization averages each source over the random effects independently, losing the fact that both data sources share the same subject-specific effect; the joint contribution averages their conditional product instead.