Read and discuss:
Dushoff, J, Plotkin, JB, Levin, SA, & Earn, DJD. (2004). Dynamical resonance can account for seasonality of influenza epidemics. Proceedings of the National Academy of Sciences of the United States of America, 101(48), 16915–6. doi:10.1073/pnas.0407293101
Code the deterministic model in Dushoff et al. 2004. If you want to do this in anything other than R, you're on your own. If you'll be using R, here are some resources to help you get started:
- Complete the R Studio tutorial and R tutorials 1-3, if you haven't already done this and are unfamiliar with either of these programs (or if you need a refresher).
- Complete the lab on ODE's in R
During the lab meeting we'll discuss issues that arose and troubleshoot each other's code, as needed. We'll also look at some of the issues we discussed today, such as what happens when you initialize at the endemic equilibrium (with and without seasonal forcing).
Your program should allow you to easily do these things:
- Reproduce Figure 1, panels A and B (minus the stochastic model results, of course!)
- Explore the effect of using different initial conditions on the dynamics observed, including what happens if you start at the disease-free equilibrium and what happens if you start at the endemic equlibrium (for the non-seasonal case)
- Explore the effect of using different parameter values on the dynamics observed
- Compare the model trajectory with and without a seasonal contact rate
- Calculate the intrinsic period of oscillation for a given set of parameters
This doesn't necessarily mean that your code already does these things, but that it would take you a negligible amount of time for you to implement them, if someone were to ask you a specific question that requires any of these things.
You should come to lab meeting with a final version of your code that meets the criteria above. Note that, if you're paying attention, you will encounter some difficulty with the first benchmark, even if your code is working properly!
Please also review the endemic equilibrium calculations if you found them difficult to follow during lab meeting.
Write code for the stochastic model in Dushoff et al. 2004. To get started with the Gillespie algorithm, and implementing it in R:
- Review page 6-11 of this presentation and the references cited therein.
- Follow this step-by-step guide to learn how to write a Gillespie implementation of the SIRS model in R; feel free to look at the example code if you get stuck
This week's lab meeting will be a troubleshooting session, guided by Tom. Benchmarks to aim for with your code will be added after the meeting, and you should make any necessary modifications to meet them by July 25.
You should come to lab meeting with a final version of your code that meets the criteria below. Note that, if you have written your code in R, running a simulation with a population of 500,000 inidividuals will take an extremely long time, even if your code is working properly! You are not expected to reproduce the blue lines in Figure 1 using your R code!
Your program should allow you to easily do these things:
- Output the time until extinction for a single run of the Gillespie simulation
- Return the amount of time it takes to run a single realization of the Gillespie algorithm
- Explore the effect of using different parameter values on the dynamics observed (for an arbitrary population size)
- Visually compare the ODE model trajectory to a single realization of the Gillespie implementation, with the same parameter values and initial conditions
Again, this doesn't necessarily mean that your code already does these things, but that it would take you a negligible amount of time for you to implement them, if someone were to ask you a specific question that requires any of these things.
Reproduce Figure 1 (both panels). You may want to use the GillespieSSA R package to do this, though you should feel free to do it through other means (e.g., coding the algorithm in C or C++, using an approximate algorithm that you write yourself, or even writing a more efficient pure R implementation