Half the difference was who lives there

Two districts are not comparable as they stand

Lesson 4 put three crude attack rates side by side and the comparison looked straightforward. It is not, and the reason is in the population file.

import pandas as pd

population = pd.read_csv("district-population-2024.v1.csv")

structure = (
    population.pivot(index="district", columns="age_band", values="population")
)
shares = structure.div(structure.sum(axis=1), axis=0)
print((shares * 100).round(1))

library(dplyr)

population |>
  mutate(share = population / sum(population), .by = district) |>
  tidyr::pivot_wider(id_cols = district, names_from = age_band, values_from = share)

District	0-4	5-14	15-44	45+
Nord	22%	30%	35%	13%
Centre	15%	25%	42%	18%
Sud	11%	20%	44%	25%

Nord has twice the share of under-fives that Sud has, and lesson 4 established that under-fives have three to four times the attack rate of working-age adults.

So Nord would have a higher crude rate than Sud even if every age-specific rate in the two districts were identical. The crude comparison mixes two things: how risky each district is, and who lives in it.

Direct standardisation

The method is one line of arithmetic applied consistently: compute each district’s age-specific rates, then apply them to one common population.

cases = pd.read_csv("cholera-line-list-2024.v1.csv")

observed = (
    cases.groupby(["district", "age_band"]).size().rename("cases").to_frame()
    .join(population.set_index(["district", "age_band"])["population"])
)
observed["rate"] = observed["cases"] / observed["population"]

standard = population.groupby("age_band")["population"].sum()
standard_share = standard / standard.sum()

standardised = (
    observed["rate"].unstack()          # district x age_band
    .mul(standard_share, axis=1).sum(axis=1)
)
crude = (
    cases.groupby("district").size()
    / population.groupby("district")["population"].sum()
)

comparison = pd.DataFrame({
    "crude_per_1000": 1000 * crude,
    "standardised_per_1000": 1000 * standardised,
})
comparison["difference"] = (
    comparison["standardised_per_1000"] - comparison["crude_per_1000"]
)
print(comparison.round(2))

observed <- cases |> count(district, age_band, name = "cases") |>
  left_join(population, by = c("district", "age_band")) |>
  mutate(rate = cases / population)

standard <- population |> summarise(pop = sum(population), .by = age_band) |>
  mutate(share = pop / sum(pop))

observed |>
  left_join(standard, by = "age_band") |>
  summarise(standardised = sum(rate * share), .by = district)

District	Crude	Standardised	Change
Nord	7.94	7.25	−0.69
Centre	6.47	6.60	+0.13
Sud	5.51	5.95	+0.44

Nord falls, Sud rises, and the gap between them narrows from 2.43 to 1.30 per 1,000. About half the apparent difference between the worst and the best district was age structure rather than risk.

The remaining 1.30 is real. Nord is genuinely worse — its age-specific rates are higher in three of four bands — and standardisation is what lets you say that rather than assert it.

What the standard population is, and why it matters

The standard is whatever population you apply every district’s rates to. Three common choices:

• The combined study population, as above. Simple, defensible, and internal — the standardised rates are comparable to each other and to nothing else.
• A national population. Lets you compare against other districts standardised the same way.
• A published world standard — the WHO or Segi world standard population. Lets you compare internationally, and produces numbers that look nothing like the crude rates.

State the standard. A standardised rate is only interpretable against others standardised to the same population, and two reports using different standards produce incomparable numbers that both look official.

print("Standard: combined population of the three districts, 145,000")

# Say it in the table caption, every time.

Indirect standardisation, and when you need it

Direct standardisation needs age-specific rates for every district, which needs enough cases in every cell. Where a district has four cases in a band, its age-specific rate is unstable and the direct method propagates that instability.

The indirect method inverts the problem: apply a standard set of rates to each district’s own population, and compare observed cases to expected.

standard_rates = observed.groupby("age_band").apply(
    lambda g: g["cases"].sum() / g["population"].sum()
)

expected = (
    population.assign(rate=population["age_band"].map(standard_rates))
    .assign(expected=lambda d: d["population"] * d["rate"])
    .groupby("district")["expected"].sum()
)
smr = cases.groupby("district").size() / expected
print(smr.round(3))

standard_rates <- observed |>
  summarise(rate = sum(cases) / sum(population), .by = age_band)

population |>
  left_join(standard_rates, by = "age_band") |>
  summarise(expected = sum(population * rate), .by = district) |>
  left_join(count(cases, district, name = "observed"), by = "district") |>
  mutate(smr = observed / expected)

The result is a standardised morbidity ratio: observed over expected, where 1.0 means the district has exactly the cases its age structure predicts. Above 1 is worse than expected, below is better.

Use indirect standardisation when cells are small, and say which method you used — the two answer slightly different questions and their numbers are not interchangeable.

Age is not the only confounder

Standardisation removes the variable you standardise on and nothing else. Cholera attack rates also vary with water source, population density, distance to a treatment centre and displacement status, and none of those is in the population file.

Standardising on age and then claiming the remaining difference is programme performance is the error this lesson creates the opportunity for, and the next lesson is entirely about it.

Report the pair

Cholera attack rate by district, weeks 1-16

  District   Crude    Age-standardised   Population
  Nord       7.94     7.25               48,000
  Centre     6.47     6.60               62,000
  Sud        5.51     5.95               35,000

  Standardised directly to the combined population of the three districts.
  About half the crude Nord-Sud difference is age structure: Nord has 22% of
  its population under five against Sud's 11%, and under-fives have three to
  four times the attack rate of adults aged 15-44.

Both columns, the standard named, and one sentence saying how much moved and why. A table with only the standardised column hides the fact that anything was adjusted; a table with only the crude column invites a comparison that is half demographic.

What comes next

Standardisation removed one alternative explanation. The next lesson lists the others — and asks what a difference between two districts is allowed to be attributed to when none of them can be removed.