玄関雑学の部屋雑学コーナー感染症数理モデル

3.基本再生産数

閉鎖集団で感染症が流行するかどうかは、感染の初期段階でI(t)を経時的に観測し、そのデータを用いて、例えばポアソン回帰分析などで係数(βS0-γ)を求めればわかります。 (βS0-γ)は単位時間あたりの感染率と隔離率の差であり、単位時間あたりに隔離される割合よりも、感染する割合の方が多ければ流行すると解釈できます。 しかしこの解釈は数学的であり、普通の人にはわかりにくいと思います。 そこで次のような値を定義します。

基本再生産数(basic reproduction number)

このR0(R noughtと発音)は、単位時間あたりの感染率と隔離率の比を表します。 でもこれも普通の人にはわかりにくいので、次のように解釈するとわかりやすいと思います。

1人の感染者が感染期間中に何人の人を感染させたかを表す値

感染期間とは感染者が感染性を持ってから、隔離または回復または死亡によってそれを失うまでの期間のことです。 つまり1人の感染者が感染期間中に1人よりも多い人を感染させていれば流行が始まるわけです。 この様子は次の図を見るとわかりやすいと思います。 (図ではR0のことを基本再生産係数と表示しています)


※上図は厚労省対策推進本部クラスター対策班の押谷仁先生が作成された
「COVID-19への対策の概念」から引用させていただきました。

感染が進行すると、感染の初期段階ではS(t)≒S0として定数扱いしたS(t)の影響が無視できなくなります。 そのためR0は次のような式で表されます。 これを実効再生産数(Effective reproduction number、効果的再生産数ともいう)と呼び、ReまたはRtで表します。

I区画の変化速度:
I区画の単位時間あたりの変化率:
実効再生産数(Effective reproduction number):

:平均感染期間(感染期間が傾きγの指数分布をする時、その平均値は1/γになる)

上記のように、実効再生産数RtI区画の単位時間あたりの変化率kとγを用いてを求めることができます。 そしてγは、感染期間が近似的に指数分布する時は感染者の平均感染期間Dから推測することができます。 ただし実際の感染では、γI(t)とそれを積分したR(t)だけが直接観測可能です。 そこで次のようにして変化率kを求め、それと平均感染期間Dから実効再生産数Rtを求めます。

より

一方、より


I(t)を対数変換してtで微分すると

1 + k=⊿Rt単位時間あたりの再生産数

以上のように、k値はI区画の単位時間あたりの変化率であると同時に、時点tの近傍で感染者関数I(t)を指数関数によって近似した時の傾きに相当します。 そのためk>0の時はI(t)が指数関数的に増加し、k=0の時はI(t)がピークになり、k<0の時はI(t)が指数関数的に減少します。

そして(1+k)は単位時間あたりの再生産数つまり1人の感染者が単位時間あたりに何人の人を感染させたかを表す値になります。 したがってk値または(1+k)を実効再生産数Rtと同様の指標として用いることが可能です。 さらにk値または(1+k)は平均感染期間Dとは無関係に求められるので、Dの値があやふやな時はRtよりもむしろ有用な指標です。

ちなみに感染初期におけるI(t)の近似指数関数の両辺を対数変換すると、次のように1次式(直線)になります。

ln{I(t)} ≒ ln(I0)+(βS0-γ)t = ln(I0)+k0t

そこでI(t)を対数変換して直線回帰分析を行えば、初期k値であるk0=(βS0-γ)を近似的に求められると考えがちです。 しかしI(t)はカウントデータ(count data)ですから「0(感染者なし)」というデータがあり、これが重要な意味を持ちます。 そのためI(t)を対数変換して直線回帰分析を行うのは不適切です。 このような場合は非線形最小2乗法またはポアソン回帰分析を用いる必要があります。 (非線形最小2乗法については当館の「統計学入門・第14章 薬物動態学 (注2)」、ポアソン回帰分析については「統計学入門・第15章 ポアソン回帰分析」参照)

参考までに8章で説明する韓国のCOVID-19について、流行初期(2020年2月28日から9日間)のデータに非線形最小2乗法とポアソン回帰分析を適用し、指数関数を求めると次のようになります。

非線形最小2乗法(緑の曲線):I(t)≒17.5389×exp(0.334896×t)   k0=0.334896
ポアソン回帰分析(青い曲線):I(t)≒7.589×exp(0.460384×t)   k0=0.460384
※ln{I(t)}の直線回帰分析:ln{I(t)}≒-0.2449+0.827271×t   k0=0.827271

また全データに簡易モデルを当てはめてk値を求めると赤い破線のようになります。 この場合、流行初期のk値は約0.3であり、2つの指数関数から求めたk0値よりも少し小さくなります。 これはS(t)が減少する影響を計算に入れているからです。 そして平均感染期間Dを5日間とすると、基本再生産数R0は、非線形最小2乗法によると1+0.334896×5≒2.7、ポアソン回帰分析によると1+0.460684×5≒3.3、簡易モデルによると1+0.3×5=2.5になります。

現実の感染ではS(t)だけでなくβもγもtによって変化するのでRtは複雑に変化します。 そのため感染の動向を把握するには、発見された感染者の行動を詳細に調査し、感染者が接触した人と、接触された人が感染しているかどうかをチェックし、Rtの値を常にモニターすることが大切です。 そしてそれに基づいて感染予防対策を色々と検討することができます。 このように感染症数理モデルと疫学的調査を併用すると、合理的な感染予防対策を検討することができます。 その具体的な内容は厚労省対策推進本部クラスター対策班の押谷仁先生が作成された「COVID-19への対策の概念」をご覧ください。

以上の説明から、感染者数の推移を予想するにはβS0とγとRtを推測し、さらに感染者の行動を詳細に調べる必要があることがわかると思います。 マスコミやインターネットなどで、色々な分野のデータ解析屋が感染者数の推移を予想することがあります。 ところが日本で感染者の行動記録を持っているのは厚労省対策推進本部クラスター対策班だけであり、一般のデータ解析屋がそれらを推測することはほとんど不可能です。 そのため、どうしても感染初期におけるk0=(βS0-γ)とRtのラフな推測値から感染者数の推移を予測することになります。

したがって感染者の行動記録を持っていないデータ解析屋が予想した感染者数の推移予想は、僕のものも含めてσ(^^;)、あまり信用できないことがわかると思います。