Skip to content

Commit

Permalink
Test to demonstrate bug GH-406
Browse files Browse the repository at this point in the history
Created a test that shows the issue of the SEIR parameters not being
randomly drawn per a slot wheras the outcome parameters are. Test
currently fails.
  • Loading branch information
TimothyWillard committed Dec 2, 2024
1 parent 70f9cd7 commit 2f93ec3
Show file tree
Hide file tree
Showing 2 changed files with 246 additions and 0 deletions.
168 changes: 168 additions & 0 deletions examples/tutorials/config_sample_2pop_vaccine_scenarios.yml
Original file line number Diff line number Diff line change
@@ -0,0 +1,168 @@
name: sample_2pop
setup_name: minimal
start_date: 2020-02-01
end_date: 2020-08-31
nslots: 10

subpop_setup:
geodata: model_input/geodata_sample_2pop.csv
mobility: model_input/mobility_sample_2pop.csv

initial_conditions:
method: SetInitialConditions
initial_conditions_file: model_input/ic_2pop.csv
allow_missing_subpops: TRUE
allow_missing_compartments: TRUE

compartments:
infection_stage: ["S", "E", "I", "R", "V"]

seir:
integration:
method: rk4
dt: 1
parameters:
sigma:
value: 1 / 4
gamma:
value: 1 / 5
Ro:
value:
distribution: truncnorm
mean: 2.5
sd: 0.1
a: 1.1
b: 3
omega_pess:
value: 0.02
omega_opt:
value: 0.01
nu_pess:
value: 0.01
nu_opt:
value: 0.03
transitions:
#infections
- source: ["S"]
destination: ["E"]
rate: ["Ro * gamma"]
proportional_to: [["S"],["I"]]
proportion_exponent: ["1","1"]
# progression to infectiousness
- source: ["E"]
destination: ["I"]
rate: ["sigma"]
proportional_to: ["E"]
proportion_exponent: ["1"]
# recovery
- source: ["I"]
destination: ["R"]
rate: ["gamma"]
proportional_to: ["I"]
proportion_exponent: ["1"]
#vaccination (offers complete protection)
- source: ["S"]
destination: ["V"]
rate: ["nu_pess + nu_opt"]
proportional_to: ["S"]
proportion_exponent: ["1"]
# waning of vaccine-induced protection
- source: ["V"]
destination: ["S"]
rate: ["omega_pess + omega_opt"]
proportional_to: ["V"]
proportion_exponent: ["1"]

seir_modifiers:
scenarios:
- no_vax
- pess_vax
- opt_vax
modifiers:
pess_vax_nu: # turn off nu_opt, only nu_pess left
method: SinglePeriodModifier
parameter: nu_opt
period_start_date: 2020-02-01
period_end_date: 2020-08-31
subpop: "all"
value: 0
pess_vax_wane: # turn off omega_opt, only omega_pess left
method: SinglePeriodModifier
parameter: omega_opt
period_start_date: 2020-02-01
period_end_date: 2020-08-31
subpop: "all"
value: 0
pess_vax: # turn off all vaccination
method: StackedModifier
modifiers: ["pess_vax_nu","pess_vax_wane"]
opt_vax_nu: # turn off nu_pess, only nu_opt left
method: SinglePeriodModifier
parameter: nu_pess
period_start_date: 2020-02-01
period_end_date: 2020-08-31
subpop: "all"
value: 0
opt_vax_wane: # turn off omega_pess, only omega_opt left
method: SinglePeriodModifier
parameter: omega_pess
period_start_date: 2020-02-01
period_end_date: 2020-08-31
subpop: "all"
value: 0
opt_vax: # turn off all vaccination
method: StackedModifier
modifiers: ["opt_vax_nu","opt_vax_wane"]
no_vax: # turn off all vaccination
method: StackedModifier
modifiers: ["pess_vax","opt_vax"]



outcomes:
method: delayframe
outcomes:
incidCase: #incidence of detected cases
source:
incidence:
infection_stage: "I"
probability:
value:
distribution: truncnorm
mean: 0.5
sd: 0.05
a: 0.3
b: 0.7
delay:
value: 5
incidHosp: #incidence of hospitalizations
source:
incidence:
infection_stage: "I"
probability:
value: 0.05
delay:
value: 7
duration:
value: 10
name: currHosp # will track number of current hospitalizations (ie prevalence)
incidDeath: #incidence of deaths
source: incidHosp
probability:
value: 0.2
delay:
value: 14

# outcome_modifiers:
# scenarios:
# - test_limits
# modifiers:
# # assume that due to limitations in testing, initially the case detection probability was lower
# test_limits:
# method: SinglePeriodModifier
# parameter: incidCase::probability
# subpop: "all"
# period_start_date: 2020-02-01
# period_end_date: 2020-06-01
# value: 0.5

Original file line number Diff line number Diff line change
@@ -0,0 +1,78 @@
import os
from pathlib import Path
import shutil

from click.testing import CliRunner
import pandas as pd
import pytest

from gempyor.simulate import _click_simulate


@pytest.fixture
def setup_sample_2pop_vaccine_scenarios(tmp_path: Path) -> Path:
tutorials_path = Path(os.path.dirname(__file__) + "/../../../../examples/tutorials")
for file in (
"config_sample_2pop_vaccine_scenarios.yml",
"model_input/geodata_sample_2pop.csv",
"model_input/mobility_sample_2pop.csv",
"model_input/ic_2pop.csv",
):
source = tutorials_path / file
destination = tmp_path / file
destination.parent.mkdir(parents=True, exist_ok=True)
shutil.copy(source, destination)
return tmp_path


def test_random_seir_parameter_draw_per_slot(
monkeypatch: pytest.MonkeyPatch, setup_sample_2pop_vaccine_scenarios: Path
) -> None:
# Test setup
monkeypatch.chdir(setup_sample_2pop_vaccine_scenarios)

# Test execution of `gempyor-simulate`
runner = CliRunner()
result = runner.invoke(_click_simulate, ["config_sample_2pop_vaccine_scenarios.yml"])
assert result.exit_code == 0

# Get the contents of 'spar' and 'hpar' directories as DataFrames
spar_directory: Path | None = None
hpar_directory: Path | None = None
for p in setup_sample_2pop_vaccine_scenarios.rglob("*"):
if p.is_dir() and p.name == "spar":
spar_directory = p
elif p.is_dir() and p.name == "hpar":
hpar_directory = p
if spar_directory is not None and hpar_directory is not None:
break

def read_directory(directory: Path) -> list[pd.DataFrame]:
dfs: list[pd.DataFrame] | pd.DataFrame = []
for i, f in enumerate(sorted(list(directory.glob("*.parquet")))):
dfs.append(pd.read_parquet(f))
dfs[-1]["slot"] = i
dfs = pd.concat(dfs)
return dfs

hpar = read_directory(hpar_directory)
spar = read_directory(spar_directory)

# Test contents of 'spar'/'hpar' DataFrames
assert (
hpar[
(hpar["subpop"] == "large_province")
& (hpar["quantity"] == "probability")
& (hpar["outcome"] == "incidCase")
]["value"].nunique()
== 10
)
assert (
hpar[
(hpar["subpop"] == "small_province")
& (hpar["quantity"] == "probability")
& (hpar["outcome"] == "incidCase")
]["value"].nunique()
== 10
)
assert spar[spar["parameter"] == "Ro"]["value"].nunique() == 10

0 comments on commit 2f93ec3

Please sign in to comment.