Skip to content

Examples

Jump to bottom
Robert Manson-Sawko edited this page Jan 7, 2022 · 3 revisions

In the preprint associated with this repository, we study four policies associated with non-pharmeceutical interventions (NPIs): combined controls on within- and between-household mixing, "support bubble" exemptions to between-household measures, short-term relaxations of restrictions on between-household mixing, and out-of-household isolation (OOHI). We also simulate an unconstrained epidemic to demonstrate that the initial conditions used in our model is suitably close to the eigenvector associated with the early exponential growth phase of the dynamics. Here we summarise the scripts found in the examples folder which can be used to replicate our results.

examples/uk

The scripts in this folder are used to simulate a single epidemic in a UK-like population of households with parameters callibrated to wild-type COVID-19. The script 'run.py' runs this simulation, and is illustrative of the steps needed to implement our model. The basic steps carried out are:

  1. Load in the scripts, classes, and objects needed to run the simulation, and create a folder to store the outputs if one does not already exist (lines 3-18);
  2. Define the doubling time we will callibrate the between-household mixing to, and from this derive a growth rate (lines 19-21);
  3. Read the set of household compositions used in our model population and the associated distribution from the 'inputs' folder (lines 22-31);
  4. Create a set of model inputs callibrated to our chosen growth rate. Start by defining the specifications of the model containing the basic numerical parameters we will use (line 36). Create a model input object which converts these specifications into input for the households model (line 37) - at this point the intensity of between household mixing defined by this model input is not yet callibrated to our chosen growth rate. Next create a household population object containing a list of all the epidemiological states a household in our system can take, and a Kolmogorov matrix defining the rates of transitions between these states in the within-household epidemiological dynamics (lines 38-39). Next create a rate equations object which can be used as the right hand side of the self-consistent equations which simulate the population-level infectious disease dynamics, including both within- and between-household transmission events (lines 40-41). Use these objects to estimate a between-household transmission rate which will give us the chosen growth rate (lines 42-43). Create a new model input object with this growth rate, which will be used later on to create a new household population object and system of equations callibrated to our chosen doubling time (lines 44-45). Once this callibrated model input is calculated, we save it in the outputs folder (lines 46-47); if the model input has already been calculated then instead of performing these calculations we load it from the outputs folder (lines 33-34);
  5. Calculate a household population object based on the callibrated model input and save to the outputs folder. Again, if such an object has been calculated before, we just load it from the outputs folder (lines 49-58);
  6. Create the system of equations for our simulations (line 59);
  7. Estimate the growth rate of the infectious disease dynamics defined by these equations. Printing these to the console demonstrates that we have successfully callibrated the model (lines 61-64);
  8. Define the initial conditions H(0) for the simulations based on a chosen starting prevalence of 0.00001 and a starting population immunity of zero. We also estimate the proportion of the population of individuals in each compartment given these initial conditions, (S,E,P,I,R), which should be equal to (0.99999, 0.00001, 0, 0, 0). This is not used anywhere subsequently in the code, but if printed to the console it serves as a useful confirmation that we have indeed chosen the correct initial conditions (lines 66-77);
  9. Define the timespan for the simulations and solve using solve_ivp (lines 78-79);
  10. Calculate the simulation outputs in terms of mean number of individuals per household in each epidemiological compartment, and save to file (lines 81-99). The outputs from examples/uk/run.py are used in two plotting scripts. The first, plot_time_series.py generates plots of the proportion of children and adults in each epidemiological compartment over the course of the simulation (not included in the accompanying preprint). The second, plot_exp_growth.py generates a plot of prevalence over time plotted on a logarithmic scale on the same axis as an exponential growth curve with the same exponential growth rate which we calibrated the model input to. The closeness of these two curves demonstrates that the initial conditions we choose for our model place us near to the eigenvector associated with the early exponential growth phase of the unconstrained epidemic.

examples/mixing_sweep

In this example we model the impact of reductions in within- and between-household transmission rates on the overall infectious disease dynamics. As in the subsequent three examples, a range of scenarios is considered and parallelisation is used to simulate multiple scenarios at once.

The simulations are carried out in examples/mixing_sweep/parameter_sweep.py. The parallelised script takes as input the number of workers to be used in the parallelisation, the range of reductions in within-household transmission rate, and the range of reductions in between-household transmission rate. These parameter ranges are defined using the same syntax as the arange function from numpy (minimum-maximum-increment). By default, the number of workers is set to 8 and both parameter ranges are set to (0, 0.99, 0.05). The script generates the following outputs:

  • a 2D array of estimated growth rates for each pair of parameter values;
  • a 2D array of peak prevalences for each pair of parameter values;
  • a 2D array of cumulative prevalences for each pair of parameter values;
  • a 2D array of the proportion of households which see at least one infection for each pair of parameter values;
  • a 2D array of attack ratios in households with at least one infection for each pair of parameter values;
  • a 3D array of attack ratios in households with at least one infection for each pair of parameter values and each household size;
  • a 3D array of "first pass" attack ratios in households with at least one infection for each pair of parameter values and each household size;
  • an array of within-household transmission rate reduction levels;
  • an array of between-household transmission rate reduction levels.

The output from examples/mixing_sweep/parameter_sweep.py is used by examples/mixing_sweep/plot_sweep_results.py to generate two figures: a grid of growth rate, peak prevalence, cumulative prevalence, and proportion of households infected as a function of reduction in internal and external mixing rates, and plots of overall and "first pass" attack rate in infected households for two specific levels of within-household mixing reduction level.

examples/long_term_bubbles

In this example we model "support bubble" exemptions to NPIs which reduce between-household transmission rates. We calculate a few epidemiological statistics over a range of reductions in between-household mixing and a range of uptake levels for support bubbles among those households which are eligible. Eligible households are those which contain no more than one adult (i.e. any single-parent household or any adult living alone), and when they join a support bubble they effectively form a single household with a randomly chosen household from the general population (which may itself also meet the eligibility criteria).

The simulations are carried out in examples/long_term_bubbles/parameter_sweep.py. The parallelised script takes as input the number of workers to be used in the parallelisation, the range of uptake rates of support bubbles among eligible households, and the range of reductions in between-household transmission rate. These parameter ranges are defined using the same syntax as the arange function from numpy (minimum-maximum-increment). By default, the number of workers is set to 8, the range of uptake rates is set to (0, 1.01, 0.05), and the range of between-household mixing reductions is set to (0, 0.99, 0.05). The script generates the following outputs:

  • a 2D array of estimated growth rates for each pair of parameter values;
  • a 2D array of peak prevalences for each pair of parameter values;
  • a 2D array of cumulative prevalences for each pair of parameter values;
  • a 2D array of the proportion of households which see at least one infection for each pair of parameter values;
  • an array of support bubble uptake rates;
  • an array of between-household transmission rate reduction levels.

The output from examples/long_term_bubbles/parameter_sweep.py is used by examples/long_term_bubbles/plot_sweep_results.py to generate a grid of growth rate, peak prevalence, cumulative prevalence, and proportion of households infected as a function of reduction in external mixing rates and uptake of support bubbles among eligible households.

examples/temp_bubbles

In this example we model short-term relaxations of NPIs which allow for restricted mixing between households over specific time periods. These relaxations differ from the exemptions studied in the long_term_bubbles example in that they are dynamic and temporary; in our simulations, we start by simulating the dynamics of an "unbubbled" population of single households, and at specific times we map the dynamics into a population of bubbles of two or three households, solve the dynamics of this population over the relaxation period, then map these dynamics back into the original unbubbled population to simulate the aftermath of the relaxation. These bubbles are interpreted as consisting of a host household and one or more guest households, interacting within the host household's home. To account for uncertainty in how interactions with visitors may differ from interactions with one's own household members, we perform all of our calculations over a range of density dependence parameters for bubbled and unbubbled household transmission. This density parameter ranges from 0 (purely density-dependent mixing) to 1 (purely frequency-dependent mixing), and by expressing transmission rates in the form beta * (S * I) / (N^(alpha)) we can tune between these extremes by varying alpha. We simulate a few different relaxation policies allowing for two or three households to form bubbles over periods of different duration and timing, which are detailed in our preprint. For each policy we calculate a few epidemiological statistics over a range of bubbled and unbubbled density parameters.

The file examples/temp_bubbles/common.py contains functions specific to this example which generate the bubbled populations from the baseline UK household composition data and define mappings between the state space of the baseline and bubbled populations.

The simulations are carried out in examples/temp_bubbles/parameter_sweep.py. The parallelised script takes as input the number of workers to be used in the parallelisation, the range of density parameters for unbubbled households, and the range of density parameters for bubbled households. These parameter ranges are defined using the same syntax as the arange function from numpy (minimum-maximum-increment). By default, the number of workers is set to 8 and both parameter ranges are set to (0, 1.01, 0.05). The script generates the following outputs:

  • a 2D array of peak prevalences for each pair of parameter values;
  • a 2D array of cumulative prevalences for each pair of parameter values;
  • a 2D array of attack ratios in households with at least one infection for each pair of parameter values;
  • a 2D array of the proportion of households which see at least one infection for each pair of parameter values;
  • an array of unbubbled density parameters;
  • an array of bubbled density parameters. The output statistics are outputted for each policy.

The output from examples/temp_bubbles/parameter_sweep.py is used by 'examples/temp_bubbles/plot_sweep_results.py' to generate two figures, visualising the peak and cumulative prevalences as a function of the two density rates associated with each relaxation policy.

examples/external_isolation

In this example we model a population composed of children, less clinically vulnerable adults, and clinically vulnerable adults to analyse an out-of-household isolation (OOHI) policy. Under the OOHI policy, whenever a case in the less clinically vulnerable adult class is detected in a household containing at least one clinically vulnerable adult, that case is given the opportunity to isolate and with a given probability is removed from the household until they have recovered. There are two main dimensions of uncertainty here: the rate at which cases can be detected, and the probability that a detected case enters external isolation. We calculate a few epidemiological statistics over a range of detection rates and isolation probabilities.

The simulations are carried out in examples/external_isolation/parameter_sweep.py. The parallelised script takes as input the number of workers to be used in the parallelisation, the range of detection rates, and the range of isolation probabilities for detected cases. These parameter ranges are defined using the same syntax as the arange function from numpy (minimum-maximum-increment). By default, the number of workers is set to 8, the range of detection rates is set to (0.01, 1.01, 0.05) (corresponding to expected detection time ranging from 10 days down to 1 day following infection), and the range of isolation probabilities is set to (0.0, 1.01, 0.05). The script generates the following outputs:

  • a 2D array of peak prevalences in the vulnerable class for each pair of parameter values;
  • a 2D array of cumulative prevalences in the vulnerable class for each pair of parameter values;
  • a 2D array of peak proportion of the less vulnerable adult class in isolation for each pair of parameter values;
  • a 2D array of cumulative proportion of the less vulnerable adult class in isolation for each pair of parameter values;
  • an array of detection rates;
  • an array of isolation probabilities.

The output from examples/external_isolation/parameter_sweep.py is used by examples/external_isolation/plot_sweep_results.py to generate plots of each output array.

Clone this wiki locally