Skip to content

Commit

Permalink
Merge pull request #244 from CITCOM-project/update-covasim
Browse files Browse the repository at this point in the history 
Update: edited example_vaccine.py
  • Loading branch information
@f-allian
f-allian authored Nov 10, 2023
2 parents 3b00417 + 4bb55d9 commit 282a225
Show file tree
Hide file tree
Showing 3 changed files with 60 additions and 157 deletions.
4 changes: 4 additions & 0 deletions examples/covasim_/vaccinating_elderly/README.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -13,6 +13,10 @@ four test cases: one focusing on each of the four previously mentioned outputs.

Further details are provided in Section 5.3 (Prioritising the elderly for vaccination) of the paper.

**Note**: this version of the CTF utilises the observational data collector in order to separate the software execution
and testing. Older versions of this framework simulate the data using the custom experimental data collector and the
`covasim` package (version 3.0.7) as outlined below.

## How to run
To run this case study:
1. Ensure all project dependencies are installed by running `pip install .` from the top
Expand Down
182 changes: 25 additions & 157 deletions examples/covasim_/vaccinating_elderly/example_vaccine.py
View file
Original file line number Diff line number Diff line change
@@ -1,33 +1,27 @@
# -*- coding: utf-8 -*-
import os
import logging
import pandas as pd
import numpy as np
import covasim as cv # Version used in our study is 3.07
import random
from causal_testing.specification.causal_dag import CausalDAG
from causal_testing.specification.scenario import Scenario
from causal_testing.specification.variable import Input, Output
from causal_testing.specification.causal_specification import CausalSpecification
from causal_testing.data_collection.data_collector import ExperimentalDataCollector
from causal_testing.data_collection.data_collector import ObservationalDataCollector
from causal_testing.testing.causal_test_case import CausalTestCase
from causal_testing.testing.causal_test_outcome import Positive, Negative, NoEffect
from causal_testing.testing.estimators import LinearRegressionEstimator
from causal_testing.testing.base_test_case import BaseTestCase

import os
import logging

logger = logging.getLogger(__name__)
logging.basicConfig(level=logging.DEBUG, format="%(message)s")

ROOT = os.path.realpath(os.path.dirname(__file__))


def test_experimental_vaccinate_elderly(runs_per_test_per_config: int = 30, verbose: bool = False):
"""Run the causal test case for the effect of changing vaccine to prioritise elderly. This uses the experimental
data collector.
def setup_test_case(verbose: bool = False):
"""Run the causal test case for the effect of changing vaccine to prioritise elderly from an observational
data collector that was previously simulated.
:param runs_per_test_per_config: Number of times to run each input configuration (control and treatment) per test.
Hence, the total number of runs per test will be twice this value.
:param verbose: Whether to print verbose details (causal test results).
:return results_dict: A dictionary containing ATE, 95% CIs, and Test Pass/Fail
"""
Expand Down Expand Up @@ -63,13 +57,10 @@ def test_experimental_vaccinate_elderly(runs_per_test_per_config: int = 30, verb
# 4. Construct a causal specification from the scenario and causal DAG
causal_specification = CausalSpecification(scenario, causal_dag)

# 5. Instantiate the experimental data collector for Covasim
covasim_parameters_dict = {"pop_size": 50000, "pop_type": "hybrid", "pop_infected": 1000, "n_days": 50}
control_input_configuration = {"covasim_parameters_dict": covasim_parameters_dict, "target_elderly": False}
treatment_input_configuration = {"covasim_parameters_dict": covasim_parameters_dict, "target_elderly": True}
data_collector = CovasimVaccineDataCollector(
scenario, control_input_configuration, treatment_input_configuration, runs_per_test_per_config
)
# 5. Instantiate the observational data collector using the previously simulated data
obs_df = pd.read_csv("simulated_data.csv")

data_collector = ObservationalDataCollector(scenario, obs_df)

# 6. Express expected outcomes
expected_outcome_effects = {
Expand All @@ -85,161 +76,38 @@ def test_experimental_vaccinate_elderly(runs_per_test_per_config: int = 30, verb
causal_test_case = CausalTestCase(
base_test_case=base_test_case, expected_causal_effect=expected_effect, control_value=0, treatment_value=1
)

# 7. Obtain the minimal adjustment set for the causal test case from the causal DAG
minimal_adjustment_set = causal_dag.identification(base_test_case)

# 8. Build statistical model
# 8. Build statistical model using the Linear Regression estimator
linear_regression_estimator = LinearRegressionEstimator(
vaccine.name, 1, 0, minimal_adjustment_set, outcome_variable.name
treatment=vaccine.name,
treatment_value=1,
control_value=0,
adjustment_set=minimal_adjustment_set,
outcome=outcome_variable.name,
df=obs_df,
)

# 9. Execute test and save results in dict
causal_test_result = causal_test_case.execute_test(linear_regression_estimator, data_collector)

if verbose:
logging.info("Causation:\n%s", causal_test_result)

results_dict[outcome_variable.name]["ate"] = causal_test_result.test_value.value

results_dict[outcome_variable.name]["cis"] = causal_test_result.confidence_intervals

results_dict[outcome_variable.name]["test_passes"] = causal_test_case.expected_causal_effect.apply(
causal_test_result
)
return results_dict


class CovasimVaccineDataCollector(ExperimentalDataCollector):
"""A custom experimental data collector for the elderly vaccination Covasim case study.

This experimental data collector runs covasim with a normal Pfizer vaccine and then again with the same vaccine but
this time prioritising the elderly for vaccination.
"""
return results_dict

def run_system_with_input_configuration(self, input_configuration: dict) -> pd.DataFrame:
"""Run the system with a given input configuration.

:param input_configuration: A nested dictionary containing Covasim parameters, desired number of repeats, and
a bool to determine whether elderly should be prioritised for vaccination.
:return: A dataframe containing results for this input configuration.
"""
results_df = self.simulate_vaccine(
input_configuration["covasim_parameters_dict"], self.n_repeats, input_configuration["target_elderly"]
)
return results_df

def simulate_vaccine(self, pars_dict: dict, n_simulations: int = 100, target_elderly: bool = False):
"""Simulate observational data that contains a vaccine that is optionally given preferentially to the elderly.
:param pars_dict: A dictionary containing simulation parameters.
:param n_simulations: Number of simulations to run.
:param target_elderly: Whether to prioritise vaccination for the elderly.
:return: A pandas dataframe containing results for each run.
"""
simulations_results_dfs = []
for sim_n in range(n_simulations):
logging.info("Simulation %s/%s.", sim_n + 1, n_simulations)

# Update simulation parameters with vaccine and optionally sub-target
if target_elderly:
logger.info("Prioritising the elderly for vaccination")
vaccine = cv.vaccinate_prob(
vaccine="Pfizer",
label="prioritise_elderly",
subtarget=self.vaccinate_by_age,
days=list(range(7, pars_dict["n_days"])),
)
else:
logger.info("Using standard vaccination protocol")
vaccine = cv.vaccinate_prob(vaccine="Pfizer", label="regular", days=list(range(7, pars_dict["n_days"])))

pars_dict["interventions"] = vaccine
pars_dict["use_waning"] = True # Must be set to true for vaccination
sim_results_df = self.run_sim_with_pars(
pars_dict=pars_dict,
desired_outputs=[
"cum_infections",
"cum_deaths",
"cum_recoveries",
"cum_vaccinations",
"cum_vaccinated",
],
n_runs=1,
)

sim_results_df["interventions"] = vaccine.label # Store label in results instead of vaccine object
sim_results_df["target_elderly"] = target_elderly
sim_results_df["vaccine"] = int(target_elderly) # 0 if standard vaccine, 1 if target elderly vaccine
sim_results_df["max_doses"] = vaccine.p["doses"] # Get max doses for the vaccine
simulations_results_dfs.append(sim_results_df)

# Create a single dataframe containing a row for every execution
obs_df = pd.concat(simulations_results_dfs, ignore_index=True)
obs_df.rename(columns={"interventions": "vaccine_type"}, inplace=True)
return obs_df

@staticmethod
def run_sim_with_pars(pars_dict: dict, desired_outputs: [str], n_runs: int = 1, verbose: int = -1):
"""Runs a Covasim COVID-19 simulation with a given dict of parameters and collects the desired outputs,
which are given as a list of output names.
:param pars_dict: A dictionary containing the parameters and their values for the run.
:param desired_outputs: A list of outputs which should be collected.
:param n_runs: Number of times to run the simulation with a different seed.
:param verbose: Covasim verbose setting (0 for no output, 1 for output).
:return results_df: A pandas df containing the results for each run
"""
results_dict = {k: [] for k in list(pars_dict.keys()) + desired_outputs + ["rand_seed"]}
for _ in range(n_runs):
# For every run, generate and use a new a random seed.
# This is to avoid using Covasim's sequential random seeds.
random.seed()
rand_seed = random.randint(0, 10000)
pars_dict["rand_seed"] = rand_seed
logger.info("Rand Seed: %s", rand_seed)
sim = cv.Sim(pars=pars_dict)
m_sim = cv.MultiSim(sim)
m_sim.run(n_runs=1, verbose=False, n_cpus=1)

for run in m_sim.sims:
results = run.results
# Append inputs to results
for param in pars_dict.keys():
results_dict[param].append(run.pars[param])

# Append outputs to results
for output in desired_outputs:
if output not in results:
raise IndexError(f"{output} is not in the Covasim outputs. Are you using v3.0.7?")
results_dict[output].append(
results[output][-1]
) # Append the final recorded value for each variable

# Any parameters without results are assigned np.nan for each execution
for param, results in results_dict.items():
if not results:
results_dict[param] = [np.nan] * len(results_dict["rand_seed"])
return pd.DataFrame(results_dict)

@staticmethod
def vaccinate_by_age(simulation):
"""A custom method to prioritise vaccination of the elderly. This method is taken from Covasim Tutorial 5:
https://github.com/InstituteforDiseaseModeling/covasim/blob/7bdf2ddf743f8798fcada28a61a03135d106f2ee/
examples/t05_vaccine_subtargeting.py
:param simulation: A covasim simulation for which the elderly will be prioritised for vaccination.
:return output: A dictionary mapping individuals to vaccine probabilities.
"""
young = cv.true(simulation.people.age < 50)
middle = cv.true((simulation.people.age >= 50) * (simulation.people.age < 75))
old = cv.true(simulation.people.age > 75)
inds = simulation.people.uid
vals = np.ones(len(simulation.people))
vals[young] = 0.1
vals[middle] = 0.5
vals[old] = 0.9
output = dict(inds=inds, vals=vals)
return output
if __name__ == "__main__":

test_results = setup_test_case(verbose=True)

if __name__ == "__main__":
test_results = test_experimental_vaccinate_elderly(runs_per_test_per_config=30, verbose=True)
logging.info("%s", test_results)
logging.info("%s", test_results)
31 changes: 31 additions & 0 deletions examples/covasim_/vaccinating_elderly/simulated_data.csv
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,31 @@
pop_size,pop_type,pop_infected,n_days,vaccine_type,use_waning,rand_seed,cum_infections,cum_deaths,cum_recoveries,cum_vaccinations,cum_vaccinated,target_elderly,vaccine,max_doses
50000,hybrid,1000,50,pfizer,True,1169,6277.0,15.0,6175.0,629466.0,530715.0,True,1,2
50000,hybrid,1000,50,pfizer,True,8888,6381.0,18.0,6274.0,630796.0,532010.0,True,1,2
50000,hybrid,1000,50,pfizer,True,370,6738.0,15.0,6621.0,631705.0,532864.0,True,1,2
50000,hybrid,1000,50,pfizer,True,9981,6784.0,18.0,6682.0,634582.0,535795.0,True,1,2
50000,hybrid,1000,50,pfizer,True,6305,6757.0,20.0,6659.0,631292.0,532464.0,True,1,2
50000,hybrid,1000,50,pfizer,True,1993,5844.0,17.0,5755.0,633314.0,534478.0,True,1,2
50000,hybrid,1000,50,pfizer,True,1938,6465.0,19.0,6353.0,627724.0,528993.0,True,1,2
50000,hybrid,1000,50,pfizer,True,4797,7044.0,15.0,6919.0,631246.0,532433.0,True,1,2
50000,hybrid,1000,50,pfizer,True,2308,6878.0,6.0,6801.0,628865.0,530038.0,True,1,2
50000,hybrid,1000,50,pfizer,True,4420,6429.0,11.0,6348.0,633803.0,535030.0,True,1,2
50000,hybrid,1000,50,pfizer,True,2314,6566.0,15.0,6477.0,629288.0,530550.0,True,1,2
50000,hybrid,1000,50,pfizer,True,7813,6913.0,17.0,6818.0,629290.0,530512.0,True,1,2
50000,hybrid,1000,50,pfizer,True,1050,6963.0,14.0,6860.0,627981.0,529212.0,True,1,2
50000,hybrid,1000,50,pfizer,True,3215,6671.0,17.0,6577.0,628802.0,530038.0,True,1,2
50000,hybrid,1000,50,pfizer,True,2286,6597.0,13.0,6505.0,628986.0,530195.0,True,1,2
50000,hybrid,1000,50,pfizer,True,3080,6926.0,16.0,6834.0,633636.0,534904.0,True,1,2
50000,hybrid,1000,50,pfizer,True,7405,6438.0,15.0,6347.0,630353.0,531540.0,True,1,2
50000,hybrid,1000,50,pfizer,True,9668,6577.0,15.0,6485.0,631257.0,532409.0,True,1,2
50000,hybrid,1000,50,pfizer,True,8211,6197.0,13.0,6103.0,633827.0,535056.0,True,1,2
50000,hybrid,1000,50,pfizer,True,4686,6761.0,16.0,6653.0,630557.0,531737.0,True,1,2
50000,hybrid,1000,50,pfizer,True,3591,7328.0,24.0,7214.0,629949.0,531124.0,True,1,2
50000,hybrid,1000,50,pfizer,True,4834,6617.0,22.0,6512.0,632609.0,533705.0,True,1,2
50000,hybrid,1000,50,pfizer,True,6142,7017.0,17.0,6902.0,635965.0,537252.0,True,1,2
50000,hybrid,1000,50,pfizer,True,6877,6845.0,15.0,6753.0,635678.0,536925.0,True,1,2
50000,hybrid,1000,50,pfizer,True,1878,6480.0,20.0,6390.0,630807.0,531999.0,True,1,2
50000,hybrid,1000,50,pfizer,True,3761,6972.0,16.0,6890.0,631100.0,532329.0,True,1,2
50000,hybrid,1000,50,pfizer,True,1741,6581.0,20.0,6491.0,632835.0,534088.0,True,1,2
50000,hybrid,1000,50,pfizer,True,5592,6561.0,19.0,6461.0,636799.0,537959.0,True,1,2
50000,hybrid,1000,50,pfizer,True,7979,7075.0,17.0,6966.0,632902.0,534140.0,True,1,2
50000,hybrid,1000,50,pfizer,True,71,6291.0,13.0,6203.0,631694.0,532901.0,True,1,2

0 comments on commit 282a225

Please sign in to comment.