---
title: "Computing prevalence with WHO growth standards"
subtitle: "SMART nutrition survey, 2024 · R"
author: "Cassion · data-analysis.cassion.dev"
format:
  html:
    toc: true
    code-fold: false
---

## What this produces

Weight-for-height z-scores against the WHO 2006 growth standards, and GAM and SAM
prevalence with a design effect for the cluster sample.

Z-scores are deliberately **not** shipped in this dataset. Computing them is the
exercise, and reading a precomputed column teaches nothing about the three
decisions that move the answer: the range check, the measurement position, and
which flagging rule you apply.

Every dataset on this platform is synthetic. No real child is described, and
these results must not be cited as a real nutrition situation.

## Use the official package, not your own LMS code

The WHO standards are an LMS table, and reimplementing the interpolation is a
well-known source of quiet error. `anthro` is maintained by WHO and applies the
length/height adjustment itself.

```{r}
#| message: false
# install.packages("anthro")
library(readr)
library(dplyr)
library(anthro)

URL <- paste0(
  "https://data-analysis.cassion.dev/datasets/files/",
  "smart-nutrition-survey-2024.v1.csv"
)

smart <- read_csv(URL, col_types = cols(
  child_id = col_character(),
  .default = col_guess()
))

glimpse(smart)
```

## Range-check before flagging

Fourteen records hold impossible measurements — weights out by a factor of ten in
both directions, and heights entered in metres. These are not outliers to be
flagged statistically; they are data entry errors, and they must be excluded
before any flagging rule is applied, because a single 114 kg child moves the
survey mean that the SMART flag is computed against.

```{r}
impossible <- smart |>
  filter(weight_kg < 2 | weight_kg > 30 | height_cm < 45 | height_cm > 130)

impossible |> select(child_id, team, weight_kg, height_cm)
```

```{r}
plausible <- smart |>
  filter(
    is.na(weight_kg) | is.na(height_cm) |
      !(weight_kg < 2 | weight_kg > 30 | height_cm < 45 | height_cm > 130)
  )

cat("kept:", nrow(plausible), "of", nrow(smart), "\n")
cat("missing age:", sum(is.na(plausible$age_months)),
    " missing weight:", sum(is.na(plausible$weight_kg)), "\n")
```

## Measurement position is not optional

Children under two years are measured lying — recumbent length — and older
children standing. Length reads about 0.7 cm greater than height for the same
child, so mixing the two without adjustment biases every z-score for the younger
half of the sample.

Do **not** adjust the column by hand. Pass the position to `anthro_zscores` via
`measure` and let it apply the WHO rule, which depends on age as well as
position.

```{r}
z <- anthro_zscores(
  sex             = ifelse(plausible$sex == "m", 1, 2),
  age             = plausible$age_months,
  is_age_in_month = TRUE,
  weight          = plausible$weight_kg,
  lenhei          = plausible$height_cm,
  measure         = ifelse(plausible$measured_lying, "l", "h")
)

scored <- plausible |>
  mutate(whz = z$zwfl, who_flag = z$fwfl)

table(scored$who_flag, useNA = "ifany")
```

## Two flagging rules, two slightly different surveys

WHO flags are fixed bounds — a weight-for-height z-score outside -5 to +5 is
biologically implausible. SMART flags are relative: more than 3 SD from the
*survey* mean. They exclude different children, and a plausibility report states
which was used.

```{r}
mean_z <- mean(scored$whz, na.rm = TRUE)
sd_z   <- sd(scored$whz, na.rm = TRUE)

scored <- scored |>
  mutate(smart_flag = abs(whz - mean_z) > 3 * sd_z)

scored |>
  summarise(
    who_flagged   = sum(who_flag == 1, na.rm = TRUE),
    smart_flagged = sum(smart_flag, na.rm = TRUE),
    both          = sum(who_flag == 1 & smart_flag, na.rm = TRUE)
  )
```

The SMART rule is relative to a mean that the flagged observations themselves
influence, which is why the range check has to come first.

## Prevalence

**Oedema overrides anthropometry.** A child with bilateral pitting oedema is
severely acutely malnourished whatever their weight-for-height, so the SAM
numerator is not simply the count below -3 z-scores.

```{r}
analysable <- scored |> filter(who_flag == 0, !is.na(whz))

prevalence <- analysable |>
  summarise(
    children = n(),
    gam = mean(whz < -2 | oedema),
    sam = mean(whz < -3 | oedema),
    mean_z = mean(whz),
    sd_z   = sd(whz)
  )

prevalence |> mutate(across(c(gam, sam), ~ round(100 * .x, 1)),
                     across(c(mean_z, sd_z), ~ round(.x, 2)))
```

Global acute malnutrition near 14.9% and severe near 3.9%. That sits just under
the 15% WHO emergency threshold — which is exactly the position where the
analytical choices above stop being academic, because a different flagging rule
or a skipped position adjustment moves the figure across the line.

Note the standard deviation of the z-score. SMART expects it between about 0.8
and 1.2; a value above that suggests measurement error inflating the spread, and
this survey sits at the top of the acceptable range for a reason the next section
identifies.

## The team effect, which is not a nutrition finding

```{r}
analysable |>
  group_by(team) |>
  summarise(
    children = n(),
    mean_z   = round(mean(whz), 2),
    gam      = round(100 * mean(whz < -2 | oedema), 1),
    .groups = "drop"
  )
```

Team 3 reports GAM near 22% against 10 to 16% for the others, with a mean z-score
of -1.09 against -0.43 to -0.69. A real difference in nutrition status between
randomly assigned clusters of that size would be extraordinary. This is a
measurement artefact — a team measuring height long or weight light — and
reporting it as a geographic finding would send resources to the wrong clusters.

## The design effect

This is a cluster sample. Children within a cluster resemble each other, so the
effective sample size is smaller than the count of children and a confidence
interval computed as though the sample were simple random understates itself.

```{r}
clusters <- analysable |>
  mutate(case = whz < -2 | oedema) |>
  group_by(cluster) |>
  summarise(m = n(), y = sum(case), .groups = "drop")

k     <- nrow(clusters)
M     <- sum(clusters$m)
p_bar <- sum(clusters$y) / M

# Ultimate-cluster variance of a ratio estimator. Note it uses the deviation of
# each cluster's case count from what the overall rate predicts for its size —
# not the variance of the cluster rates, which ignores that clusters differ in
# size.
var_cluster <- (k / ((k - 1) * M^2)) * sum((clusters$y - p_bar * clusters$m)^2)
var_srs     <- p_bar * (1 - p_bar) / M

deff <- var_cluster / var_srs
icc  <- (deff - 1) / (mean(clusters$m) - 1)

cat(sprintf("clusters: %d   children: %d   mean cluster size: %.1f\n",
            k, M, mean(clusters$m)))
cat(sprintf("design effect: %.2f   ICC: %.3f\n", deff, icc))
cat(sprintf("effective sample size: %.0f of %d\n", M / deff, M))
```

A design effect near 2.3 and an ICC around 0.045 are ordinary for a nutrition
cluster survey. A DEFF below 1 or above about 4 usually means the calculation is
wrong rather than the survey unusual — this is a check worth running on your own
arithmetic before you report it.

```{r}
se_cluster <- sqrt(var_cluster)
se_srs     <- sqrt(var_srs)

cat(sprintf("GAM %.1f%%  (95%% CI %.1f - %.1f)  accounting for clustering\n",
            100 * p_bar, 100 * (p_bar - 1.96 * se_cluster), 100 * (p_bar + 1.96 * se_cluster)))
cat(sprintf("           (95%% CI %.1f - %.1f)  if clustering is ignored\n",
            100 * (p_bar - 1.96 * se_srs), 100 * (p_bar + 1.96 * se_srs)))
```

Two things to read there. The correct interval is wider, and a report that ignores
clustering claims precision the design cannot deliver — the most common way a
survey overstates what it knows.

And the upper bound crosses 15%. The point estimate sits below the WHO emergency
threshold; the interval does not rule out being above it. That is the sentence
the report needs, not a bare "14.9%, below the emergency threshold".

## What to report

Prevalence with its interval and the design effect used, the flagging rule named,
the exclusions counted, and the team comparison — because a survey where one team
differs from the others by half a z-score has a measurement problem that outranks
every prevalence figure in the report.
