Epidemiology: The SIR model: Difference between revisions
A WASSERMANN (talk | contribs) No edit summary |
A WASSERMANN (talk | contribs) No edit summary |
||
(55 intermediate revisions by 2 users not shown) | |||
Line 1: | Line 1: | ||
Simulation of differential equations with turtle graphics. | Simulation of differential equations with turtle graphics using [http://jsxgraph.org JSXGraph]. | ||
===SIR model without vital dynamics=== | ===SIR model without vital dynamics=== | ||
The SIR model measures the number of susceptible, infected, and recovered individuals in a host population. | |||
Given a fixed population, let <math>S(t)</math> be the fraction that is susceptible to an infectious, but not deadly, disease at time t; let <math>I(t)</math> be the fraction that is infected at time <math>t</math>; | |||
and let <math>R(t)</math> be the fraction that has recovered. Let <math>\beta</math> be the rate at which an infected person infects a susceptible person. Let <math>\gamma</math> be the rate at which infected people recover from the disease. | |||
A single epidemic outbreak is usually far more rapid than the vital dynamics of a population, thus, if the aim is to study the immediate consequences of a single epidemic, one may neglect birth-death processes. In this case the SIR system can be expressed by the following set of differential equations: | |||
:<math> \frac{dS}{dt} = - \beta I S </math> | |||
:<math> \frac{dR}{dt} = \gamma I </math> | |||
:<math> \frac{dI}{dt} = -(\frac{dS}{dt}+\frac{dR}{dt}) </math> | |||
====Example Hong Kong flu==== | |||
* initially 7.9 million people, | |||
* 10 infected, | |||
* 0 recovered. | |||
* estimated average period of infection: 3 days, so <math>\gamma = 1/3</math> | |||
* infection rate: one new person every other day, so <math>\beta = 1/2</math> | |||
Thus S(0) = 1, I(0) = 1.27E-6, R(0) = 0, see [http://www.cs.princeton.edu/introcs/94diffeq/]. | |||
The lines in the JSXGraph-simulation below have the following meaning: | |||
* <span style="color:Blue">Blue: Rate of susceptible population</span> | |||
* <span style="color:red">Red: Rate of infected population</span> | |||
* <span style="color:green">Green: Rate of recovered population (which means: immune, isolated or dead) | |||
<html> | |||
<form><input type="button" value="clear and run a simulation of 100 days" onClick="clearturtle();run()"> | |||
<input type="button" value="stop" onClick="stop()"> | |||
<input type="button" value="continue" onClick="goOn()"></form> | |||
</html> | |||
<jsxgraph box="box" width="600" height="450"> | |||
var brd = JXG.JSXGraph.initBoard('box', {axis: true, boundingbox: [-5, 1.2, 100, -1.2]}); | |||
var S = brd.create('turtle',[],{strokeColor:'blue',strokeWidth:3}); | |||
var I = brd.create('turtle',[],{strokeColor:'red',strokeWidth:3}); | |||
var R = brd.create('turtle',[],{strokeColor:'green',strokeWidth:3}); | |||
var s = brd.create('slider', [[0,-0.3], [30,-0.3],[0,1.27E-6,1]], {name:'s'}); | |||
brd.create('text', [40,-0.3, "initially infected population rate (on load: I(0)=1.27E-6)"]); | |||
var beta = brd.create('slider', [[0,-0.4], [30,-0.4],[0,0.5,1]], {name:'β'}); | |||
brd.create('text', [40,-0.4, "β: infection rate"]); | |||
var gamma = brd.create('slider', [[0,-0.5], [30,-0.5],[0,0.3,1]], {name:'γ'}); | |||
brd.create('text', [40,-0.5, "γ: recovery rate = 1/(days of infection)"]); | |||
var t = 0; // global | |||
brd.create('text', [40,-0.2, | |||
function() {return "Day "+t+": infected="+(7900000*I.Y()).toFixed(1)+" recovered="+(7900000*R.Y()).toFixed(1);}]); | |||
S.hideTurtle(); | |||
I.hideTurtle(); | |||
R.hideTurtle(); | |||
function clearturtle() { | |||
S.cs(); | |||
I.cs(); | |||
R.cs(); | |||
S.hideTurtle(); | |||
I.hideTurtle(); | |||
R.hideTurtle(); | |||
} | |||
function run() { | |||
S.setPos(0,1.0-s.Value()); | |||
R.setPos(0,0); | |||
I.setPos(0,s.Value()); | |||
delta = 1; // global | |||
t = 0; // global | |||
loop(); | |||
} | |||
function turtleMove(turtle,dx,dy) { | |||
turtle.moveTo([dx+turtle.X(),dy+turtle.Y()]); | |||
} | |||
function loop() { | |||
var dS = -beta.Value()*S.Y()*I.Y(); | |||
var dR = gamma.Value()*I.Y(); | |||
var dI = -(dS+dR); | |||
turtleMove(S,delta,dS); | |||
turtleMove(R,delta,dR); | |||
turtleMove(I,delta,dI); | |||
t += delta; | |||
if (t<100.0) { | |||
active = setTimeout(loop,10); | |||
} | |||
} | |||
function stop() { | |||
if (active) clearTimeout(active); | |||
active = null; | |||
} | |||
function goOn() { | |||
if (t>0) { | |||
if (active==null) { | |||
active = setTimeout(loop,10); | |||
} | |||
} else { | |||
run(); | |||
} | |||
} | |||
</jsxgraph> | |||
===The underlying JavaScript code=== | |||
<source lang="javascript"> | |||
var brd = JXG.JSXGraph.initBoard('box', {axis: true, boundingbox: [-5, 1.2, 100, -1.2]}); | |||
var S = brd.create('turtle',[],{strokeColor:'blue',strokeWidth:3}); | |||
var I = brd.create('turtle',[],{strokeColor:'red',strokeWidth:3}); | |||
var R = brd.create('turtle',[],{strokeColor:'green',strokeWidth:3}); | |||
var s = brd.create('slider', [[0,-0.3], [30,-0.3],[0,1.27E-6,1]], {name:'s'}); | |||
brd.create('text', [40,-0.3, "initially infected population rate (on load: I(0)=1.27E-6)"]); | |||
var beta = brd.create('slider', [[0,-0.4], [30,-0.4],[0,0.5,1]], {name:'β'}); | |||
brd.create('text', [40,-0.4, "β: infection rate"]); | |||
var gamma = brd.create('slider', [[0,-0.5], [30,-0.5],[0,0.3,1]], {name:'γ'}); | |||
brd.create('text', [40,-0.5, "γ: recovery rate = 1/(days of infection)"]); | |||
var t = 0; // global | |||
brd.create('text', [40,-0.2, | |||
function() {return "Day "+t+": infected="+(7900000*I.Y()).toFixed(1)+" recovered="+(7900000*R.Y()).toFixed(1);}]); | |||
S.hideTurtle(); | |||
I.hideTurtle(); | |||
R.hideTurtle(); | |||
function clearturtle() { | |||
S.cs(); | |||
I.cs(); | |||
R.cs(); | |||
S.hideTurtle(); | |||
I.hideTurtle(); | |||
R.hideTurtle(); | |||
} | |||
function run() { | |||
S.setPos(0,1.0-s.Value()); | |||
R.setPos(0,0); | |||
I.setPos(0,s.Value()); | |||
delta = 1; // global | |||
t = 0; // global | |||
loop(); | |||
} | |||
function turtleMove(turtle,dx,dy) { | |||
turtle.moveTo([dx+turtle.X(),dy+turtle.Y()]); | |||
} | |||
function loop() { | |||
var dS = -beta.Value()*S.Y()*I.Y(); | |||
var dR = gamma.Value()*I.Y(); | |||
var dI = -(dS+dR); | |||
turtleMove(S,delta,dS); | |||
turtleMove(R,delta,dR); | |||
turtleMove(I,delta,dI); | |||
t += delta; | |||
if (t<100.0) { | |||
active = setTimeout(loop,10); | |||
} | |||
} | |||
function stop() { | |||
if (active) clearTimeout(active); | |||
active = null; | |||
} | |||
function goOn() { | |||
if (t>0) { | |||
if (active==null) { | |||
active = setTimeout(loop,10); | |||
} | |||
} else { | |||
run(); | |||
} | |||
} | |||
</source> | |||
===See also=== | |||
* [[Epidemiology: The SEIR model]] | |||
* [[Population growth models]] | |||
* [[Lotka-Volterra equations]] | |||
===References=== | |||
* [http://en.wikipedia.org/wiki/Compartmental_models_in_epidemiology http://en.wikipedia.org/wiki/Compartmental_models_in_epidemiology] | |||
* [http://mathworld.wolfram.com/SIRModel.html http://mathworld.wolfram.com/SIRModel.html] | |||
* [http://www.cs.princeton.edu/introcs/94diffeq http://www.cs.princeton.edu/introcs/94diffeq] | |||
[[Category:Examples]] | |||
[[Category:Turtle Graphics]] | |||
[[Category:Calculus]] |
Latest revision as of 14:59, 20 February 2013
Simulation of differential equations with turtle graphics using JSXGraph.
SIR model without vital dynamics
The SIR model measures the number of susceptible, infected, and recovered individuals in a host population. Given a fixed population, let [math]\displaystyle{ S(t) }[/math] be the fraction that is susceptible to an infectious, but not deadly, disease at time t; let [math]\displaystyle{ I(t) }[/math] be the fraction that is infected at time [math]\displaystyle{ t }[/math]; and let [math]\displaystyle{ R(t) }[/math] be the fraction that has recovered. Let [math]\displaystyle{ \beta }[/math] be the rate at which an infected person infects a susceptible person. Let [math]\displaystyle{ \gamma }[/math] be the rate at which infected people recover from the disease.
A single epidemic outbreak is usually far more rapid than the vital dynamics of a population, thus, if the aim is to study the immediate consequences of a single epidemic, one may neglect birth-death processes. In this case the SIR system can be expressed by the following set of differential equations:
- [math]\displaystyle{ \frac{dS}{dt} = - \beta I S }[/math]
- [math]\displaystyle{ \frac{dR}{dt} = \gamma I }[/math]
- [math]\displaystyle{ \frac{dI}{dt} = -(\frac{dS}{dt}+\frac{dR}{dt}) }[/math]
Example Hong Kong flu
- initially 7.9 million people,
- 10 infected,
- 0 recovered.
- estimated average period of infection: 3 days, so [math]\displaystyle{ \gamma = 1/3 }[/math]
- infection rate: one new person every other day, so [math]\displaystyle{ \beta = 1/2 }[/math]
Thus S(0) = 1, I(0) = 1.27E-6, R(0) = 0, see [1].
The lines in the JSXGraph-simulation below have the following meaning:
* Blue: Rate of susceptible population * Red: Rate of infected population * Green: Rate of recovered population (which means: immune, isolated or dead)
The underlying JavaScript code
var brd = JXG.JSXGraph.initBoard('box', {axis: true, boundingbox: [-5, 1.2, 100, -1.2]});
var S = brd.create('turtle',[],{strokeColor:'blue',strokeWidth:3});
var I = brd.create('turtle',[],{strokeColor:'red',strokeWidth:3});
var R = brd.create('turtle',[],{strokeColor:'green',strokeWidth:3});
var s = brd.create('slider', [[0,-0.3], [30,-0.3],[0,1.27E-6,1]], {name:'s'});
brd.create('text', [40,-0.3, "initially infected population rate (on load: I(0)=1.27E-6)"]);
var beta = brd.create('slider', [[0,-0.4], [30,-0.4],[0,0.5,1]], {name:'β'});
brd.create('text', [40,-0.4, "β: infection rate"]);
var gamma = brd.create('slider', [[0,-0.5], [30,-0.5],[0,0.3,1]], {name:'γ'});
brd.create('text', [40,-0.5, "γ: recovery rate = 1/(days of infection)"]);
var t = 0; // global
brd.create('text', [40,-0.2,
function() {return "Day "+t+": infected="+(7900000*I.Y()).toFixed(1)+" recovered="+(7900000*R.Y()).toFixed(1);}]);
S.hideTurtle();
I.hideTurtle();
R.hideTurtle();
function clearturtle() {
S.cs();
I.cs();
R.cs();
S.hideTurtle();
I.hideTurtle();
R.hideTurtle();
}
function run() {
S.setPos(0,1.0-s.Value());
R.setPos(0,0);
I.setPos(0,s.Value());
delta = 1; // global
t = 0; // global
loop();
}
function turtleMove(turtle,dx,dy) {
turtle.moveTo([dx+turtle.X(),dy+turtle.Y()]);
}
function loop() {
var dS = -beta.Value()*S.Y()*I.Y();
var dR = gamma.Value()*I.Y();
var dI = -(dS+dR);
turtleMove(S,delta,dS);
turtleMove(R,delta,dR);
turtleMove(I,delta,dI);
t += delta;
if (t<100.0) {
active = setTimeout(loop,10);
}
}
function stop() {
if (active) clearTimeout(active);
active = null;
}
function goOn() {
if (t>0) {
if (active==null) {
active = setTimeout(loop,10);
}
} else {
run();
}
}