Why the Reference Category Matters
A categorical predictor does not have a single numerical slope. Its coefficient table depends on how its categories are represented in the model matrix. With treatment coding, one category is selected as the reference; each displayed coefficient for that predictor compares one other category with the reference. Choosing and documenting the reference makes those comparisons deliberate and the output easier to interpret.
The choice changes the names and meanings of individual coefficients, not the underlying data or, when the recoding is full rank, the fitted model. A reference category is therefore a parameterization decision. It is not a declaration that the reference group is a control, a standard, or scientifically more important.
Reference Levels and GLMM Coefficients
For illustration, suppose a binary outcome \(Y_{ij}\) is measured for observation \(i\) in cluster \(j\), and a GLMM includes a categorical treatment predictor and a cluster-specific random intercept. A logistic random-intercept model can be written as
Here \(c_1\) is the reference level, \(I(X_{ij}=c_k)\) is 1 when the observation is in category \(c_k\) and 0 otherwise, and \(b_j\) is the cluster random intercept. This example assumes the specified conditional outcome distribution and link are appropriate, observations are conditionally independent given the random effects and covariates, and the random effects follow the stated distribution. The fixed-effect design must also be identifiable: for example, the categories being compared cannot be perfectly confounded with other model terms. Reference selection does not establish or repair these assumptions.
Under this coding, \(\beta_0\) is the linear predictor for the reference category when the other fixed-effect predictors are at their zero values and \(b_j=0\). Each \(\beta_k\) is the difference in conditional linear predictors between \(c_k\) and \(c_1\), at the same values of the other predictors and random effect. For a logistic model that scale is log odds. In a model with interactions, the meaning of a main-effect coefficient depends on the values of the interacting predictors; it is not automatically an overall comparison.
Worked Example: Reading a Treatment-Coded Coefficient Table
Suppose a hypothetical GLMM compares three treatments, with Placebo as the reference. Its fixed-effect estimates are: intercept \(-1.20\), Drug A coefficient \(0.55\), and Drug B coefficient \(0.90\). The model includes a cluster random intercept, and no other fixed predictors are included in this illustration.
For a cluster with \(b_j=0\), the estimated linear predictor for Placebo is \(-1.20\). For Drug A it is \(-1.20+0.55=-0.65\), and for Drug B it is \(-1.20+0.90=-0.30\).
Thus, the Drug A row reports the estimated conditional log-odds difference between Drug A and Placebo, and the Drug B row reports the estimated conditional log-odds difference between Drug B and Placebo. The coefficients do not directly compare Drug A with Drug B: that log-odds difference is \(0.55-0.90=-0.35\). These arithmetic comparisons describe the fitted parameterization; assessing uncertainty requires the relevant standard error or interval.
Setting and Verifying Levels in R
In R, a factor stores both category values and their level order. For an unordered factor, the first level is the reference under treatment contrasts. If levels are not explicitly supplied, R generally sorts the observed values when constructing a factor; the resulting order may not match the intended scientific comparison and can depend on how labels sort.
For a hypothetical dataset named dat with a character variable arm, an explicit declaration is: arm levels are Placebo, Drug A, Drug B; then use factor(arm, levels = c("Placebo", "Drug A", "Drug B")). If arm is already an unordered factor, use relevel(arm, ref = "Placebo") to make Placebo the reference. Assigning the factor back to dat, rather than changing a temporary copy, ensures the model uses the intended coding.
Check the result before fitting. The R expression levels(dat$arm) displays the order, and table(dat$arm, useNA = "ifany") displays the observed counts, including missing values. For an unordered factor, contrasts(dat$arm) shows its contrast matrix. With treatment contrasts and Placebo first, the rows for Drug A and Drug B have a 1 in their respective columns, while the Placebo row has zeros. In a fitted model, model.matrix(fit) can be inspected to confirm the columns actually used in the fixed-effects design.
Worked Example: Preventing an Accidental Alphabetical Reference
A hypothetical dataset contains 12 observations: four labelled “Usual care,” four “Dose 1,” and four “Dose 2.” The analysis intends to compare each dose with Usual care. If the labels are turned into a factor without specifying levels, a sorted order may put Dose 1 first, making it the treatment reference rather than Usual care.
Set the intended order explicitly with factor(group, levels = c("Usual care", "Dose 1", "Dose 2")). The expected level order is then Usual care, Dose 1, Dose 2, and the expected treatment-coded coefficient columns correspond to Dose 1 and Dose 2. Each column compares that dose with Usual care.
Check the counts as well as the order: table(dat$group) should show four observations for each of the three categories. If the model-matrix columns instead show Dose 2 as the omitted category, inspect the factor levels and contrasts used in the actual model frame. A printed coefficient table is evidence about the coding that was fitted, not about the coding the analyst intended.
What Changes When the Reference Changes?
Changing the reference under treatment coding changes the coefficient coordinates. It does not change fitted values if the new coding represents the same model space and the model remains identifiable. The result follows directly by expressing each category’s linear predictor in terms of the new baseline.
Worked Example: Re-Expressing Coefficients With Drug B as Reference
Use the hypothetical estimates from the first example: Placebo is the reference, the intercept is \(-1.20\), Drug A has coefficient \(0.55\), and Drug B has coefficient \(0.90\). The category-specific fixed linear predictors are therefore Placebo \(-1.20\), Drug A \(-0.65\), and Drug B \(-0.30\).
With Drug B as the new reference, the intercept becomes Drug B’s linear predictor, \(-0.30\). The Placebo coefficient is Placebo minus Drug B: \(-1.20-(-0.30)=-0.90\). The Drug A coefficient is Drug A minus Drug B: \(-0.65-(-0.30)=-0.35\). Equivalently, these are \(0-0.90=-0.90\) for Placebo and \(0.55-0.90=-0.35\) for Drug A.
For example, the new Drug A linear predictor is \(-0.30+(-0.35)=-0.65\), exactly the old value. The random intercept and its distribution have not changed. Only the fixed-effect reference and resulting coefficient representation have changed.
Checks, Alternatives, and Common Pitfalls
A short coding audit can catch errors before the coefficient table is interpreted. Confirm that the variable is a factor rather than a numeric code, that the observed labels match the declared levels, and that the intended reference is first for treatment coding. Then inspect the contrast matrix or model matrix and compare its columns with the coefficient names. Also check for empty or unused levels: a category with no observations cannot supply an empirical comparison, and may produce a rank-deficient design or an aliased coefficient. Removing unused levels with droplevels() can prevent stale categories from confusing the audit, but it cannot create information about an unobserved category.
The standard treatment contrasts are not the only valid parameterization. Sum contrasts express category effects relative to a constrained overall reference point; other contrasts can encode planned comparisons. With a full-rank design, these alternatives can describe the same fitted means while producing different coefficient values and interpretations. Use treatment coding when comparisons against one named category are the desired coefficient-table summaries. Choose other contrasts when their comparisons answer the scientific question more directly, and document that choice.
A frequent mistake is to assume the first category shown in a raw data table is the reference, or to assume R always uses the desired category. Neither is guaranteed: the fitted factor levels and contrasts determine the design matrix. Another common mistake is to apply relevel() to an ordered factor and expect ordinary treatment coding. Ordered factors can use polynomial contrasts by default; relevel() is intended for unordered factors. If nominal categories are intended, represent them as an unordered factor and verify the contrasts rather than relying on defaults.
A further limitation is that reference coding changes interpretation, not model adequacy. It cannot fix sparse data, separation, confounding, an unsuitable link, or a misspecified random-effects structure. Nor does choosing a clinically familiar baseline by itself justify a causal conclusion. It only defines how the fitted model’s parameters are expressed.
Use the intended labels and an explicit level order when creating or recoding the factor.
Check counts, missing values, and whether any declared levels have no observations.
Inspect the factor contrasts or fitted model matrix, then confirm the coefficient names represent the planned comparisons.
Document the reference level and contrast system alongside the analysis code so the comparisons can be reproduced.
Check Your Understanding
Answer each question, then compare with the solutions below.
- An unordered factor has levels “Standard,” “Low,” and “High,” in that order, and uses treatment contrasts. Which category is the reference, and what does the coefficient for High compare?
- A model with category A as reference has an intercept of \(0.40\), a B coefficient of \(0.30\), and a C coefficient of \(-0.20\). What are the new intercept and the coefficient for A if C becomes the reference?
- Why is it not enough to inspect the raw values in a column to confirm the reference category used by a GLMM?
- A factor has a declared “Pilot” level but no observations in that category. Can its coefficient be interpreted as an estimated comparison with the reference? Explain.
Solutions. (1) Standard is the reference. The High coefficient is the High-minus-Standard difference on the linear-predictor scale, holding other terms fixed. (2) The category-specific predictors are A \(0.40\), B \(0.70\), and C \(0.20\). The new intercept is \(0.20\); the A coefficient is \(0.40-0.20=0.20\). (3) The raw column does not show the factor’s level order or contrast system used to construct the model matrix; inspect the factor and fitted design. (4) No. Without observations in Pilot, the data provide no comparison for that category; the coefficient is not an estimable data-supported effect.