Interpreting regression coefficients in R

A regression coefficient is the difference in the average outcome that goes with a one-unit difference in one predictor, comparing cases that have the same values on every other predictor in the model. The intercept is the model’s prediction when all predictors are zero, which is often not a case anyone observed.

This post reads one small regression coefficient by coefficient: the slope of a numeric predictor, the intercept, what centring does to it, and a categorical predictor. The data are R’s swiss data, and every example number is either a deterministic engine output or comes from a code block shown in the post, so you can run it yourself.

A coefficient is a comparison, not a rate of change

The sentence “for each one-unit increase in X, Y goes up by b” is the usual textbook reading, and it hides two things. First, the data are a set of different cases, not one case that was changed, so the coefficient describes how averages differ between cases that differ by one unit. Second, “the other predictors held fixed” is a statement about the model: it compares cases with equal values of the other predictors, which may be rare or absent in your data.

A sound reading of any coefficient therefore has five parts:

1. The outcome and its units.
2. The predictor and the size of one unit.
3. The direction and size of the difference, with its confidence interval.
4. The adjustment set: which other predictors are in the model.
5. The reach of the claim: association in this sample, not a causal effect, unless the design supports one.

The rest of the post fills those five parts in for real numbers.

The data and the model

R’s swiss data hold one row for each of 47 French-speaking Swiss provinces around 1888. The outcome is a standardised fertility measure (Fertility). Two of the predictors are Education, the percentage of army draftees with education beyond primary school, and Agriculture, the percentage of males working in agriculture. The rows are provinces, so everything below is about provinces, not about the people in them. The definitions are in the dataset documentation.

The question here is how fertility relates to education and agriculture across the 47 provinces. The outcome is continuous and each province is treated as one independent row, so a multiple linear regression fits; the post on multiple linear regression and when to switch explains how that choice is made, and the linear regression method reference lists what the engine computes. The engine fitted the model on all 47 rows, with 44 residual degrees of freedom.

Together the two predictors give R-squared = 0.4492. The Shapiro-Wilk check on the residuals found no evidence against normality, which is one of the assumptions that the post on linear regression assumptions walks through. The engine also noted one province above the 4/n Cook’s distance screening threshold (largest 0.092) and five Education values beyond Tukey’s fences; these are screening notes, not reasons to drop rows.

The table holds three coefficient tests with no multiplicity correction. Only the education slope is the focus of the run; the other two rows are read here for what they teach about interpretation, not as separate findings.

The same fit in R is one call. Below is an excerpt of what summary() and confint() print (call, residuals and significance codes trimmed), so you can match each column of the table above:

fit <- lm(Fertility ~ Education + Agriculture, data = swiss)
summary(fit)
confint(fit)

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 84.08005    5.78180  14.542  < 2e-16 ***
Education   -0.96276    0.18906  -5.092  7.1e-06 ***
Agriculture -0.06648    0.08005  -0.830    0.411    

Residual standard error: 9.479 on 44 degrees of freedom
Multiple R-squared:  0.4492,	Adjusted R-squared:  0.4242 

                 2.5 %     97.5 %
(Intercept) 72.4276058 95.7325022
Education   -1.3437949 -0.5817303
Agriculture -0.2277953  0.0948453

Reading a slope

Take the Education row. The estimate is -0.9628: across provinces, each additional percentage point of draftees with education beyond primary school goes with a 0.9628-unit lower average fertility measure, among provinces with the same share of males in agriculture. The 95% interval, -1.344 to -0.5817, is the range of slopes the data are compatible with, and it excludes zero. Its t statistic is -5.092 on 44 degrees of freedom, with p < 0.001. Read the interval first; the p-value only says that a t statistic this far from zero would be unusual if the true slope were zero and the model’s assumptions held (Greenland and colleagues 2016 list the common ways these quantities are misread).

Now the Agriculture row. The estimate is -0.06648, but the interval, -0.2278 to 0.09485, runs from a modest negative to a small positive slope. The data do not settle the direction. That is not evidence that agriculture and fertility are unrelated: a wide interval after adjusting for education says the data cannot pin this slope down. The same coefficient looks different in the next model, which is the point of the section on adjustment below.

Two habits make slopes easier to read.

Rescale to a useful unit. One percentage point of education is a small step. Ten points is a more natural comparison, and the slope and its interval simply scale by ten:

10 * coef(fit)["Education"]
10 * confint(fit)["Education", ]

Education 
-9.627626 

     2.5 %     97.5 % 
-13.437949  -5.817303 

So a province ten points higher in education has an average fertility measure about 9.628 units lower, with an interval for the ten-point slope from -13.44 to -5.817, at the same agriculture share. The conclusion has not changed; only the unit has.

Compare predictors on a common scale only with care. A slope per standard deviation puts predictors that were measured in different units side by side. Gelman (2008) argues for dividing by two standard deviations so that a numeric predictor is comparable with a binary one. Either way, the standardised slope depends on the spread of the sample as well as on the relationship, so a predictor that varied little in your data will look weak. Report the raw slope in its natural unit as well.

Reading the intercept

The intercept, 84.08, is the model’s predicted average fertility measure for a province with 0 percent education beyond primary school and 0 percent of males in agriculture. Its interval is 72.43 to 95.73.

Whether that number means anything depends on whether the zero point exists in the data. Check the ranges:

sapply(swiss[c("Education", "Agriculture")], range)

     Education Agriculture
[1,]         1         1.2
[2,]        53        89.7

No province has 0 for either predictor: the smallest education value is 1 and the smallest agriculture value is 1.2. The intercept is therefore a prediction slightly outside the data on both axes. Here the extrapolation is short, so the number is not absurd, but in many models zero is nowhere near the data (height at age 0 in a model fitted to adults, a blood pressure with a weight of 0 kg) and the intercept is a construction that anchors the line, not a finding. Do not interpret it as a result unless zero is a plausible, observed case. You still keep it in the model, because the predictions need it.

Centring gives the intercept a meaning

Subtract each predictor’s mean and the zero point becomes “an average province”. The slopes keep their meaning, because shifting a predictor by a constant does not change how the outcome differs between cases one unit apart. Only the intercept moves:

sw <- swiss
sw$Education_c   <- sw$Education   - mean(sw$Education)
sw$Agriculture_c <- sw$Agriculture - mean(sw$Agriculture)
fit_c <- lm(Fertility ~ Education_c + Agriculture_c, data = sw)
summary(fit_c)$coefficients
confint(fit_c)
c(mean(swiss$Education), mean(swiss$Agriculture), mean(swiss$Fertility))

                 Estimate Std. Error    t value     Pr(>|t|)
(Intercept)   70.14255319 1.38262043 50.7316047 1.104087e-40
Education_c   -0.96276262 0.18906343 -5.0922731 7.104842e-06
Agriculture_c -0.06647502 0.08004511 -0.8304695 4.107570e-01

                   2.5 %     97.5 %
(Intercept)   67.3560648 72.9290416
Education_c   -1.3437949 -0.5817303
Agriculture_c -0.2277953  0.0948453

[1] 10.97872 50.65957 70.14255

The slopes, their standard errors, t values and intervals match the uncentred fit digit for digit. The intercept is now 70.14 with an interval from 67.36 to 72.93: the model’s average prediction for a province at the mean of both predictors, and a much narrower interval than before, because that prediction sits in the middle of the data, not beyond its edge. The last line shows why 70.14 is the number to expect: when every predictor is centred at its mean, the intercept of a least-squares fit equals the mean of the outcome.

The means here are 10.98 for education and 50.66 for agriculture. Centring at a round value you can explain, such as 10 or 50, works the same way and may read better in a paper. Schielzeth (2010) discusses centring and scaling as ways to make coefficients interpretable. In a model with no interaction or squared terms, as here, centring changes only the intercept.

A prediction for any specific case is the intercept plus each slope times the predictor’s value. For a province with education 10 and agriculture 50 the model predicts

84.08 - 0.9628 x 10 - 0.06648 x 50 = 71.13

and predict() on either fit gives it:

predict(fit,   data.frame(Education = 10, Agriculture = 50))
predict(fit_c, data.frame(Education_c   = 10 - mean(sw$Education),
                          Agriculture_c = 50 - mean(sw$Agriculture)))

       1 
71.12868 
       1 
71.12868 

Both fits are the same model written two ways. A prediction describes the average of provinces like this one, not any single province.

Reading a categorical predictor

A categorical predictor enters the model as indicator columns. For a two-level factor there is one, and its coefficient is the difference in average outcome between that level and the reference level, at equal values of the other predictors.

To have something to read, the code below makes a two-level factor from Catholic, the percentage of Catholic residents: provinces with more than 50 percent are the “yes” group. That cut-off is ours, chosen to give a simple illustration; splitting a percentage into two groups throws information away, and a real analysis would usually keep the number. It then prints the coding R uses, so you can see which level is the reference:

sw$catholic_majority <- factor(ifelse(sw$Catholic > 50, "yes", "no"),
                               levels = c("no", "yes"))
table(sw$catholic_majority)
contrasts(sw$catholic_majority)
fit_g <- lm(Fertility ~ Education + Agriculture + catholic_majority, data = sw)
summary(fit_g)$coefficients
confint(fit_g)
tapply(sw$Fertility, sw$catholic_majority, mean)

 no yes 
 29  18 

    yes
no    0
yes   1

                      Estimate Std. Error   t value     Pr(>|t|)
(Intercept)          87.224351 5.11346960 17.057763 9.356933e-21
Education            -1.061408 0.16703857 -6.354268 1.115557e-07
Agriculture          -0.187135 0.07662275 -2.442290 1.878003e-02
catholic_majorityyes 10.578348 2.75579495  3.838583 4.010123e-04

                          2.5 %      97.5 %
(Intercept)          76.9120567 97.53664524
Education            -1.3982732 -0.72454244
Agriculture          -0.3416594 -0.03261046
catholic_majorityyes  5.0207579 16.13593824

       no       yes 
 66.22069  76.46111 

The contrasts() output shows the coding: “no” is coded 0 on the yes column and “yes” is coded 1, so “no” is the reference level. R names the coefficient by pasting the variable name and the level it stands for, catholic_majorityyes. R’s contrasts help page describes the coding options. Because the levels were set explicitly above, the reference is the one we chose. If you leave the levels unspecified, check the printed coding before you read the sign.

Reading it: among provinces with the same education and agriculture, those with a Catholic majority had an average fertility measure 10.58 units higher than those without (95% CI 5.021 to 16.14). The intercept, 87.22, is now the prediction for a reference-level province with 0 on both numeric predictors, and the same caution about zero applies.

Three details are worth noticing. The adjusted gap, 10.58, differs from the raw gap in group means, 76.46 minus 66.22, which is 10.24, because the groups also differ on the other predictors. The Education slope moved from -0.9628 to -1.061 and the Agriculture slope from -0.06648 to -0.1871, with the agriculture interval now excluding zero. A coefficient is a property of the model it sits in. And with more than two levels you get one coefficient per non-reference level, each against the same reference, so differences between two non-reference levels are not read directly off the table.

Adjustment changes what each coefficient means

The agriculture slope moved when a third predictor arrived, and that is expected. “Holding the other predictors fixed” depends on which predictors are there. In the two-predictor model the agriculture slope compares provinces with equal education; in the three-predictor model it compares provinces with equal education and equal Catholic majority status.

Westreich and Greenland (2013) call the temptation to read every coefficient in a model as an effect of its predictor the Table 2 fallacy. A model built to estimate one exposure’s association adjusts for other variables so the exposure’s coefficient is cleaner; those other coefficients were not set up to answer their own questions and are usually not interpretable as effects. If you have one focal predictor, interpret that coefficient and describe the rest as adjustment. If you need several coefficients to be meaningful, the model needs to be designed for that.

When predictors overlap strongly, individual coefficients become imprecise even though the model fits well; the post on multicollinearity and VIF shows what that does to an interval.

Reporting the result

A results sentence needs the outcome, the predictor with its unit, the direction and size, the interval, the adjustment, and a statement of scope. For the education slope:

In a linear regression across 47 provinces, each additional percentage point of draftees with education beyond primary school was associated with a 0.9628-unit lower average fertility measure (95% CI -1.344 to -0.5817, t(44) = -5.09, p < 0.001), adjusting for the share of males working in agriculture. The model explained 44.92% of the variation in the fertility measure (R-squared = 0.4492).

For the intercept, only report it when zero is meaningful or after centring: “At the mean of both predictors, the predicted fertility measure was 70.14 (95% CI 67.36 to 72.93).” For a categorical predictor name the reference level: “Provinces where more than 50 percent of residents were Catholic had a 10.58-unit higher average fertility measure than provinces at or below 50 percent, the reference group (95% CI 5.021 to 16.14), adjusting for education and agriculture.”

Check these before you submit:

• The unit of the predictor is stated, and the effect is quoted for a unit that makes sense.
• The interval is reported, not only the p-value.
• The other predictors in the model are named.
• The reference level of any factor is named.
• The wording is “associated with”, unless the design supports “effect of”. For provinces, the claim is about provinces: do not restate it for individuals. This is an ecological association, and individuals within a province can differ from the province pattern in either direction.
• The wording about the interval is about the data, not about a probability for the coefficient: it is the range of slopes compatible with the data, in the sense described in the overview of p-values and confidence intervals.

If the outcome is binary, a logistic model has coefficients on the log-odds scale and the reading changes; the logistic regression guide covers it. For help deciding on the model before you read it, how to choose a statistical test starts from the outcome type.

What the example run recomputes

The reproduce.R file in the example folder reads the bundled swiss data and refits the two-predictor model with lm(), so it recomputes the estimates, standard errors, t statistics, p-values and R-squared, and prints the sensitivity power calculation. It does not print the confidence intervals (confint(fit) in the first code block gives them), and it does not recompute the Shapiro-Wilk check or the engine’s influence and outlier notes. The centred model, the categorical model and the rescaled slope above are computed by the code blocks in this post, not by the engine, and you can run them in any R session with the built-in data.

# Results

- Test: Linear regression
- n: 47
- statistic: -5.092
- p-value: < 0.001
- effect size: R-squared = 0.4492
- 95% CI (Education): [-1.344, -0.5817]
- seed: 20261024

## Assumption checks
- normality (Shapiro-Wilk): passed=TRUE (no evidence against normality)

## Power (sensitivity, not observed power)
- minimum detectable effect at n=47: Cohen's f2 = 0.2196
- to detect a conventional medium effect: 68 total to detect f2=0.15
- alpha = 0.05, target power = 0.8

# Celeus - reproduce.R  (run from this package folder; bundled data in data/)
# Test: Linear regression  |  engine 0.31.7  |  catalog c39-d3cd
# Research context: not provided for this run.
# Reproduces the recorded result deterministically.

set.seed(20261024)
.data_path <- "data/interpreting-regression-coefficients.csv"
if (!file.exists(.data_path)) stop("Run reproduce.R from the package folder (the bundled data is in data/).")
data <- readr::read_csv(.data_path, show_col_types = FALSE)  # needs the readr package
clean_cells <- function (df) 
{
    ctrl_collapse <- "[\r\n\t\v\f]+"
    zero_width_chars <- intToUtf8(c(8203L, 65279L), multiple = TRUE)
    uni_space_chars <- intToUtf8(c(160L, 8239L, 8199L), multiple = TRUE)
    clean_text <- function(x) {
        out <- x
        for (ch in zero_width_chars) out <- gsub(ch, "", out, 
            fixed = TRUE)
        out <- gsub(ctrl_collapse, " ", out)
        for (ch in uni_space_chars) out <- gsub(ch, " ", out, 
            fixed = TRUE)
        trimws(gsub(" +", " ", out))
    }
    changed_cols <- character(0)
    cells_cleaned <- 0L
    retyped <- character(0)
    for (j in seq_along(df)) {
        col <- df[[j]]
        if (!is.character(col)) 
            next
        cleaned <- clean_text(col)
        n_changed <- sum(!is.na(col) & !is.na(cleaned) & col != 
            cleaned, na.rm = TRUE)
        if (n_changed > 0L) {
            changed_cols <- c(changed_cols, names(df)[j])
            cells_cleaned <- cells_cleaned + n_changed
        }
        df[[j]] <- cleaned
        non_blank <- cleaned[!is.na(cleaned) & nzchar(cleaned)]
        if (length(non_blank) && all(grepl("^[+-]?(0|[1-9][0-9]*)([.][0-9]+)?$", 
            non_blank))) {
            parsed <- suppressWarnings(as.numeric(cleaned))
            was_blank <- is.na(cleaned) | !nzchar(cleaned)
            lost <- (is.na(parsed) | is.infinite(parsed)) & !was_blank
            if (!any(lost)) {
                finite <- parsed[!is.na(parsed)]
                all_whole <- all(finite == floor(finite))
                in_range <- all(abs(finite) <= .Machine$integer.max)
                df[[j]] <- if (all_whole && in_range) 
                  as.integer(parsed)
                else parsed
                retyped <- c(retyped, names(df)[j])
            }
        }
    }
    names(df) <- clean_text(names(df))
    if (length(changed_cols)) {
        attr(df, "ingest_clean_columns") <- changed_cols
        attr(df, "ingest_clean_cells_cleaned") <- cells_cleaned
    }
    if (length(retyped)) 
        attr(df, "ingest_clean_retyped") <- retyped
    df
}
data <- clean_cells(data)
# Complete-case analysis: drop rows missing any analysis value.
.acols <- c("Fertility", "Education", "Agriculture")
data <- data[stats::complete.cases(data[, .acols, drop = FALSE]), , drop = FALSE]
result <- summary(lm(Fertility ~ Education + Agriculture, data = data))
print(result)

# Sensitivity: minimum detectable effect at n=47 (power=0.8, alpha=0.05). NOT observed power.
power <- pwr::pwr.f2.test(u = 2, v = 44, sig.level = 0.05, power = 0.8)
print(power)