Introduction and Overview

Earlier lessons produced a clean, descriptive view of the data. This lesson takes the next step from description to inference: linear regression is the workhorse model for explaining or predicting a continuous outcome from one or more predictors, fit by ordinary least squares (Stigler, 1981). Across four content sections we walk through this in order: the simple and multivariable model and what its coefficients mean (this section), the ANOVA decomposition and how to test the model and its individual coefficients (a later section), how to handle different types of predictor variables and detect collinearity (a later section), and how to detect and model interactions and give a regression a defensible causal interpretation (a later section). Model diagnostics, the residual and influence checks that show whether the fit can be trusted, run alongside the R work throughout.

Why Linear Regression?

Up to this point, most examples of relating an outcome to an exposure have been based on qualitative outcome variables, that is, variables that are categorical or dichotomous. Linear regression is suitable for modelling the outcome when it is measured on a continuous or near-continuous scale. Examples include birth weight, blood pressure, body mass index, and disease frequency at a regional level.

The Simple Regression Model

When only one predictor variable is used, the model is called a simple regression model. The term “model” denotes the formal statistical formula that describes the relationship between the predictor and the outcome.

In this equation, β₀ is the intercept (or constant), β₁ is the regression coefficient, and ε is the error term. The errors are assumed to be normally and independently distributed (ε ~ N(0, σ²)). We estimate these errors by residuals, the difference between the observed value and the value predicted by the model.

The Multivariable Model

Almost without exception, the regression models used by epidemiologists will contain more than one predictor variable. These are known as multiple regression or multivariable models.

A major difference from simple regression is that in the multivariable model, β₁ is an estimate of the effect of X₁ on Y after controlling for the effects of X₂. This is the key advantage of multivariable analysis: it accounts for confounding by extraneous variables.

# 0. Load the cleaned data + packages we will use ---------------------------
library(corrplot);  library(regclass);  library(caret)
phaa <- read.csv("phaa_survey_clean.csv", stringsAsFactors = FALSE)

# 1. Bivariate correlation between two numeric variables --------------------
cor.test(phaa$age, phaa$systolic_bp,
         method = "pearson")

# 2. Correlation matrix + visual --------------------------------------------
keep_num <- c("age", "bmi", "systolic_bp", "diastolic_bp",
              "phys_act_min", "discrimination_score",
              "social_support_score", "dep_score", "anx_score")
cor_mat <- cor(phaa[, keep_num], use = "complete.obs")
round(cor_mat, 2)
corrplot(cor_mat, method = "color", type = "upper",
         addCoef.col = "black", tl.col = "black", tl.srt = 45)

# 3. Set the reference level on a factor before fitting lm() ----------------
phaa$gender <- as.factor(phaa$gender)
phaa$gender <- relevel(phaa$gender, ref = "Woman")

# 4. Multivariable linear model for systolic BP -----------------------------
model_3 <- lm(systolic_bp ~ age + gender + smoker + bmi
                          + dep_score + phys_act_min,
              data = phaa)
summary(model_3)
confint(model_3)

# 5. Diagnostics: linearity, equal variance, normal residuals, outliers -----
par(mfrow = c(2, 2));  plot(model_3);  par(mfrow = c(1, 1))
VIF(model_3)            # multicollinearity
varImp(model_3)         # variable importance

Introduction and Overview

An earlier section set up the regression model. This section turns to the question of whether the model is doing useful work: how much of the variation in the outcome does it actually explain, and which individual coefficients are meaningfully different from zero? The ANOVA decomposition and the formal tests of model significance are how those questions get answered.

The ANOVA Table

The idea behind regression is that information in the X-variables can be used to predict the value of Y. The formal way this is approached is to ascertain how much of the sums of squares (SS) of Y we can explain with knowledge of the X-variable(s).

Here, k is the number of predictor variables in the model (not counting the intercept). When the SS are divided by their degrees of freedom (df), the result is a mean square, denoted MSM (model), MSE (error), and MST (total). The MSE is our estimate of the error variance σ², and the square root of σ² is called the root MSE or the standard error of prediction.

Assessing the Significance of a Linear Regression Model

We use the F-test from the ANOVA table to assess whether the predictors collectively have a statistically significant relationship with the outcome. The null hypothesis is H₀: β₁ = β₂ = … = β_k = 0.

In plain terms, the overall F-test asks a single yes-or-no question: taken as a set, do the predictors track the outcome better than simply predicting its overall mean for everyone? If the answer is no, the F ratio sits near 1; a large F with a small p-value says the predictors are doing real work.

Testing Individual Regression Coefficients

A t-test with n−(k+1) degrees of freedom is used to evaluate the significance of any individual regression coefficient. The usual null hypothesis is H₀: β_j = 0.

Introduction and Overview

An earlier section evaluated the model's overall fit. This section turns to a practical question that often determines whether your model gives sensible answers: are your predictors entered correctly? Continuous, categorical, indicator, and polynomial predictors all need different handling, and highly correlated predictors (multicollinearity) can destabilize coefficient estimates without obvious warning signs.

Types of Predictor Variables

The X-variables can be continuous or categorical. Categorical variables can be either nominal (levels with no meaningful numerical representation, e.g., race or city of residence) or ordinal (ordered levels, e.g., severity: low, medium, high). Nominal and ordinal variables with more than 2 levels must be converted to indicator variables before entering the regression.

Detecting Highly Correlated (Collinear) Variables

If the predictor variables are too highly correlated, a number of problems arise. The estimated effect of each variable depends on the other predictors in the model. With highly correlated predictors, the βs will be highly and negatively correlated, and in extreme cases none of the individual coefficients will be significantly different from zero despite a significant overall F-test.

Introduction and Overview

Earlier sections set up a model with main effects only. This section takes two final design steps: testing whether the effect of one predictor depends on another (interaction) and giving the resulting coefficients a defensible causal interpretation. Both push linear regression beyond a curve-fitting exercise into a tool for answering causal questions, anchored in the DAG-based framework you met in an earlier lesson.

Detecting and Modelling Interaction

Given the component cause model, we might expect to see interaction when 2 factors act synergistically or antagonistically. In previous sections, models contained only main effects, assuming the association of X₁ to Y is the same at all levels of X₂. An interaction term tests whether the effect of one variable depends on the level of another.

We assess interaction by testing whether β₃ = 0. If the interaction is absent (i.e., β₃ is not significantly different from 0), the main effects (additive) model is deemed adequate. If the interaction is needed, centring becomes useful because it allows us to interpret β₁ and β₂ as linear effects when the centred version of the other variable is zero.

Causal Interpretation of a Multivariable Linear Model

So far, we have focused on the technical interpretation of regression coefficients. When making causal inferences, extra care is needed to ensure that only the appropriate variables are included in the analysis. A causal diagram is very helpful in this regard.

Introduction and Overview

Earlier lessons worked through linear regression for continuous outcomes. This lesson takes the most common alternative in epidemiology: the binary outcome (disease present/absent, vaccinated/unvaccinated, alive/dead). The logistic model was formalised by Cox (1958) and has since become the workhorse logistic regression for dichotomous outcomes in epidemiology. The four content sections walk through this in order: why ordinary linear regression breaks down for binary outcomes and what the logistic model replaces it with (this section), the assumptions and how to test the overall model and individual coefficients (a later section), goodness-of-fit and predictive ability (a later section), and finally extensions including generalised linear models and exact logistic regression for sparse data (a later section). The odds ratios you computed by hand in an earlier course reappear here as exponentiated regression coefficients.

Why We Cannot Use Linear Regression for Dichotomous Outcomes

When the outcome variable is dichotomous (e.g., disease present/absent), ordinary linear regression is inappropriate for three fundamental reasons:

1. Non-normal errors: The residuals from a linear model with a binary outcome follow a binomial distribution, not a normal distribution, violating a key assumption of linear regression.
2. Heteroscedasticity: The variance of the residuals depends on the predicted probability, so the assumption of constant variance is violated.
3. Predictions outside 0–1: A linear model can produce predicted values less than 0 or greater than 1, which are nonsensical for probabilities.

The Logistic Model

The logistic regression model expresses the log odds of the outcome as a linear combination of predictors:

The inverse logit (or logistic function) converts back to the probability scale:

Note that, unlike linear regression, the logistic model has no error term because it models on the logit scale. The randomness enters through the binomial distribution of the outcome.

Odds and Odds Ratios

The odds of the outcome are p / (1 − p). The odds ratio for the kth predictor is obtained by exponentiating its coefficient:

For a dichotomous predictor, this is the odds ratio comparing the group coded 1 to the group coded 0, adjusted for all other variables in the model. Care is needed in interpretation: odds ratios are not risk ratios and can exaggerate the apparent strength of association when the outcome is common (Norton, Dowd, & Maciejewski, 2018).

Maximum Likelihood Estimation (MLE)

Unlike linear regression, which uses least squares, logistic regression uses maximum likelihood estimation (MLE). MLE is an iterative process that finds the parameter values most likely to have produced the observed data. The algorithm starts with initial estimates and refines them until convergence, the point at which the change in the log-likelihood between iterations falls below a specified criterion.

Introduction and Overview

An earlier section set up the model and showed how its coefficients become odds ratios on the exponentiated scale. This section turns to the same questions you asked of linear regression in an earlier lesson: are the assumptions met, is the overall model significant, what do the individual coefficients mean, and how do confounding and interaction enter the picture? Most of the framework is identical; the differences are mostly in how we compute and interpret coefficients on the log-odds scale.

Assumptions of Logistic Regression

Logistic regression requires two key assumptions: (1) independence of observations, and (2) linearity on the logit scale, that is, the relationship between each continuous predictor and the log odds of the outcome is linear. Note that the relationship on the probability scale will be non-linear (S-shaped).

Testing the Overall Model

The likelihood ratio test (LRT) compares the fitted model to the null model (intercept only). The test statistic is:

This statistic follows an approximate chi-squared distribution with degrees of freedom equal to the number of predictors. It can also be used to compare any two nested models (Eq 16.10), a full model versus a reduced model, to test whether the excluded variables contribute significantly.

The Wald test divides the coefficient by its standard error (following a Z distribution) and is more commonly reported by software. However, Wald tests can be unreliable when the true probability is near 0 or 1, or when the sample size is small (Vittinghoff & McCulloch, 2007).

Interpreting Coefficients

# 0. Load + ensure outcome has the REFERENCE level FIRST ("No")
phaa <- read.csv("phaa_survey_clean.csv", stringsAsFactors = FALSE)
phaa$hypertension <- factor(phaa$hypertension, levels = c("No", "Yes"))
phaa$smoker       <- factor(phaa$smoker,       levels = c("No", "Yes"))

# 1. Crude (unadjusted) model -----------------------------------------------
glm_1 <- glm(hypertension ~ smoker,
             data   = phaa,
             family = binomial(link = "logit"))
summary(glm_1)
exp(coef(glm_1))      # crude odds ratio
exp(confint(glm_1))   # 95% CI for the OR

# 2. Adjusted (multivariable) model -----------------------------------------
glm_2 <- glm(hypertension ~ smoker + age + gender + bmi + dep_score,
             data   = phaa,
             family = binomial)
summary(glm_2)
or  <- exp(coef(glm_2))
ci  <- exp(confint(glm_2))
round(cbind(OR = or, ci, p = summary(glm_2)$coef[,"Pr(>|z|)"]), 3)

# 3. Likelihood-ratio test for nested models ---------------------------------
anova(glm_1, glm_2, test = "Chisq")

# 4. Goodness-of-fit and discrimination --------------------------------------
library(generalhoslem);  library(DescTools);  library(pROC)
logitgof(obs = glm_2$y, fitted(glm_2))   # Hosmer-Lemeshow
PseudoR2(glm_2, which = "all")             # McFadden / Nagelkerke
phaa$pred_htn <- predict(glm_2, type = "response")
auc(roc(phaa$hypertension, phaa$pred_htn,
        levels = c("No", "Yes")))

Assessing Confounding and Interaction

Introduction and Overview

An earlier section covered model construction and interpretation. This section turns to model evaluation: residuals, formal goodness-of-fit tests (Hosmer & Lemeshow, 1980), predictive ability via discrimination (ROC curves; Hanley & McNeil, 1982) and calibration, the question of overdispersion, pseudo-R² statistics, and influential observations. Each gives a different angle on whether the model you've built is actually fit for the question you're asking (Steyerberg et al., 2010).

Model-Building Process

The model-building process for logistic regression follows the same general principles as for linear regression (Chapter 15): develop a causal diagram, perform unconditional (univariable) analyses, evaluate linearity of continuous predictors on the logit scale, and use automated selection methods with caution. Subject matter knowledge should guide decisions at every step.

Covariate Patterns and Data Structure

A covariate pattern is a unique combination of predictor values. Whether the data are treated as binary (one observation per row) or binomial/grouped (multiple observations per covariate pattern) has implications for how residuals and goodness-of-fit statistics are computed and interpreted.

Residuals

Pearson residuals and deviance residuals are used to assess model fit at the level of individual covariate patterns (Eq 16.16). Both types compare observed outcomes to predicted probabilities, but they differ in how discrepancies are scaled. These residuals are the building blocks of several goodness-of-fit tests.

Goodness-of-Fit Tests

Predictive Ability

Discrimination is most often summarised by the area under the receiver operating characteristic (ROC) curve (Hanley & McNeil, 1982), while calibration assesses agreement between predicted and observed risks (Steyerberg et al., 2010).

Selecting a cutpoint involves a trade-off between sensitivity and specificity. A lower cutpoint increases sensitivity but decreases specificity, and vice versa. The ROC curve provides a visual summary of this trade-off across all possible cutpoints.

Overdispersion

Pseudo-R² and Influential Observations

Pseudo-R² measures (e.g., McFadden’s, Cox-Snell, Nagelkerke) provide an indication of how much of the variation in the outcome is explained by the model. They are analogues of R² in linear regression but are not directly comparable. Values tend to be lower for logistic regression than for linear regression.

Influential observations can be identified using several diagnostic measures: outliers (large residuals), leverage (unusual covariate patterns), delta-betas (influence on individual coefficients), delta-χ², and delta-deviance (influence on overall fit). These diagnostics help identify observations that disproportionately affect the model.

Introduction and Overview

Earlier sections covered standard logistic regression. This section places it in a wider context. Logistic regression is one example of the generalised linear model (GLM) family, which also includes the linear, Poisson, and other regressions you'll meet later in this course. The section closes with exact logistic regression, the small-sample alternative when standard maximum-likelihood methods fail.

Generalised Linear Models (GLMs)

Logistic regression is a member of the broader family of Generalised Linear Models (GLMs). A GLM is defined by two key components: (1) a link function that relates the expected value of the outcome to the linear combination of predictors, and (2) the distribution of the outcome variable.

The canonical link is the “natural” link function for each distribution. For binary data, the canonical link is the logit. Non-canonical links (e.g., probit, complementary log-log for binary data) can also be used. GLMs are estimated using maximum likelihood, often with iterative algorithms such as Newton-Raphson or iteratively reweighted least squares. Quasi-likelihood estimation can be used when the full distribution is not specified, requiring only the mean-variance relationship.

Exact Logistic Regression

Standard logistic regression relies on large-sample approximations. When the dataset is very small or severely unbalanced (e.g., very few events), these approximations may be poor, and ML estimates can be biased or fail to converge. Exact logistic regression uses conditional maximum likelihood to produce exact P-values and confidence intervals without relying on large-sample theory. A widely used alternative is the penalised-likelihood approach of Firth (1993), which reduces small-sample bias and handles separation gracefully.

Introduction and Overview

An earlier lesson walked through linear regression with a fixed set of predictors. Real data rarely arrive with a clean “here are the four predictors you should use” instruction. This lesson turns to the broader question of how to build a defensible model when you have many candidate predictors and need to decide which to include. The four content sections move from goals and frameworks (this section), to reducing predictors and handling missing values (a later section), to modelling continuous predictor–outcome relationships flexibly with categorisation, polynomials, and splines (a later section), and finally to interactions, moderation, and selection criteria (a later section). Throughout, the central tension is between data-driven flexibility and theory-driven discipline, and the answer is usually closer to theory than the data alone would suggest (Babyak, 2004; Heinze, Wallisch, & Dunkler, 2018).

Why Model-Building Strategies Matter

When building a regression model, we need to decide on the goals of the analysis, incorporate both statistical considerations and subject matter knowledge, and balance the desire for parsimony (simplicity) with the desire for a model that “best fits” the data. The definition of “best fit” depends on the goal of the analysis, and the principles discussed in this chapter apply to all types of regression models.

Goals of the Analysis

If the goal is prediction, we want to keep any variables whose relationship with the dependent variable is questionable, because excluding them might lead to inaccurate predictions when future observations have extreme values for those variables. The details of specific predictors are of little consequence; we just want overall accuracy. Reporting guidelines such as TRIPOD (Collins, Reitsma, Altman, & Moons, 2015) provide a checklist for transparently documenting prediction models.

If the goal is understanding biological relationships, we want precise estimates of coefficients for the variables of interest. Careful attention must be paid to interaction and confounding effects. Factors likely to be confounders should be retained in the model regardless of statistical significance, while factors that are almost certainly not confounders should generally be excluded, especially if they are intervening variables, as their inclusion may bias results.

Steps in Building a Regression Model

The process of building a regression model follows a systematic set of steps. While statistical software handles the computation, the researcher must make many decisions along the way that require both subject matter expertise and statistical reasoning.

Building a Causal Model

Before beginning the model-building process, it is imperative to have a causal model in place, usually presented as a causal diagram. The diagram identifies potential causal relationships among the predictors and the outcome of interest.

Introduction and Overview

An earlier section settled the conceptual goals and the workflow. This section turns to two practical issues that determine which predictors actually make it into the final model. First: when you have more candidate predictors than your sample can support, how do you reduce the set without introducing bias? Second: missing values reduce your effective sample further. How should you handle them?

Reducing the Number of Predictors

It is sometimes necessary to reduce the number of predictors in the model-building process. Before undertaking any reduction, it is essential to identify the primary variables of interest and any variables that might be confounders or interacting variables; these should always be retained for consideration.

The Problem of Missing Values

Missing data are common in observational studies. Statistical programs use complete case analysis by default: only observations with no missing values for any variable are included. Even a relatively low overall percentage of missing values can result in a substantial reduction of the sample if missing data are spread across observations.

Dealing with Missing Data: Imputation

The two main alternatives to complete case analysis are imputation and analysis methods where missing data are ignorable. Imputation involves replacing missing data points with values predicted from available data.

Introduction and Overview

An earlier section reduced your candidate-predictor list. This section zooms into how to model the predictors that survived, specifically the continuous ones. Linear regression assumes a linear relationship between predictor and outcome by default, but real relationships are often curved. Categorisation, polynomials, and splines are three different ways to allow curvature; each has its own trade-offs.

Evaluating Continuous Predictor–Outcome Relationships

It is important to evaluate the structure of the relationship between a continuous predictor and the outcome before starting model building. The underlying assumption of linearity can be evaluated through diagnostics after fitting the model, but it is useful to explore the nature of the relationship beforehand.

Scatterplots & Smoothed Lines

Scatterplots are 2-way plots of the outcome (Y-axis) versus the continuous predictor (X-axis). They are most useful for continuous outcomes; scatterplots of dichotomous outcomes present as two lines of dots and are rarely informative by themselves.

Scatterplots can be greatly improved by adding a smoothed line through the centre of the data. All smoothed lines have a local-influence property: the position of the line at any value of x is influenced by nearby points but not by distant points.

Categorising Continuous Predictors

The assumption of linearity can be avoided by converting the continuous predictor into categories. However, this is generally not advisable for three reasons:

1. Categorisation involves the loss of information
2. It is unlikely that biological processes have a step-function relationship (i.e., sudden changes at specific cutpoints)
3. The choice of cutpoints is arbitrary and, if data-driven, may lead to biased results

That said, about 5 categories will usually suffice to control for confounding effects. A model with a categorised variable can be compared to one with a continuous (linear) variable using AIC or BIC.

Polynomial Models

Polynomials allow the regression line to follow a curve rather than a straight line. Power terms (e.g., x² or x³) are added to the model. Unlike smoothed lines, polynomial models have a global-influence property: the shape of the entire line is influenced by all the data.

Fractional Polynomials

Fractional polynomials (FPs) extend the idea of polynomial models by allowing power terms that are not restricted to positive integers. The most common set of powers to consider is: −3, −2, −1, −0.5, 0 (= ln), 0.5, 1, 2, 3. A 2-degree FP can fit a wide range of non-linear shapes and may be the most parsimonious way to model non-linearity.

Splines

An alternative to polynomial models is to fit a piecewise linear function. Points where the slope changes are called knot points. In the absence of prior evidence, knots may be chosen at percentiles of the predictor (e.g., 25th, 50th, 75th). Cubic splines allow for smoother transitions across knots compared to linear splines, producing more biologically plausible curves; the same logic extends to generalised additive models (Hastie & Tibshirani, 1986).

Introduction and Overview

Earlier sections settled which predictors enter the model and how each one is shaped. This section closes the lesson with the most consequential remaining design choices: which interaction terms (if any) to include, how to think about moderation as a substantive question, and which selection criterion to use when several plausible models compete.

Identifying Interaction Terms

It is important to consider including interaction terms when specifying the maximum model. There are five general strategies for creating and evaluating 2-way interactions:

Moderation: Interactions With a Substantive Story

In the regression literature, an interaction term is often called a moderation when it captures a substantive claim that the effect of one variable depends on another. A moderator W changes the slope of X → Y. Mathematically it is identical to an interaction (Y ~ X * W); the difference is conceptual.

# install.packages(c("MASS", "ggplot2", "interactions"))
library(MASS);  library(ggplot2);  library(interactions)
data("birthwt", package = "MASS")
bw <- birthwt
bw$smoke <- factor(bw$smoke, levels = c(0, 1),
                   labels = c("non-smoker", "smoker"))

# 1. Main-effects model (no moderation)
m_main <- lm(bwt ~ smoke + age, data = bw)

# 2. Moderation model: smoke * age
m_mod  <- lm(bwt ~ smoke * age, data = bw)

# 3. Is the moderation needed? Compare nested models with a likelihood-ratio
#    test (here, the partial F-test from anova() because models are linear).
anova(m_main, m_mod)

summary(m_mod)

# 4. Plot the moderation: smoking slopes at low / mean / high age
interact_plot(m_mod, pred = age, modx = smoke,
              interval = TRUE) +
  labs(y = "Birth weight (g)",
       title = "Effect of maternal age on birth weight, by smoking status")

# 5. Simple-slopes / Johnson-Neyman: across the range of maternal age,
#    where is the effect of smoking on birth weight distinguishable from zero?
sim_slopes(m_mod, pred = smoke, modx = age, johnson_neyman = TRUE)

Building the Model: Selection Criteria

Once the maximum model has been specified, you need to decide how to determine which predictors to retain. Both non-statistical and statistical criteria should be considered.

Specifying the Selection Strategy

Once criteria are established, there are several strategies for selecting which variables to include in the final model.