Cumulative Link Likelihoods For Ordinal Responses

From Ordered Categories to a Likelihood

An ordinal response records which of several ordered categories applies, such as “none,” “mild,” “moderate,” or “severe.” The category labels have an order, but the distances between them need not be equal. A cumulative link model represents the probability of being at or below each category boundary. It then obtains the probability of an exact response category by taking a difference between neighboring cumulative probabilities.

For a GLMM, these probabilities are usually specified conditional on the random effects. The likelihood contribution for an observation is the conditional probability of its observed category. To obtain a marginal likelihood for a cluster, the product of its conditional contributions is integrated over the random-effects distribution.

Definition: Let $Y_{ij}$ be an ordinal response for observation $i$ in cluster $j$, with ordered categories $1,\ldots,K$. Define the cumulative probability $C_{ijk}=P(Y_{ij}\leq k\mid b_j)$, for $k=1,\ldots,K-1$, conditional on random effects $b_j$. A cumulative link model specifies $$ g(C_{ijk})=\theta_k-\eta_{ij}, $$ where $g$ is a link function, $\theta_1<\cdots<\theta_{K-1}$ are thresholds, and $\eta_{ij}$ is the linear predictor, which may include fixed effects and random effects. For the cumulative logit link, $g(p)=\log\{p/(1-p)\}$ and its inverse is the logistic CDF $F(t)=1/(1+e^{-t})$. For the cumulative probit link, $g(p)=\Phi^{-1}(p)$ and its inverse is the standard normal CDF $\Phi(t)$.

Write $F$ for the inverse link CDF, so $C_{ijk}=F(\theta_k-\eta_{ij})$. The first and last category probabilities require the endpoints $C_{ij0}=0$ and $C_{ijK}=1$. Thus, for every $r=1,\ldots,K$,

P(Y_{ij}=r\mid b_j) =C_{ijr}-C_{ij,r-1} =F(\theta_r-\eta_{ij})-F(\theta_{r-1}-\eta_{ij}),

where the endpoint values are used when $r=1$ or $r=K$. This difference, rather than a cumulative probability by itself, is the likelihood contribution when the observed response is category $r$.

Assumptions and Probability Coherence

The cumulative-link likelihood is well-defined when the response categories are ordered and the model assigns them the intended order. The link CDF must be increasing, the thresholds must satisfy $\theta_1<\cdots<\theta_{K-1}$, and the linear predictor must be defined on the scale used in the cumulative equation. These conditions ensure that the cumulative probabilities are ordered. In particular, for each observation,

0\leq C_{ij1}\leq C_{ij2}\leq\cdots\leq C_{ij,K-1}\leq1.

The category probabilities are then nonnegative and sum to one: summing $C_{ijr}-C_{ij,r-1}$ over $r=1,\ldots,K$ cancels all intermediate terms, leaving $C_{ijK}-C_{ij0}=1$. This telescoping identity is the basic justification for turning cumulative probabilities into a categorical likelihood.

In the GLMM formulation, responses are assumed independent conditional on the random effects and the modeled covariates. A common specification takes $b_j$ to have a mean-zero multivariate normal distribution with covariance matrix $D$, independently across clusters. The conditional probabilities are evaluated at the realized $b_j$; marginally, observations from the same cluster need not be independent because they share that random effect.

The displayed equation uses a shared linear predictor across thresholds. It is a useful cumulative-link likelihood specification, but the restriction it imposes on how covariates relate to cumulative probabilities is a distinct modeling assumption. The likelihood formulas below do not by themselves establish that the restriction is suitable for a particular application.

Worked Example: Turning Cumulative Logits Into Category Probabilities

Suppose a hypothetical four-category response has thresholds $\theta_1=-0.8$, $\theta_2=0.4$, and $\theta_3=1.5$. For an observation with $\eta=0.2$, calculate the cumulative probabilities and the probability of each category under the cumulative logit model. The thresholds are strictly ordered, and the logistic CDF is increasing, so the resulting cumulative probabilities must be ordered.

The three arguments to the logistic CDF are $-0.8-0.2=-1.0$, $0.4-0.2=0.2$, and $1.5-0.2=1.3$. Hence $C_1=F(-1.0)\approx0.268941$, $C_2=F(0.2)\approx0.549834$, and $C_3=F(1.3)\approx0.785835$.

The category probabilities are the successive differences, including the endpoint probabilities: $p_1=C_1\approx0.268941$; $p_2=C_2-C_1\approx0.549834-0.268941=0.280893$; $p_3=C_3-C_2\approx0.785835-0.549834=0.236001$; and $p_4=1-C_3\approx1-0.785835=0.214165$.

All four values are nonnegative, and their sum is approximately $0.268941+0.280893+0.236001+0.214165=1.000000$, up to rounding. For example, the model assigns about 0.2360 probability to category 3 at this predictor value. These are model probabilities, not observed frequencies.

Likelihood Contributions and Their Derivation

Let $y_{ij}$ denote the observed category. Conditional on $b_j$, define $p_{ijr}=P(Y_{ij}=r\mid b_j)$ using the threshold-difference formula. Because $Y_{ij}$ takes exactly one category, its conditional likelihood contribution is $p_{ij,y_{ij}}$. Equivalently, using an indicator that equals one when $y_{ij}=r$ and zero otherwise, this contribution is

L_{ij}(b_j)=\prod_{r=1}^{K}p_{ijr}^{\,I(y_{ij}=r)}.

Only the factor corresponding to the observed category remains, because all other indicators are zero. Under conditional independence, the likelihood for all observations in cluster $j$ is the product of these contributions:

L_j(b_j)=\prod_i L_{ij}(b_j) =\prod_i\prod_{r=1}^{K}p_{ijr}^{\,I(y_{ij}=r)}.

For a GLMM, this is a conditional likelihood: it treats $b_j$ as given. If the goal is a likelihood for the observed responses without conditioning on the unobserved random effects, integrate the cluster likelihood against the random-effects density $f(b_j;D)$:

L_j=\int L_j(b_j)f(b_j;D)\,db_j.

The product over observations is inside the integral because all observations in that cluster share $b_j$. In general, this integral has no simple closed form and is evaluated numerically. The marginal likelihood across independent clusters is the product of their cluster marginal likelihoods.

Worked Example: A Conditional Likelihood for Observed Categories

Consider five hypothetical conditionally independent responses with the same covariates and random effects held fixed. Suppose the model gives category probabilities $(p_1,p_2,p_3,p_4)=(0.20,0.30,0.35,0.15)$, and the observed categories are $1,2,4,2,3$. These probabilities are nonnegative and sum to one, so they define a valid categorical distribution for each response.

The likelihood is the product of the probabilities assigned to the observed categories: $L=0.20\times0.30\times0.15\times0.30\times0.35=0.000945$. The category-2 probability appears twice because category 2 was observed twice.

The log likelihood is $\ell=\log(0.20)+\log(0.30)+\log(0.15)+\log(0.30)+\log(0.35)\approx-6.96432$. The likelihood itself is small because it is a product of five probabilities; its magnitude alone is not evidence that the model fits poorly. For model fitting, log likelihoods are used to avoid numerical underflow and to turn products into sums.

Integrating a Random Intercept

A conditional category probability at $b_j=0$ is not generally the marginal probability for a randomly selected cluster. The latter averages the conditional probability over the random-intercept distribution. This distinction matters for both likelihood construction and interpretation.

Worked Example: A Marginal Category Probability by Quadrature

For a hypothetical single observation, suppose a cumulative logit model has random intercept $b\sim N(0,0.25)$, so its standard deviation is $0.5$. Let the two thresholds of interest satisfy $\theta_1-x\beta=0.3$ and $\theta_2-x\beta=1.1$, where $x\beta$ is the fixed part of the predictor. With $\eta=x\beta+b$, the conditional probability of category 2 is $H(b)=F(1.1-b)-F(0.3-b)$.

The assumptions here are an ordered pair of thresholds, the logistic CDF, and the stated normal distribution for the random intercept. The marginal probability is the normal average $E\{H(b)\}$, not simply $H(0)$. A three-node Gauss-Hermite approximation for a standard normal uses nodes $0,\pm\sqrt{3}$ and weights $2/3,1/6,1/6$. Scaling the nodes by the random-intercept standard deviation gives $0,\pm0.5\sqrt{3}$. Therefore,

E\{H(b)\}\approx \frac{1}{6}H(-\sqrt{3}/2) +\frac{2}{3}H(0) +\frac{1}{6}H(\sqrt{3}/2).

Equivalently, the three conditional category-2 probabilities at the quadrature nodes, ordered by the arguments’ common shift, are approximately $0.11475773$, $0.17581759$, and $0.19607380$. Substituting these values and the weights gives $(1/6)(0.11475773)+(2/3)(0.17581759)+(1/6)(0.19607380)\approx0.16901698$, or approximately $0.1690$ to four decimal places.

At $b=0$, the conditional probability is $F(1.1)-F(0.3)\approx0.1758$, which differs from the integrated estimate. The difference illustrates why a conditional prediction at the mean random effect should not automatically be reported as a population-averaged probability. In practice, numerical integration accuracy should be checked by increasing the number of quadrature points or comparing with another suitable integration method.

Diagnostics, Alternatives, and Limitations

A likelihood can be computed exactly from the stated model and still be a poor description of the response process. Before interpreting fitted probabilities, inspect whether category coding preserves the real order, whether all categories are represented as intended, and whether fitted cumulative probabilities increase with the threshold. The corresponding category probabilities should be nonnegative and sum to one; violations indicate a coding or model-implementation problem. Also examine observed and fitted category frequencies across meaningful covariate ranges and clusters. Marked systematic discrepancies can signal an inadequate predictor structure, an unsuitable link, an inappropriate random-effects distribution, or a mismatch between the model and response process.

The cumulative logit and cumulative probit differ in their inverse link CDFs: one uses the logistic distribution and the other the standard normal distribution. Both construct category probabilities from adjacent cumulative probabilities, but their tail behavior and coefficient scales differ. A cumulative link uses the ordering by modeling probabilities below successive cutpoints. A nominal multinomial model does not use the category order in the same way. Adjacent-category likelihoods and continuation-ratio likelihoods express other comparisons among categories or stages; they are alternatives when the response-generating process calls for those comparisons, not algebraic substitutes that can be interchanged without changing the model.

Treating ordinal scores as continuous instead imposes a numerical scale and spacing interpretation that may not be justified. Conversely, a cumulative link retains the order but does not use a claim that adjacent categories are equally spaced. A further limitation is that the likelihood depends on threshold and predictor specifications, and GLMM marginal likelihoods require integration over random effects. Sparse or empty categories can also make threshold estimation unstable. These issues call for careful model checking and sensitivity analysis rather than treating the likelihood formula as automatic evidence of adequacy.

Common pitfall: It is incorrect to use $P(Y\leq k\mid b_j)$ as the likelihood contribution when the observed response is exactly category $k$. That is the probability of category $k$ or anything below it. The correct contribution is $P(Y=k\mid b_j)=C_{ijk}-C_{ij,k-1}$. In a GLMM, another common error is to multiply probabilities evaluated at $b_j=0$ and call the result the marginal likelihood. The marginal cluster likelihood integrates the product of conditional probabilities over $b_j$.

Check Your Understanding

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

For three ordered categories with cumulative probabilities $C_1=0.25$ and $C_2=0.70$, find all category probabilities and verify they sum to one.
Why must thresholds be ordered when the inverse link CDF is increasing?
In a cluster with a shared random intercept, why is the marginal likelihood not generally the product of probabilities evaluated at the mean random intercept?
For the cumulative probit model, what function converts $\theta_k-\eta_{ij}$ into $P(Y_{ij}\leq k\mid b_j)$?

Solutions. (1) The probabilities are $p_1=0.25$, $p_2=0.70-0.25=0.45$, and $p_3=1-0.70=0.30$; their sum is $0.25+0.45+0.30=1$. (2) Ordered thresholds give ordered CDF arguments at a fixed predictor, and an increasing CDF therefore gives ordered cumulative probabilities; their adjacent differences are then nonnegative. (3) The random intercept is shared, so conditional observation probabilities are multiplied at a common value of that intercept and the product is then averaged over its distribution. Averaging each factor first and multiplying generally gives a different result. (4) The standard normal CDF $\Phi$ gives $P(Y_{ij}\leq k\mid b_j)=\Phi(\theta_k-\eta_{ij})$.