<jupyter_start><jupyter_text>SIR modelling<jupyter_code>import pandas as pd
import numpy as np

from datetime import datetime
import pandas as pd 

from scipy import optimize
from scipy.optimize import curve_fit
from scipy import integrate

%matplotlib inline
import matplotlib as mpl
import matplotlib.pyplot as plt

import seaborn as sns


sns.set(style="darkgrid")

mpl.rcParams['figure.figsize'] = (16, 9)
pd.set_option('display.max_rows', 500)


df_analyse=pd.read_csv('../data/processed/COVID_small_flat_table.csv',sep=';')  
df_analyse.sort_values('date',ascending=True).head()
country_list=['Italy','US','India','Germany','China']
N0=1000000 #max susceptible population
beta=0.4   # infection spread dynamics
gamma=0.1  # recovery rate


# condition I0+S0+R0=N0
I0=df_analyse.Germany[35]
S0=N0-I0
R0=0
def SIR_model(SIR,beta,gamma):
    ''' Simple SIR model
        S: susceptible population
        I: infected people
        R: recovered people
        beta: 
        
        overall condition is that the sum of changes (differnces) sum up to 0
        dS+dI+dR=0
        S+I+R= N (constant size of population)
    
    '''
    
    S,I,R=SIR
    dS_dt=-beta*S*I/N0          #S*I is the 
    dI_dt=beta*S*I/N0-gamma*I
    dR_dt=gamma*I
    return([dS_dt,dI_dt,dR_dt])<jupyter_output><empty_output><jupyter_text>Simulative approach to curves<jupyter_code>SIR=np.array([S0,I0,R0])
propagation_rates=pd.DataFrame(columns={'susceptible':S0,
                                        'infected':I0,
                                        'recoverd':R0})



for each_t in np.arange(100):
   
    delta_vec=SIR_model(SIR,beta,gamma)
   
    SIR=SIR+delta_vec
    
    propagation_rates=propagation_rates.append({'susceptible':SIR[0],
                                                'infected':SIR[1],
                                                'recovered':SIR[2]}, ignore_index=True)
fig, ax1 = plt.subplots(1, 1)

ax1.plot(propagation_rates.index,propagation_rates.infected,label='infected',color='k')
ax1.plot(propagation_rates.index,propagation_rates.recovered,label='recovered',color='r')
ax1.plot(propagation_rates.index,propagation_rates.susceptible,label='susceptible')

ax1.set_ylim(10, 1000000)
ax1.set_yscale('linear')
ax1.set_title('Szenario SIR simulations  (demonstration purposes only)',size=16)
ax1.set_xlabel('time in days',size=16)
ax1.legend(loc='best',
           prop={'size': 16});<jupyter_output><empty_output><jupyter_text>Fitting the Suitable parameters of SIR model<jupyter_code>for each in df_analyse:
    fitdata=np.array(df_analyse[each]) #individual country
fitdata
# N0=5000000
# t=np.arange(len(fitdata))

# for each in countrylist['India']:

#     #df_plot=df_analyseSIR[df_analyseSIR['country']==each]
#     #df_plot=df_plot[['state','country','confirmed','date']].groupby(['country','date']).agg(np.mean).reset_index()
    
#     fitdata=np.array(df_analyseSIR[each])#.reshape(-1,1)  
# #t=np.arange(len(fitdata))

# #N0=1000000
#     I0=fitdata[0]
#     S0=N0-I0
#     R0=0

#     popt=[0.4,0.1]
    
#     popt,pcov = curve_fit(fit_odeint,t,fitdata,maxfev=15000) #,bounds=(1,len(t)))

#     perr = np.sqrt(np.diag(pcov))
#     #print(popt)
#     fitted=fit_odeint(t, *popt)#.reshape(-1,1)
#     print('fitted', fitted)
#     df_analyseSIR[each + '_SIR']=fitted

#     #fitted=fit_odeint(t, *popt)
#     #df_analyseSIR[each + '_sir']
    
t=np.arange(len(fitdata))
I0=fitdata[0]
S0=1000000
R0=0
beta
def SIR_model_t(SIR,t,beta,gamma):
    ''' Simple SIR model
        S: susceptible population
        t: time step, mandatory for integral.odeint
        I: infected people
        R: recovered people
        beta: 
        
        overall condition is that the sum of changes (differnces) sum up to 0
        dS+dI+dR=0
        S+I+R= N (constant size of population)
    
    '''
    
    S,I,R=SIR
    dS_dt=-beta*S*I/N0          #S*I is the 
    dI_dt=beta*S*I/N0-gamma*I
    dR_dt=gamma*I
    return dS_dt,dI_dt,dR_dt
#defining the fit function

def fit_odeint(x, beta, gamma):
    return integrate.odeint(SIR_model_t,(S0,I0,R0),t,args=(beta,gamma)) [:,1] #this is only for infection(dI)
#example for the differential equation
popt=[0.4,0.1]
fit_odeint(t, *popt)
popt, pcov = optimize.curve_fit(fit_odeint,t,fitdata)
perr = np.sqrt(np.diag(pcov))
    
print('standard deviation errors : ',str(perr), ' start infect:',fitdata[0])
print("Optimal parameters: beta =", popt[0], " and gamma = ", popt[1])
fitted=fit_odeint(t, *popt)
plt.semilogy(t, fitdata, 'o')
plt.semilogy(t, fitted)
plt.title("Fit of SIR model for Germany cases")
plt.ylabel("Population infected")
plt.xlabel("Days")
plt.show()
print("Optimal parameters: beta =", popt[0], " and gamma = ", popt[1])
print("Basic Reproduction Number R0 " , popt[0]/ popt[1])
print("This ratio is derived as the expected number of new infections (these new infections are sometimes called secondary infections from a single infection in a population where all subjects are susceptible. @wiki")<jupyter_output><empty_output><jupyter_text>Dynamic beta in SIR(Infection rate)<jupyter_code>df_analyseSIR=pd.read_csv('../data/processed/COVID_full_flat_table.csv',sep=';') 
df_analyseSIR.sort_values('date',ascending=True).iloc[35:].head()

t_initial=24
t_intro =14
t_hold  =21
t_relax =21

beta_max=0.4
beta_min=0.11
gamma=0.1
pd_beta=np.concatenate((np.array(t_initial*[beta_max]),
                       np.linspace(beta_max,beta_min,t_intro),
                       np.array(t_hold*[beta_min]),
                       np.linspace(beta_min,beta_max,t_relax),
                       ))
pd_beta
SIR=np.array([S0,I0,R0])
propagation_rates=pd.DataFrame(columns={'susceptible':S0,
                                        'infected':I0,
                                        'recoverd':R0})



for each_beta in pd_beta:
   
    new_delta_vec=SIR_model(SIR,each_beta,gamma)
   
    SIR=SIR+new_delta_vec
    
    propagation_rates=propagation_rates.append({'susceptible':SIR[0],
                                                'infected':SIR[1],
                                                'recovered':SIR[2]}, ignore_index=True)
fig, ax1 = plt.subplots(1, 1)

ax1.plot(propagation_rates.index,propagation_rates.infected,label='infected',linewidth=4)

t_phases=np.array([t_initial,t_intro,t_hold,t_relax]).cumsum()
ax1.bar(np.arange(len(fitdata)),fitdata, width=0.8,label=' current infected Germany',color='r')
ax1.axvspan(0,t_phases[0], facecolor='b', alpha=0.2,label='no measures')
ax1.axvspan(t_phases[0],t_phases[1], facecolor='b', alpha=0.3,label='hard measures introduced')
ax1.axvspan(t_phases[1],t_phases[2], facecolor='b', alpha=0.4,label='hold measures')
ax1.axvspan(t_phases[2],t_phases[3], facecolor='b', alpha=0.5,label='relax measures')
ax1.axvspan(t_phases[3],len(propagation_rates.infected), facecolor='b', alpha=0.6,label='repead hard measures')

ax1.set_ylim(10, 1.5*max(propagation_rates.infected))
ax1.set_yscale('log')
ax1.set_title('Scenario SIR simulations  (demonstration purposes only)',size=16)
ax1.set_xlabel('time in days',size=16)
ax1.legend(loc='upper right', prop={'size': 16});<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.14401907165895386, 'fraction_numerical': 0.03000981629504978, 'mean_word_length': 3.8026936026936027, 'pattern_counts': {'":': 0, '<': 17, '<?xml version=': 0, '>': 17, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 12, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '4', 'dir': 'jupyter-structured-clean-dedup', 'id': '650209', 'license': 'FTL', 'n_tokens_mistral': 2816, 'n_tokens_neox': 2643, 'n_words': 436, 'path': 'notebooks/Step9_ SIR Modelling.ipynb', 'repo_name': 'sindunarra/COVID19-prototype', 'size': '7131'}