The SIRS model is
Where \(S=S(t)\) is the population of susceptible individuals, \(I=I(t)\), the infected population, and \(R=R(t)\) the recovered population. This diagram shows the model where now some of the recovered population can become susceptible again and become infected. The parameter \(\mu \) indicates how much of the recovered population could become susceptible again.
The units of \(S,R,I\) are population measured in person. These values can not be negative since they are population amount. The units of \(\mu \) is \(\frac {1}{\text {time}}\) where time can be day or week or any other unit of time. The units of \(\upsilon \) is also \(\frac {1}{\text {time}}\). The units of \(\beta \) is \(\frac {1}{\left ( \text {time}\right ) \left ( \text {person}\right ) }\)
The following program allows one to do analysis on this model and it also displays the critical points.
The program is written using Mathematica. Here is the notebook
After downloading the notebook and opening it inside Mathematica, you can now use the sliders and analyze the system behaviour for different parameters.
Source code
Manipulate[ tick; Module[{ic, eq1, eq2, eq3, S, R, I0, t, initialR, sol}, If[initialS + initialI > 1000, initialI = 1000 - initialS]; initialR = 1000 - (initialS + initialI); ic = {S[0] == initialS, I0[0] == initialI, R[0] == initialR}; eq1 = S'[t] == -beta S[t] *I0[t] + mu*R[t]; eq2 = I0'[t] == beta*S[t]*I0[t] - v I0[t]; eq3 = R'[t] == v*I0[t] - mu*R[t]; sol = NDSolve[{eq1, eq2, eq3, ic}, {S, I0, R}, {t, 0, 60}]; Grid[{ {Grid[{ {"Transmission rate (beta) ", padIt2[beta, {5, 4}]}, {"Recovery rate (v) ", padIt2[v, {3, 3}]}, {"reinfection rate (mu) ", padIt2[mu, {3, 3}]}, {"v/beta ", padIt2[v/beta, {6, 3}]}, {}, {If[1000 >= v/beta, Row[{Style["case (a)", Bold], " population >= v/beta"}], Row[{Style["case (b)", Bold], " population <= v/beta"}] ]}}, Frame -> True ] }, {Plot[{ Evaluate[S[t] /. sol], Evaluate[I0[t] /. sol], Evaluate[R[t] /. sol], If[addSI, Evaluate[S[t] /. sol] + Evaluate[I0[t] /. sol], Nothing], If[addSI, 1000] }, {t, 0, 60}, PlotLegends -> If[addSI, {"susceptible (S)", "Infected (I)", "Recovered (R)", "S+I", "Total population"}, {"susceptible (S)", "Infected (I)", "Recovered (R)"}], PlotStyle -> If[addSI, {Blue, Red, Black, Dashed, Green}, {Blue, Red, Black}], GridLines -> Automatic, GridLinesStyle -> LightGray, AxesLabel -> {"Time (days)", "Population size"}, ImageSize -> 500, ImagePadding -> 30, PlotRange -> {{0, 60}, {0, 1000}} ] }}] ], {{beta, 0.0005, "Transmission rate (beta)"}, 0.0001, 0.003, 0.0001, Appearance -> "Labeled"}, {{v, 0.05, "Recovery rate (v)"}, 0.001, 0.4, 0.001, Appearance -> "Labeled"}, {{mu, 0.05, "Recovery to Susceptible rate (mu)"}, 0.001, 0.4, 0.001, Appearance -> "Labeled"}, {{initialI, 5, "Initial infections size I(0)"}, 0, 1000, 1, Appearance -> "Labeled"}, {{initialS, 995, "Initial susceptible size S(0)"}, 0, 1000, 1, Appearance -> "Labeled"}, Grid[{{ Button[ Text[Style["First critical point", 11]], {initialI = 0, initialS = 1000, tick = Not[tick]}, ImageSize -> {100, 35}], Button[Text[Style["second critical point", 11]], tick = Not[tick]; If[1000 >= v/beta, {initialI = mu/(v + mu)*(1000 - v/beta), initialS = v/beta} , Nothing ], ImageSize -> {130, 35}], Button[ Text[Style["Add S+I curve", 11]], {addSI = True, tick = Not[tick]}, ImageSize -> {100, 35}], Button[ Text[Style["Remove S+I curve", 11]], {addSI = False, tick = Not[tick]}, ImageSize -> {110, 35}] } } ], {{tick, False}, None}, {{addSI, True}, None}, TrackedSymbols :> {beta, v, mu, initialI, initialS, tick}, SaveDefinitions -> False, Initialization :> (padIt2[v_, f_List] := AccountingForm[v, f, NumberSigns -> {"", ""}, NumberPadding -> {"0", "0"}, SignPadding -> True]) ]