{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# T5 - Using interventions\n", "\n", "Interventions are one of the most critical parts of Covasim. This tutorial shows how to implement standard interventions, as well as how to define your own custom interventions." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "
\n", " \n", "Click [here](https://mybinder.org/v2/gh/institutefordiseasemodeling/covasim/HEAD?urlpath=lab%2Ftree%2Fdocs%2Ftutorials%2Ftut_interventions.ipynb) to open an interactive version of this notebook.\n", " \n", "
" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Change beta and clip edges\n", "\n", "The most basic intervention in Covasim is \"change beta\": this changes overall transmission rates on the specified day(s), no questions asked. For example:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import covasim as cv\n", "cv.options(jupyter=True, verbose=0)\n", "\n", "# Define baseline parameters and sim\n", "pars = dict(\n", " start_day = '2020-03-01',\n", " end_day = '2020-06-01',\n", " pop_type = 'hybrid',\n", ")\n", "orig_sim = cv.Sim(pars, label='Baseline')\n", "\n", "# Define sim with change_beta\n", "cb = cv.change_beta(days=['2020-04-15', '2020-05-01', '2020-05-15'], changes=[0.2, 1.5, 0.7])\n", "sim = cv.Sim(pars, interventions=cb, label='With beta changes')\n", "\n", "# Run and plot\n", "msim = cv.parallel(orig_sim, sim)\n", "msim.plot()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "A few things to note here:\n", "\n", "* By default, interventions are shown with vertical dashed lines. You can turn this off by passing `do_plot=False` to the intervention.\n", "* Look at the new infections plot (middle) for the clearest illustration of the impact the intervention has.\n", "* Note that like other \"parameters\", you can either pass interventions to the sim directly or as part of the pars dictionary; the examples below illustrate these options.\n", "* All interventions handle dates automatically, so days can be specified as strings, integer numbers of days from the start of the sim, or datetime objects. You can even mix and match (e.g. `days=[10, '2020-05-05']`), although this is not recommended for obvious reasons!\n", "\n", "If your simulation has different layers (discussed in Tutorial 1), you can also define the beta change to only apply to a single layer. For example, an intervention that closed 90% of schools (`'s'`) on September 1st would be written as:\n", "\n", "```python\n", "close_schools = cv.change_beta(days='2020-09-01', changes=0.1, layers='s')\n", "```\n", "\n", "The \"clip edges\" intervention works similarly to the \"change beta\" intervention, with two important distinctions:\n", "\n", "1. Whereas beta changes affect the per-contact transmission probability, clip edges change the number of contacts. For example, in the `close_schools` example above, each school-aged child has an unchanged number of contacts (about 20), but 10 times lower probability of infection per contact (about 0.1% per contact per day). If clip edges were used instead, then each child would have unchanged probability of infection per contact (about 1% per contact per day), but 10 times fewer contacts (about 2).\n", "2. While beta changes can be any value (greater or less than 1), as the name implies, clip edges can only *clip* edges, not create new ones. Edges that are clipped are saved and can later be restored, but you can never have a change of greater than 1.0 with clip edges." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "
\n", "\n", "Overall, it's best to use `change_beta` for interventions that reduce per-contact risk (e.g., masks), and `clip_edges` for interventions that reduce the number of contacts per person. However, in practice, it tends to make relatively little difference to the simulation results. If in doubt, use `change_beta`.\n", "\n", "
" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can use `change_beta` to do some gnarly things, like create sinusoidally varying infection probability:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import sciris as sc\n", "import covasim as cv\n", "\n", "pars = sc.objdict(\n", " beta = 0.008,\n", " n_agents = 50e3,\n", " n_days = 180,\n", " verbose = 0,\n", ")\n", "\n", "beta_days = np.arange(pars.n_days)\n", "beta_vals = np.cos(2*np.pi*beta_days/90)**2+0.5\n", "beta = cv.change_beta(beta_days, beta_vals, do_plot=False)\n", "\n", "s1 = cv.Sim(pars, label='Normal')\n", "s2 = cv.Sim(pars, interventions=beta, label='Waves')\n", "\n", "msim = cv.parallel(s1, s2)\n", "msim.plot(['r_eff'])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Testing interventions\n", "\n", "There are two types of testing interventions in Covasim:\n", "\n", "1. `test_num` performs the specified *number* of tests per day;\n", "2. `test_prob` performs testing with the specified *probability*.\n", "\n", "Typically, if you have actual data on number of tests, you would use `test_num`; otherwise, you would usually use `test_prob`. For `test_num`, other than the number of tests performed, the main parameter is `symp_prob`, which determines how much more likely someone with (true positive) COVID symptoms is to test than a person who is asymptomatic or uninfected. Higher values of `symp_prob` give more diagnoses per number of tests, since symptomatic people preferentially test. The default value is 100, meaning that on any given day, a person experiencing COVID symptoms is 100 times more likely to test than someone who isn't.\n", "\n", "Likewise, `test_prob` also lets you specify the (daily) probability of symptomatic vs. asymptomatic testing. Unlike `test_num`, however, the total number of tests is determined by these probabilities. For example, since symptomatic testing probabilities are usually much higher than asymptomatic, if the epidemic is growing, more people will be symptomatic and thus go to be tested." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "
\n", "\n", "Note that test probabilities for `test_prob` are daily probabilities, not total. If you set `symp_prob = 0.3`, this means a person has a 30% chance of testing *per day of symptoms*, meaning that on average they will go to be tested after 2-3 days. Given that symptoms last for up to 8 days, the total probability that they will test at some point during their illness is much higher than the daily probability.\n", "\n", "
" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The example below illustrates (a) `test_num` with data, (b) `test_num` with a fixed number of tests per day, and (c) `test_prob`." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import covasim as cv\n", "\n", "# Define the testing interventions\n", "tn_data = cv.test_num('data')\n", "tn_fixed = cv.test_num(daily_tests=500, start_day='2020-03-01')\n", "tp = cv.test_prob(symp_prob=0.2, asymp_prob=0.001, start_day='2020-03-01')\n", "\n", "# Define the default parameters\n", "pars = dict(\n", " pop_size = 50e3,\n", " pop_infected = 100,\n", " start_day = '2020-02-01',\n", " end_day = '2020-03-30',\n", ")\n", "\n", "# Define the simulations\n", "sim1 = cv.Sim(pars, datafile='example_data.csv', interventions=tn_data, label='Number of tests from data')\n", "sim2 = cv.Sim(pars, interventions=tn_fixed, label='Constant daily number of tests')\n", "sim3 = cv.Sim(pars, interventions=tp, label='Probability-based testing')\n", "\n", "# Run and plot results\n", "msim = cv.parallel(sim1, sim2, sim3)\n", "msim.plot(['new_infections', 'new_tests', 'new_diagnoses', 'cum_diagnoses'])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Contact tracing\n", "\n", "Contact tracing is the process of notifying the contacts of a diagnosed positive COVID case and letting them know that they have been exposed, so they can be tested, go into quarantine, or both. Since contact tracing uses diagnosed people to start with, and since people are only diagnosed if they are tested, you **must** include at least one testing intervention in your sim before including contact tracing. The main option for contact tracing is `trace_prob`, the probability that a person in each layer will be traced (usually high for households and lower for other layers). For example:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import covasim as cv\n", "\n", "# Define the testing and contact tracing interventions\n", "tp = cv.test_prob(symp_prob=0.2, asymp_prob=0.001, symp_quar_prob=1.0, asymp_quar_prob=1.0, do_plot=False)\n", "ct = cv.contact_tracing(trace_probs=dict(h=1.0, s=0.5, w=0.5, c=0.3), do_plot=False)\n", "\n", "# Define the default parameters\n", "pars = dict(\n", " pop_type = 'hybrid',\n", " pop_size = 50e3,\n", " pop_infected = 100,\n", " start_day = '2020-02-01',\n", " end_day = '2020-05-01',\n", " interventions = [tp, ct],\n", ")\n", "\n", "# Create and run the simulation\n", "sim = cv.Sim(pars)\n", "sim.run()\n", "\n", "# Plot the simulation\n", "cv.options.set(fontsize=8) # Make the font size a bit smaller here so the labels don't overlap\n", "sim.plot(['new_infections', 'new_tests', 'new_diagnoses', 'cum_diagnoses', 'new_quarantined', 'test_yield'])\n", "cv.options.set(fontsize='default') # Reset to the default value" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Since it's assumed that known contacts are placed into quarantine (with efficacy `sim['quar_factor']`), the number of contacts who are successfully traced each day is equal to the number of people newly quarantined (bottom left panel). As is commonly seen using testing and tracing as the only means of epidemic control, these programs stop the epidemic from growing exponentially, but do not bring it to zero.\n", "\n", "Since these interventions happen at `t=0`, it's not very useful to plot them. Note that we have turned off plotting by passing `do_plot=False` to each intervention." ] }, { "cell_type": "markdown", "metadata": { "tags": [] }, "source": [ "## Simple vaccination\n", "\n", "Vaccines can do one of two things: they can stop you from becoming infected in the first place (acquisition blocking), or they can stop you from developing symptoms or severe disease once infected (symptom blocking). The Covasim `simple_vaccine` intervention lets you control both of these options. In its simplest form, a vaccine is like a change beta intervention. For example, this vaccine:\n", "\n", "```python\n", "vaccine = cv.simple_vaccine(days=30, prob=1.0, rel_sus=0.3, rel_symp=1.0)\n", "```\n", "\n", "is equivalent to this beta change:\n", "\n", "```python\n", "cb = cv.change_beta(days=30, changes=0.3)\n", "```\n", "\n", "But that's not very realistic. A vaccine given on days 30 and 44 (two weeks later), with efficacy of 50% per dose which accumulates, given to 60% of the population, and which blocks 50% of acquisition and (among those who get infected even so) 90% of symptoms, would look like this:\n", "\n", "```python\n", "vaccine = cv.simple_vaccine(days=[30, 44], cumulative=[0.5, 0.5], prob=0.6, rel_sus=0.5, rel_symp=0.1)\n", "```" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "
\n", "\n", "For more advanced vaccine usage, see *T8 - Vaccines and variants*.\n", "\n", "
" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Age targeting and other subtargeting\n", "\n", "For some interventions, you want them to preferentially apply to some people in the population, e.g. the elderly. For testing and vaccine interventions, this can be done using the `subtarget` argument. Subtargeting is quite flexible in terms of what it accepts, but the easiest approach is usually to define a function to do it. This function must return a dictionary with keys `inds` (which people to subtarget) and `vals` (the values they are subtargeted with). This example shows how to define a 20% transmission-blocking, 94% symptom-blocking vaccine with varying coverage levels by age (for the vaccine intervention, `subtarget` modifies `prob`, the probability of a person being vaccinated) applied on day 20:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import covasim as cv\n", "import numpy as np\n", "\n", "# Define the vaccine subtargeting\n", "def vaccinate_by_age(sim):\n", " young = cv.true(sim.people.age < 50) # cv.true() returns indices of people matching this condition, i.e. people under 50\n", " middle = cv.true((sim.people.age >= 50) * (sim.people.age < 75)) # Multiplication means \"and\" here\n", " old = cv.true(sim.people.age >= 75)\n", " inds = sim.people.uid # Everyone in the population -- equivalent to np.arange(len(sim.people))\n", " vals = np.ones(len(sim.people)) # Create the array\n", " vals[young] = 0.1 # 10% probability for people <50\n", " vals[middle] = 0.5 # 50% probability for people 50-75\n", " vals[old] = 0.9 # 90% probability for people >75\n", " output = dict(inds=inds, vals=vals)\n", " return output\n", "\n", "# Define the vaccine\n", "vaccine = cv.simple_vaccine(days=20, rel_sus=0.8, rel_symp=0.06, subtarget=vaccinate_by_age)\n", "\n", "# Create, run, and plot the simulations\n", "sim1 = cv.Sim(label='Baseline')\n", "sim2 = cv.Sim(interventions=vaccine, label='With age-targeted vaccine')\n", "msim = cv.parallel(sim1, sim2)\n", "msim.plot()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "If you have a simple conditional, you can also define subtargeting using a lambda function, e.g. this is a vaccine with 90% probability of being given to people over age 75, and 10% probability of being applied to everyone else (i.e. people under age 75%):\n", "```python\n", "vaccine = cv.simple_vaccine(days=20, prob=0.1, subtarget=dict(inds=lambda sim: cv.true(sim.people.age>50), vals=0.9))\n", "```" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "However, lambda functions are not known for their easy readability, so their use is discouraged." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Dynamic parameters\n", "\n", "A slightly more advanced intervention is the dynamic parameters intervention. This is not an intervention the same way the others are; instead, it provides a consistent way to modify parameters dynamically. For example, to change the number of imported infections on day 15 from 0 to 100 (representing, say, resuming international travel), and then change back to 0 on day 30:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import covasim as cv\n", "\n", "# Define the dynamic parameters\n", "imports = cv.dynamic_pars(n_imports=dict(days=[15, 30], vals=[100, 0]))\n", "\n", "# Create, run, and plot the simulations\n", "sim1 = cv.Sim(label='Baseline')\n", "sim2 = cv.Sim(interventions=imports, label='With imported infections')\n", "msim = cv.parallel(sim1, sim2)\n", "msim.plot()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can see the sudden jump in new infections when importations are turned on." ] }, { "cell_type": "markdown", "metadata": { "tags": [] }, "source": [ "## Dynamic triggering\n", "\n", "Another option is to replace the ``days`` arguments with custom functions defining additional criteria. For example, perhaps you only want your beta change intervention to take effect once infections get to a sufficiently high level. Here's a fairly complex example (feel free to skip the details) that toggles the intervention on and off depending on the current number of people who are infectious.\n", "\n", "This example illustrates a few different features:\n", "\n", "* The simplest change is just that we're supplying `days=inf_thresh` instead of a number or list. If we had `def inf_thresh(interv, sim): return [20,30]` this would be the same as just setting `days=[20,30]`.\n", "* Because the first argument this function gets is the intervention itself, we can stretch the rules a little bit and name this variable `self` -- as if we're defining a new method for the intervention, even though it's actually just a function.\n", "* We want to keep track of a few things with this intervention -- namely when it toggles on and off, and whether or not it's active. Since the intervention is just an object, we can add these attributes directly to it.\n", "* Finally, this example shows some of the flexibility in how interventions are plotted -- i.e. shown in the legend with a label and with a custom color." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "def inf_thresh(self, sim, thresh=500):\n", " ''' Dynamically define on and off days for a beta change -- use self like it's a method '''\n", "\n", " # Meets threshold, activate\n", " if sim.people.infectious.sum() > thresh:\n", " if not self.active:\n", " self.active = True\n", " self.t_on = sim.t\n", " self.plot_days.append(self.t_on)\n", "\n", " # Does not meet threshold, deactivate\n", " else:\n", " if self.active:\n", " self.active = False\n", " self.t_off = sim.t\n", " self.plot_days.append(self.t_off)\n", "\n", " return [self.t_on, self.t_off]\n", "\n", "# Set up the intervention\n", "on = 0.2 # Beta less than 1 -- intervention is on\n", "off = 1.0 # Beta is 1, i.e. normal -- intervention is off\n", "changes = [on, off]\n", "plot_args = dict(label='Dynamic beta', show_label=True, line_args={'c':'blue'})\n", "cb = cv.change_beta(days=inf_thresh, changes=changes, **plot_args)\n", "\n", "# Set custom properties\n", "cb.t_on = np.nan\n", "cb.t_off = np.nan\n", "cb.active = False\n", "cb.plot_days = []\n", "\n", "# Run the simulation and plot\n", "sim = cv.Sim(interventions=cb)\n", "sim.run().plot()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Custom interventions\n", "\n", "If you're still standing after the previous example, Covasim also lets you do things that are even *more* complicated, namely define an arbitrary function or class to act as an intervention instead. If a custom intervention is supplied as a function (as it was in Tutorial 1 as well), then it receives the `sim` object as its only argument, and is called on each timestep. It can perform arbitrary manipulations to the sim object, such as changing parameters, modifying state, etc.\n", "\n", "This example reimplements the dynamic parameters example above, except using a hard-coded custom intervention:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "def dynamic_imports(sim):\n", " if sim.t == 15:\n", " sim['n_imports'] = 100\n", " elif sim.t == 30:\n", " sim['n_imports'] = 0\n", " return\n", "\n", "sim = cv.Sim(interventions=dynamic_imports)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "However, function-based interventions only take you so far.\n", "\n", "We saw in Tutorial 1 how you could define a simple \"protect the elderly\" intervention with just a few lines of code:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "def protect_elderly(sim):\n", " if sim.t == sim.day('2020-04-01'):\n", " elderly = sim.people.age>70\n", " sim.people.rel_sus[elderly] = 0.0" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This example explains how to create an intervention object that does much the same thing, but is more fully-featured because it uses the `Intervention` class." ] }, { "cell_type": "markdown", "metadata": { "tags": [] }, "source": [ "
\n", "\n", "You must include the line `super().__init__(**kwargs)` in the `self.__init__()` method, or else the intervention won't work. You must also include `super().initialize()` in the `self.initialize()` method.\n", "\n", "
" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import pylab as pl\n", "import covasim as cv\n", "\n", "class protect_elderly(cv.Intervention):\n", "\n", " def __init__(self, start_day=None, end_day=None, age_cutoff=70, rel_sus=0.0, *args, **kwargs):\n", " super().__init__(**kwargs) # NB: This line must be included\n", " self.start_day = start_day\n", " self.end_day = end_day\n", " self.age_cutoff = age_cutoff\n", " self.rel_sus = rel_sus\n", " return\n", "\n", " def initialize(self, sim):\n", " super().initialize() # NB: This line must also be included\n", " self.start_day = sim.day(self.start_day) # Convert string or dateobject dates into an integer number of days\n", " self.end_day = sim.day(self.end_day)\n", " self.days = [self.start_day, self.end_day]\n", " self.elderly = sim.people.age > self.age_cutoff # Find the elderly people here\n", " self.exposed = np.zeros(sim.npts) # Initialize results\n", " self.tvec = sim.tvec # Copy the time vector into this intervention\n", " return\n", "\n", " def apply(self, sim):\n", " self.exposed[sim.t] = sim.people.exposed[self.elderly].sum()\n", "\n", " # Start the intervention\n", " if sim.t == self.start_day:\n", " sim.people.rel_sus[self.elderly] = self.rel_sus\n", "\n", " # End the intervention\n", " elif sim.t == self.end_day:\n", " sim.people.rel_sus[self.elderly] = 1.0\n", "\n", " return\n", "\n", " def plot(self):\n", " fig = pl.figure()\n", " pl.plot(self.tvec, self.exposed)\n", " pl.xlabel('Day')\n", " pl.ylabel('Number infected')\n", " pl.title('Number of elderly people with active COVID')\n", " return fig" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "While this example is fairly long, hopefully it's fairly straightforward:\n", "\n", "- `__init__()` just does what init always does; it's needed to create the class instance. For interventions, it's usually used to store keyword arguments and perform some basic initialization (the first two lines).\n", "- `initialize()` is similar to init, with the difference that it's invoked when the *sim* itself is initialized. It receives the sim as an input argument. This means you can use it to do a fairly powerful things. Here, since `sim.people` already exists, we can calculate up-front who the elderly are so we don't have to do it on every timestep.\n", "- `apply()` is the crux of the intervention. It's run on every timestep of the model, and also receives `sim` as an input. You almost always use `sim.t` to get the current timestep, here to turn the intervention on and off. But as this example shows, its real power is that it can make direct modifications to the sim itself (`sim.people.rel_sus[self.elderly] = self.rel_sus`). It can also perform calculations and store data in itself, as shown here with `self.exposed` (although in general, analyzers are better for this, since they happen at the end of the timestep, while interventions happen in the middle).\n", "- `plot()` is a custom method that shows a plot of the data gathered during the sim. Again, it's usually better to use analyzers for this, but for something simple like this it's fine to double-dip and use an intervention.\n", "\n", "Here is what this custom intervention looks like in action. Note how it automatically shows when the intervention starts and stops (with vertical dashed lines)." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Define and run the baseline simulation\n", "pars = dict(\n", " pop_size = 50e3,\n", " pop_infected = 100,\n", " n_days = 90,\n", " verbose = 0,\n", ")\n", "orig_sim = cv.Sim(pars, label='Default')\n", "\n", "# Define the intervention and the scenario sim\n", "protect = protect_elderly(start_day='2020-04-01', end_day='2020-05-01', rel_sus=0.1) # Create intervention\n", "sim = cv.Sim(pars, interventions=protect, label='Protect the elderly')\n", "\n", "# Run and plot\n", "msim = cv.parallel(orig_sim, sim)\n", "msim.plot()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Because we collected information in the intervention itself, we can also plot that:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Plot intervention\n", "protect = msim.sims[1].get_intervention(protect_elderly) # Find intervention by type\n", "protect.plot()" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.11.7" }, "pycharm": { "stem_cell": { "cell_type": "raw", "metadata": { "collapsed": false }, "source": [] } }, "toc": { "base_numbering": 1, "nav_menu": {}, "number_sections": true, "sideBar": true, "skip_h1_title": false, "title_cell": "Table of Contents", "title_sidebar": "Contents", "toc_cell": false, "toc_position": {}, "toc_section_display": true, "toc_window_display": false } }, "nbformat": 4, "nbformat_minor": 4 }