Python で SIR モデルを理解

感染症の予測を行うための基本となるモデル（微分方程式）として，SIR モデルというものが知られている。

Python を用いて，SIR モデルでグラフを作成してみる。

SIR モデルとは

SIR モデルとは，全人口，すなわち感受性保持者 (Susceptible)，感染者 (Infected)，免疫保持者 (Recovered) または隔離者 (Removed) それぞれの頭文字である S，I，R によって命名されたものである。

この SIR モデルによって，人口全体 (S) が感染者 (I)，回復・隔離者 (R) に移行していきます。

免疫保持者または隔離者の増分

感染者 (I) のうち一定数が免疫保持者/隔離者 (R) になるとして計算できます。

定数 gamma を用いて，次式で計算できる。

R = I * gamma

感染者の増分

次に，感染者 (I) は，人口全体 (S) から減った分だけ増える。

そのうち回復・隔離者 (R) に変化した分だけ減るとして計算できます。

I = S - R

人口全体の増分

最後に「人口全体 (S) のうち感染者 (I) に変化するものがどれくらいか」を考える。

人口全体 (S) が多ければ多いほど，また，感染者 (I) が多ければ多いほど，それらに比例して増える。

定数 beta を用いて，次式で計算できる。

S = S * I * beta

数値シミュレーション

SIR モデルによる数値シミュレーションを行います。

具体的には，全人口が 200,000 人として，初期値 2 人の感染者が広がり，徐々に回復（あるいは隔離）される様子をシミュレーションします。

定数として beta = 0.000003，gamma = 0.1 を設定します。

SIR モデルの Python コード

ここで，SIR モデルの Python コードを示します。

import numpy as np
import matplotlib.pyplot as plt
import japanize_matplotlib

# パラメータ設定
S, I, R = 200000, 2, 0
dt = 1
beta = 0.000003
gamma = 0.1
num = 100

# 初期化
sus = np.zeros(num)
inf = np.zeros(num)
rec = np.zeros(num)

sus[0], inf[0], rec[0] = S, I, R

day = np.arange(num)  # 日数の配列

# SIRモデルの時間発展
for t in range(1, num):
    S, I, R = sus[t-1], inf[t-1], rec[t-1]
    alpha = I / (S + I + R)

    # 変化分の計算
    delta_R = I * gamma
    delta_S = -beta * S * I
    delta_S = min(delta_S, 0)  # Sは増加しないよう制限
    delta_I = -delta_S - delta_R

    # 更新
    sus[t] = max(S + delta_S * dt, 0)  # Sの値が負にならないよう制限
    inf[t] = I + delta_I * dt
    rec[t] = R + delta_R * dt

    print(f"日 {t}: S = {sus[t]:.2f}, I = {inf[t]:.2f}, R = {rec[t]:.2f}")

# 結果の可視化
plt.figure(figsize=(10, 6))
plt.plot(day, sus, label="S (Susceptible)", linestyle='-', marker='o')
plt.plot(day, inf, label="I (Infected)", linestyle='-', marker='s')
plt.plot(day, rec, label="R (Recovered)", linestyle='-', marker='d')

plt.legend()
plt.xlabel('日数')
plt.ylabel('個体数')
plt.grid(True)
plt.savefig("./images/SIRmodel.png", dpi=300)
plt.show()