Vytautas Jancauskas 2020
In this tutorial we will show you how you can use EasyVVUQ to investigate the properties of a simple epidemiological model. The model is very simplistic and is not intended to realisticly portray a real epidemic such as COVID-19. However it is also fun enough to experiment with and we can use it as an example to show how you can use EasyVVUQ to answer questions about your scientific models. EasyVVUQ is successfully used by researchers in various different fields, such as plasma physics, weather and air pollution modelling and materials science.
Before we do anything else we need to install EasyVVUQ for use in this notebook. Please skip this if you already have EasyVVUQ installed in your system. This is meant for situations where this notebook is hosted externally to your local computer. We also want to clone the git repository because it contains some files used in this tutorial. Please note that you might need to restart the runtime after the installation is complete. If the following code examples don't work please try that first. I have added a command that kills the runtime (which causes it to be restarted). But I'm not sure if it will always work.
#!pip install git+https://github.com/UCL-CCS/EasyVVUQ
#!git clone https://github.com/UCL-CCS/EasyVVUQ
#import os
#os.kill(os.getpid(), 9)
In our model individuals are placed on a two dimensional square grid. During each turn, the individual moves at random one square in one of eight possible directions (N, NE, E, SE, S, SW, W or NW). If the individual encounters another individual and that person is ill, the individual will also become sick for a specified number of turns (that number is a parameter to the model). Once ill, the individual cannot become sick anymore and after getting over the disease they will have immunity. Immunity means that this person cannot get infected anymore. The simulation continues until no one is sick. At that point the disease counts as erradicated.
Let us load the model and see how it operates. We create a population on a 10x10 grid containing 20 individuals. Once infected the individual can transmit the disease for 28 turns. After that they have a 20% chance of dying.
import matplotlib.pyplot as plt
import numpy as np
import epidemic.epidemic as epidemic
population = epidemic.Population(grid_size=10, n=20, duration=28, mortality=0.2)
The code cell below is supposed to be run multiple times. Each time you run it the image below will update to show the state of the population. The black squares are empty, the white squares represent individuals who are not immune and are not sick. The red squares represent ill individuals and the intensity of the red shows how many turns they will stay ill for. Green squares represent individuals who are immune. Once the disease is eradicated, the graphic will change to a plot showing the evolution of the disease over time.
population.move()
if np.count_nonzero(population.ill):
plt.imshow(population.get_im())
plt.show()
else:
plt.plot(population.ill_history, label="Sick Individuals")
plt.plot(population.immune_history, label="Immune Individuals")
plt.plot(population.n_history, label="Population")
plt.xlabel("Time")
plt.ylabel("Count")
plt.legend()
plt.show()
Since EasyVVUQ is meant to be a general framework (non-Python specific) we don't call Python functions directly to get results of the simulation. After all, many simulations are still written in Fortran and operate by taking an input file and producing a file with outputs of the simulation. To do statistical analysis we need to be able to provide an appropriate input file to the simulation and be able to parse the outputs of the simulation. You can run our simulation as in the following example.
!cd epidemic/; time python3 epidemic.py example.json output.csv
1.887u 0.711s 0:01.25 207.2% 0+0k 0+0io 525pf+0w
This will have produced a file called output.csv that consists of four columns labeled "iteration", "ill", "immune" and "population". These should be fairly self explanatory.
!cat epidemic/output.csv
iteration,ill,immune,population 0,2,1,20 1,2,1,20 2,2,1,20 3,2,1,20 4,2,1,20 5,2,1,20 6,3,2,20 7,3,2,20 8,3,2,20 9,3,2,20 10,3,2,20 11,3,2,20 12,3,2,20 13,4,3,20 14,5,4,20 15,6,5,20 16,6,5,20 17,7,6,20 18,7,6,20 19,7,6,20 20,8,7,20 21,9,8,20 22,10,9,20 23,11,10,20 24,13,12,20 25,14,13,20 26,15,14,20 27,14,15,20 28,14,15,20 29,14,15,20 30,16,17,20 31,16,17,20 32,16,17,20 33,16,18,20 34,16,18,20 35,16,18,20 36,16,18,20 37,16,18,20 38,17,19,20 39,17,19,20 40,16,19,20 41,16,20,20 42,15,20,20 43,15,20,20 44,14,20,20 45,14,20,20 46,14,20,20 47,13,20,20 48,12,20,20 49,11,20,20 50,10,20,20 51,8,20,20 52,7,20,20 53,6,20,20 54,5,19,19 55,5,19,19 56,5,19,19 57,3,19,19 58,3,19,19 59,3,19,19 60,2,19,19 61,2,19,19 62,2,19,19 63,2,19,19 64,2,19,19 65,1,18,18 66,1,18,18 67,1,18,18 68,0,18,18
Let us plot this data and see what the evolution of the disease looks like in our population.
import matplotlib.pyplot as plt
import pandas as pd
df = pd.read_csv("epidemic/output.csv")
plt.plot(df['ill'], label="Sick Individuals")
plt.plot(df['immune'], label="Immune Individuals")
plt.plot(df['population'], label="Population")
plt.xlabel("Time")
plt.ylabel("Count")
plt.legend()
plt.show()
A short summary is in order so that we can start exploring the sample space of our model:
We will use EasyVVUQ to help us answer some questions about the model. Here are some simple ones that arise from toying with the model:
But first we need to set EasyVVUQ up to produce the configuration files in the suitable format and read in the output of the simulation. We also need to give a description of the parameter space. We also need to specify how we will execute our simulation. The next sections is concerened with these tasks.
For the examples in this tutorial we import some libraries that will be used throughout. EasyVVUQ will be referred to as 'uq' in the code. We also need Chaospy because we use it for the probability distribution classes. We use numpy for certain small tasks and we use pandas DataFrame as the standard data exchange format as is customary in the Python data science infrastucture.
import easyvvuq as uq
import chaospy as cp
import easyvvuq.collate
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import os
from pathlib import Path
from easyvvuq.actions import CreateRunDirectory, Encode, Decode, CleanUp, ExecuteLocal, Actions
While we are at it we also want to describe the arguments to our model. This takes a form of a Python dictionary. The dictionary is based on the Cerberus validator dictionary format. For more information refer to Cerberus documentation. This dictionary is used in both validation of the results and for the default values when we don't want to vary a certain parameter.
params = {
"grid_size": {"type": "float", "default": 10},
"n": {"type": "float", "min": 1, "max": 100, "default": 20},
"duration": {"type": "float", "min": 0, "max": 365, "default": 28},
"mortality": {"type": "float", "min": 0.0, "max": 1.0, "default": 0.2},
"ensemble" : {"type": "integer", "default": 0},
"ensemble_id" : {"type": "integer", "default": 0}
}
We will also want to set-up some elements that will stay the same for all the examples. These components are the encoder - which is responsible for creating input files for our simulation and the decoder - which is responsible for parsing the output of the simulation.
For the Encoder we use the GenericEncoder class. It is a very simple template based encoder. It finds a specified delimiter, and replaces the variable name that follows that delimiter with the corresponding value. In our case the template file looks like follows:
{
"grid_size" : $grid_size,
"n" : $n,
"duration" : $duration,
"mortality" : $mortality
}
From this template, a JSON file will be created and then passed to the simulation script as an argument. EasyVVUQ has other encoders as well. For example the Jinja encoder.
encoder = uq.encoders.GenericEncoder(
template_fname='epidemic/epidemic.template',
delimiter='$',
target_filename='epidemic_in.json')
Since the quantity of interest (number of turns until the disease is erradicated) is a function of the simulation output (it is the iteration number of the last row) we need to extend the Decoder class to take this in to account. To this end we inherit from SimpleCSV decoder and redefine the parse_sim_output
method to take the last value of the iteration
column in the file produced by the simulation. This gives us the length in turns for which the simulation ran or in other words before the disease disappeared in our simulation.
class EpidemicDecoder(uq.decoders.SimpleCSV): # , decoder_name='epidemic_decoder'
def parse_sim_output(self, run_info={}):
result = super().parse_sim_output(run_info)
return {'iteration': result['iteration'][-1]}
decoder = EpidemicDecoder(
target_filename="output.csv", output_columns=["iteration"])
We will also define a helper function that will execute the simulation with the provided input files. This function takes a campaign object, creates the directories with input files and then runs our script in them with those files as inputs. The exact details of this process can be found here and here.
We start by creating an EasyVVUQ Campaign. Here we call it 'epidemic_basic'. :
execute = ExecuteLocal('{}/epidemic/epidemic.py epidemic_in.json output.csv'.format(os.getcwd()))
actions = Actions(CreateRunDirectory('/tmp'), Encode(encoder), execute, Decode(decoder))
campaign = uq.Campaign(name='epidemic_basic', params=params, actions=actions)
print(campaign)
db_location = sqlite:///./epidemic_basicw4x3b3o2/campaign.db active_sampler_id = None campaign_name = epidemic_basic campaign_dir = /Volumes/Samsung8TB/dpc/GIT/EasyVVUQ/tutorials/./epidemic_basicw4x3b3o2 campaign_id = 1
We then want to describe our application. This means passing parameter dictionary, enoder, decoder and collater to the campaign object.
For this particular task we are not interested in the relationship between input parameters and the behavior of the simulation. All we want is to see how much the result varies between runs that are identical but for the random seed.
from easyvvuq.sampling import EmptySampler
campaign.set_sampler(EmptySampler())
EmptySampler is a convenience class for such cases. However, another option is, if your simulation provides the option to specify different seeds, to draw the seeds from a probability distribution. In cases like these you could specify the sampler like this:
from easyvvuq.sampling import RandomSampler
vary = {'seed' : cp.DiscreteUniform(0, MAX_SEED)}
sampler = RandomSampler(vary)
In the above example MAX_SEED
is the maximum value the seed can take plus one.
We now want to execute the simulations and collate the results.
campaign.execute(nsamples=20).collate()
We can now see what the result is. It will be a DataFrame containing the number of iterations before the disease is erradicated. Calling collate on the campaign object will find and parse the output files of our simulation and then produce a single DataFrame with the results.
df = campaign.get_collation_result()
print(df)
run_id iteration grid_size n duration mortality ensemble 0 0 0 0 0 0 0 0 1 73.0 10 20 28 0.2 0 1 2 80.0 10 20 28 0.2 0 2 3 89.0 10 20 28 0.2 0 3 4 73.0 10 20 28 0.2 0 4 5 57.0 10 20 28 0.2 0 5 6 88.0 10 20 28 0.2 0 6 7 72.0 10 20 28 0.2 0 7 8 81.0 10 20 28 0.2 0 8 9 81.0 10 20 28 0.2 0 9 10 67.0 10 20 28 0.2 0 10 11 53.0 10 20 28 0.2 0 11 12 59.0 10 20 28 0.2 0 12 13 73.0 10 20 28 0.2 0 13 14 55.0 10 20 28 0.2 0 14 15 56.0 10 20 28 0.2 0 15 16 61.0 10 20 28 0.2 0 16 17 65.0 10 20 28 0.2 0 17 18 70.0 10 20 28 0.2 0 18 19 61.0 10 20 28 0.2 0 19 20 75.0 10 20 28 0.2 0 20 21 77.0 10 20 28 0.2 0 21 22 75.0 10 20 28 0.2 0 22 23 71.0 10 20 28 0.2 0 23 24 76.0 10 20 28 0.2 0 24 25 58.0 10 20 28 0.2 0 25 26 71.0 10 20 28 0.2 0 26 27 89.0 10 20 28 0.2 0 27 28 89.0 10 20 28 0.2 0 28 29 77.0 10 20 28 0.2 0 29 30 65.0 10 20 28 0.2 0 30 31 76.0 10 20 28 0.2 0 31 32 71.0 10 20 28 0.2 0 32 33 85.0 10 20 28 0.2 0 33 34 101.0 10 20 28 0.2 0 34 35 65.0 10 20 28 0.2 0 35 36 88.0 10 20 28 0.2 0 36 37 64.0 10 20 28 0.2 0 37 38 75.0 10 20 28 0.2 0 38 39 27.0 10 20 28 0.2 0 39 40 56.0 10 20 28 0.2 0
This collated data is stored in the campaign database. An analysis element, here EnsembleBoot (which performs bootstraping, can then be applied to the campaign's collation result. :
analysis = uq.analysis.EnsembleBoot(qoi_cols=[("iteration", 0)], stat_func=np.mean, stat_name='mean', alpha=0.05)
campaign.apply_analysis(analysis)
The output of this is dependent on the type of analysis element. :
# Get Descriptive Statistics
results = campaign.get_last_analysis()
print(results)
iteration 0 mean low high True 71.35 66.974375 75.075625
The above gives the mean value of the number of iterations before the disease eradication and the 95% confidence region for that estimator.
Suppose we want to examine the behaviour of the model by doing a sweep across a range of parameter values. Let's say that given a fixed population size and fixed mortality rate we want to see how does the model behaves if we vary the disease duration parameter from 4 days to, for example, 28 days. We will do this in 2 day increments. Also, in order to perform some kind of statistical analysis of the results we will want to do sample the same parameter set several times. We will do this using the replicas mechanism in the EasyVVUQ.
First we define a new campaign.
campaign = uq.Campaign(name='epidemic_sweep', params=params, actions=actions)
We need to define a different kind of sampler. In our case we want BasicSweep which allows one to sample an n-dimensional grid of values. In this case we specify a single parameter and give it a list of values. If you add more parameters the sampler will then sample a Cartesian product of those lists. For examples {'a' : [1, 2], 'b' : [3, 4]}
will result in {'a': 1, 'b': 3}
, {'a': 1, 'b': 4}
, {'a': 2, 'b': 3}
and {'a': 2, 'b': 4}
being sampled.
sweep = {
"duration" : list(range(2, 30, 2))
}
sweep_sampler = uq.sampling.BasicSweep(sweep=sweep)
We will also want to sample each point multiple times in order to construct confidence intervals. To this end we use a ReplicaSampler element which wraps around our sweep_sampler
. In essence it will simply run the simulation with the same inputs multiple times with different random number generator seeds.
from easyvvuq.sampling import ReplicaSampler
sampler = ReplicaSampler(sweep_sampler)
campaign.set_sampler(sampler)
The following steps are the same as in our more basic example. Note that execution might take some time. That is because we are running 280 simulations on the machine that is hosting this notebook.
campaign.execute(nsamples=20 * 14).collate(progress_bar=True)
100%|██████████| 744/744 [03:20<00:00, 3.70it/s]
df = campaign.get_collation_result()
print(df)
run_id iteration duration ensemble_id grid_size n mortality ensemble 0 0 0 0 0 0 0 0 0 1 2.0 2 0 10 20 0.2 0 1 2 3.0 4 1 10 20 0.2 0 2 3 8.0 6 2 10 20 0.2 0 3 4 7.0 8 3 10 20 0.2 0 4 5 9.0 10 4 10 20 0.2 0 .. ... ... ... ... ... .. ... ... 835 836 85.0 20 9 10 20 0.2 0 836 837 70.0 22 10 10 20 0.2 0 837 838 74.0 24 11 10 20 0.2 0 838 839 72.0 26 12 10 20 0.2 0 839 840 71.0 28 13 10 20 0.2 0 [840 rows x 8 columns]
One important difference in the analysis stage below is the addition of the groupby
keyword argument. This means, essentially, that bootstrapping will be done on rows with the same ensemble_id
value separately. So in the end we will get (30 - 2) / 2 = 14
triples of the mean value and lower/upper confidence bounds.
analysis = uq.analysis.EnsembleBoot(qoi_cols=[("iteration", 0)], groupby=('ensemble_id', 0), stat_func=np.mean, stat_name='mean')
campaign.apply_analysis(analysis)
We can now view the resulting DataFrame
.
results = campaign.get_last_analysis()
print(results)
iteration 0 mean low high (ensemble_id, 0) 0 1.333333 1.183333 1.500417 1 3.950000 3.600000 4.350000 2 9.283333 8.082917 10.767500 3 14.666667 12.850000 16.866667 4 28.800000 25.733333 32.384167 5 33.850000 30.048750 37.917083 6 43.283333 38.098750 47.784167 7 51.150000 46.316250 55.767500 8 52.583333 48.648750 56.183333 9 64.525000 60.316250 68.450833 10 64.975000 61.815417 68.134583 11 65.808333 62.082917 69.600417 12 71.500000 68.265417 74.800417 13 75.658333 72.150000 79.435000
#plt.errorbar(list(range(2, 30, 2)),
# results[('iteration', 0)]['mean'].values,
# np.array([results['iteration']['mean'].values - results['iteration']['low'].values,
# results['iteration']['high'].values - results['iteration']['mean'].values]))
#plt.xlabel('duration (days)')
#plt.ylabel('days until eradication')
#plt.show()
This part of the tutorial assumes that you have a Kubernetes cluster access configured. In other words the ~/.kube/config
file needs to be populated with the data that Kubernetes API can use to connect to a cluster. The exact details for how to do this will depend on your cloud service provider. For example, in order to start a cluster on GKE I would need to run the following command (or similar):
!gcloud container clusters create easyvvuq
To use this functionality you would first need to create a Docker image for your simulation. There are many resource on how to do this. See for example here. For your convenience below is the Dockerfile we have used when creating an EasyVVUQ image that will be used in his example. You also need to publish your Docker image to where it can be downloaded by the Kubernetes cluster. An obvious place is Docker Hub.
FROM ubuntu:latest
RUN apt-get update && \
apt-get install -y python3-pip && \
apt-get install -y git && \
apt-get install -y tini && \
pip3 install easyvvuq && \
git clone https://github.com/UCL-CCS/EasyVVUQ.git
ENTRYPOINT ["tini", "--"]
After this, the only other bit of boring admin you need to do is to create a Kubernetes pod configuration. Again here is the one that was used in this tutorial:
apiVersion: v1
kind: Pod
metadata:
name: epidemic
spec:
restartPolicy: Never
containers:
- name: epidemic
image: orbitfold/easyvvuq:latest
command: ["/bin/sh", "-c"]
args: ["python3 /EasyVVUQ/docs/epidemic/epidemic.py /config/example.json out.csv && cat out.csv"]
It will likely be the same for your simulation. Please note the image name, and the command that will be used to execute the simulation. In the current implementation your simulation needs to put all output to the standard output. So after creating the output CSV file we use cat
to print it to the screen. Your way of doing this is likely to be different.
We have, so far, drawn 20 samples in our parameter sweep which has resulted in fairly wide error bars. We can shrink them by adding more samples. To this end will draw 80 more samples and create the corresponding input files.
campaign.draw_samples(60 * 14)
campaign.populate_runs_dir()
Now we call apply_for_each_run_dir
with the ExecuteKubernetes
action. This will submit the jobs to the Kubernetes cluster and the execution will automatically start.
statuses = campaign.apply_for_each_run_dir(
uq.actions.ExecuteKubernetes('../tutorials/kubernetes/epidemic.yaml',
['epidemic_in.json'], 'output.csv'), batch_size=8)
statuses = campaign.sample_and_apply(
60 * 14,
uq.actions.ExecuteKubernetes(
'../tutorials/kubernetes/epidemic.yaml',
['epidemic_in.json'],
'output.csv'), 8)
statuses.start()
print(statuses.progress())
statuses.actions
campaign.collate()
df = campaign.get_collation_result()
analysis = uq.analysis.EnsembleBoot(qoi_cols=["iteration"], groupby='ensemble', stat_func=np.mean, stat_name='mean')
campaign.apply_analysis(analysis)
results = campaign.get_last_analysis()
print(results)
plt.errorbar(list(range(2, 30, 2)),
results['iteration']['mean'].values,
np.array([results['iteration']['mean'].values - results['iteration']['low'].values,
results['iteration']['high'].values - results['iteration']['mean'].values]))
plt.xlabel('duration (days)')
plt.ylabel('days until eradication')
plt.show()
campaign = uq.Campaign(name='sobol_method', params=params, actions=actions)
vary = {
"duration": cp.DiscreteUniform(7, 14),
"mortality": cp.Uniform(0.1, 0.3),
}
qmc_sampler = uq.sampling.QMCSampler(vary, 128)
campaign.set_sampler(qmc_sampler)
campaign.execute().collate(progress_bar=True)
100%|██████████| 512/512 [01:40<00:00, 5.09it/s]
df = campaign.get_collation_result()
print(df)
run_id iteration duration mortality grid_size n ensemble ensemble_id 0 0 0 0 0 0 0 0 0 1 13.0 6.511719 0.175293 10 20 0 0 1 2 14.0 10.082031 0.175293 10 20 0 0 2 3 12.0 6.511719 0.197363 10 20 0 0 3 4 9.0 10.082031 0.197363 10 20 0 0 4 5 24.0 10.511719 0.275293 10 20 0 0 .. ... ... ... ... ... .. ... ... 507 508 13.0 10.394531 0.173926 10 20 0 0 508 509 12.0 6.574219 0.258105 10 20 0 0 509 510 55.0 14.394531 0.258105 10 20 0 0 510 511 5.0 6.574219 0.273926 10 20 0 0 511 512 24.0 14.394531 0.273926 10 20 0 0 [512 rows x 8 columns]
from easyvvuq.analysis import QMCAnalysis
analysis = QMCAnalysis(qmc_sampler, ['n', 'duration', 'mortality'])
campaign.apply_analysis(analysis)
results = campaign.get_last_analysis()
results.sobols_first()
{'n': {'duration': array([0.]), 'mortality': array([0.])}, 'duration': {'duration': array([1.01554484]), 'mortality': array([0.])}, 'mortality': {'duration': array([0.]), 'mortality': array([1.01224875])}}
results.sobols_first('mortality')
{'duration': array([0.]), 'mortality': array([1.01224875])}
results.sobols_first('duration')
{'duration': array([1.01554484]), 'mortality': array([0.])}
results.describe()
n | duration | mortality | |
---|---|---|---|
0 | 0 | 0 | |
mean | 20.0 | 10.484375 | 0.199609 |
var | 0.0 | 5.333023 | 0.003333 |
std | 0.0 | 2.309334 | 0.057733 |
campaign = uq.Campaign(name='epidemic_sc', params=params, actions=actions)
vary = {
"n": cp.DiscreteUniform(10, 100),
"duration": cp.DiscreteUniform(7, 56),
"mortality": cp.Uniform(0.0, 1.0),
}
sampler = uq.sampling.SCSampler(vary=vary, polynomial_order=9)
campaign.set_sampler(sampler)
campaign.execute().collate(progress_bar=True)
100%|██████████| 1000/1000 [03:52<00:00, 4.30it/s]
df = campaign.get_collation_result()
print(df)
run_id iteration n duration mortality grid_size ensemble ensemble_id 0 0 0 0 0 0 0 0 0 1 12.0 14.0 9.0 0.013047 10 0 0 1 2 10.0 14.0 9.0 0.067468 10 0 0 2 3 22.0 14.0 9.0 0.160295 10 0 0 3 4 13.0 14.0 9.0 0.283302 10 0 0 4 5 8.0 14.0 9.0 0.425563 10 0 0 .. ... ... ... ... ... ... ... ... 995 996 69.0 96.0 54.0 0.574437 10 0 0 996 997 71.0 96.0 54.0 0.716698 10 0 0 997 998 74.0 96.0 54.0 0.839705 10 0 0 998 999 65.0 96.0 54.0 0.932532 10 0 0 999 1000 69.0 96.0 54.0 0.986953 10 0 0 [1000 rows x 8 columns]
And the analysis can be done with: :
analysis = uq.analysis.SCAnalysis(sampler=sampler, qoi_cols=["iteration"])
campaign.apply_analysis(analysis)
results_sc = campaign.get_last_analysis()
results_sc.describe()
iteration | |
---|---|
0 | |
mean | 62.304993 |
std | 24.803309 |
var | 615.204157 |
analysis.surrogate('iteration', np.array([5, 7, 0.7]))
array([133196.93383935])
n = np.linspace(5, 100, 20)
duration = np.linspace(7, 50, 20)
mortality = 0.7
grid = np.meshgrid(n, duration, mortality)
z = np.array([analysis.surrogate('iteration', np.array([grid[0].flatten()[i], grid[1].flatten()[i], mortality])) for i in range(n.shape[0] * duration.shape[0])]).flatten()
fig = plt.figure(figsize=(10, 10), dpi= 80, facecolor='w', edgecolor='k')
ax = fig.gca(projection='3d')
ax.set_zlim(0, 200)
ax.plot_trisurf(grid[0].flatten(), grid[1].flatten(), z)
ax.set_xlabel('n')
ax.set_ylabel('duration')
ax.set_zlabel('time')
plt.show()
Let us quickly check the correspondence of the surrogate model to the data we have collected during our parameter sweep experiment. What we're looking for is the curve should be within the confidence interval we have calculated.
n = [20] * 14
duration = list(range(2, 30, 2))
mortality = [0.2] * 14
plt.plot(([analysis.surrogate('iteration', np.array([n_, duration_, mortality_]))[0] for n_, duration_, mortality_ in zip(n, duration, mortality)]))
plt.ylim(0, 100)
plt.show()