deephub翻譯組:劉欣然
當今世界正在與一個新的敵人作鬥爭,那就是Covid-19病毒。
該病毒自首次在中國出現以來,在世界範圍內迅速傳播。不幸的是,義大利的Covid-19感染人數是歐洲最高的,為19人。我們是西方世界第一個面對這個新敵人的國家,我們每天都在與這種病毒帶來的經濟和社會影響作鬥爭。
在本文中,我將用Python向您展示感染增長的簡單數學分析和兩個模型,以更好地理解感染的演變。
義大利民防部門每天都會更新感染者的累積資料。這些資料在GitHub上作為開放資料公開在Github這裡:
https://raw.githubusercontent.com/pcm-dpc/COVID-19/master/dati-andamento-nazionale/dpc-covid19-ita-andamento-nazionale.csv我的目標是建立迄今為止受感染人數(即實際感染人數加上已感染人數)的時間序列模型。這些模型具有引數,這些引數將透過曲線擬合進行估算。
我們用Python來做。
首先,讓我們匯入一些庫。
importnumpy as np
from datetime import datetime,timedelta
from sklearn.metrics import mean_squared_error
from scipy.optimize import curve_fit
from scipy.optimize import fsolve
import matplotlib.pyplot as plt
%matplotlib inline
現在,讓我們看一下原始資料。
url = https://raw.githubusercontent.com/pcm-dpc/COVID-19/master/dati-andamento-nazionale/dpc-covid19-ita-andamento-nazionale.csv
df =pd.read_csv(url)
我們需要的列是" totale_casi ",它包含到目前為止的累計感染人數。
這是原始資料。現在,讓我們為分析做準備。
首先,我們需要將日期改為數字。我們將從一月一日起開始算。
df =df.loc[:,["data","totale_casi"]]
FMT ="%Y-%m-%d %H:%M:%S"
date =df["data"]
df["data"]= date.map(lambda x : (datetime.strptime(x, FMT) -datetime.strptime("2020-01-01 00:00:00", FMT)).days )
現在,我們可以分析要參加測試的兩個模型,分別是邏輯函式(logistic function)和指數函式(exponential function)。
每個模型都有三個引數,這些引數將透過對歷史資料進行曲線擬合計算來估計。
logistic模型被廣泛用於描述人口的增長。感染可以被描述為病原體數量的增長,因此使用logistic模型似乎是合理的。
這個公式在資料科學家中非常有名,因為它被用於邏輯迴歸分類器,並且是神經網路的一個啟用函式。
logistic函式最一般的表示式為:
在這個公式中,我們有變數x(它是時間)和三個引數:a,b,c。
•a為感染速度•b為感染髮生最多的一天•c是在感染結束時記錄的感染者總數
在高時間值時,被感染的人數越來越接近c值,也就是我們說感染已經結束的時間點。這個函式在b點也有一個拐點,也就是一階導數開始下降的點(即感染開始減弱並下降的峰值)。
讓我們在Python中定義模型:
def logistic_model(x,a,b,c):
return c/(1+np.exp(-(x-b)/a))
我們可以使用scipy庫中的curve_fit函式從原始資料開始估計引數值和錯誤。
x =list(df.iloc[:,0])
y =list(df.iloc[:,1])fit = curve_fit(logistic_model,x,y,p0=[2,100,20000])
這裡是一些值:
· a: 3.54
· b: 68.00
· c: 15968.38
該函式也返回協方差矩陣,其對角值是引數的方差。取它們的平方根,我們就能計算出標準誤差。
errors= [np.sqrt(fit[1][i][i]) for i in [0,1,2]]
· a的標準誤差:0.24
· b的標準誤差:1.53
· c的標準誤差:4174.69
這些數字給了我們許多有用的見解。
預計感染人數在感染結束時為15968+/-4174。
感染高峰預計在2020年3月9日左右。
預期的感染結束日期可以計算為受感染者累計計數四捨五入約等於到最接近整數的c引數的那一天。
我們可以使用scipy的fsolve函式來計算出定義感染結束日的方程的根。
求解出來時間是2020年4月15日。
logistic模型描述了未來將會停止的感染增長,而指數模型描述了不可阻擋的感染增長。例如,如果一個病人每天感染2個病人,1天后我們會有2個感染,2天后4個,3天后8個,等等。
最通用的指數函式是:
變數x是時間,我們仍然有引數a, b, c,但是它的意義不同於logistic函式引數。
讓我們在Python中定義這個函式,並執行與logistic增長相同的曲線擬合過程。
def exponential_model(x,a,b,c):
return a*np.exp(b*(x-c))exp_fit =curve_fit(exponential_model,x,y,p0=[1,1,1])
引數及其標準差為:
· a: 0.0019 +/- 64.6796
· b: 0.2278 +/- 0.0073
· c: 0.50 +/- 144254.77
我們現在有了所有必要的資料來視覺化我們的結果。
pred_x= list(range(max(x),sol))
plt.rcParams["figure.figsize"]= [7, 7]
plt.rc("font",size=14)
## Realdata
plt.scatter(x,y,label="Real data",color="red")
#Predicted logistic curve
plt.plot(x+pred_x,[logistic_model(i,fit[0][0],fit[0][1],fit[0][2]) for i inx+pred_x], label="Logistic model" )
#Predicted exponential curve
plt.plot(x+pred_x,[exponential_model(i,exp_fit[0][0],exp_fit[0][1],exp_fit[0][2])for i in x+pred_x], label="Exponential model" )
plt.legend()
plt.xlabel("Days since 1 January 2020")
plt.ylabel("Total number of infected people")
plt.ylim((min(y)*0.9,c*1.1))plt.show()
這兩條理論曲線似乎都很接近實驗趨勢。哪一個更好?讓我們看一下殘差(residuals.)。
殘差是指各實驗點與相應理論點的差值。我們可以透過分析兩種模型的殘差來驗證最佳擬合曲線。在第一次近似中,理論和實驗資料的均方誤差越小,擬合越好。
y_pred_logistic=[logistic_model(i,fit[0][0],fit[0][1],fit[0][2])
for iin x]y_pred_exp = [exponential_model(i,exp_fit[0][0], exp_fit[0][1], exp_fit[0][2]) for iin x]
mean_squared_error(y,y_pred_logistic)
mean_squared_error(y,y_pred_exp)
Logistic模型MSE(均方誤差):8254.07
指數模型MSE: 16219.82
殘差分析似乎指向邏輯模型。很可能是因為感染應該會在將來的某一天結束;即使每個人都會被感染,他們也會適當地發展出免疫防禦措施以避免再次感染。只要病毒沒有發生太多變異(例如,流感病毒),這就是正確的模型。
但是有些事情仍然讓我擔心。自感染開始以來,我每天都在擬合logistic曲線,而且每天都有不同的引數值。感染的人數最終會增加,最大感染日通常是當天或第二天(與該引數的1天標準誤差是一致的)。
這就是為什麼我認為,儘管邏輯模型似乎是最合理的模型,但是曲線的形狀可能會由於新的感染熱點,政府約束感染的行動措施等外在影響而發生變化。
因此,我認為這個模型的預測只有在感染高峰期之後的幾周內才會開始有用。
有。中國疫情發生後也有國外專家做過數學模型,到中國透過努力,成功的大臉了外國專家,數學模型更使用於各國的疫情模型計算。
deephub翻譯組:劉欣然
當今世界正在與一個新的敵人作鬥爭,那就是Covid-19病毒。
該病毒自首次在中國出現以來,在世界範圍內迅速傳播。不幸的是,義大利的Covid-19感染人數是歐洲最高的,為19人。我們是西方世界第一個面對這個新敵人的國家,我們每天都在與這種病毒帶來的經濟和社會影響作鬥爭。
在本文中,我將用Python向您展示感染增長的簡單數學分析和兩個模型,以更好地理解感染的演變。
資料收集(Data collection)義大利民防部門每天都會更新感染者的累積資料。這些資料在GitHub上作為開放資料公開在Github這裡:
https://raw.githubusercontent.com/pcm-dpc/COVID-19/master/dati-andamento-nazionale/dpc-covid19-ita-andamento-nazionale.csv我的目標是建立迄今為止受感染人數(即實際感染人數加上已感染人數)的時間序列模型。這些模型具有引數,這些引數將透過曲線擬合進行估算。
我們用Python來做。
首先,讓我們匯入一些庫。
importpandas as pdimportnumpy as np
from datetime import datetime,timedelta
from sklearn.metrics import mean_squared_error
from scipy.optimize import curve_fit
from scipy.optimize import fsolve
import matplotlib.pyplot as plt
%matplotlib inline
現在,讓我們看一下原始資料。
url = https://raw.githubusercontent.com/pcm-dpc/COVID-19/master/dati-andamento-nazionale/dpc-covid19-ita-andamento-nazionale.csv
df =pd.read_csv(url)
我們需要的列是" totale_casi ",它包含到目前為止的累計感染人數。
這是原始資料。現在,讓我們為分析做準備。
資料準備(Data preparation)首先,我們需要將日期改為數字。我們將從一月一日起開始算。
df =df.loc[:,["data","totale_casi"]]
FMT ="%Y-%m-%d %H:%M:%S"
date =df["data"]
df["data"]= date.map(lambda x : (datetime.strptime(x, FMT) -datetime.strptime("2020-01-01 00:00:00", FMT)).days )
現在,我們可以分析要參加測試的兩個模型,分別是邏輯函式(logistic function)和指數函式(exponential function)。
每個模型都有三個引數,這些引數將透過對歷史資料進行曲線擬合計算來估計。
logistic模型(The logistic model)logistic模型被廣泛用於描述人口的增長。感染可以被描述為病原體數量的增長,因此使用logistic模型似乎是合理的。
這個公式在資料科學家中非常有名,因為它被用於邏輯迴歸分類器,並且是神經網路的一個啟用函式。
logistic函式最一般的表示式為:
在這個公式中,我們有變數x(它是時間)和三個引數:a,b,c。
•a為感染速度•b為感染髮生最多的一天•c是在感染結束時記錄的感染者總數
在高時間值時,被感染的人數越來越接近c值,也就是我們說感染已經結束的時間點。這個函式在b點也有一個拐點,也就是一階導數開始下降的點(即感染開始減弱並下降的峰值)。
讓我們在Python中定義模型:
def logistic_model(x,a,b,c):
return c/(1+np.exp(-(x-b)/a))
我們可以使用scipy庫中的curve_fit函式從原始資料開始估計引數值和錯誤。
x =list(df.iloc[:,0])
y =list(df.iloc[:,1])fit = curve_fit(logistic_model,x,y,p0=[2,100,20000])
這裡是一些值:
· a: 3.54
· b: 68.00
· c: 15968.38
該函式也返回協方差矩陣,其對角值是引數的方差。取它們的平方根,我們就能計算出標準誤差。
errors= [np.sqrt(fit[1][i][i]) for i in [0,1,2]]
· a的標準誤差:0.24
· b的標準誤差:1.53
· c的標準誤差:4174.69
這些數字給了我們許多有用的見解。
預計感染人數在感染結束時為15968+/-4174。
感染高峰預計在2020年3月9日左右。
預期的感染結束日期可以計算為受感染者累計計數四捨五入約等於到最接近整數的c引數的那一天。
我們可以使用scipy的fsolve函式來計算出定義感染結束日的方程的根。
sol =int(fsolve(lambda x : logistic_model(x,a,b,c) - int(c),b))求解出來時間是2020年4月15日。
指數模型(Exponential model)logistic模型描述了未來將會停止的感染增長,而指數模型描述了不可阻擋的感染增長。例如,如果一個病人每天感染2個病人,1天后我們會有2個感染,2天后4個,3天后8個,等等。
最通用的指數函式是:
變數x是時間,我們仍然有引數a, b, c,但是它的意義不同於logistic函式引數。
讓我們在Python中定義這個函式,並執行與logistic增長相同的曲線擬合過程。
def exponential_model(x,a,b,c):
return a*np.exp(b*(x-c))exp_fit =curve_fit(exponential_model,x,y,p0=[1,1,1])
引數及其標準差為:
· a: 0.0019 +/- 64.6796
· b: 0.2278 +/- 0.0073
· c: 0.50 +/- 144254.77
畫圖我們現在有了所有必要的資料來視覺化我們的結果。
pred_x= list(range(max(x),sol))
plt.rcParams["figure.figsize"]= [7, 7]
plt.rc("font",size=14)
## Realdata
plt.scatter(x,y,label="Real data",color="red")
#Predicted logistic curve
plt.plot(x+pred_x,[logistic_model(i,fit[0][0],fit[0][1],fit[0][2]) for i inx+pred_x], label="Logistic model" )
#Predicted exponential curve
plt.plot(x+pred_x,[exponential_model(i,exp_fit[0][0],exp_fit[0][1],exp_fit[0][2])for i in x+pred_x], label="Exponential model" )
plt.legend()
plt.xlabel("Days since 1 January 2020")
plt.ylabel("Total number of infected people")
plt.ylim((min(y)*0.9,c*1.1))plt.show()
這兩條理論曲線似乎都很接近實驗趨勢。哪一個更好?讓我們看一下殘差(residuals.)。
殘差分析殘差是指各實驗點與相應理論點的差值。我們可以透過分析兩種模型的殘差來驗證最佳擬合曲線。在第一次近似中,理論和實驗資料的均方誤差越小,擬合越好。
y_pred_logistic=[logistic_model(i,fit[0][0],fit[0][1],fit[0][2])
for iin x]y_pred_exp = [exponential_model(i,exp_fit[0][0], exp_fit[0][1], exp_fit[0][2]) for iin x]
mean_squared_error(y,y_pred_logistic)
mean_squared_error(y,y_pred_exp)
Logistic模型MSE(均方誤差):8254.07
指數模型MSE: 16219.82
哪個是正確的模型?殘差分析似乎指向邏輯模型。很可能是因為感染應該會在將來的某一天結束;即使每個人都會被感染,他們也會適當地發展出免疫防禦措施以避免再次感染。只要病毒沒有發生太多變異(例如,流感病毒),這就是正確的模型。
但是有些事情仍然讓我擔心。自感染開始以來,我每天都在擬合logistic曲線,而且每天都有不同的引數值。感染的人數最終會增加,最大感染日通常是當天或第二天(與該引數的1天標準誤差是一致的)。
這就是為什麼我認為,儘管邏輯模型似乎是最合理的模型,但是曲線的形狀可能會由於新的感染熱點,政府約束感染的行動措施等外在影響而發生變化。
因此,我認為這個模型的預測只有在感染高峰期之後的幾周內才會開始有用。