Workshop 3 - Exploring Agent-Based models
In this workshop we will explore agent-based models using NetLogo.
1 Setting up NetLogo
If you’re on a university computer you can run NetLogo from AppsAnywhere. If you’re using your own computer, you can download NetLogo (for free) here.
Once you have NetLogo installed and open, you can open a model via the ‘File’ menu. There is a large selection of models in the library, many of which have useful and interesting information with them and demonstrate well the sorts of dynamics ABMs can simulate.
To open a netlogo model file, first download and save the file somewhere, and then navigate to it from within NetLogo. In this workshop we will be working with the file ‘Epidemic_ABM_RHO.nlogo’, which is available in Ultra under ‘Resources’ and ‘Rfiles’, and also on the module webpage. Save this file somewhere on your computer, and open it in NetLogo.
2 The mass action version
We will first briefly explore our recreation of a mass action model, by keeping all the turtles still and having all turtles equally likely to be infected by any infectious turtle. Remember that this was to (imperfectly) mimic the mass action of the compartmental models. In the model Epidemic_ABM.nlogo, there are various ways to do this - either the turtles move completely at random, or they infect other people at random (rather than nearby) or both.
First of all, it will be good to play around with the values of infectious-period and prob-transmission-pp and see what effect different settings have on the system.
Vary the input values and run the model at these different settings. As you’re experimenting with different input values, make a note of things like the value of \(R_0\), the final epidemic size (or spread), how long (in time steps) until the epidemic dies out and so on.
- What is the effect of changing each parameter?
- What happens if you keep the parameters the same?
2.1 Repeated experiments with BehaviourSpace
The explorations above gave us an idea of how much the outcome varied as the input parameters are changed (or as we move around the input space, in Lecture 3 parlance). To explore the uncertainty coming from the randomness that is built into the model, we perform repeated experiments using BehaviourSpace.
BehaviourSpace can be found in the Tools menu, and the window is mostly intuitive - if you have problems or questions about it, speak to a tutor. The ‘table’ output is stored as a .csv (comma separated value) file, which can be easily read into R. It is also in a fairly convenient form for plotting in R, since it is in ‘long’ format (more rows, fewer columns). The ‘spreadsheet’ file is more human friendly, as it is in wide format (a column for each repeat).
2.1.1 Using BehaviorSpace
To use BehaviorSpace, go to ‘Tools -> BehaviorSpace’. This will open the BehaviorSpace window, shown in Figure 2.1.
Figure 2.1: The BehaviourSpace window.
To edit an existing experiment, select that experiment and click ‘edit’. To start a new experiment, click ‘new’. Either of these actions will open the Experiment window, shown in Figure 2.2.
Figure 2.2: The Experiment window.
The Experiment window is where you specify the parameter values and output variables for your experiment. In the second box (‘Vary variables as follows’), give the name of each variable followed by the value you wish that variable to take, for each variable:
["<variable name>" "<variable value>" ]
In the fourth box (‘Measure runs using these reporters’), give the values you would like BehaviorSpace to monitor at each tick. In Figure 2.2, the outputs will be
- The count of infectious turtles
- The count of removed/recovered turtles
- The estimate of \(R_n\).
After pressing ‘OK’, you will return to the BehaviorSpace window. To run an experiment, select that experiment and click ‘run’. This will open the ‘Run options’ window, show in Figure 2.3.
Figure 2.3: The Run options window.
Enter filenames for the output files you wish to create, and tick ‘OK’ to run the experiment.
2.1.2 BehaviorSpace output in R
The following few lines of code should give you a framework for processing and visualising the output files from your BehaviourSpace experiments, but if you run into problems please speak to a tutor. It may be helpful to open the file in Excel (or similar) if the code below doesn’t work, as some of the numbers may need to be changed.
library(ggplot2)
## Read the file into R - the filepath will need to be changed.
## Setting your working directory to the folder you're saving the output
## files in will make this simpler.
## Naming your files and objects well can also help as you try to distinguish
## between many similar plots and files
exp_anywhere_ip6_prob06 = read.csv("Data/ABM_MA1_p06-table", skip=6, header=T)
## A quick way to check that the data has imported correctly (and to see
## what the columns are) is to use 'summary'
summary(exp_anywhere_ip6_prob06) ## X.run.number. movement infectious.period transmit.to
## Min. : 1.00 Length:2003 Min. :6 Length:2003
## 1st Qu.:14.00 Class :character 1st Qu.:6 Class :character
## Median :26.00 Mode :character Median :6 Mode :character
## Mean :26.24 Mean :6
## 3rd Qu.:40.00 3rd Qu.:6
## Max. :50.00 Max. :6
## initial.I prob.transmission.pp prob.move num.people X.step.
## Min. :1 Min. :0.6 Min. :1 Min. :4000 Min. : 0.00
## 1st Qu.:1 1st Qu.:0.6 1st Qu.:1 1st Qu.:4000 1st Qu.:14.00
## Median :1 Median :0.6 Median :1 Median :4000 Median :30.00
## Mean :1 Mean :0.6 Mean :1 Mean :4000 Mean :31.16
## 3rd Qu.:1 3rd Qu.:0.6 3rd Qu.:1 3rd Qu.:4000 3rd Qu.:47.00
## Max. :1 Max. :0.6 Max. :1 Max. :4000 Max. :86.00
## count.turtles.with...status....I.. count.turtles.with...status....R..
## Min. : 0.0 Min. : 0
## 1st Qu.: 3.0 1st Qu.:1606
## Median : 36.0 Median :3872
## Mean : 373.7 Mean :2848
## 3rd Qu.: 378.5 3rd Qu.:3993
## Max. :2591.0 Max. :4000## Rename the columns of interest to be easier to work with
## These will likely be:
## - X.run.number
## - X.step
## - count.turtles.with...Status..I
## - count.turtles.with...Status..R
## but may include others #
## The following line selects the columns of interest and renames them
names(exp_anywhere_ip6_prob06)[c(1,9,10,11)] = c("runID", "TimeStep", "CountInfected", "CountRecovered")
summary(exp_anywhere_ip6_prob06)## runID movement infectious.period transmit.to
## Min. : 1.00 Length:2003 Min. :6 Length:2003
## 1st Qu.:14.00 Class :character 1st Qu.:6 Class :character
## Median :26.00 Mode :character Median :6 Mode :character
## Mean :26.24 Mean :6
## 3rd Qu.:40.00 3rd Qu.:6
## Max. :50.00 Max. :6
## initial.I prob.transmission.pp prob.move num.people TimeStep
## Min. :1 Min. :0.6 Min. :1 Min. :4000 Min. : 0.00
## 1st Qu.:1 1st Qu.:0.6 1st Qu.:1 1st Qu.:4000 1st Qu.:14.00
## Median :1 Median :0.6 Median :1 Median :4000 Median :30.00
## Mean :1 Mean :0.6 Mean :1 Mean :4000 Mean :31.16
## 3rd Qu.:1 3rd Qu.:0.6 3rd Qu.:1 3rd Qu.:4000 3rd Qu.:47.00
## Max. :1 Max. :0.6 Max. :1 Max. :4000 Max. :86.00
## CountInfected CountRecovered
## Min. : 0.0 Min. : 0
## 1st Qu.: 3.0 1st Qu.:1606
## Median : 36.0 Median :3872
## Mean : 373.7 Mean :2848
## 3rd Qu.: 378.5 3rd Qu.:3993
## Max. :2591.0 Max. :4000## Plot [in this case] the CountInfected data against time
## the 'group' argument means that the plot shows each runID as a separate line
## but doesn't colour them or make a legend
ggplot(data=exp_anywhere_ip6_prob06, aes(x=TimeStep, y=CountInfected, group = runID)) + geom_path()
Figure 2.4: The count of infected turtles over time for 50 ABM runs with the same input values.
Figure 2.4 shows the count of Infected turtles for 50 runs of the same model. The model is a mass action style one, where there is no movement and infectious turtles can infect other turtles in any part of the space. The probability of transmission per person (in effective contact) is set to 0
You can also tabulate and plot particular characteristics of each repeat. For example, here we find the time step of the peak (where the number of infected turtles is at its maximum) for each run:
# create a vector of runIDs in the dataset
runIDs = unique(exp_anywhere_ip6_prob06$runID)
# Create a vector the same length as runIDs
# This vector will hold the time of maximum for each run
times_maxI = rep(NA, length(runIDs))
# For each i in runIDs
for (i in runIDs){
## create a dataframe containing that run's data
exp_i = exp_anywhere_ip6_prob06[exp_anywhere_ip6_prob06$runID == i,]
## find the maximum of the infected values
max_inf_i = max(exp_i$CountInfected)
## find the time step at which that occurs
time_maxinf_i = exp_i$TimeStep[exp_i$CountInfected == max_inf_i]
## record the time step of maximum infection in the times_maxI vector
## Sometimes there will be more than one timestep for which the infectious value
## is at its highest, so we use 'min' to choose the first of these
times_maxI[i] = min(time_maxinf_i)
}
# Plot this as a histogram
hist(times_maxI, breaks=20)
Figure 2.5: Histogram of the time-step of the peak of infection for each of the 50 runs.
For the model in Figures 2.4 and 2.5 the number of initial Infectious was set to one, and you can see that this means the epidemic dies out immediately for almost 20 of the 50 runs. Setting the initial number of infectious to a higher number will make this less likely.
Examples of quantities that could be interesting to plot are:
- The infected count over time
- The cured count over time
- The time of the maximum (as a histogram)
- The time of the epidemic dying out (as a historgram)
- The final epidemic size (as a histogram)
- The maximum number of infected turtles at any one time (as a histogram)
and you may well think of more. If there is something you’d like to plot and you’re not sure how, please speak to your tutor.
For some different settings of input values, run a repeated experiment in BehaviourSpace. Do these as separate experiments for each input parameter combination. Make sure to save the output files with different names. You could team up with a friend to save replicating the same (or very nearby) points.
- Does the amount of uncertainty in things like location and height of peak, final epidemic size or time to end of epidemic vary for the different points?
3 Local movement
We will now move onto the [slightly more realistic] version in which turtles can only move a short distance at each timestep, and can only infect turtles on their own patch. Before you begin to experiment with this model, please think about the following questions (feel free to discuss them with others).
What effect do you think the change in movement will have (and why) on:
- The speed and extent of the epidemic?
- The way the epidemic spreads?
To visually compare the different versions of the model, it helps to start with just one initially infectious person. Remember that this is the starting point for SIR compartmental models too.
Setting initial-I to one, try various configurations of the model, watching the progression of the epidemic in the main window. Which settings do you think are more realistic?
We can now compare the different versions more systematically.
Replicate your model exploration from Exercise 2.1, and compare your findings. Do model runs with the same parameter values behave similarly? If not, how?
We can go one step further and compare the uncertainty at each point, using BehaviourSpace.
For the same input points that you used in Exercise 2.2, run experiments in BehaviourSpace with the local movement model. How do the results compare? Were your answers in Exercise 3.1 correct?
Finally, we can think about how this links to the compartmental SIR models we studied in week 1.
Looking back at some of the example plots for the SIR model from Lecture 1, or using the code from the workshop:
- Can you find values of the ABM parameters that correspond to the SIR example?
- Run the ABM model at these parameters, for different versions of the settings. How do the results compare?
- Which settings in the ABM result in greater differences between ABM and SIR output?