感染症の流行を表すSIRモデルについて

新型コロナウイルス感染症の拡大を防止するために、政府は市民に人との接触の機会を8割削減するよう求めました。その要請の背後にあるとされる数理モデルを理解し、リスクの定量評価と対策の妥当性を評価できるようになっておくことが市民側にも求められていると思います。すでに多くの解説記事が出されていますが、多くの方にできる限りシンプルに理解できるように、このメモにまとめなおしてみました。

モデルの前提

接触機会の削減要請の背後にあるのがSIRモデルというものです。ある閉じられたコミュニティー内に、新しい感染症が発生したときに、その感染者数がどのように推移していくのかを表すシンプルな数理モデルです。このモデルによって、再生産数というパラメータや集団免疫などの基本的な概念とあわせて、感染症という現象の基本的なふるまいを理解することができると思います。

　まず感染という現象は、病原体を持っている人（感染者）が他者と接触したときに、ある確率で起こるものと仮定します。そして感染者はある一定の期間が過ぎれば回復して感染力を失うとともに、回復後は抗体を獲得しているために二度と感染はしない（病原体への感度を失う）とします。このような特徴を持つ感染症があるコミュニティーに発生したとき、はじめはみなS（Susceptible/感受性者=未感染者）であった人々が、I(Infectious/感染者)となり、最終的にはR(Recovered／回復者)の状態へと一方向の変化を経験することになります。

　上記のように想定した感染という現象を数式で表すために、以下のパラメータで感染を特徴づけられるとします。
B：一人の感染者が一日の間に人と接触する量［接触／日］
P：接触した人に感染させる確率［１／接触］
A：感染した人が回復に要する日数［日］
I：ある時点での感染者数［人］
ここで、接触量という曖昧な量を持ち出しましたが、感染事象に寄与する主には人数・時間・距離などが関係してくる量だと考えておけばよいと思います。

　そうしたとき一日の新規感染者数は、BPI［人／日］と表すことが出来ます。一方、一日ごとに感染者のうち(１／A) I［人／日］は回復していくとみなせるので、感染者数の1日ごとの増加人数は、BPI－(1／A)I＝（BP－1／A）I［人／日］であると表せます。そして、（BP－1／A）＞０である限り感染者は増え続けることになります。逆に感染者数を減少に転じさせるためには（BP－1／A）＜０にならなければならないことがわかります。

再生産数

感染者数が減少する条件（BP－1／A）＜０を書き換えると、BP＜1／A、さらに、ABP＜１となりますが、ここで出てきたABPは、「ひとりの感染者が感染してから回復するまでの間に感染させる人の数」を表していて、このパラメータには「再生産数」という名前がついています。感染者数の増減の条件を再生産数で表現すると、再生産数が１よりも大きければ感染者数は増加し続けるし、再生産数が１よりも小さければ感染者数は減少しつづけ収束に向かいます。

集団免疫という考え

上で仮定したように、一度かかった感染症に対しては抗体が作られて、もう感染することはないとした場合、おもしろいことに再生産数ABPは特別な対策を取らずとも徐々に小さくなっていくことが起こります。一人の感染者が一日の間に人と接触する量Bは変化しなくても、接触した相手のうち、すでに感染を経験した回復者の割合は時間とともに増えていくことになるでしょう。回復者は感染症への感度を持っていませんから、回復者の割合の増加は１接触当たりの感染のおこる確率Pを減少させることになります。定量的には、確率Pは全人口のうちの未感染者数の割合に比例することになると予想できます。
　たとえば、感染症が最初に発生したときの１接触当たりの感染の起こる確率をPo、最初の未感染者数をSoとすると、時間がたって未感染者数がSとなったときの感染の起こる確率は、P=Po (S/So)のように書けます。そのとき、再生産数が１よりも小さいという条件は、ABP<1すなわちABPo(S/So)<1であることから、S/So < 1/(ABPo)が得られます。ここででてきたABPoは、再生産数の初期値で基本再生産数と呼ばれています。つまり感染拡大がすすんで未感染者数割合がABPo(基本再生産数)の逆数よりも小さくなれば、再生産数は１を下回りいずれ感染の流行は終息するであろうことが予想できます。
　新型コロナウイルス感染症の場合、基本再生産数として得られている2.5という値を用いると1/2.5=0.4。つまり未感染者の割合が40％以下になれば、あるいは逆にすでに感染を経験して免疫を獲得した回復者の割合が60％以上になれば、感染は終息に向かうことになります。これが集団免疫という考え方で、イギリスのジョンソン首相が当初打ち出しすぐに取り下げた、“国民の6割に感染させる”という政策の数字の根拠はここにあると思われます。

S,I,Rの時間発展グラフ

さて次に感染者数が時間とともにどのように推移していくのか、実際に数値計算をしてその様子をみてみることにします。S,I,Rの時間変化は、これまで用いてきたパラメータを用いてそれぞれ以下のように書くことができます。

S: 未感染者Sは、一日ごとにBPI人減少する
    dS/dt = - B P I = - B Po (S/So) I
I : 感染者Iは、一日ごとに(BP-1/A)I人増加する
    dI/dt = ( B P – 1 / A ) I = ( B Po (S/So) – 1 / A ) I
R : 回復者Rは、一日ごとに(1/A)I人増加する
    dR/dt = ( 1 / A ) I　　

　それほど計算精度を求めなければ、オイラー法という微分方程式の解法を用いて比較的簡単にS,I,Rそれぞれの時間発展を計算することができます。初期値として、基本再生産数 ABPo=2.5（Imai et al. によると1.5～3.5）、回復に要する日数 A=10[日] （ある感染者の感染日と感染源となった人の感染日との間隔の平均が4.8日という研究結果Nishiura et al.から、平均の2倍を回復時間と仮定してみました) の特徴をもつ感染症が、総人口1億人の国に発生した場合を想定してみます。

第0日目に1人の感染者が発生したとすると、SIRそれぞれの初期値は
So=99999999, Io=1, Ro=0
とおけます。
翌日の第1日目は、
S1 = So – BPoIo = 99999999 - 2.5/10*1 = 99999998.75
I1 = Io + (BPo(So/So) – 1 / A )Io = 1 + (2.5/10-1/10)*1 = 1.15
R1 = Ro + (1/A)Io = 0 + (1/10)*1 = 0.1
（ABPo=2.5, A=10より、BPo=2.5/10が代入されていることに注意）

さらに第2日目は、
S2 = S1 – BPo(S1/So) I1
     = 99999998.75 - 2.5/10* (99999998.75/99999999)*1.15 = 99999998.46
I2 = I1 + (BPo(S1/So) - 1/A)I1 
     = 1.15 +(2.5/10*(99999998.75/99999999)-1/10)*1.15 = 1.32
R2 = R1 + (1/A)I1
     = 0.1 + (1/10)*1.15 = 0.215

のように順番に数字をいれて計算することができます。そのような方法で作ったグラフが下の図です。

　図中の赤い線で示したのは、再生産数ABP=ABPo(S/So)です。未感染者数が減少するのに呼応して再生産数も減少し、再生産数が１となるときに、感染者数がピークを迎え、それ以降減少に転じています。そして最終的には感染者数はゼロに、回復者数と未感染者数はそれぞれおよそ9000万人と1000万人とに収束していっています。再生産数が１のときの未感染者数は4000万人で、集団免疫のことを考えたときに求めた40％にちょうど一致しています。

最終規模

この図を見てよくわかることは、再生産数が１を下回るところは中間地点に過ぎないということです。感染者数はピークを迎えますが、その後も感染者数は積み上がり、最終的には全人口の9割が感染してようやく流行は終息することになります。基本再生産数2.5の感染症の場合、何の対策も取らずに自然に任せた場合、終息時点で感染を免れることができるのは全人口のわずか1割ほどになる(S∞/So=0.1)という予想が立てられるということです。
　感染が終息したときの未感染者あるいは回復者の割合は「最終規模」と呼ばれています。これも基本再生産数の値に応じて決まっていて、両者の間には次のような関係式が成り立つことが分かっています。
　ABPo（基本再生産数）= ln(S∞/So)/(S∞/So – 1)
この式は最終規模方程式と呼ばれていますが、以下のように微分方程式を解くことで導くことが出来ます。
dS/dt = - B P I ー①
dI/dt = (B P – 1 / A ) I　ー②
②／①より、
dI/dS = - (B P – 1 / A ) / B P = - 1 + 1 / ( A B P )
∫dI = ∫ ( - 1 + 1 / ( A B P ) ) dS
I – Io = - ( S – So ) + So / (A B Po ) ∫1/S dS
I – Io = - ( S – So ) + So / ( A B Po ) [lnS – lnSo]
I – Io = - ( S – So ) + So / ( A B Po ) ln(S/So)
のようにIとSの間の関係式が求まります。ここで、t→∞のときのことを考えると、
I∞ - Io = - ( S∞ - So ) + So / (A B Po ) ln (S∞/So)
A B Po / So (I∞ - Io + S∞ - So ) = ln (S∞/So)
A B Po (I∞/So – Io/So + S∞/So - 1 ) = ln (S∞/So) ここで、I∞=0, Io/So～0, なので、
A B Po (S∞/So – 1 ) = ln (S∞/So) 、従って、
A B Po = ln (S∞/So) / (S∞/So – 1 )　 のようにして最終規模方程式を求めることができました。

　グラフにしてみると、下のようになります。この図からも基本再生産数A B Po = 2.5 のときの最終規模は S∞/So=0.1となることがわかります。

流行を早く収めるには、再生産数を小さくする

以上SIRモデルを用いて、感染症が発生してから流行がどのように拡大して終息するのか、何の対策も取らずに自然に任せたときにどうなるのか、その様子が予想できました。それでは次に、流行を早く収めるためにどのような対策を取ればよいのかということですが、これは再生産数を人為的に小さくすることが必要だということになります。
　再生産数は、ABP=ABPo(S/So)のように表せます。ここで、それぞれのパラメータの意味は、
A：感染した人が回復に要する日数［日］
B：一人の感染者が一日の間に人と接触する量［接触／日］
Po：感染者が未感染者に1回接触したときに感染が起こる確率［１／接触］
S ：未感染者数［人］
So：最初の未感染者数（全人口）［人］
だったことから、これまで言われてきた感染防止対策は、それぞれ再生産数のどの部分を小さくすることに対応してきたのかについて以下のようにまとめられると思います。

①治療薬の開発⇒早く回復できるようにしてAを小さくする。
②人と人との接触の機会を減らす⇒Bを小さくする。
③マスクと手洗い⇒Poを小さくする。
④ワクチン開発⇒Sを小さくする。

②については、感染が起こるのはあくまでも感染者から未感染者へですから、感染者だけ隔離して出歩かないようにしておけばよいのかもしれません。しかしながら新型コロナウイルス感染症では感染していても全く症状の出ない人も大勢報告されていますし、症状がなくても唾液などにウイルスが含まれていて感染力を有しているとも言われています。だれが感染者なのかは外から見ただけではわかりませんし、自分が知らない間に感染者となってしまっていることも起こり得ます。従って人との接触量を減らす対策はあらゆる人が取り組まなければならないと理解できます。

参考資料

・西浦・稲葉　統計数理2006,第54巻第2号461-480
・稲葉寿　数理科学No.538, 2008
・「計算ウイルス学・免疫学の展開①」岩見真吾（九州大学大学院）集中講義＠広島大学．１２．７．２３( http://www.math.sci.hiroshima-u.ac.jp/~ryo/Lectures/Iwami/Note_1.pdf )
・「新型コロナウイルス感染症の流行予測」大橋順（東大ヒトゲノム多様性研究室）( http://www.bs.s.u-tokyo.ac.jp/content/files/SEIRモデルによるCOVID-19流行予測%20200331ver3.1.pdf )
・國谷紀良　「応用数学特論IIIb」 2018年12月17日講義メモ( http://www2.kobe-u.ac.jp/~tkuniya/18am3b )
・霊媒師 飯島章「SIRモデル」( https://note.com/mentalist007/n/n38acc092a99b )
・Natsuko Imai, Anne Cori, Ilaria Dorigatti, Marc Baguelin, Christl A. Donnelly, Steven Riley, Neil M. Ferguson. Transmissibility of 2019-nCoV. Imperial College London (25-01-2020)
・Nishiura, Linton, Akhmetzhanov, International Journal of Infectious Diseases, 93(2020)284–286