Scenario Modelling

Dr. Arun Mitra Peddireddy
JHAPSMCON 2026
4^th Annual Jharkhand State Conference of IAPSM

The Counterfactual

What would have happened without lockdown / masks / vaccines?

We can never observe both as there is only one history.

Models are the only way to answer counterfactual questions.

Build a baseline. Build an alternative.

The comparison is the answer.

Three Intervention Types

A Hypothetical COVID-Like Wave

A pathogen seeds an age-structured population (POLYMOD survey: 0–20, 20–40, 40+). COVID-like timescales, with one-in-a-million initially infectious in each age band.

Parameter	Value	Meaning
`transmission_rate`	0.3	per-contact infection rate
`infectiousness_rate`	1/5	rate of E → I (5-day latent period)
`recovery_rate`	1/7	rate of I → R (7-day infectious period)
Population	POLYMOD demography (3 age groups)	from `socialmixr::polymod`
`time_end`	365 days	one-year horizon

The question: how do three interventions compare against doing nothing?

Building the Population

library(epidemics); library(socialmixr); library(tidyverse)

contact_data <- contact_matrix(
  socialmixr::polymod,
  age_limits = c(0, 20, 40),
  symmetric  = TRUE
)

contact_matrix    <- t(contact_data[["matrix"]])
demography_vector <- contact_data[["demography"]][["population"]]
names(demography_vector) <- rownames(contact_matrix)

initial_i <- 1e-6
init_row  <- c(S = 1 - initial_i, E = 0, I = initial_i, R = 0, V = 0)
initial_conditions <- rbind(init_row, init_row, init_row)
rownames(initial_conditions) <- rownames(contact_matrix)

uk_population <- population(
  name               = "Population",
  contact_matrix     = contact_matrix,
  demography_vector  = demography_vector,
  initial_conditions = initial_conditions
)

1: POLYMOD survey of who-meets-whom — used to build a contact matrix per day, per age group.
2: Three age groups → a 3 × 3 contact matrix C, plus a demography vector with one population count per group.
3: Initial conditions as proportions per age group (epidemics expects rows summing to 1).
4: The population object bundles structure + demography + starting state.

A Helper To Run Scenarios

run_scenario <- function(npi = NULL, vac = NULL, t_end = 365) {

  args <- list(
    population          = uk_population,
    transmission_rate   = 0.3,                # COVID-like
    infectiousness_rate = 1 / 5,
    recovery_rate       = 1 / 7,
    time_end            = t_end,
    increment           = 1.0
  )

  if (!is.null(npi)) args$intervention <- list(contacts = npi)
  if (!is.null(vac)) args$vaccination  <- vac

  do.call(model_default, args)
}

1: Baseline arguments — same for every scenario.
2: Optional levers — pass an intervention() object for an NPI, a vaccination() object for a campaign, or both. No arguments → counterfactual.

Scenario 0 — Counterfactual (Do Nothing)

counterfactual <- run_scenario()

Scenario A — 60-Day Lockdown (NPI)

lockdown_60 <- intervention(
  name       = "Lockdown",
  type       = "contacts",
  reduction  = matrix(c(0.5, 0.5, 0.3)),
  time_begin = 30,
  time_end   = 90
)

scenario_npi <- run_scenario(npi = lockdown_60)

1: Reduction per age group: kids 50%, working-age 50%, older 30%. Older population is already isolating, so a smaller marginal change.
2: The intervention is on for 60 days; β returns to baseline on day 91.

Scenario A — Plot

Scenario B — Vaccination Campaign

campaign <- vaccination(
  name        = "Campaign",
  time_begin  = matrix(60,    nrow = 3),
  time_end    = matrix(360,   nrow = 3),
  nu          = matrix(0.005, nrow = 3)
)

scenario_vac <- run_scenario(vac = campaign)

1: Start day per age group. Matrix shape n_age × n_doses (here 3 × 1). Same start day for everyone here.
2: End day per age group. A 300-day campaign window.
3: Daily vaccination rate — fraction of remaining S moved to V per day. 0.5%/day.

Scenario B — Plot

Scenario C — Combined (NPI + Vaccination)

scenario_combined <- run_scenario(npi = lockdown_60, vac = campaign)

All Four on One Plot

Outcomes

scenario	Peak day	Peak infectious	Final size	Vaccinated	Cases averted
0. Counterfactual	182	25,165,113	215,180,509	0	0
1. NPI only	264	25,166,324	214,756,837	0	423,672
2. Vaccination only	162	879,869	10,265,527	261,494,441	204,914,981
3. Combined	168	2,346	37,843	271,722,125	215,142,666

NPIs buy time.

Vaccination prevents.

Combined wins.

Hands-on

Open scripts/scenario_explorer.R. Setup loads the population + helpers; pick one task per pair.

Pair	Task
A	Weak lockdown — change `reduction` to `c(0.3, 0.3, 0.2)`. Does the wave still flatten?
B	Long lockdown — extend to days 30–150. Smaller rebound?
C	Late lockdown — start day 60 instead of day 30. How does peak day shift?
D	Earlier vaccination — start day 30 instead of day 60. Final size change?
E	Slow vaccination — `nu = 0.001` instead of 0.005. Still useful?
F	Vaccinate only — drop the NPI, but double `nu` to 0.01. Trade-off?

For your tweak: predict first, then run, then write down peak day / peak I / final size.

Reveal — Pair Results

Closing

A trajectory tells you what could happen.
A counterfactual tells you what to do.