How can an observational study imitate randomisation?

Start here

Why can't you compare treated and untreated patients directly?

In routine practice doctors choose who gets a treatment, and they choose for reasons. If the same reasons also change the outcome, the treated and untreated groups differ before any treatment is given.

A confounder is linked to the treatment and, separately, to the outcome. The bias and validity chapter introduces confounding and the Mantel-Haenszel estimate. This chapter shows how an observational study can compare like with like: first by stratification, then with a propensity score that handles many confounders at once. Austin (2011) describes the propensity score as a way to design and analyse an observational study “so that it mimics some of the particular characteristics of a randomized controlled trial”.

Stratification

How does comparing like with like remove confounding?

Split the patients by the confounder, so that within each group (stratum) the treated and untreated patients have the same severity. Compare inside each stratum, then combine the stratum results into one pooled estimate.

OR_MH = Σ (a_i × d_i ÷ n_i) ÷ Σ (b_i × c_i ÷ n_i)

a, b

treated patients who died, and who survived, in stratum i

c, d

untreated patients who died, and who survived, in stratum i

n

all patients in stratum i

Σ

add over all strata

Each stratum contributes in proportion to its size, so a small stratum does not dominate. Before pooling, check that the stratum odds ratios are similar (here 0.64 and 0.70, Breslow-Day test p = 0.88). If they differ clearly, the confounder is modifying the effect, and you report each stratum instead of one pooled number. The full working is in the worked example.

Where stratification stops working

• Each extra confounder multiplies the number of strata. Sex × 3 severity grades × 4 age bands × diabetes gives 48 strata, and many contain few or no treated or untreated patients.
• A continuous confounder must be cut into categories, and patients inside one category still differ.
• Only the confounders you stratify on are handled.

The propensity score

How can many confounders become one number?

Instead of stratifying on every characteristic, estimate for each patient the probability that a patient with those characteristics receives the treatment. Patients with the same probability can then be compared, whatever combination of characteristics produced it.

The score is a balancing score: among patients with the same propensity score, the distribution of the measured baseline characteristics is similar in treated and untreated patients (Austin, 2011). A treated and an untreated patient with scores of 0.62 may differ in age or diabetes, but across many such pairs the measured characteristics balance out, as they would after randomisation.

Two conditions must hold

1. No unmeasured confounders: every variable that affects both treatment and outcome has been measured and put into the model. This is the main assumption, and the data cannot prove it.
2. Positivity: every patient has a probability of receiving each treatment that is above 0 and below 1 (see below).

Build the model from characteristics measured before treatment was decided. A variable changed by the treatment is not a confounder and does not belong in the score.

Methods

What can you do with the propensity score?

Austin (2011) describes four methods: matching on the score, stratification on the score, inverse probability of treatment weighting, and covariate adjustment using the score.

1. Matching

Each treated patient is paired with an untreated patient whose score is close. Austin (2011) describes the choices: nearest-neighbour matching picks the untreated patient with the closest score, and nearest-neighbour matching within a caliper also requires the difference to be below a set maximum. Matching without replacement uses each untreated patient once. A treated patient with no untreated patient inside the caliper is excluded. Austin (2011) recommends matching on the logit of the score with a caliper of 0.2 of its standard deviation. Matching estimates the effect in the treated patients (the ATT).

2. Stratification on the score

Patients are ranked by the score and divided into five equal groups (quintiles). Rosenbaum and Rubin (1984) stated that stratifying on the quintiles of the propensity score removes approximately 90% of the bias due to measured confounders when estimating a linear treatment effect, extending Cochran's (1968) result for a single confounder (Austin, 2011).

3. Inverse probability of treatment weighting (IPTW)

Each patient is weighted by the inverse of the probability of the treatment they actually received: 1 ÷ e for a treated patient and 1 ÷ (1 − e) for an untreated patient, where e is the score. A treated patient with a score of 0.2 was unusual for being treated and gets weight 5. An untreated patient with a score of 0.2 gets 1.25. The weighted sample estimates the effect in the whole population (the ATE). Weights of 1 for treated and e ÷ (1 − e) for untreated patients estimate the ATT instead (Austin and Stuart, 2015). Very large weights make the estimate unstable, and stabilised or truncated weights are options.

4. Covariate adjustment using the score

The outcome is regressed on treatment and the estimated score, with a linear model for a continuous outcome or a logistic model for a yes/no outcome (Austin, 2011). Unlike the other three methods, it does not separate the design of the comparison from the analysis of the outcome.

Interactive

Change how strongly severity decides treatment, and compare the crude, stratified, matched and weighted estimates

All four estimates use the same patients. They differ in how much of the difference between treated and untreated patients they remove.

Balance

How do you know the propensity score worked?

Compare the baseline characteristics of treated and untreated patients after matching or weighting. If the score did its job, they look alike. The measure used is the standardised mean difference, not a p-value.

SMD = (mean₁ − mean₀) ÷ √[(SD₁² + SD₀²) ÷ 2]

SMD

standardised mean difference: the gap between groups in units of the pooled standard deviation

1, 0

treated and untreated patients

Yes/no variable

(p₁ − p₀) ÷ √{[p₁(1 − p₁) + p₀(1 − p₀)] ÷ 2}, where p is the proportion

Threshold

there is no universally agreed cut-off. An SMD below 0.1 has been taken to indicate a negligible difference (Austin, 2011)

• Report the SMD for every covariate before and after matching or weighting, in a table or a dot plot (love.plot() in the cobalt package draws one).
• For continuous covariates also compare variances and distributions, not only means (Austin, 2009).
• The c-statistic of the propensity score model “provides no information as to whether the propensity score model has been correctly specified” (Austin, 2011). A score that predicts treatment very well is not the aim. Balance is.
• If balance is poor, change the score model (add interactions or non-linear terms) and check again, before looking at the outcome.

Positivity

What if some patients could only ever have received one treatment?

A comparison needs patients in both groups who resemble each other. If every patient with a certain profile was treated, there is no untreated patient to compare them with, and no method can create one.

Rosenbaum and Rubin's conditions include that every subject has a non-zero probability of receiving either treatment, written 0 < P(treatment | characteristics) < 1 (Austin, 2011). This is positivity. In practice, look at the distribution of the score in the two groups. Where the distributions do not overlap:

• matching leaves many treated patients without a partner, and the matched result no longer describes all treated patients;
• weighting gives very large weights to the few patients in the sparse region, and the estimate depends heavily on them (Austin and Stuart, 2015);

When overlap is sparse, Austin (2011) describes two choices: restrict the analysis to the treated and untreated patients with similar covariate patterns, or stop the analysis. If you restrict, state clearly which patients the result applies to.

Interactive

Make treatment depend more strongly on severity and watch the two score distributions separate

The stronger the link between patient characteristics and treatment, the fewer patients have comparable partners in the other group.

Unmeasured confounding

What about confounders nobody measured?

Randomisation balances characteristics you measured and characteristics you did not. A propensity score balances only what went into it. The remaining question is how strong an unmeasured confounder would have to be to change the conclusion.

E-value = RR + √[RR × (RR − 1)]

RR

the observed risk ratio. If it is below 1, use 1 ÷ RR first

Rare outcome

an odds or hazard ratio with an outcome under 15% is used as a risk ratio (EValue package)

Common OR

RR ≈ √OR

Common HR

RR ≈ (1 − 0.5^√HR) ÷ (1 − 0.5^√(1/HR))

CI limit

repeat for the limit closer to 1. If the interval already includes 1, its E-value is 1

Example: an observed RR of 2 gives an E-value of 2 + √(2 × 1) = 3.41 (VanderWeele, 2022). An unmeasured confounder associated with both treatment and outcome by a risk ratio of 3.41 each could explain the association away. Weaker confounding could not. VanderWeele and Ding propose that observational studies intended to support causal conclusions report the E-value for the estimate and for the confidence limit closest to the null.

The E-value comes from a bounding factor (Ding and VanderWeele, 2016): a confounder with risk ratio RR_EU with the exposure and RR_UD with the outcome can shift the risk ratio by at most RR_EU × RR_UD ÷ (RR_EU + RR_UD − 1). For example, (2, 2) gives 2 × 2 ÷ 3 = 1.33. The E-value is the value that both would need, if equal, to make the bounding factor as large as the observed risk ratio.

A large E-value means that considerable unmeasured confounding would be needed to explain the result. It does not show that such confounding is absent. Judge it against what is known about the confounders that were not measured.

Interactive

Enter a published estimate and its confidence limit, and see how strong an unmeasured confounder would need to be

The E-value for the confidence limit tells you how strong confounding would need to be before the result stops excluding no effect.

Decide

Which adjustment method fits your study?

Answer the questions in order.

Comparison

Which of these ideas do residents mix up?

Randomisation or a propensity score?

Matching or weighting?

Standardised difference or p-value for balance?

The formula

Worked example: stratified odds ratio, balance, weights and E-value

Made-up teaching data: the pneumonia thesis from the problem at the top. 1000 patients, 500 given the drug and 500 not. Rows: treated, untreated. Columns: died, survived.

Step 1: crude and stratum odds ratios

1. Crude OR = (125 × 432) ÷ (375 × 68) = 54000 ÷ 25500 = 2.12 (95% CI 1.53 to 2.93).
2. Severe: (120 × 60) ÷ (280 × 40) = 7200 ÷ 11200 = 0.64. Mild: (5 × 372) ÷ (95 × 28) = 1860 ÷ 2660 = 0.70.
3. Breslow-Day test of equal odds ratios: χ² = 0.02, p = 0.88. The stratum odds ratios are similar, so pooling is reasonable.

Step 2: Mantel-Haenszel pooled odds ratio

1. Severe: a × d ÷ n = 120 × 60 ÷ 500 = 14.40 and b × c ÷ n = 280 × 40 ÷ 500 = 22.40.
2. Mild: 5 × 372 ÷ 500 = 3.72 and 95 × 28 ÷ 500 = 5.32.
3. OR_MH = (14.40 + 3.72) ÷ (22.40 + 5.32) = 18.12 ÷ 27.72 = 0.654.
4. 95% CI 0.43 to 0.99. Cochran-Mantel-Haenszel test χ² = 4.12, p = 0.042.

Crude OR 2.12 → OR_MH 0.65

Step 3: how unbalanced was severity?

Severe illness: 80% of treated and 20% of untreated patients. SMD = (0.80 − 0.20) ÷ √[(0.80 × 0.20 + 0.20 × 0.80) ÷ 2] = 0.60 ÷ 0.40 = 1.50, far above 0.1.

Step 4: inverse probability weights

A treated patient with a propensity score of 0.8 gets weight 1 ÷ 0.8 = 1.25. An untreated patient with the same score gets 1 ÷ 0.2 = 5. The untreated patient resembles the treated majority, so counts for more in the untreated group.

Step 5: E-value

Suppose a matched analysis gave RR 0.70 (95% CI 0.55 to 0.89).

1. RR below 1, so invert: 1 ÷ 0.70 = 1.429. E-value = 1.429 + √(1.429 × 0.429) = 2.21.
2. Confidence limit closer to 1 is 0.89: 1 ÷ 0.89 = 1.124. E-value = 1.124 + √(1.124 × 0.124) = 1.50.

Odds ratios, confidence intervals and tests computed in Python with statsmodels (Table2x2, StratifiedTable). The drug, the patients and the numbers are invented for teaching.

In a real paper

How are propensity score methods reported in real studies?

How to run it

How do you run propensity score analyses?

Estimating the score is ordinary logistic regression. Matching, weighting, balance tables and E-values are easiest in R, with packages written for the purpose.

R

# 1. Matching: nearest neighbour on the logit of the score, caliper 0.2 SD
library(MatchIt)
m.out <- matchit(drug ~ severity + age + diabetes + creatinine, data = icu,
                 method = "nearest", distance = "glm", link = "linear.logit",
                 caliper = 0.2, std.caliper = TRUE)
summary(m.out, un = TRUE)
matched <- match_data(m.out)

# 2. Balance: standardised differences before and after, threshold 0.1
library(cobalt)
bal.tab(m.out, un = TRUE, thresholds = c(m = 0.1))

# 3. Weighting for the average treatment effect
library(WeightIt)
W <- weightit(drug ~ severity + age + diabetes + creatinine, data = icu,
              method = "glm", estimand = "ATE")
summary(W)
ESS(W$weights[icu$drug == 1])

# 4. E-values
library(EValue)
evalues.RR(0.70, lo = 0.55, hi = 0.89)
evalues.HR(0.73, lo = 0.66, hi = 0.81, rare = TRUE)

• matchit defaults: method "nearest" (greedy nearest-neighbour), distance "glm" with a logit link, estimand "ATT", matching without replacement, and when a score is estimated, matching in descending order of the score. link = "linear.logit" matches on the logit. With std.caliper = TRUE the caliper is in standard deviation units.
• bal.tab: un = TRUE prints the unadjusted sample as well. thresholds = c(m = 0.1) flags mean differences above 0.1. cobalt computes the standard deviation for the standardised difference in the full, unadjusted sample.
• ESS gives the effective sample size of a weighted group, (Σw)² ÷ Σw².
• evalues.OR and evalues.HR need rare = TRUE if the outcome is under 15% at the end of follow-up and FALSE otherwise.

SPSS

1. Analyze > Regression > Binary Logistic. Select the treatment (coded 0/1) as the dependent variable and the baseline characteristics as covariates.
2. Click Save and tick Probabilities under Predicted Values. SPSS saves each patient's predicted probability of the higher-coded category, which is the estimated propensity score.
3. Use the saved score to form strata or weights, or take it to R for matching and balance tables.

How to report it

What should you write in the thesis or paper?

Say which confounders went into the score and why, which method you used and with what settings, how balance looked before and after, how many patients were lost, and how sensitive the result is to unmeasured confounding.

A worked sentence with the teaching numbers

“The propensity score for receiving the drug was estimated by logistic regression on severity, age, diabetes and creatinine, measured at admission. Treated patients were matched 1:1 without replacement to untreated patients by nearest-neighbour matching on the logit of the score, with a caliper of 0.2 standard deviations. After matching, all standardised mean differences were below 0.1. Death within 30 days had a risk ratio of 0.70 (95% CI 0.55 to 0.89). The E-values for the estimate and the upper confidence limit were 2.21 and 1.50.”

Also give the numbers matched and unmatched, the balance table before and after matching, and the crude estimate alongside the adjusted one (STROBE 16a). For weighting, give the weight type (ATE, ATT, stabilised or truncated) and the largest weights.

By specialty

Confounding by indication by specialty

Common myths

Common mistaken beliefs

• ✗Myth: A propensity score study is as good as a randomised trial.

  ✓Randomisation balances measured and unmeasured characteristics. A propensity score balances only the characteristics in the model, and assumes there are no unmeasured confounders (Austin, 2011).

• ✗Myth: Balance after matching is shown by p-values above 0.05.

  ✓Significance tests for balance are discouraged: p-values depend on sample size, and the matched sample is smaller. Use standardised differences (Austin, 2011).

• ✗Myth: A propensity score model with a high c-statistic is a good model.

  ✓The c-statistic provides no information about whether the model is correctly specified (Austin, 2011). Judge the model by the balance it produces.

• ✗Myth: Matching uses all the treated patients.

  ✓Treated patients with no untreated patient inside the caliper are excluded from the matched sample (Austin, 2011). Report how many.

• ✗Myth: Mantel-Haenszel stratification can handle any number of confounders.

  ✓Each added confounder multiplies the number of strata, and many strata end up with few or no patients in one group.

• ✗Myth: A large E-value proves there is no unmeasured confounding.

  ✓It shows that considerable unmeasured confounding would be needed to explain the result away (VanderWeele and Ding, 2017), not that none exists.

Why this matters to you

Five consequences of ignoring confounding by indication

1. A conclusion in the wrong direction

   In the teaching example the crude odds ratio of 2.12 made a beneficial drug look harmful.

2. A thesis that cannot answer its question

   If the variables that decided treatment were never recorded, no analysis can adjust for them. Decide what to collect before data collection starts.

3. Questions at the viva

   An examiner will ask how you checked balance. "The p-values were not significant" is the answer they expect to correct.

4. Misreading papers

   A matched analysis that dropped half the treated patients, or a weighted analysis driven by a few huge weights, can be presented as if it described everyone.

5. Harm to patients

   Observational results that are really confounding by indication can lead to useful treatments being withheld or useless ones being adopted.

History

History of stratification and propensity score methods

1. 1959

   Nathan Mantel and William Haenszel

   Published "Statistical aspects of the analysis of data from retrospective studies of disease" in the Journal of the National Cancer Institute, the basis of the Mantel-Haenszel pooled estimate.

2. 1968

   William G. Cochran

   Published "The effectiveness of adjustment by subclassification in removing bias in observational studies". Stratifying on the quintiles of a continuous confounder removed approximately 90% of the bias due to it.

3. 1983

   Paul Rosenbaum and Donald Rubin

   Published "The central role of the propensity score in observational studies for causal effects" in Biometrika, defining the propensity score.

4. 1984

   Rosenbaum and Rubin

   Published "Reducing bias in observational studies using subclassification on the propensity score", extending Cochran's quintile result to the propensity score.

5. 1985

   Rosenbaum and Rubin

   Published "Constructing a control group using multivariate matched sampling methods that incorporate the propensity score". Austin (2011) summarises it as suggesting matching on the logit of the score with calipers.

6. 2009

   Peter Austin

   Described balance diagnostics for propensity-score matched samples, including standardised differences and variance ratios.

7. 2011

   Peter Austin

   Recommended, from simulations, matching on the logit of the propensity score with a caliper of 0.2 of its standard deviation.

8. 2016

   Peng Ding and Tyler VanderWeele

   Derived a bounding factor for unmeasured confounding that makes no assumptions about the type of confounder.

9. 2017

   Tyler VanderWeele and Peng Ding

   Introduced the E-value in Annals of Internal Medicine and proposed that it be reported in observational studies intended to show causality.

Guidelines

Propensity scores and confounding in guidance

Reading a paper

Five checks on any propensity score study

• Which variables went into the score, were they measured before treatment, and is any important reason for choosing treatment missing?
• Is balance shown with standardised differences before and after, not with p-values?
• How many treated patients were left unmatched, or how large were the largest weights, and does the result still describe the population of interest?
• Does the stated method match the estimand (ATT for matching, ATE or ATT for weighting), and were matched pairs analysed as matched?
• Is there a sensitivity analysis for unmeasured confounding, such as an E-value for the estimate and the confidence limit?

Why is it called that?

Word surgery

Exam answers

Short-note outlines

Further reading

Related Stateazy articles