{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# SMART plausibility report\n",
    "\n",
    "*SMART nutrition survey, 2024*\n",
    "\n",
    "Cassion · data-analysis.cassion.dev\n",
    "\n",
    "[![Open In Colab](https://colab.research.google.com/assets/colab-badge.svg)](https://colab.research.google.com/github/thecassion/cassion-learning-platform/blob/main/apps/data-analysis/public/datasets/examples/smart-nutrition-survey-2024/plausibility.python.en.ipynb)"
   ],
   "id": "cell-000"
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## What a plausibility report is for\n",
    "\n",
    "Before a SMART survey's prevalence figure is used, the survey itself is\n",
    "assessed: were the measurements taken properly, were the children sampled\n",
    "properly, and is the resulting distribution shaped like a real population? A\n",
    "survey that fails these checks does not get a caveat — it gets rejected, or\n",
    "re-run.\n",
    "\n",
    "This report reproduces the standard checks and ends with a verdict.\n",
    "\n",
    "Every dataset on this platform is synthetic. No real child is described.\n",
    "\n",
    "## Setup"
   ],
   "id": "cell-001"
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "\n",
    "BASE = \"https://data-analysis.cassion.dev/datasets/\"\n",
    "\n",
    "smart = pd.read_csv(BASE + \"files/smart-nutrition-survey-2024.v1.csv\",\n",
    "                    dtype={\"child_id\": \"string\", \"sex\": \"string\"})\n",
    "reference = pd.read_csv(BASE + \"reference/who-2006-weight-for-lenhei.csv\")\n",
    "\n",
    "print(f\"{len(smart)} children, {smart['cluster'].nunique()} clusters, \"\n",
    "      f\"{smart['team'].nunique()} teams\")"
   ],
   "id": "cell-002"
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Check 1: completeness and impossible values"
   ],
   "id": "cell-003"
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "impossible = (\n",
    "    (smart[\"weight_kg\"] < 2) | (smart[\"weight_kg\"] > 30)\n",
    "    | (smart[\"height_cm\"] < 45) | (smart[\"height_cm\"] > 130)\n",
    ")\n",
    "\n",
    "pd.DataFrame({\n",
    "    \"count\": [\n",
    "        len(smart),\n",
    "        int(smart[\"age_months\"].isna().sum()),\n",
    "        int(smart[\"weight_kg\"].isna().sum()),\n",
    "        int(impossible.sum()),\n",
    "    ]\n",
    "}, index=[\"children\", \"missing age\", \"missing weight\", \"impossible measurement\"])"
   ],
   "id": "cell-004"
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Fourteen impossible measurements out of 930 is 1.5% — acceptable for a field\n",
    "survey, and they are excluded rather than corrected because you cannot know what\n",
    "the enumerator meant.\n",
    "\n",
    "## Check 2: z-scores, so the later checks have something to work on"
   ],
   "id": "cell-005"
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plausible = smart[~impossible.fillna(False)].copy()\n",
    "plausible[\"standard\"] = np.where(plausible[\"age_months\"] < 24, \"L\", \"H\")\n",
    "\n",
    "lying_should_stand = (plausible[\"standard\"] == \"H\") & plausible[\"measured_lying\"]\n",
    "stand_should_lie = (plausible[\"standard\"] == \"L\") & ~plausible[\"measured_lying\"]\n",
    "\n",
    "plausible[\"lenhei\"] = plausible[\"height_cm\"]\n",
    "plausible.loc[lying_should_stand, \"lenhei\"] -= 0.7\n",
    "plausible.loc[stand_should_lie, \"lenhei\"] += 0.7\n",
    "plausible[\"lenhei_key\"] = (plausible[\"lenhei\"] * 10).round() / 10\n",
    "\n",
    "lms = reference.set_index([\"sex\", \"lorh\", \"lenhei\"])[[\"l\", \"m\", \"s\"]]\n",
    "scored = plausible.join(lms, on=[\"sex\", \"standard\", \"lenhei_key\"])\n",
    "\n",
    "raw_z = ((scored[\"weight_kg\"] / scored[\"m\"]) ** scored[\"l\"] - 1) / (scored[\"l\"] * scored[\"s\"])\n",
    "\n",
    "def sd_at(row, n):\n",
    "    return row[\"m\"] * (1 + row[\"l\"] * row[\"s\"] * n) ** (1 / row[\"l\"])\n",
    "\n",
    "def who_adjust(row, z):\n",
    "    if pd.isna(z):\n",
    "        return z\n",
    "    if z > 3:\n",
    "        sd3, sd2 = sd_at(row, 3), sd_at(row, 2)\n",
    "        return 3 + (row[\"weight_kg\"] - sd3) / (sd3 - sd2)\n",
    "    if z < -3:\n",
    "        sd3, sd2 = sd_at(row, -3), sd_at(row, -2)\n",
    "        return -3 + (row[\"weight_kg\"] - sd3) / (sd2 - sd3)\n",
    "    return z\n",
    "\n",
    "scored[\"whz\"] = [who_adjust(r, z) for r, z in zip(scored.to_dict(\"records\"), raw_z)]"
   ],
   "id": "cell-006"
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Check 3: flagged records"
   ],
   "id": "cell-007"
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mean_z, sd_z = scored[\"whz\"].mean(), scored[\"whz\"].std()\n",
    "\n",
    "flags = pd.DataFrame({\n",
    "    \"flagged\": [\n",
    "        int((scored[\"whz\"].abs() > 5).sum()),\n",
    "        int(((scored[\"whz\"] - mean_z).abs() > 3 * sd_z).sum()),\n",
    "    ],\n",
    "    \"share\": [\n",
    "        (scored[\"whz\"].abs() > 5).mean(),\n",
    "        ((scored[\"whz\"] - mean_z).abs() > 3 * sd_z).mean(),\n",
    "    ],\n",
    "}, index=[\"WHO (fixed -5 to +5)\", \"SMART (3 SD from survey mean)\"]).round(4)\n",
    "flags"
   ],
   "id": "cell-008"
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "SMART treats above 2.5% flagged as a problem and above 5% as grounds for\n",
    "rejection. Both rules clear that comfortably here.\n",
    "\n",
    "## Check 4: the standard deviation of the z-score\n",
    "\n",
    "This is the single most informative number in the report. SMART expects the SD of\n",
    "weight-for-height z between 0.8 and 1.2. A real population has a spread close to\n",
    "1; measurement error widens it."
   ],
   "id": "cell-009"
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "analysable = scored[(scored[\"whz\"].abs() <= 5) & scored[\"whz\"].notna()]\n",
    "\n",
    "print(f\"mean weight-for-height z: {analysable['whz'].mean():.3f}\")\n",
    "print(f\"SD                      : {analysable['whz'].std():.3f}\")"
   ],
   "id": "cell-010"
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "At the top of the acceptable band. That is a warning, not a failure, and the next\n",
    "two checks locate its source.\n",
    "\n",
    "## Check 5: digit preference\n",
    "\n",
    "An enumerator reading a height board under pressure rounds. A team whose\n",
    "measurements pile up on `.0` and `.5` is not measuring to the millimetre they are\n",
    "recording."
   ],
   "id": "cell-011"
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "scored[\"last_digit\"] = ((scored[\"height_cm\"] * 10).round() % 10).astype(\"Int64\")\n",
    "\n",
    "digits = pd.crosstab(scored[\"team\"], scored[\"last_digit\"], normalize=\"index\") * 100\n",
    "digits.round(1)"
   ],
   "id": "cell-012"
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rounded = digits[[0, 5]].sum(axis=1).round(1)\n",
    "rounded.name = \"% ending .0 or .5\"\n",
    "rounded.to_frame().assign(expected=20.0)"
   ],
   "id": "cell-013"
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Team 2 records about 69% of its heights on a whole or half centimetre, against\n",
    "17 to 23% for the other teams. With ten possible last digits, 20% is what an\n",
    "unbiased team produces. This is the source of the wide standard deviation, and it\n",
    "is a training issue with a name attached.\n",
    "\n",
    "## Check 6: age heaping\n",
    "\n",
    "Ages reported by carers rather than documents pile up on whole years. It matters\n",
    "because age determines which growth standard applies — the length/height rule\n",
    "switches at exactly 24 months."
   ],
   "id": "cell-014"
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ages = scored[\"age_months\"].dropna()\n",
    "whole_years = ages.isin([12, 24, 36, 48, 60])\n",
    "\n",
    "print(f\"children at an exact whole year: {int(whole_years.sum())} ({whole_years.mean():.1%})\")\n",
    "\n",
    "ages.value_counts().reindex([23, 24, 25, 35, 36, 37, 47, 48, 49]).to_frame(\"children\")"
   ],
   "id": "cell-015"
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Sixty-six children recorded at exactly 24 months against 15 and 19 either side;\n",
    "eighty at 36 months against 28 and 17. That is not a birth pattern, it is\n",
    "rounding. Since 24 months is the boundary between the length and height\n",
    "standards, some of those children are being scored against the wrong reference.\n",
    "\n",
    "## Check 7: sex ratio"
   ],
   "id": "cell-016"
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "counts = scored[\"sex\"].value_counts()\n",
    "ratio = counts.get(\"m\", 0) / counts.get(\"f\", 1)\n",
    "\n",
    "chi_square = (counts.get(\"m\", 0) - len(scored) / 2) ** 2 / (len(scored) / 4)\n",
    "\n",
    "print(f\"boys: {counts.get('m', 0)}   girls: {counts.get('f', 0)}\")\n",
    "print(f\"ratio: {ratio:.3f}   chi-square against 1:1: {chi_square:.2f}\")"
   ],
   "id": "cell-017"
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A ratio near 1.0 and a chi-square well under 3.84 — no evidence that one sex was\n",
    "preferentially sampled or preferentially skipped.\n",
    "\n",
    "## Check 8: bias between teams"
   ],
   "id": "cell-018"
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "by_team = analysable.groupby(\"team\").agg(\n",
    "    children=(\"whz\", \"size\"),\n",
    "    mean_z=(\"whz\", \"mean\"),\n",
    "    sd_z=(\"whz\", \"std\"),\n",
    "    mean_height=(\"height_cm\", \"mean\"),\n",
    "    mean_weight=(\"weight_kg\", \"mean\"),\n",
    ")\n",
    "by_team[\"gam\"] = analysable.groupby(\"team\").apply(\n",
    "    lambda g: ((g[\"whz\"] < -2) | g[\"oedema\"]).mean(), include_groups=False\n",
    ")\n",
    "by_team.round(3)"
   ],
   "id": "cell-019"
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Team 3's mean z-score is -1.10 against -0.43 to -0.69 for the others, and its GAM\n",
    "comes out near 22% against 10 to 16%. Clusters were assigned to teams\n",
    "independently of nutrition status, so a real difference of half a z-score between\n",
    "teams is not a plausible reading. Team 3 is measuring long, light, or both."
   ],
   "id": "cell-020"
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "spread = by_team[\"mean_z\"].max() - by_team[\"mean_z\"].min()\n",
    "print(f\"spread in mean z-score across teams: {spread:.2f}\")\n",
    "print(\"SMART treats a between-team spread above ~0.3 z as a supervision problem.\")"
   ],
   "id": "cell-021"
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## The verdict"
   ],
   "id": "cell-022"
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "verdict = pd.DataFrame([\n",
    "    (\"Impossible measurements\", \"1.5% excluded\", \"pass\"),\n",
    "    (\"Flagged records\", f\"{flags.loc['SMART (3 SD from survey mean)', 'share']:.1%} SMART\", \"pass\"),\n",
    "    (\"SD of weight-for-height z\", f\"{analysable['whz'].std():.2f}\", \"warning\"),\n",
    "    (\"Digit preference\", \"team 2 at 69% on .0/.5\", \"fail\"),\n",
    "    (\"Age heaping\", f\"{whole_years.mean():.0%} on whole years\", \"warning\"),\n",
    "    (\"Sex ratio\", f\"{ratio:.2f}\", \"pass\"),\n",
    "    (\"Between-team bias\", f\"{spread:.2f} z spread\", \"fail\"),\n",
    "], columns=[\"check\", \"value\", \"result\"])\n",
    "verdict"
   ],
   "id": "cell-023"
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Would this survey be accepted?** Not as it stands. The prevalence figure is\n",
    "computable and the sampling looks sound, but two checks fail on the same root\n",
    "cause: one team measured differently from the others, and a second rounded its\n",
    "heights. Both are supervision and training problems, and both inflate the spread\n",
    "that the prevalence estimate rests on.\n",
    "\n",
    "The defensible action is not to publish 14.9% with a footnote. It is to re-measure\n",
    "team 3's clusters if the survey is still in the field, or to publish with the team\n",
    "comparison in the body of the report rather than an annex — so the reader can see\n",
    "that a sixth of the sample was measured by someone whose results do not match\n",
    "anyone else's.\n",
    "\n",
    "## What to report\n",
    "\n",
    "Every check with its value and its threshold, the verdict, and the action. A\n",
    "plausibility annex that reports only the checks that passed is not a plausibility\n",
    "report."
   ],
   "id": "cell-024"
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "name": "python",
   "file_extension": ".py"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}