快捷導(dǎo)航

Python實現(xiàn)新型冠狀病毒傳播模型及預(yù)測代碼實例

更新時間：2020年02月05日 14:39:58 作者：jiazhen.

在本篇文章里小編給大家整理的是關(guān)于Python實現(xiàn)新型冠狀病毒傳播模型及預(yù)測代碼內(nèi)容，有興趣的朋友們可以學(xué)習(xí)下。

1.傳染及發(fā)病過程

一個健康人感染病毒后進(jìn)入潛伏期（時間長度為Q天），潛伏期之后進(jìn)入發(fā)病期（時間長度為D天），發(fā)病期之后該患者有三個可能去向，分別是自愈、接收隔離、死亡。

2.模型假設(shè)

潛伏期Q=7天，根據(jù)報道潛伏期為2~14天，取中間值；發(fā)病期D=10天，根據(jù)文獻(xiàn)報告，WHO認(rèn)定SARS發(fā)病期為10天，假設(shè)武漢肺炎與此相同；潛伏期的患者不具有將病毒傳染給他人的能力；發(fā)病期的患者具有將病毒傳染給他人的能力；患者在發(fā)病期之后不再具有將病毒傳染他人的能力；假設(shè)處于發(fā)病期的患者平均每天密切接觸1人，致使該人患病的概率為γ最初只有一個人類感染者；病情自然發(fā)展，沒有外部干擾。

3.模型公式

$P_{N} = P_{N-1} + P_{N-Q} \cdot \gamma + P_{N-Q-D} \cdot \gamma$

：人類感染該病毒的天數(shù)

$P_{N}$ : 第N天感染該病毒并且處于發(fā)病期的患者數(shù)量

$P_{N-Q} \cdot \gamma$ ： Q天前新被感染患者，當(dāng)日進(jìn)入發(fā)病期的數(shù)量

$P_{N-Q-D} \cdot \gamma$ ：當(dāng)日發(fā)病期滿，不再具有傳染能力的患者數(shù)量

4.模型初始值

根據(jù)假設(shè)，最初只有一個人類感染者，所以：

$P_{1} = P_{2} = P_{3} ... = P_{7} = 1$

5.實際疫情數(shù)據(jù)

人類感染病毒且發(fā)病的初始日期：根據(jù)財新網(wǎng)的報道，官方通報首例不明原因肺炎是在12月8日，考慮到確診之前肯定已經(jīng)嘗試過各種治療方案無效后認(rèn)定為不明原因肺炎，所以有理由認(rèn)為該名患者在12月8日已經(jīng)處于發(fā)病期末端，根據(jù)假設(shè)發(fā)病期為10天，所以可以假設(shè)該名患者在11月29號發(fā)病，即N=1對應(yīng)11月29日。

日新增發(fā)病數(shù) $\Delta P_{n}$

近期疫情防控大事記：

從以上信息可以判斷核酸檢測試劑是在1月16日、17日大幅使用的，18日、19日確診大量病例，因此1月20日之前的確診病例數(shù)對模型參考意義不大。1月20日之后，可以認(rèn)為新發(fā)病例，發(fā)病即檢測。

截至1月20日24時，國家衛(wèi)健委公告累計確診病例291

截至1月21日24時，國家衛(wèi)健委公告累計確診病例440

截至1月22日24時，國家衛(wèi)健委公告累計確診病例571

截至1月23日24時，國家衛(wèi)健委公告累計確診病例830

截至1月24日24時，國家衛(wèi)健委公告累計確診病例1287

所以：

1月21日新增確診病例：440-291=1491月22日新增確診病例：571-440=1311月23日新增確診病例：830-571=2591月24日新增確診病例：1287-830=457

考慮到：

1月20日之后發(fā)病即檢測確診檢測用時2日國外新增病例在個位數(shù)，且不能保證發(fā)病即檢測等因素，暫不考慮新增病例應(yīng)該遞增，所以1月22日新增數(shù)據(jù)異常，舍去

使用上溯每日新增數(shù)據(jù),同時考慮到確診需要2天，可以得到:

$\Delta P_{52} = 149$ , $\Delta P_{54} = 259$ , $\Delta P_{55} = 457$

（ $P_{1}$ 對應(yīng)11月29日）

6.擬合確定 $\gamma$

根據(jù)近期每日新增數(shù)據(jù)、模型初始值及模型公式，用最小二乘擬合得到 $\gamma = 0.55$

7.預(yù)測

患者數(shù)量：根據(jù)上文確定的模型及參數(shù)，從11月29日（N=1）至1月27日（N=60）人群中累計處于發(fā)病期的人數(shù)如下圖所示：

根據(jù)模型，近期人群中患者數(shù)量計算如下：

人群中感染了病毒并處于發(fā)病期的患者數(shù)量，注意Pn一般需要延期2日才能確診

每日新增患者數(shù)量：

根據(jù)模型，近期人群中每日新增患者數(shù)量計算如下：

注意：圖中是人群中新增發(fā)病患者數(shù)量，可與晚2日的政府發(fā)布新增數(shù)量進(jìn)行對比。即22日新增患者數(shù)量可與24日政府發(fā)布的新增病例進(jìn)行對比。截至目前模型計算22日新增為369人，政府公布的24日新增病例457人

根據(jù)以上模型預(yù)計未來幾日的情況如下：

1月25日將新確診432例，人群中發(fā)病患者為4068人；

1月26日將新確診505例，人群中發(fā)病患者為4759人；

1月27日將新確診590例，人群中發(fā)病患者為5568人；

1月28日將新確診691例，人群中發(fā)病患者為6514人；

1月29日將新確診809例，人群中發(fā)病患者為7621人。

由于1月20日之后采取的各種措施，將導(dǎo)致發(fā)病期D下降，感染概率 $\gamma$ 下降，1月29日之后日均增長勢頭會減弱。

8.模型代碼

import numpy as np
import matplotlib.pyplot as plt
 
gamma = 0.55
Q = 7
D = 10
P = np.zeros(300, dtype=np.float)
Psum = np.zeros(300, dtype=np.float)
for i in range(Q):
  P[i] = 1
for j in range(300-Q):
  P[j+Q] = P[j+Q-1]+P[j]*gamma
  if j+Q-D-Q >= 0:
    P[j+Q] -= P[j+Q-D-Q]*gamma
  if j+Q == D:
    P[j+Q] -= 1
plt.xlabel("N")
plt.ylabel("PN")
plt.rcParams['font.sans-serif'] = ['SimHei']
plt.rcParams['axes.unicode_minus'] = False
plt.plot(range(1, 61), P[0:60])
plt.grid()
plt.show()

以上就是腳本之家小編給大家整理的全部內(nèi)容，如果大家有任何補(bǔ)充可以聯(lián)系小編。

您可能感興趣的文章: