%% Ein einfaches SIR-Modell
%

%% Kontext
% Das SIR-Modell (susceptible-infected-removed model) beschreibt in sehr
% vereinfachter Form die Ausbreitung von ansteckenden Krankheiten mit
% Immunitätsbildung.
%  
% Das Modell ist formuliert als System von drei Differentialgleichungen 1. Ordnung:
%
% $$\frac{dS}{dt}=-\frac{\beta IS}{N} $$
%
% $$\frac{dI}{dt}= \frac{\beta IS}{N}-\gamma I $$
% 
% $$\frac{dR}{dt}=\gamma I $$
%
% Darin steht $S(t)$ für die Anzahl der suszeptiblen (anfälligen) Personen,
% $I(t)$ ist die Zahl der Infizierten und $R(t)$ die aus dem Modell
% rausgenommenen Personen (alle, die nicht mehr infiziert werden können,
% weil sie immun oder tot sind). $N = S(t) + I(t) + R(t)$ ist die konstante
% Gesamtgröße der Population. (Dieses Modell berücksichtigt also nicht die
% allgemeine demografische Entwicklung, also Geburten und nicht durch die
% Epidemie bedingte Todesfälle; das ließe sich aber mit wenigen
% Zusatz-Termen mit berücksichtigen).
%
% Der Parameter $\beta$ gibt an, wie viele enge Kontakte
% ein Individuum pro Tag durchschnittlich hat. (wobei nur ein Kontakt
% zwischen einem S- und einem I-Individuum zu einer Ansteckund führt; Kontakte vom Typ S-S,
% S-R, I-I und I-R bleiben folgenlos).
%
% Der Parameter $\gamma$ lässt sich einfacher in Form seines Reziprokwerts interpretieren
% als
%
% $\displaystyle \frac{1}{\gamma}$ Typische Zeitdauer der Infektion
%
% Maßgeblich für die Dynamik der Epidemie ist die Basisreproduktionszahl 
%
% $$R_0 = \frac{\beta}{\gamma}\qquad . $$
%
% Dieser Wert gibt an, wie viele weitere Personen eine infizierte Person in
% der Anfangsphase der Epidemie, (wenn in der allgemeinen Population noch
% kaum jemand immun ist), im Durchscnitt ansteckt.
%
% Die Urheber dieses Modells, der Biochemiker William Ogilvy Kermack und
% der Militärarzt Anderson Gray McKendrick, konnten trotz der Einfachheit
% des Modells gut die Daten einer Pestepidemie in Bombay 1905/06
% modellieren. Nur etwas über 100 Jahre später, 2014, publizierten Harko
% and Koautoren eine exakte analytische Lösung dieses
% Differentialgleichungssystems. In aller Regel wird dieses
% Differentialgleichungssystem aber  numerisch gelöst.
 
%% Modellparameter

beta = 0.4 ;
gamma = 0.04;
N = 1000000;

R0=beta/gamma

%% Die Funktion f definiert die Differentialgleichung
% Bezeichnung: y = [S;I;R], also y(1) steht für S, y(2) für I, y(3) für R.

f = @(t,y)[ 
    - beta*y(2)*y(1)/N;
      beta*y(2)*y(1)/N - gamma*y(2);
      gamma*y(2) ];
%% Aufruf des Lösers mit passenden Anfangsbedingungen und Plot
% Der DG-Löser ode45 ist ein Allzweck-Standard-Verfahren.

ode45(f,[0 100],[N-10000, 10000,0])

%%
% Beschriftung und Kosmetik
title 'SIR-Modell'
xlabel 'Zeit t (Tage)'
ylabel 'Individuen in S, I, R'