You are what you ATE: Choosing an effect size measure for binary outcomes

Introduction

When summarising the effect of a treatment or intervention on an outcome measured on a continuous scale, it’s almost ubiquitous to represent the Average Treatment Effect (ATE) in the form of a difference of means between the treatment and control groups. When the outcome is binary (yes or no, event occurred or event did not occur), there are three commonly-used measures. The risk difference, also known as absolute risk reduction or difference of proportions, is an absolute measure of effect size: the proportion in one group minus the proportion in the other group. The relative risk, also known as the risk ratio, is a relative measure of effect size: the proportion in one group divided by the proportion in the other group. Finally, there is the odds ratio, another relative measure of effect: the odds in one group divided by the odds in the other group. The odds ratio is the odd one out in a few ways, and a bit more controversial.

In the first half of this post, I’ll dive into the controversy. In the middle, I’ll summarise the arguments using a meme. In the second half, I’ll work through an example of estimating odds ratios, risk differences, and relative risks in a simulated example, showing R code examples and hand-calculations. I’m taking a particular interest in the application to randomised controlled trials (RCTs) but most of the issues discussed here apply more broadly.

The great odds ratio debate

A lot of noise has been made about why an odds ratio may not be a desirable summary for a treatment effect on a binary variable². One of the biggest practical problems is that normal humans are unable to correctly interpret an odds ratio (Altman et al., 1998). But for statisticians, there are also deeper concerns to worry about.

In recent years, people³ have become more aware of the mathematical fact of noncollapsibility (Daniel et al., 2021; Greenland et al., 1999; Morris et al., 2022): the marginal⁴ odds ratio (population-level odds ratio) is not any kind of average of conditional odds ratios (individual-level odds ratios). Marginal odds ratios are the ones you might calculate directly from summary tables. Where do conditional odds ratios come from? They are artefacts of statistical models — logistic regression models, specifically. In the context of a randomised trial, these models predict the probability of an event occurring for each participant based on their treatment assignment and other variables measured prior to treatment (usually demographic characteristics, but could be any variable that is informative of the outcome)⁵. The second half of this post shows a worked example where every participant has the same conditional odds ratio of treatment, but the conditional odds ratio is not equal to the marginal odds ratio.

The problem of noncollapsibility is relevant to the analysis and reporting of randomised controlled trials. The ICH E9(R1) estimand framework for clinical trials requires a “population-level summary” be defined for each outcome without reference to method used to estimate it. In the language of causal inference, that’s describing an Average Treatment Effect (ATE), and the odds ratio isn’t one⁶.

Risk differences (the difference in probability attributable to a particular treatment) and relative risks (the ratio of probabilities due to a particular treatment) avoid both the interpretability and noncollapsibility problems of odds ratios. However, statisticians are usually taught to analyse binary outcomes using logistic regression, the direct output⁷ of which is odds ratios, so odds ratios are ubiquitous in the research literature. Someone on Stats Twitter rediscovers all of this every few months and starts a heated argument where nobody goes away happy.

Personally, I’m a fan of risk differences for communicating potential risks to an individual, ideally presented alongside the baseline level of risk. There’s some evidence that patients and clinicians find these measures easier to interpret than relative measures (Zipkin et al., 2014), especially when presented in the form of a natural frequency, like “6 in 1000 people” rather than “0.6%”.

So, with all those disadvantages, will anyone speak up in support of odds ratios? Should we all hang up our logistic regression hats for good?

In defense of logistic regression

Statistical analyses are often adjusted for covariates, either to reduce confounding (in observational studies) or improve power (in randomised trials). Adjusting for covariates requires specifying some kind of model. If we’re using a generalised linear model⁸, the mathematical form of covariate adjustment will depend on whether you are modelling log-odds (logistic regression model, where model parameters correspond to an odds ratio), probability (linear probability model, where model parameters correspond to a risk difference) or log-probability (quasi-Poisson or log-link binomial model, where model parameters correspond to a relative risk). When modelling risk differences or relative risks directly, it’s possible to end up with impossible predicted probabilities: “probabilities” which are less than zero or greater than one. Using logistic regression avoids this problem, because any real number on the log-odds scale translates to a probability between zero and one. Any odds ratio can be applied to any level of baseline risk without making mathematicians sad⁹.

There are also compelling — but not universally accepted — arguments that despite the difficulty in interpreting odds ratios, conditional odds ratios are more likely to be transportable between different levels of baseline risk than risk differences or relative risks (Doi et al., 2022; Senn, 2011). This is an empirical matter, not a mathematical one, and the evidence is not clear-cut. If we accept this, though, it is another reason to prefer logistic regression for statistical modelling. Effect size measures being “transportable” is another way of saying that the effects are closer to being additive on the scale that effect measure lives on (probability, log-probability, or log-odds). In the context of regression models, using a scale where the effect size is more transportable reduces the need for interaction terms, which can only be estimated well in large samples.

It’s possible to model the data using logistic regression, and then use that model to produce other quantities of interest, such as average risk differences. This approach is described in Permutt (2020), which has some of the best writing I’ve encountered in a statistics paper — being written in the form of a dialogue between a randomiser (the “causal inference” perspective) and a parameteriser (the “statistical modelling” perspective) walking through the Garden of Eden, planning to conduct and analyse a clinical trial. Does this method give us the best of both worlds: the modelling advantages of logistic regression and the interpretability advantages of risk differences and relative risks?

Permutt (2020) also considered what information different audiences might want from the results of a clinical trial: regulator, patient, and scientist. There’s a fourth audience which I think is worth considering, only very briefly mentioned by Permutt: the meta-analyst, trying to aggregate information from multiple trials.

Permutt argues in favour of the ATE being the main quantity of regulatory interest:

The average treatment effect should be of regulatory interest, however. The primary analysis of outcome should be of a variable that is reasonably linear in utility to the patient. Then, if and only if the average effect is positive, the ensemble of patients can be said to be better off under the test condition than under the control. This is perhaps the weakest imaginable statistical condition for approval of a drug product, but it is surely a necessary condition.

Permutt notes that studies designed to be able to detect average treatment effects are unlikely to be adequate for patient-specific decision making or providing a more detailed scientific understanding.

The most clearly expressed counter-argument to this comes from a blog post by Frank Harrell arguing against single-number summaries for treatment effects on binary outcomes:

Marginal adjusted estimates may be robust, but may not accurately estimate RD for either any patient in the RCT or for the clinical population to which RCT results are to be applied, because in effect they assume that the RCT sample is a random sample from the clinical population, something not required and never realized for RCTs.

This refers to the distinction between the population average treatment effect (PATE) and the sample average treatment effect (SATE). RCT participants are not random samples from any population, not even from the eligible pool of participants for a particular study. But the regulator’s decision that Permutt described earlier is motivated by generalising to a broader population, implicitly relying on properties of the PATE. It’s not clear to me whether there are likely to be any real-world scenarios where both (1) the practical conclusions drawn from the SATE and PATE would be different; and (2) the conditional odds ratio derived from the RCT sample is in agreement with the PATE but the risk difference SATE is not.

I am a coward and not (yet?) willing to take a strong position in this fight, but am always sympathetic to the idea that a single number is rarely sufficient to describe scientific evidence (see also: p-values). At some point I might write another blog about Harrell’s idea of plotting the distribution of estimated patient-level treatment effects, which is intriguing, although I struggle to see the practical purpose of it. Harrell’s other writing on this topic is also worth reading:

• Incorrect Covariate Adjustment May Be More Correct than Adjusted Marginal Estimates makes a similar point to White et al. (2021) with a detailed simulated example;
• Assessing Heterogeneity of Treatment Effect, Estimating Patient-Specific Efficacy, and Studying Variation in Odds ratios, Risk Ratios, and Risk Differences argues in favour of the use of odds ratios for modelling, and points out the extreme difficulty of identifying heterogeneous treatment effects even in very large studies; and
• Unadjusted Odds Ratios are Conditional demonstrates noncollapsibility and argues that conditional odds ratios are more useful than marginal odds ratios.

One issue which I remain unclear about is whether randomised trials powered to detect main treatment effects are likely to provide reasonable estimates of patient-specific baseline risk — a simpler task than patient-specific treatment effects, but still outside of the usual design remit for an efficacy trial. Common analytical approaches for clinical trials have good properties for estimating average treatment effects when the covariates are regarded as nuisance parameters (White et al., 2021), but are not guaranteed to perform so well if the effects of the covariates are themselves of interest.

Finally, there is the meta-analytic perspective to consider. An effect size which is less heterogeneous between studies is once again desirable. If arguments about the transportability of odds ratios by the Harrell, Senn, Doi, and others are to be believed, we should report conditional odds ratios, as those are likely to be the most useful for meta-analysts. The CONSORT guidelines for reporting randomised trials (Moher et al., 2010) state that both absolute and relative effect sizes should be reported for binary outcomes. I think there is merit in reporting all commonly-used effect size measures: risk difference, relative risk, and odds ratio. This provides the most flexibility for future meta-analysts.

Summary, in meme format

If all of that was a bit much to take in, maybe this will help:

In the rest of this post I’ll demonstrate how odds ratios are noncollapsible and show how to calculate risk differences and relative risks from logistic regression models using R and the marginaleffects package.

An example

This is the story of some hypothetical researchers who did a randomised controlled trial (RCT) where 560 patients were randomly assigned to either a treatment or control condition. The primary outcome of the trial was a binary measure, indicating whether a patient’s condition had worsened after 6 weeks. The scenario is based on one from Frank Harrell’s blog, with the details changed and embellished.

The hidden R code block below loads some R packages and sets up some simulated data for this trial.

library(tidyverse)
library(gt)
library(gtsummary)
library(broom)
library(cowplot)
library(marginaleffects)

theme_set(
  theme_cowplot(font_size = 11, rel_large = 1, rel_small = 1, rel_tiny = 1)
)

rct_data <- tribble(
  ~treatment, ~risk, ~outcome_prob,
  "Control", "Low risk", 1/7,
  "Control", "High risk", 3/5,
  "Treatment", "Low risk", 1/10,
  "Treatment", "High risk", 1/2
) %>%
  mutate(across(c(treatment, risk), fct_inorder)) %>%
  rowwise(everything()) %>%
  reframe(
    outcome = rep(c(1, 0), times = round(140*c(outcome_prob, 1 - outcome_prob)))
  )

Looking at a 2-way table of our hypothetical trial, we can see that the worse outcome (“1”) is more common in the control group (37%) than the treatment group (30%). This is an absolute risk reduction of 7%, a promising sign that our treatment may be beneficial.

tbl_cross(rct_data, outcome, treatment, percent = "column")

Just looking at a table won’t convince anybody, though. At this point, our hypothetical researchers asked a statistician to help out¹⁰. The hypothetical statistician immediately realised the need to estimate the treatment effect in some way that also quantified its uncertainty. A binary logistic regression model seemed like a suitable way to examine the effect of the treatment on the outcome.

lrm_unadj <- glm(
  outcome ~ treatment,
  family = binomial,
  data = rct_data
)

The table below shows the results of this logistic regression: an odds ratio of 0.73 (95% CI: 0.51 to 1.03). The odds ratio being less than 1 indicates that the treatment was beneficial (since the outcome occurring was bad, in this case) but the p-value of 0.074 indicates this effect is not statistically significant. The researchers were sad, their dreams of publishing in the British Medical Journal scuppered.

tbl_regression(lrm_unadj, exponentiate = TRUE)

Adding a covariate

Just as the statistician was finishing writing up the results, the researchers mentioned that the patients they recruited came from two different groups, one of which was known to have much worse outcomes than the other.

“Does that matter?” they asked.

“Well, using this information could improve your statistical power” said the statistician.

The statistician made another 2-way table, stratified by the risk group. Through some fluke, exactly half of the trial sample was low risk and the other half high risk, perfectly balanced across treatment arms. It turned out that while only 14% of the low risk patients in the control arm experienced the worse outcome, 60% of the high risk patients in the control arm experienced that outcome.

tbl_strata(rct_data, risk, tbl_cross, outcome, treatment, percent = "column")

The statistician decided to fit another logistic regression model, this time including risk as an additional variable in the model.

lrm_adj <- glm(
  outcome ~ treatment + risk,
  family = binomial,
  data = rct_data
)

Looking at the coefficients in this model, we can see that the estimated effect of treatment has increased: the odds ratio is further away from 1, now being 0.67 (95% CI: 0.45 to 0.99) instead of 0.73 (95% CI: 0.51 to 1.03). The p-value has decreased from 0.074 to 0.045, meaning that the treatment effect is now statistically significant.

tbl_regression(lrm_adj, exponentiate = TRUE)

The researchers went home happy, dreaming once again of publishing in a high impact factor journal, but the statistician was left with a nagging feeling that something was not quite right. The covariate was balanced between treatment groups, with equal numbers of low risk and high risk assigned to each treatment. Undergraduate linear models courses taught that adding a covariate to the model which is balanced between treatment groups shouldn’t change the coefficient for treatment, only the standard error. Why, then, did the estimates change like that? Is logistic regression different from linear regression in this regard?

Hand calculations: back to the 2-way table

Let’s take another look at the 2-way table and do some hand calculations.

The marginal risk difference is the difference between treatment and control groups in the proportion of patients experiencing the outcome event. From the table above, ignoring the covariate, we can calculate a risk reduction of 0.071:

104/280 - 84/280

[1] 0.07142857

The marginal odds ratio is the ratio of the oddses for treatment and the control group, where the odds is itself the number who experience the outcome event divided by the number who do not¹¹. We can calculate the marginal odds ratio here as 0.73, the same as the unadjusted logistic regression:

(84/196) / (104/176)

[1] 0.7252747

Now let’s look at the table stratified by risk, and see how that affects our calculations.

Using the same formula as before, we can calculate a risk difference of 0.043 in the low risk group, 0.100 in the high risk group, and the equally-weighted average of those two (because the two risk groups are equally common in this example) is 0.071, same as the marginal risk difference calculated above.

# conditional risk difference in the low risk group
20/140 - 14/140

[1] 0.04285714

# conditional risk difference in the high risk group
84/140 - 70/140

[1] 0.1

# average of conditional risk differences is the marginal risk difference
0.5*((20/140 - 14/140) + (84/140 - 70/140))

[1] 0.07142857

Using the formula for odds ratios, we can calculate an odds ratio of 0.67 in the low risk group and 0.67 in the high risk group. The conditional odds ratios in both groups are equal, and equal to the conditional odds ratio from the adjusted logistic regression, but are different from the marginal odds ratio. This is noncollapsibility in action.

# conditional odds ratio in the low risk group
(14/126) / (20/120)

[1] 0.6666667

# conditional odds ratio in the high risk group
(70/70) / (84/56)

[1] 0.6666667

In this example, we observed a treatment-by-risk-group interaction on the risk difference scale (absolute risk reduction of 0.041 in the low-risk group, absolute risk reduction of 0.100 in the high-risk group) but not on the odds ratio scale (conditional odds ratio of 0.67 in both groups). There is no mathematical requirement for the odds ratio scale to have a smaller interaction than other scales, but it is mathematically impossible for there to be no interaction on any scale unless the effect size is zero. This issue is discussed in detail from the perspective of psychological research in Rohrer & Arslan (2021).

Calculating risk difference and relative risks using the marginaleffects package

The researchers returned to the statistician, requesting adjusted and unadjusted risk differences and relative risks to report in their paper. Fortunately, the statistician was familiar with the marginaleffects package, which makes it fairly easy to calculate risk differences from logistic regression models¹².

avg_comparisons(lrm_unadj, variables = "treatment")

      Term Contrast Estimate Std. Error     z Pr(>|z|)   S  2.5 %  97.5 %
 treatment    T - C  -0.0714     0.0398 -1.79   0.0727 3.8 -0.149 0.00657

Columns: term, contrast, estimate, std.error, statistic, p.value, s.value, conf.low, conf.high 

avg_comparisons(lrm_adj, variables = "treatment")

      Term Contrast Estimate Std. Error     z Pr(>|z|)   S  2.5 %   97.5 %
 treatment    T - C  -0.0714     0.0354 -2.02   0.0437 4.5 -0.141 -0.00201

Columns: term, contrast, estimate, std.error, statistic, p.value, s.value, conf.low, conf.high 

avg_comparisons(lrm_unadj, variables = "treatment",
                comparison = "lnratioavg", transform = "exp")

      Term              Contrast Estimate Pr(>|z|)   S 2.5 % 97.5 %
 treatment ln(mean(T) / mean(C))    0.808   0.0749 3.7 0.639   1.02

Columns: term, contrast, estimate, p.value, s.value, conf.low, conf.high, predicted, predicted_hi, predicted_lo 

avg_comparisons(lrm_adj, variables = "treatment",
                comparison = "lnratioavg", transform = "exp")

      Term              Contrast Estimate Pr(>|z|)   S 2.5 % 97.5 %
 treatment ln(mean(T) / mean(C))    0.808   0.0458 4.4 0.655  0.996

Columns: term, contrast, estimate, p.value, s.value, conf.low, conf.high, predicted, predicted_hi, predicted_lo 

In pictures

One last way of showing the same thing: below is a plot of unadjusted and adjusted odds ratios, risk differences, and relative risks. You can see that for the odds ratio, the adjustment has moved the point estimate but not made the confidence interval narrower, whereas for the other measures, the point estimate has not changed but the confidence intervals are narrower for the adjusted estimate.

# odds ratio estimates
or_estimates <- bind_rows(
  "Unadjusted" = 
    tidy(lrm_unadj, conf.int = TRUE, exponentiate = TRUE) %>%
    filter(term == "treatmentT"),
  "Adjusted" = 
    tidy(lrm_adj, conf.int = TRUE, exponentiate = TRUE) %>%
    filter(term == "treatmentT"),
  .id = "model"
) %>%
  select(model, estimate, conf.low, conf.high)

# risk difference estimates
rd_estimates <- bind_rows(
  "Unadjusted" =
    avg_comparisons(lrm_unadj, variables = "treatment"),
  "Adjusted" =
    avg_comparisons(lrm_adj, variables = "treatment"),
  .id = "model"
) %>%
  select(model, estimate, conf.low, conf.high)

# relative risk estimates
rr_estimates <- bind_rows(
  "Unadjusted" =
    avg_comparisons(lrm_unadj, variables = "treatment",
                    comparison = "lnratioavg", transform = "exp"),
  "Adjusted" =
    avg_comparisons(lrm_adj, variables = "treatment",
                    comparison = "lnratioavg", transform = "exp"),
  .id = "model"
) %>%
  select(model, estimate, conf.low, conf.high)

# combine all estimands into a single data frame
all_estimates <- bind_rows(
  "Odds ratio" = or_estimates,
  "Risk difference" = rd_estimates,
  "Relative risk" = rr_estimates,
  .id = "estimand"
)

# plot!
all_estimates %>%
  mutate(model = fct_inorder(model), estimand = fct_inorder(estimand)) %>%
  ggplot(aes(x = estimate, xmin = conf.low, xmax = conf.high,
             y = model)) +
  geom_pointrange() +
  scale_y_discrete(limits = rev) +
  facet_wrap(vars(estimand), nrow = 3, scales = "free") +
  panel_border() +
  labs(y = NULL, x = "Estimate (95% CI)")

Acknowledgements

Thanks to those who provided feedback on earlier versions of this post, some which resulted in substantial revisions and improvements:

• Isabella Ghement
• Lachlan Cribb
• Solomon Kurz

This does not represent an endorsement of this post by any of the above. Mistakes and opinions are entirely my own.

Feedback welcome, especially any corrections in cases where I may have misunderstood or misattributed arguments, or made clear factual errors.