書記が数学やるだけ#904 SIRモデル
感染症モデルの基本であるSIRモデルについて扱っていく。
問題
説明
SIRモデルは、感染症の短期的な流行過程を決定論的に記述する古典的なモデル方程式である。名称は、このモデルが対象とする三つの集団――感受性保持者（Susceptible）、感染者（Infected）、免疫保持者（Recovered、あるいは隔離者 Removed）――の頭文字に由来する。原型となるモデルは、W・O・カーマックとA・G・マッケンドリックによって1927年に提案された。SIRモデルでは、全人口を感受性保持者、感染者、免疫保持者の三群に分ける。感受性保持者 S は、感受性保持者 S と感染者 I の積に比例して定率で感染者 I に移行し、感染者 I は一定の割合で免疫保持者 R に移行すると仮定する。すなわち、感染期間は指数分布に従うとされる。
解答
感染初期では、感受性保持者 S(t) がほぼ全人口 N に等しいと近似できる。このときの感染者数の変化から、感染が拡大する条件は基本再生産数 R0 = >1 の場合であり、逆に R0 <1 のとき感染は減少する。またこのときの一般解は指数関数形で与えられ、感染者数が時間とともに増加または減少する様子が示される。また式変形により、感染者数が減少するのに必要な免疫人口の割合は1-(1/R0)であることがわかる。
各方程式を加えて積分することで，任意の時刻において全人口が一定であることが言える。また，R(t) を消去すると、I(t) と S(t) の関係式が導かれる。これを相平面上に描くことで、感染のピークや I(t) の最大値を解析的に求めることができる。
S(t)の正値性について背理法で示す。同様の方法からI(t)とR(t)の正値性も確認できる。
定義式と各項の正値性から、S(t) は単調減少関数であり、R(t) は単調増加関数であることが示される。また、I(t) については、基本再生産数 R0 が 1 より小さい場合には単調減少し、R0 が 1 より大きい場合には、ある時点まで単調増加したのち単調減少に転じる。
SIRモデルの時間発展と相図を以下に示す。
import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint
# SIR model differential equations
def SIR(y, t, sigma, rho, N):
S, I, R = y
dSdt = -sigma * S * I
dIdt = sigma * S * I - rho * I
dRdt = rho * I
return [dSdt, dIdt, dRdt]
# Parameters
N = 1.0 # total population (normalized)
I0 = 0.01 # initial infected
R0 = 0.0 # initial recovered
S0 = N - I0 - R0
t = np.linspace(0, 100, 400)
# Two cases for comparison (R0 < 1 and R0 > 1)
params = [
(0.05, 0.1, "R0 < 1 (decay)"), # sigma*N/rho = 0.5
(0.2, 0.1, "R0 > 1 (outbreak)") # sigma*N/rho = 2
]
fig, axes = plt.subplots(2, 2, figsize=(10, 8))
axes = axes.ravel()
for i, (sigma, rho, label) in enumerate(params):
# Numerical integration
sol = odeint(SIR, [S0, I0, R0], t, args=(sigma, rho, N))
S, I, R = sol.T
# --- Time evolution ---
axes[2*i].plot(t, S, label="S(t): Susceptible", linewidth=2)
axes[2*i].plot(t, I, label="I(t): Infected", linewidth=2)
axes[2*i].plot(t, R, label="R(t): Recovered", linewidth=2)
axes[2*i].set_title(f"Time evolution: {label}")
axes[2*i].set_xlabel("Time")
axes[2*i].set_ylabel("Population fraction")
axes[2*i].legend()
axes[2*i].grid(True)
# --- Phase plane (S-I plot) ---
axes[2*i+1].plot(S, I, color="purple", linewidth=2)
axes[2*i+1].set_title(f"Phase plane (S–I): {label}")
axes[2*i+1].set_xlabel("S(t)")
axes[2*i+1].set_ylabel("I(t)")
axes[2*i+1].grid(True)
plt.tight_layout()
plt.show()
基本再生産数 R0 =ρσN を低下させるには、感染率 σ を減少させるか、回復率 ρ を増加させる必要がある。具体的には、
感染率の低下：マスク着用、接触機会の削減、ワクチン接種などによる感染確率の抑制。
回復率の上昇：早期治療や医療体制の改善による感染期間の短縮。
これらの対策により R0 <1 が実現すれば、感染症の流行は収束に向かう。
本記事のもくじはこちら：
