感染病の数学予測モデル　(SIRモデル）：事例紹介(1)

16380 ワード
scipy コロナウイルス COVID-19 Python SIRモデル Python テキストリンク

1. 初めに

前回の記事に、感染病の数学予測モデルのSIRモデルを紹介した。

タイトル：感染病の数学予測モデルの紹介 (SIRモデル)
https://qiita.com/kotai2003/items/3078f4095c3e94e5325c

今回は、実際に発生した感染事例を挙げ、このSIRモデルの予測性能を確認することにする。

2. イギリスのある寄宿学校でのインフルエンザ感染　(1978年)

1978年1月、イギリス北部にある全寮制の寄宿学校に、インフルエンザが発生した。全校生の10歳から18歳の763人の少年の中、30人を除く全員にインフルエンザ感染が確認された。その中,512人は1月22日から2月4日まで続いた流行の間に、病状にあった。また、感染は1人の少年から始まった。


https://www.ncbi.nlm.nih.gov/pmc/articles/PMC1603269/pdf/brmedj00115-0064.pdf

表に、図から推定された、病状にあった(Infectious感染者)人数の値を示す。

Day Number of Infectious
1 3
2 8
3 28
4 75
5 221
6 291
7 255
8 235
9 190
10 125
11 70
12 28
13 12
14 5

図にこれらのデータをプロットする。感染者(Infectious)の数は、6日目にピークに達し、その値は255名であった。

3. パラメータ推定

ここでSIRモデルの微分方程式を再確認する。このモデルを解くために、各変数の初期条件を知る必要がある。今回の感染は、763人の全校生の中、感染者が1人から始まったことなので、S0 = 762, I0, R0=0 である。


\begin{align}

\frac{dS}{dt} &= -\beta SI \\
\frac{dI}{dt} &=  \beta SI -\gamma I \\
\frac{dR}{dt} &=  \gamma I \\

\end{align} 

\begin{align}
S &: 感染可能者 \quad \text{(Susceptible)}　\\
I &: 感染者 \quad \text{(Infectious)}　\\
R &: 感染後死亡者、もしくは免疫を獲得した者 \quad \text{(Removed)}　\\
\beta &: 感染率\quad \text{(The infectious rate)} \quad [1/day] \\
\gamma &:除去率\quad \text{(The Recovery rate)} \quad [1/day]　\\
\end{align}

問題は、感染率と除去率である。一般的に、SIRモデルが必要となる時期は、感染が発生した初期である。初期の感染者のデータから、素早く感染率と除去率のパラメータを推定し、SIRモデルの計算で、この感染が拡大するか収束するかを確認する。もし、SIRモデルで、感染が拡大することが分かったら（R0>1)、この拡大を防ぐために感染者の隔離や人の移動制限などの保健政策の優先順位を決めるプロセスに入る。

初期の感染の様子から、1人の病気の少年から、1日後にさらに2人の少年に感染することが分かった。大まかな概算で、感染可能者は、1日当り2人のペースで減少すると考えられる。

\frac{dS}{dt} \simeq -2 \\

初期条件として、S_0 = 762, I_0=1がわかってるので、下記の式により感染率を計算する。

\beta = 
\frac{-dS/dt}{S_0I_0}\simeq\frac{2}{762\times  1}=0.0026

次は除去率の推定である。除去率は、感染者が死亡及び免疫により、感染者の数から減る割合である。当時の報告によると、入院した少年たちは、1~2日以内に回復することが分かった。そのため、2日で回復すると仮定し、感染した人口の約1/2が毎日除去(Removed)されることになる。

\gamma = \frac{1}{2\text{days}}=0.5[\frac{1}{day} ]

参考として、除去率と感染率との比を計算する。その値は192である。

\rho = \frac{\gamma}{\beta}=\frac{0.5}{0.0026}= 192

4. モデルの予測と実測データとの比較

上記の条件を入れて、SIRモデルを数値積分で解く。

import numpy as np
from scipy.integrate import odeint
from scipy.optimize import minimize
import matplotlib.pyplot as plt
%matplotlib inline
#Define differential equation of SIR model

'''
dS/dt = -beta * S * I
dI/dt = beta * S * I - gamma * I
dR/dt = gamma * I

[v[0], v[1], v[2]]=[S, I, R]

dv[0]/dt = -beta * v[0] * v[1]
dv[1]/dt = beta * v[0] * v[1] - gamma * v[1]
dv[2]/dt = gamma * v[1]

'''
def SIR_EQ(v, t, beta, gamma):
    return [-beta*v[0]*v[1], beta * v[0] * v[1] - gamma * v[1], gamma * v[1]]

ここ部分に初期条件 S_0, I_0, R_0と感染率及び除去率を入力する。

#parameters
t_max = 14
dt = 0.01

beta_const = 0.0026
gamma_const = 0.5


#initial_state
S_0=762
I_0=1
R_0=0
ini_state = [S_0,I_0,R_0] #[S[0], I[0], R[0]]

#numerical integration
times =np.arange(1,t_max, dt)
args  =(beta_const, gamma_const)

#R0
N_total = S_0+I_0+R_0
R0 = N_total*beta_const *(1/gamma_const)
print(R0)

#Numerical Solution using scipy.integrate
#Solver SIR model
result = odeint(SIR_EQ, ini_state, times, args)
#plot
plt.plot(times,result)
plt.legend(['Susceptible','Infectious', 'Removed'])
plt.title("Influenza 1978: a=0.0026, b=0.5, $S_0$=762, $I_0$=1")
plt.xlabel('time(days)')
plt.ylabel('population')

plt.show()

実際の感染者のデータをSIRモデルと一緒にプロットしてみる。

# Real Data
data_day = [1,2,3,4,5,6,7,8,9,10,11,12,13,14]
data_infectious = [3,8,28,75,221,291,255,235,190,125,70,28,12,5]

#SIR model Plot
plt.plot(times,result)

# Real Data Plot

plt.scatter(data_day, data_infectious,s=10, c="pink", alpha=0.5, linewidths="2",
            edgecolors="red")
plt.title("Influenza 1978: a=0.0026, b=0.5, $S_0$=762, $I_0$=1")
plt.legend(['Susceptible','Infectious', 'Removed','data'])
plt.xlabel('time(days)')
plt.ylabel('population')
plt.show()

図の中のSIRモデルInfectious（感染者）と　実際の感染者のデータ(data)を比較すると概ねトレンドが一致することが確認できる。
このSIRモデルが予測した最大感染者数(Infections_Max)は、306名で、実際の最大の感染者は291名であって、その誤差は5%であった。
そして、SIRモデルが予測した最後の感染可能者（＝最後まで感染されなかった人数）は、16名であった。実際には、30名であった。

5.　まとめ

今回は実際の感染症事例を用いて、SIRモデルの有効性を確認した。
初期段階での感染者の情報より、初期条件および感染率・除去率のパラメータを推定した。
そして、実際の感染者データと比較し、SIRモデルの予測結果と概ね一致していることが確認できた。

6.参考資料

  1. https://qiita.com/kotai2003/items/3078f4095c3e94e5325c
  2. https://www.researchgate.net/publication/336701551_On_parameter_estimation_approaches_for_predicting_disease_transmission_through_optimization_deep_learning_and_statistical_inference_methods
  3. https://www.ncbi.nlm.nih.gov/pmc/articles/PMC1603269/pdf/brmedj00115-0064.pdf