Come creare un grafico previsionale sulla diffusione del Covid-19 con Linux, Bash e Python

All’inizio dell’emergenza sanitaria da Coronavirus (nCoV-19), durante i primi giorni di lockdown, mi sono dato subito da fare mettendo in piedi una semplice webapp www.coronavirus-italy.it, per il monitoraggio, tramite API nazionali e internazionali, della pandemia.

n

Poi arricchendola di alcune funzionalità più specifiche, una di queste: le previsioni sulla diffusione in Italia del virus, in base ai dati forniti giorno per giorno dalla Protezione Civile.

n

In questo articolo voglio condividere la mia esperienza per mettere in piedi un sistema previsionale automatico sul nostro server Linux.

n

Iniziamo con la preparazione dell’ambiente.

n

n

Nel nostro server Linux dobbiamo accertarci di avere installato Python3 e questi moduli necessari: scipy, venv, pandas, sklearn, numpy, matplotlib.

n

n

Ora possiamo iniziare con la creazione del file Python che genererà i grafici previsionali con i dati della Protezione Civile.

n

Teniamo presente che i nuovi dati vengono diffusi ogni giorno intorno alle 18.00 circa, e comunque entro le 18.30.

n

n

Per fare i calcoli previsionali tramite modello matematico logistico, mi sono affidato a questo progetto: https://github.com/marcello-dev/coronavirus-forecast, facendone le opportune modifiche che condivido qui.

n

n

Possiamo inserire questo codice all’interno di un nuovo file di testo che chiameremo covid-forecast.py

n

import pandas as pdnimport numpy as npnfrom datetime import datetime, timedeltanfrom sklearn.metrics import mean_squared_errornfrom scipy.optimize import curve_fitnfrom scipy.optimize import fsolvenimport matplotlib.pyplot as pltnnurl = "https://raw.githubusercontent.com/pcm-dpc/COVID-19/master/dati-andamento-nazionale/dpc-covid19-ita-andamento-nazionale.csv"ndf = pd.read_csv(url)nndf = df.loc[:, ['data', 'totale_casi']]nnFMT = '{a99a1f178160d0ccabf421e3b6d1240d1fe1df588bebbd34e3de16e354822ce9}Y-{a99a1f178160d0ccabf421e3b6d1240d1fe1df588bebbd34e3de16e354822ce9}m-{a99a1f178160d0ccabf421e3b6d1240d1fe1df588bebbd34e3de16e354822ce9}dT{a99a1f178160d0ccabf421e3b6d1240d1fe1df588bebbd34e3de16e354822ce9}H:{a99a1f178160d0ccabf421e3b6d1240d1fe1df588bebbd34e3de16e354822ce9}M:{a99a1f178160d0ccabf421e3b6d1240d1fe1df588bebbd34e3de16e354822ce9}S'ndate = df['data']ndf['data'] = date.map(lambda x2: (datetime.strptime(x2, FMT) - datetime.strptime("2020-01-01T00:00:00", FMT)).days)nndef logistic_model(x3, a3, b3, c3):n    return c3/(1+np.exp(-(x3-b3)/a3))nnx = list(df.iloc[:, 0])ny = list(df.iloc[:, 1])nnfit = curve_fit(logistic_model, x, y, p0=[2, 100, 20000])nna = fit[0][0]nb = fit[0][1]nc = fit[0][2]nnfirst_january = datetime.strptime('2020/01/01', "{a99a1f178160d0ccabf421e3b6d1240d1fe1df588bebbd34e3de16e354822ce9}Y/{a99a1f178160d0ccabf421e3b6d1240d1fe1df588bebbd34e3de16e354822ce9}m/{a99a1f178160d0ccabf421e3b6d1240d1fe1df588bebbd34e3de16e354822ce9}d")na = fit[0][0]nb = fit[0][1]nc = fit[0][2]nnfirst_january = datetime.strptime('2020/01/01', "{a99a1f178160d0ccabf421e3b6d1240d1fe1df588bebbd34e3de16e354822ce9}Y/{a99a1f178160d0ccabf421e3b6d1240d1fe1df588bebbd34e3de16e354822ce9}m/{a99a1f178160d0ccabf421e3b6d1240d1fe1df588bebbd34e3de16e354822ce9}d")ninfection_peak_date = first_january + timedelta(days=int(b))nprint('{"status": "ok","previsioni": [{ "piccoContagio": "',datetime.strftime(infection_peak_date,"{a99a1f178160d0ccabf421e3b6d1240d1fe1df588bebbd34e3de16e354822ce9}d/{a99a1f178160d0ccabf421e3b6d1240d1fe1df588bebbd34e3de16e354822ce9}m/{a99a1f178160d0ccabf421e3b6d1240d1fe1df588bebbd34e3de16e354822ce9}Y"),'",')nerrors = [np.sqrt(fit[1][i][i]) for i in [0, 1, 2]]nprint('"totaleInfetti": "{}", "minInfetti": "{}", "maxInfetti": "{}",'.format(int(c), int(c-errors[2]),int(c+errors[2])))nnsol = int(fsolve(lambda x : logistic_model(x, a, b, c) - int(c),b))nninfection_end_date = first_january + timedelta(days=int(sol))nprint('"fineContagio": "',datetime.strftime(infection_end_date,"{a99a1f178160d0ccabf421e3b6d1240d1fe1df588bebbd34e3de16e354822ce9}d/{a99a1f178160d0ccabf421e3b6d1240d1fe1df588bebbd34e3de16e354822ce9}m/{a99a1f178160d0ccabf421e3b6d1240d1fe1df588bebbd34e3de16e354822ce9}Y"),'",')nndef add_real_data(df, label, color=None):n    x = df['data'].tolist()n    y = df['totale_casi'].tolist()n    plt.scatter(x, y, label="Dati reali (" + label + ")", c=color)nnpred_x = list(range(max(x), sol))nplt.rcParams['figure.figsize'] = [7, 7]nplt.rc('font', size=14)n# Real datanadd_real_data(df[-1:], "oggi")nadd_real_data(df[-2:-1], "ieri")nadd_real_data(df[:-2], "2 giorni fa")n# Predicted curve of todaynplt.plot(x+pred_x, [logistic_model(i,fit[0][0],fit[0][1],fit[0][2]) for i in x+pred_x], label="Previsione dati oggi")nn# Predicted curve of yesterdaynx = list(df[:-1].iloc[:, 0])ny = list(df[:-1].iloc[:, 1])npred_x = list(range(max(x), sol))nfit = curve_fit(logistic_model, x, y, p0=[2, 100, 20000])nplt.plot(x+pred_x, [logistic_model(i,fit[0][0],fit[0][1],fit[0][2]) for i in x+pred_x],n         label="Previsione dati ieri", dashes=[4, 4])nn# Predicted curve of 2 days ago curvenx = list(df[:-2].iloc[:, 0])ny = list(df[:-2].iloc[:, 1])npred_x = list(range(max(x), sol))nfit = curve_fit(logistic_model, x, y, p0=[2, 100, 20000])nplt.plot(x+pred_x, [logistic_model(i,fit[0][0],fit[0][1],fit[0][2]) for i in x+pred_x],n         label="Previsione dati 2 giorni fa",dashes=[8, 8])nntoday_date = datetime.today().strftime('{a99a1f178160d0ccabf421e3b6d1240d1fe1df588bebbd34e3de16e354822ce9}Y-{a99a1f178160d0ccabf421e3b6d1240d1fe1df588bebbd34e3de16e354822ce9}m-{a99a1f178160d0ccabf421e3b6d1240d1fe1df588bebbd34e3de16e354822ce9}d')nplt.title("Previsioni per casi accertati in Italia del " + today_date)nnplt.legend()nplt.xlabel("Giorni da 1 gennaio 2020")nplt.ylabel("Numero totale persone infette")nplt.ylim((min(y)*0.9,c*1.1))nnfilename = 'plot-' + today_date + '.png'nplt.savefig('plots/'+filename, bbox_inches="tight")ny_pred_logistic = [logistic_model(i,fit[0][0],fit[0][1],fit[0][2]) for i in x]nprint('"erroreModello": "', mean_squared_error(y,y_pred_logistic),'"}]}')n

n

n

Creiamo allo stesso livello del nostro file Python una cartella che chiameremo plots/.

n

n

A questo punto, prima di procedere facciamo il test, per verificare che tutti gli import e le istruzioni Python vadano a buon fine, digitiamo questi 3 comandi uno dopo l’altro, il tutto si dovrebbe eseguire senza errori.

n

python3 -m venv venvnnsource venv/bin/activatennpython3 main.py > plots/previsione.jsonn

n

n

Tutta la parte di codice di generazione del file Json assieme al grafico previsionale è una modifica che ho implementato io e non troverete nel progetto originale, ad ogni modo la condivido in quanto ritenuta utile a future implementazioni di ogni genere, rimane pur sempre un JSON aggiornato quotidianamente contenente dei dati.

n

n

Possiamo ora procedere all’automatizzazione del tutto semplicemente creando uno script Bash che chiameremo covidscript:

n

#!/bin/bashncd /home/Covid/coronavirus-forecastnpython3 -m venv venvnsource venv/bin/activatenpython3 main.py > plots/previsione.jsonn

n

Ovviamente ricordate di personalizzare l’indirizzo della cartella dove risiede il vostro applicativo Python per questo progetto, dichiarato nella seconda riga di questo script appena creato.

n

n

Se vogliamo inserire lo script in un Cron automatico ricordiamoci di non inserire alcuna estensione al file Bash e di renderlo eseguibile:

n

chmod +x covidscript

n

n

Poi effettuiamo un link simbolico dello script alla cartella di Cron:

n

ln -s /home/covid/covidscript /etc/cron.daily

n

n

A questo punto ogni giorno avremo sempre aggiornati il grafico previsionale e il file Json pronto per qualsiasi utilizzo.

n

n

Se vogliamo modificare l’orario dei Cron giornalieri possiamo inserire una riga all’interno del file /etc/crontab ricordandoci di inserire un orario non inferiore alle 18.30, perché altrimenti rischiamo ogni giorno di essere troppo in anticipo rispetto alla diffusione pubblica dei dati di quel giorno:

n

35 18   * * *   root    test -x /usr/sbin/anacron || ( cd / && run-parts --report /etc/cron.daily )