Estimating parameters of SIR model and problem with real-life data

73 Views Asked by At

I tried to make an SIR model based on real-world data. But, I ran into a snag when I'm trying to estimate the parameters of $\beta$ and $\gamma$. With equations:

$$ \begin{cases} \frac{dS(t)}{dt}=-\beta\frac{S(t)I(t)}{N}\\[6pt] \frac{dI(t)}{dt}=\beta\frac{S(t)I(t)}{N}- \gamma I(t) \\[6pt] \frac{dR(t)}{dt}=\gamma I(t) \end{cases} $$

So my data that I am using is divided into two the number of infected individuals and recovered individuals. I already tried to get the value of $(\beta - \gamma)$ by integrating the infected differential equation and yield:

$$I = I_0 e^{(\beta - \gamma)t}$$

I then linearize the equation using $ln()$ and do linear regression from my data of daily cases of infected individuals which yields an $m$ = $(\beta - \gamma)$ value of

$$ln(I(t))=(\beta - \gamma)t+ln(I_0)$$

$$m = (\beta - \gamma) = 2.808$$

Then, I estimated the value of $\gamma$ by using:

$$\gamma = \frac{R(t+1) - R(t)}{I(t)}$$

By evaluating the average $\gamma$ from my data I calculated a value of:

$$\gamma = 0.01224$$

Thus, $\beta$ should be

$$\beta = 2.808 + 0.01224 = 2.8202$$

Though when I model in MATLAB with my parameters of initial population $N = 100$ and initial infected $I = 1$, the value is very abnormal, with a simple ODE code:

clear all, close all, clc

A = 2.808 + 0.01224965083;
B = 0.01224965083;

tspan = [0:1:300];
y0 = [100 1 0];
[t,y] = ode45(@(t,y)odefcn(y,A,B),tspan,y0);

figure(1)
hold on
plot(t,y(:,1),'-o')
plot(t,y(:,2),'-*')
plot(t,y(:,3),'-+')

function dydt = odefcn(y,A,B)

    dydt = zeros(3,1);
    dydt(1) = -A*y(2)*y(1);
    dydt(2) = A*y(2)*y(1)-B*y(2);
    dydt(3) = B*y(2);
end

Question Why does the value of my susceptible and infected individuals immediately jump to max and min on the very first day, then it exponentially decreases? Is this due to the parameters or is it something else? You can refer to image below.

My reasoning at first would most likely be because of the $\beta$ value is too high, as when I change the $\beta$ value to be smaller than the value of $\gamma$ the graph works it does not jump from min to max but increases and decreases exponentially. But, I don't know why this occur, any ideas?

enter image description here

Additionally, @user10354138 mentioned that I miss dividing $\beta / N$ but still though there isn't much of a difference. If only jumps to min and max at x = 3 not the first.

enter image description here