Skip to content

Instantly share code, notes, and snippets.

@arivero
Last active April 28, 2020 15:10
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save arivero/c629d41e95717d9eab98d68225880ecd to your computer and use it in GitHub Desktop.
Save arivero/c629d41e95717d9eab98d68225880ecd to your computer and use it in GitHub Desktop.
descarga y ajusta los datos de fallecimientos oficiales, para estimar los datos finales
#!/usr/bin/env python
# coding: utf-8
# In[1]:
from scipy.optimize import leastsq,curve_fit
import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
def logcorrected(x,s0,k0,dk,mu,i0,v):
k = k0 - dk * (x > v)
k1 = k0 - dk
log0=np.log(s0-i0)
log1= np.log (s0-v) - k1/k0 * (np.log (s0-v) - log0 )
l = log0 * ( x <= v) + log1 * (x > v)
return (s0-x) * ( k * x /s0 + mu * (np.log (s0-x) - l ))
dataframe=pd.read_csv("https://covid19.isciii.es/resources/serie_historica_acumulados.csv", skipfooter=6, engine="python", encoding='iso8859')
dataframe['date']=pd.to_datetime(dataframe['FECHA'],format='%d/%m/%Y')
assert(dataframe['CCAA'][len(dataframe)-1]=='RI')
data=dataframe.drop(columns='FECHA').set_index(['CCAA','date']).sort_index().unstack(level=0).fillna(0)
f=data['Fallecidos'].sum(axis=1)
# In[2]:
fit,cov=curve_fit(logcorrected, f[~np.isnan(f.diff())],f.diff()[~np.isnan(f.diff())],p0=[600000,1,0.1,0,100,8000.5],
bounds=np.array([[0,np.inf],[0,np.inf],[0,10],[0,10],[1,10000],[2000,30001]]).T,
method='trf')
s0,k0,k1,mu,i0,corte = fit
print(fit)
print(s0,s0/46000000*100)
print(k0,(k0-k1))
print(k0/mu,(k0-k1)/mu)
print(i0,corte)
print(cov)
xp=np.linspace(0, 300000, 300000)
plt.xlim(0,35000)
plt.ylim(0,1200)
plt.plot(f[~np.isnan(f.diff())],f.diff()[~np.isnan(f.diff())],'.',xp,logcorrected(xp,*fit))
# In[3]:
fit2=fit.copy()
fit2[2]=0
xp=np.linspace(0, 300000, 3000)
plt.xlim(0,80000)
plt.ylim(0,2500)
plt.xlabel("total de fallecidos")
plt.ylabel("fallecidos diarios")
plt.plot(f[~np.isnan(f.diff())],f.diff()[~np.isnan(f.diff())],'.',xp,logcorrected(xp,*fit2),xp,logcorrected(xp,*fit),'-')
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment