いろいろなことを考察する

西浦氏の Scientific Reports 論文について(3)活動度は報告率によって変化すべきではない


私がネット上でしていることの まとめ は、こちらに。
https://sarkov28.hatenablog.com/entry/2022/03/29/160915



目次


(1)序論

西浦氏らのグループは、2023-10 に Scientific Reports にコロナワクチンに関する論文を発表しました。
後に「コロナワクチンで死者9割以上減 京都大チームが推計」と共同が報じた論文です。(以下、西浦論文、と書きます。)
Evaluating the COVID-19 vaccination program in Japan, 2021 using the counterfactual reproduction number
(論文には補足情報があります。補足情報の所在については 本稿末尾 に。)

図3.1 西浦論文の Figure 1
図3.1 西浦論文の Figure 1

本稿で Figure 1 について直接言及することはありませんが、補足的事項で、計算した活動度が Figure 1 に与える影響を見るために、重ね合わせたグラフを示しています。

西浦論文は、幾つかの要素から構成されています。
論文 式(4)で計算される活動度は一つの要素です。 また、reporting coverage 報告率も重要な要素であり、西浦論文はこの値を複数用意し、それぞれについて計算を示しています。 他にも論文を構成する要素はあるのですが、本稿が注目するのは今挙げた「活動度」と「報告率」の2つです。

本稿はまず、西浦論文が、

  • 活動度と報告率を独立したものとして扱うべきなのに、
  • 独立したものとして扱っていない、

ということを示します。この点は、明確に西浦論文に示されています。

本稿は次に、西浦論文において活動度と報告率が独立しているべき理由を示します。

ただし、活動度と報告率が独立していないことが、西浦論文の結果に定量的に大きな影響を与えているかどうかは不明です。本稿が述べるのは、(定量的に大きな問題の指摘というより、)定性的に、モデル構成のあるべき姿として問題ではないか、という話です。

また、報告率が論文 式(4)で計算している活動度  h_t に与えている影響をグラフで示します。

西浦論文については、複数のエントリで書いています。

(2)活動度と報告率の関係 ─ あるべき構成と西浦論文の構成

西浦論文のモデルは、多数の要素から構成されています。

論文で、活動度は  h_t と表記されています。式(4)で計算され、式(3)で用いられています。  h_t は、「google mobility の6指標のうち3指標と、その3指標を調整するパラメータから計算される数値」です。

本来、西浦論文のモデルは図3.2の構成になっているべきです。図3.2のポイントは、活動度と報告率とが独立していること、です。ところが西浦論文のモデルは、図3.3の構成になっています。西浦論文のモデルは、図3.3に示したように、活動度が報告率から独立していません。

図3.2 西浦論文のモデルのあるべき構成
図3.2 西浦論文のモデルのあるべき構成

図3.3 西浦論文のモデル構成
図3.3 西浦論文のモデル構成

(3)西浦論文の活動度が報告率から独立していないことの根拠

西浦論文の活動度が報告率から独立していないことは、論文に明記されています。
補足情報Table S1 には、4つの報告率毎の  ω_h,\ ω_w が記載されています。(西浦論文の活動度に関する表記には不統一があるので、別項で整理しています。

表3.1 西浦論文の補足情報の Table S1
表3.1 西浦論文の補足情報の Table S1。(クリックすると大きく表示される。)

表3.1 にあるように、活動度(を決めるパラメータ  ω_h,\ ω_w )は報告率毎に異なります。これは、活動度が報告率から独立していないことを示しています。

 ω_h,\ ω_w が報告率毎に異なっているのは、西浦論文がこれらを報告率毎に計算したからです。西浦氏ら論文執筆グループがそのような計算方法を選んだのであって、これは必然ではありません。
報告率によって変化しない、(4つの)報告率に共通な  ω_h,\ ω_w を計算するような尤度関数を作成し、計算することも選択できます。

西浦氏ら論文執筆グループは、(意識して行ったかは不明ですが、)図3.2 の構成ではなく、図3.3 の構成を選択しています。
この選択は不適切だと思います。

(4)活動度が報告率から独立しているべき理由

活動度と報告率は、以下のものを表しています。

  • 活動度は、人々の活動がコロナに与える影響をモデル化したものです。
  • 報告率は、行政がコロナに感染した人を捕捉する割合を示したものです。

実社会における「人々の活動がコロナに与える影響(活動度)」は、「行政がコロナに感染した人を捕捉する割合(報告率)」と独立していると考えられます。両者の間に特段の関係があるとは考えられません。ですから、西浦論文のモデルが良いモデルであるなら、実社会におけるコロナの感染態様を適切にモデル化したものであるなら、活動度は報告率から独立しているべきです。

本稿の主たる考察は以上です。次節は、以上に関連する参考事項です。

(5)活動度は報告率からどのような影響を受けているのか

論文 式(4) で計算される活動度  h_t が、どのように報告率の影響を受けるのかを実際に計算します。

西浦論文が用いたのは、google mobility の、6個ある指標のうち、"retail_and_recreation", "workplaces", "residential" の3個です。(ただし西浦論文は、それぞれ "community", "workplace", "household" と表記しています。この付近の表記には他にも問題があるので、別に整理しています。
西浦論文は、式(4) の計算のために google mobility の指標を加工し、Figure S4 に示しています。

図3.4 西浦論文の補足情報 Figure S4
図3.4 西浦論文の補足情報 Figure S4

私は、図3.4 を再現しました。

図3.5 図3.4の再現
図3.5 google mobility データを加工して、図3.4を私が再現したもの
西浦論文 式(4) の  α_t^{community},\ α_t^{house},\ α_t^{work}

図3.5は、google mobility データを加工して、図3.4を私が再現したものです。よく再現できていると思います。加工の詳細は(10)に。 図3.5を描くのに用いたデータはこちら。

再現した活動度データと、西浦論文の補足情報 Table S1 を使うと、論文 式(4)を計算できます。

式(4)は、こちら。

 h_t = ω^{community}\ α_t^{community} + ω^{house}\ α_t^{house} + ω^{work}\ α_t^{work} ...... 論文 式(4)

再現した図3.5 が 式(4)の  α_t^{community},\ α_t^{house},\ α_t^{work} です。  ω^{community} は定義から  1 であると本文に書いてあり、 ω^{house},\ ω^{work}Table S1 に( ω_h,\ ω_w の表記で)記載されています。

reporting coverage  ω_h  ω_w
1 −0.490 −0.083
0.5 −0.483 −0.088
0.25 −0.467 −0.096
0.125 −0.427 −0.118

表2 西浦論文の補足情報 Table S1 のうち、式(4) に関係する部分

これらを用いると  h_t が計算できます。

図3.6 論文 式(4)のグラフ
図3.6 論文 式(4)を計算してグラフにしたもの
西浦論文 式(3)、式(4) の  h_t

図3.6 は論文 式(4) の  h_t です。描くのに用いたデータはこちら。4本の曲線が描かれているのは、表2 に示されている4つの reporting coverage 報告率毎の  ω_h,\ ω_w を使って、論文 式(4) を計算したからです。

これが、報告率 reporting coverage が、活動度  h_t に与えている影響です。

ただし、論文 式(3)の左辺  R_{ab,t} に対して、図3.6 が示している報告率によるパラメータの変化が与える影響を考える際は、図3.6 以外の要因にも注意する必要があります。  R_{ab,t} が報告率から受ける影響は、(報告率によって変化する)活動度  h_t 経由のもの(=図3.6)だけではなく、他の(報告率によって変化する)パラメータ経由のものもあるからです。

ややこしく書きましたが要するに、図3.6 からすぐに何かを言えるということはないと思います。

次節以降は、補足的な事項です。

(6)(参考グラフ)図3.7

図3.6 に示した論文 式(4) の活動度 が、西浦論文 Figure 1 にどのように影響したのかは気になるところです。論文 Figure 1 と 図3.6 を重ね合わせた図を作成しましたので、ご参考までに示します。こちらです。

図3.7 西浦論文 Figure 1 と図3.6 の一部を重ね合わせたグラフ
図3.7 西浦論文 Figure 1 と図3.6 の一部を重ね合わせたグラフ

図3.7 の重ね合わせに使ったのは、図3.6 の一部、報告率  1 の曲線だけです。西浦論文 Figure 1 が報告率  1 における曲線なので、この曲線だけを使いました。

ご説明するまでもなく、Figure 1 の形状に影響を与えるのは、図3.6 の活動度だけではありません。この他、感染やワクチンによる免疫獲得や、市中におけるデルタ株の比率なども影響を与えています。
ただしこれらは、活動度ほど「細かく」上下には動きません。(例えばワクチン接種による免疫保持者の割合を示した、補足情報 Figure S2 を参照。上下しますが、ゆったりとした動きです。)したがって、Figure 1 曲線の細かな上下動に大きく寄与しているのは、活動度ではないかと思われます。

(7)(参考グラフ)図3.8

本稿の内容を把握して頂いたなら、図3.6 に示した論文 式(4) の(報告率毎の)活動度  h_t に、論文 式(3) の(報告率毎の) p×r を掛算したグラフにご興味があるかも知れません。このグラフから読み取れる特記すべき事項はないと考えますが、ご参考として示しておきますと、こちらです。

図3.8 図3.6(論文 式(4)の h_t)に、報告率毎の p×r を掛けたもの
図3.8 (参考グラフ)図3.6(論文 式(4)の h_t)に、報告率毎の  p×r を掛けたもの

図3.8 のグラフの形状は、図3.6 と似ていますが、違いも少しあり、報告率 0.125 の曲線と、その他の曲線との相対的位置がやや近づいたように見えます。グラフ縦軸が変化していますが、これは本質的問題ではありません。

(8)論文の活動度に関する表記の不統一

西浦論文の活動度  h_t に関する表記には不統一があるので、整理します。
 h_t は論文 式(4)で計算され、式(3)で用いられています。式(4) の計算では、「google mobility の6指標のうち3指標("retail_and_recreation", "workplaces", "residential")と、それらを調整するパラメータ」を使います。

google mobility の3指標を指し示すはずの論文の表記は不統一で、以下の表のような対応関係があります。

google mobility での表記 retail_and_recreation residential workplaces
論文本文での表記 community household workplace
式(4)での活動度の表記  α^{community}  α^{house}  α^{work}
式(4)でのパラメータの表記  ω^{community}  ω^{house}  ω^{work}
補足情報 Table S1 での
パラメータの表記
(表記なし *1)  ω_{h}  ω_{w}

表3 論文の活動度に関する表記の対応
(*1 Table S1 に community に関する項目はありません。この値は常に 1 と論文本文で説明されているので、記載を省略したと思われます。)

西浦論文の表記は不統一ですが、本稿では論文の表記を踏襲しています。

(9)論文の補足情報の所在

論文の補足情報の所在の詳細は、論文
https://www.nature.com/articles/s41598-023-44942-6
をご確認頂いた方がいいですが、以下に概要をご説明します。

論文の下部には、補足情報 Supplementary Information として幾つかの url が示されていて、そのうちの1つに本稿で何度か言及した Table S1 が入っています。

Table S1 が入っているファイルの url は
https://static-content.springer.com/esm/art:10.1038/s41598-023-44942-6/MediaObjects/41598_2023_44942_MOESM1_ESM.docx
になりますが、私の環境でこの url をクリックすると、Microsoft Word のファイルである 41598_2023_44942_MOESM1_ESM.docx が、ブラウザのデフォルトのダウンロードフォルダにダウンロードされます。そのファイルを開くと、ほぼ最下部に Table S1 があります。

(10)google mobility の加工してパラメータαを得る方法の詳細

西浦論文は、google mobility を加工してパラメータ  α^{community},\ α^{house},\ α^{work} を得ていますが、加工の細部は示されていません。私はこれを以下のように推測しました。図3.5が図3.4を非常によく再現していることから、私の推測は正しいと思われます。

 α = (google\ mobility\ 指標の7日間移動平均) / 100 + 1

ここで、移動平均を計算する際には、python で、pandas.DataFrame.rolling(7, center=True).mean() を使いました。デフォルトではない center=True を指定している点を特に注意します。

(99)更新履歴

  • 2023-12-30
    公開。
  • 2024-01-04
    別稿「西浦氏の Scientific Reports 論文について(2)」に、「更新履歴 2024-01-04」の大幅な修正をしました。これに合わせて本稿(7)の関連する箇所を修正しました。
  • 2024-01-10
    図3.7 を鮮明なものに差し替えました。グラフの内容は同じです。
  • 2024-01-11
    図3.5 の凡例の順序が 図3.4 の凡例と不一致だったのを修正しました。
  • 2024-01-23
    図3.5 と 図3.6 のデータを github に置き、リンクを貼りました。
  • 2024-01-25
    (10)を加筆し、ここへの言及を追記しました。


西浦氏の Scientific Reports 論文について(2)論文 式(3)、式(8)の無駄なパラメータ


私がネット上でしていることの まとめ は、こちらに。
https://sarkov28.hatenablog.com/entry/2022/03/29/160915



目次


(1)序論

西浦氏らのグループは、2023-10 に Scientific Reports にコロナワクチンに関する論文を発表しました。
後に「コロナワクチンで死者9割以上減 京都大チームが推計」と共同が報じた論文です。(以下、西浦論文、と書きます。)
Evaluating the COVID-19 vaccination program in Japan, 2021 using the counterfactual reproduction number
(論文には補足情報があります。補足情報の所在については 本稿末尾 に。)

図2.1 西浦論文の Figure 1
図2.1 西浦論文の Figure 1(本稿で Figure 1 について直接言及することはない。)

西浦論文の式(3)と式(8)には、無駄なパラメータがあるのでご説明します。

この無駄なパラメータの存在は、西浦論文の結果に実質的な影響を与えることはないと思いますが、論文の見通しを低下させています。

途中で出てくる数式は、「多少」ややこしいですが、奇妙な点の核心は非常に単純で明確です。

西浦論文については、複数のエントリで書いています。

(2)西浦論文の式(3)には無駄なパラメータがある

本項は、西浦論文 式(3) に無駄なパラメータがあることを示します。と言っても単に、式(3) に論文に定義されてる式を代入し、式変形するだけです。この「代入」「式変形」は中学校の数学の範囲の話だと思います。

数式に  Σ 記号が出てきますが、本項は  Σ が関係する部分は操作しませんので、仮に  Σ のところが分からなくても、本項が指摘する問題はご理解頂けると思います。

西浦論文 式(3)は、以下です。

 \displaystyle R_{ab, t} = \left(1 - l_{a, t} - \frac{\sum_{k=1}^{t-1}{i_{a, k}}}{n_a}\right) K_{ab}\ p\ h_t\ d_t\ c_t ...... 論文 式(3)

論文によると、 d_t d_t = r u_t です。
これを 式(3) に代入すると、以下になります。

 \displaystyle R_{ab, t} = \left(1 - l_{a, t} - \frac{\sum_{k=1}^{t-1}{i_{a, k}}}{n_a}\right) K_{ab}\ p\ h_t\ r u_t\ c_t ...... (3)'

式(3)' のうち、 p,\ r は単一の値を取るパラメータです。
 p,\ r が単一の値を取るのは、西浦論文の 補足情報 にある Table S1 からも明らかです。Table S1 には、各条件において、 p,\ r の値が「1つ」だけ示されています。

表2.1 西浦論文の補足情報の Table S1
表2.1 西浦論文の 補足情報 の Table S1。(クリックすると大きく表示される。)

 p,\ r は単一の値のパラメータなので、前に移動します。

 \displaystyle R_{ab, t} = p\ r \left(1 - l_{a, t} - \frac{\sum_{k=1}^{t-1}{i_{a, k}}}{n_a}\right) K_{ab}\ h_t\ u_t\ c_t ...... (3)' '

 p\ r は単一の値のパラメータの掛算なので、まとめることができます。これを  q と書くことにします。この  q = p×r を使って、式(3)' ' を書き換えます。

 \displaystyle R_{ab, t} = q \left(1 - l_{a, t} - \frac{\sum_{k=1}^{t-1}{i_{a, k}}}{n_a}\right) K_{ab}\ h_t\ u_t\ c_t ...... (3)' ' '

2つのパラメータが1つになりました。

(なお、 c_t = e β_t は特殊な動きをするパラメータなので、 d_t = r u_t と同様に扱うことはできません。変更履歴の 2024-01-04 の項目に関連事項があります。)

(3)西浦論文の式(8)は無駄なパラメータを最尤推定している

西浦論文は、尤度関数(論文式(8))という関数を定義し、 θ を推定し、先ほど参照した西浦論文 補足情報 の Table S1 を得ています。

 θ = \{p,\ ω_{house},\ ω_{work},\ r,\ e\}

この  θ には、先ほど  q に( q = p×r と)まとめた  p,\ r が入っています。

しかし(2)で示したように、 p,\ r のそれぞれを推定するのではなく、掛算された  q の1つを推定すれば足りています。

 q = p×r の関係が成立していますから、

  •  p,\ r をそれぞれ最尤推定し、(式の中で実質的に) p×r を計算すること(西浦論文の方法)と、
  •  q最尤推定すること(私が本稿で示している方法)は、

論理的には等価です。しかし西浦論文の方法は無駄です。論文のように計算する理由はありませんので、一般にこのようなことはしません。

(4)西浦論文が計算したパラメータは、本稿が導入したパラメータよりバラつく

西浦論文が計算した  p,\ r と、上で導入した  q = p×r は、4つの報告率 reporting coverage において以下になります。

報告率  p  r  q = p×r
1 2.749 0.994 2.733
0.5 2.722 1.005 2.736
0.25 2.660 1.029 2.737
0.125 2.509 1.080 2.710

表2.2: 西浦論文 補足情報 Table S1 のパラメータ  p,\ r と、上で導入した  q = p×r。4つの報告率 reporting coverage における値。

値のバラつきを見やすくするため、 p,\ r,\ q のそれぞれを最小値  min(\ ) で割ると以下になります。

報告率  p/min(p)  r/min(r)  q/min(q)
1 1.096 1.000 1.008
0.5 1.085 1.011 1.010
0.25 1.060 1.035 1.010
0.125 1.000 1.087 1.000

表2.3: 表2.2の  p,\ r,\ q を、それぞれの最小値  min(\ ) で割った値。


表2.3をグラフにすると以下になります。
 p/min(p) r/min(r) はバラついていて、掛算した  q/min(q) ではバラつきが相殺されていることが分かります。

図2.1 パラメータ p, r とその掛算 q の、報告率によるバラつき(それぞれを最小値で割って正規化した)
図2.1 パラメータ  p,\ r とその掛算  q の、報告率によるバラつき(それぞれを最小値で割って正規化した)

 q/min(q) の、報告率=0.125 以外の3つの値がかなり近い値になっているのが目を引きます。(これは偶然なのかも知れませんが、)このモデルの作成者ならここは気になるところではないかと思います。

(5)まとめ

西浦論文がパラメータ  q p,\ r の3つに分けているのは無駄です。

これは単なる無駄であり、論文の結果に実質的な影響を与えることはないと思います。
しかしこうした無駄は、数値的な計算を行う際には極力避けるものです。
計算誤差を拡大する上、計算の(したがって論文の)見通しを低下させるからです。

(6)論文の補足情報の所在

論文の補足情報の所在の詳細は、論文
https://www.nature.com/articles/s41598-023-44942-6
をご確認頂いた方がいいですが、以下に概要をご説明します。

論文の下部には、補足情報 Supplementary Information として幾つかの url が示されていて、そのうちの1つに本稿で何度か言及した Table S1 が入っています。

Table S1 が入っているファイルの url は
https://static-content.springer.com/esm/art:10.1038/s41598-023-44942-6/MediaObjects/41598_2023_44942_MOESM1_ESM.docx
になりますが、私の環境でこの url をクリックすると、Microsoft Word のファイルである 41598_2023_44942_MOESM1_ESM.docx が、ブラウザのデフォルトのダウンロードフォルダにダウンロードされます。そのファイルを開くと、ほぼ最下部に Table S1 があります。

(99)更新履歴

  • 2023-12-28
    公開。
  • 2023-12-30
    図2.1 の配色を変更。
  • 2024-01-04
    大幅に修正しました。修正前は、パラメータ  p,\ r の他に  e も同様の議論の対象になると考え、そう記述していました。したがって  q = p×r ではなく、 q = p×r×e としていました。表2.2、表2.3、図2.1や関連箇所も同様です。しかし  e は対象にできないことに気付いたので、 e に関する部分を削除しました。


西浦氏の Scientific Reports 論文について(1)論文の最大の特徴は、類のない非現実的な想定


私がネット上でしていることの まとめ は、こちらに。
https://sarkov28.hatenablog.com/entry/2022/03/29/160915



目次


(1)序論

西浦氏らのグループは、2023-10 に Scientific Reports にコロナワクチンに関する論文を発表しました。
後に「コロナワクチンで死者9割以上減 京都大チームが推計」と共同が報じた論文です。(以下、西浦論文、と書きます。)
Evaluating the COVID-19 vaccination program in Japan, 2021 using the counterfactual reproduction number

図1.1 西浦論文の Figure 1
図1.1 西浦論文の Figure 1(本稿では Figure 1 に直接言及することはない。)

私は、

  • 西浦論文は、「何も感染症対策を変えなければ」との非現実的な状況を想定したものである。
  • これは極めて特殊な論文である。
  • この論文を紹介した記事には不備があった。
  • 西浦氏は不備のある記事を「ぜひ多くの方に読んでいただきたい」と紹介したことは、不適切である。

と考えますので、以下に、その根拠を書きます。

この論文については、複数のエントリで書いています。

西浦論文は、ワクチン接種の4つのシナリオでの感染被害を計算していますが、本稿ではこのうち、「ワクチン接種なしシナリオ」に注目します。
このシナリオに注目するのは、このシナリオが一番極端な結果を示しているから、この論文を報じた記事がこのシナリオでの計算を特に示したから、また他のシナリオに関しても同様の議論が適用できるから、です。

(2)西浦氏は西浦論文の紹介を適切に行うべき

西浦氏は一般の研究者ではないし、この研究は一般的な研究ではないのだから、西浦論文を多くの人向けに紹介するのであれば、適切な紹介をするべきだと、私は考えます。
この点を本項で説明します。

一般の専門家が自分の論文を(あるいは自分の論文を報じた記事を)紹介した際、その周辺に多少の誤りがあっても一々批判する必要があるとは思いません。

一般の専門家は、しばしば人々に名前を記憶されていませんし、「何を研究しているのか分からない」と思われがちです。だから、「少しでも多くの人に自分(の研究)について知ってもらいたい。そこに多少の誤りがあってもいい」と考えることにも一定の正当性があります。また、一般的な研究は、多くの人々の生活に影響を与えるとしても、ほとんどが間接的なものです。したがって、紹介に多少の誤りがあってもその影響は大抵は小さなものになります。

しかし西浦氏は、この意味での「一般の専門家」ではありません。また今回の紹介は、この意味での「一般的な研究の紹介」には該当しません。

まず西浦氏は、日本における理論疫学の代表者であることを既に多くの人に知られています。 「8割おじさん」との(適切なのかが不明な)呼称もあります。 一般の研究者とは違うのであり、今回のような研究についての人々の関心はどちらかと言うと高いと言えます。 この研究やこの記事や西浦氏の紹介に関して、twitter.com 上で多くの投稿がなされたことも、人々の関心の高さを示しています。

また この研究を報じた記事 は、「ワクチン接種政策が我々の社会に与えた影響を、この研究が示した」かのように記述していました。 したがってこの記事は、社会への影響が間接的だったり、小さかったりするものとは異なります。 もしかすると、この記事を読んでワクチンへの評価やワクチン接種の判断を変えた読者がいたかも知れません。この意味で、むしろ大きな直接的な影響を与えうる記事だったと言えます。

私が、「西浦氏は一般の研究者ではないし、この研究は一般的な研究ではないのだから、西浦論文(や関連の記事)を多くの人向けに紹介するのであれば、適切な紹介をするべきだ」と、考る理由は以上です。

(3)西浦論文は、非現実的な状況を想定している

西浦論文は「何も感染症対策を変えなければ」という非現実的な状況を想定して計算しています。各項では以下を説明します。

  • (3-1): 西浦論文が前提とする「何も感染症対策を変えなければ」について。
  • (3-2) と (3-3): 西浦氏が twitter.com への投稿で、「何も感染症対策を変えなければ」が非現実的であると述べていることを。
  • (3-4): 西浦論文を「知的体操」と表現した西浦氏の投稿について。

(3-1)「ワクチン接種なしシナリオ」 ─ 「何も感染症対策を変えなければ」とは

まず、西浦論文が想定している「何も感染症対策を変えなければ」を説明します。
「ワクチン接種なしシナリオ」では、ワクチン接種をしないこと以外の全ての感染症対策が、現実のものと全く同じだったと想定しています。
つまりこの論文は、以下のことを想定しています。

  • 政府の対応について
    西浦論文が計算対象とした 2021-02-17 から 2021-11-30 の間、緊急事態宣言が2度出されました。この宣言の期間は(「ワクチン接種なし」でどんなに感染が拡大しても)宣言期間を延長するなどせず、全く同じと想定します。
  • 病院・施設の対応について
    病院・施設での面会制限も、(「ワクチン接種なし」でどんなに感染が拡大しても)強化などせず、全く同じと想定します。
  • 個人の対策について
    一般の人の対策や行動も同様であり、(「ワクチン接種なし」でどんなに感染が拡大しても)全く同じと想定します。
    例えば「ワクチン接種に出掛ける」との行動についても同じと考えるしかありません。
    つまり、現実にこの期間にワクチンを接種した人の「ワクチン接種なしシナリオ」での行動は、「同じ日に接種会場に出向き、そこで何かを注射してもらった」とでも考えるべきなのです。これは「接種会場に出掛ける」ことの感染に対する影響を、現実と全く同じにするための想定です。現実にワクチン接種の副反応が出た人は、何かを注射してもらった事により、同じ何らかの副反応が出たと想定するしかありません。
    (この段落は削除します。「ワクチン接種なし」では、「ワクチンの効果だけが無くなる」と考えてこのように書いたのですが、「ワクチン接種に伴う一連の行動が無くなる」と考えるべきでした。失礼しました。)

西浦論文が「そこまで同じと考えなくてもいい」と考えているはずはありません。もし「その感染症対策は同じでなくてもいい」というなら、その感染症対策を無駄だと考えていることになります。もし「その行動は同じでなくてもいい」というなら、その行動は感染に影響しないと考えていることになります。西浦氏は感染症対策を強く主張し、行動に注意せよと呼びかけていたので、このように考えているはずはないのです。

政府や、病院・施設や、一般の人は、「どんなに感染状況が酷くなっても(逆に軽くなっても)、実際の 2021-02-17 から 2021-11-30 と同じ様に行動した」と考えるのが、西浦論文が想定する「何も感染症対策を変えなければ」です。
私が極端なことを言っているのではありません。極端なのは西浦論文の想定です。

なお上の説明の内容は、論文のロジックにも見て取れます。
西浦論文は以下のように計算しています。

  • 「現実の感染被害」を Figure 1 で「再現」した後、
  • ワクチンの接種率「だけ」を変更し、
  • 「ワクチン接種なしシナリオの感染被害」を計算しています。

この操作は、上で説明した「何も感染症対策を変えなければ」を、計算モデル上で表現したものです。

西浦論文が「何も感染症対策を変えなければ」を採用した事情については、西浦論文の中にもある程度の説明があります。この件は(5) に書きました。

(3-2) 西浦論文は、非現実的状況を想定している(その1)

そもそも西浦氏は、西浦論文が非現実的状況を想定としていると端的に説明しています。

非現実的だとしても仮想的に作っている
https://twitter.com/nishiurah/status/1719494672786321896 (2023-11-01)

(3-3) 西浦論文は、非現実的状況を想定している(その2)

「ワクチン接種なしシナリオ」での計算は、「10カ月で36万人が死亡」でした。西浦論文は、こんなことになっても社会は全ての行動を変えないと想定しています。

しかし西浦氏は、そんなことは現実ではあり得ないと 2021 年に述べています。

「十万単位死亡なんかあるかい」はおっしゃる通り。仮に対策なくとも、日本で情緒的にそういう流行に耐えられません。皆が接触を止めます。
https://twitter.com/nishiurah/status/1437551244436914177 (2021-09-14)

西浦氏は、
「数十万人が死亡するような状況では、皆が接触を止めます」
と述べていて、これは、
「数十万人が死亡するような状況では、社会は行動を変える」
という意味です。
この説明は、西浦論文の想定が非現実的であると示しています。

もしも感染状況が現実よりも数倍の規模に拡大したら、政府は緊急事態宣言を延長するに違いありません。(こうした対策が、社会に与える影響に見合う適切なものかは別の議論です。)何も対策を変えないということは現実には起きないのです。

(3-4) 西浦氏は、西浦論文を「知的体操」だと説明している

西浦氏は、西浦論文を「知的体操」だと説明しています。
https://twitter.com/nishiurah/status/1719494672786321896 (2023-11-01)

これらは単なる知的体操の問題なんですが、要するに「他の条件が全部一緒だとして、予防接種だけなかった場合に」というのを非現実的だとしても仮想的に作っているということですよね。

西浦氏が「知的体操」とわざわざ書いたのは、この論文が、「非現実的想定と、現実とを比較する」という特殊な計算であることを、西浦氏が良く知っているからだと思われます。

なお、ここに示した、西浦氏による twitter.com への 投稿 は、「FAQにお答えしましょう」と書き出されたスレッド の一部ですから、通常の書き込み以上に西浦氏の真意を表しているはずです。

(4)西浦論文の最大の特徴は、類のない非現実的な想定

(4-1) 他の分野の専門家による、被害想定

専門家が何らかの被害想定を計算することは、他の分野でも行われています。よく目にするのは、巨大地震による被害想定です。死者数、負傷者数、経済的被害などが計算されます。

これは低い確率で起きる状況ですが、「決して実現する事はない」とその専門家本人が述べるような状況が計算されることはありません。確率は低くても、起こる可能性があるからこそ、被害を計算する価値があるのですし、社会が注意を向けるのです。しばしば、その確率も示されます。「こういう地震が何十年以内に何%の確率で起きる」「何十年の一度の大雨が降るとこうなる」など。地震災害の専門家が「決して実現する事はない地震の被害を計算しました」「ぜひ多くの方に読んでいただきたい」と述べたら、どちらかと言うと不思議がられるのではないでしょうか。

西浦論文や、西浦氏が2020年4月に発表した「42万死亡推計」は、「決して実現する事はない感染状況での被害を計算したもの」です。ただしこれらの発表は、あまり不思議がられることはありませんでした。

(4-2) 西浦論文が非現実的な想定で計算したのは、「理論疫学の都合」

地震や水害では「ある程度は起こるかも知れない状況」を計算するのに、理論疫学が「決して実現する事はない状況」を計算するのは、なぜなのか。単に西浦氏が甚大な被害を計算して見せたいからではないと思います。

現在の理論疫学においては、「何も感染症対策を変えなければ」「何も感染症対策をしなければ」といった非現実的な被害想定ならば計算できることになっていますが、「実際と少し違う特定の状況で」などのより現実的な被害想定は計算できないと思われます。これは理論疫学の都合です。

「実効再生産数が〇〇だけ変化したら」は計算できますが、ある感染症対策が実効再生産数をどれだけ変えるのかは計算できないのです。西浦論文や「42万死亡推計」が非現実的な想定で計算されたのはこれが理由でしょう。

(上で「計算できることになっています」と回りくどく書いたのは、西浦氏はそう考えているはずですが、私は異論があるからです。この件は今は措きます。)

私は、非現実的な状況を想定して計算すること自体を批判しているのではありません。しかしそれは、専門家の「知的体操」であり、専門家の間で議論の対象にするような特殊な計算です。一般の人が見聞して内容を理解することが要請されたり、奨励されたりするものではありません。

(4-3) 西浦論文の想定は「実現したかも知れない非現実」ではない

「誰もワクチンを接種しない」との想定だけなら、(これもかなり特殊ですが、)実現するかも知れません。しかしこの論文はそれと同時に、(3-3) などで述べたように、「10万人単位が死亡するような酷い感染状況下でも、政府や人々の感染対策が、実際と全く変わらない」との状況を想定しています。この状況は非現実的です。この「何も感染症対策を変えなければ」が非現実的で特殊な想定なのです。

ここでいう非現実が、「現実的な状況から少し外れた、実現したかも知れないけどたまたま実現しなかった非現実」とは違うことを特に指摘します。西浦論文の非現実は、「どんなことがあっても実現することのない非現実」です。

(4-4) 西浦論文の最大の特徴は、他に類を見ない非現実的で特殊な想定

論文は、実現することがないとはっきり分かっている状況を計算機の中に仮想的に作り、そこでの感染者数や死者数を計算しています。

西浦論文の最大の特徴は、他の分野で類を見ない非現実的で特殊な想定にあります。

(5)非現実的な想定で計算した事情に関する、西浦論文の記述

(5-1) 非現実的な想定で計算した事情に関する、西浦論文の記述

「理論疫学の都合」について述べた (4-2) で、

現在の理論疫学においては、「何も感染症対策を変えなければ」「何も感染症対策をしなければ」といった非現実的な被害想定ならば計算できることになっていますが、「実際と少し違う特定の状況で」などのより現実的な被害想定は計算できないと思われます。

と書きました。関連する記述が西浦論文にもあるのでご説明します。

Our study involved several technical limitations. ... during the research period, Japan experienced three state of emergency declarations ... . Rather than incorporating the specific variable of a state of emergency into the model (e.g., quantified effectiveness of PHSM), we tried to indirectly capture its impact via estimating the effective reproduction number using several explanatory variables, including mobility.

本研究にはいくつかの技術的制約があった。(略)研究期間中に日本では(略)3回、非常事態宣言が発令された。緊急事態宣言の具体的な変数(例えば、定量化した公衆衛生的・社会的対策の効果)をモデルに組み込むのではなく、活動度を含むいくつかの説明変数を用いて実効再生産数を推定することで、間接的にその影響を捉えようとした。

この一節は、西浦論文が「何も感染症対策を変えなければ」という非現実的な想定で計算した事情を説明した重要な箇所です。しかしこの記述は、曖昧な表現を使っていると私には見えます。以下の (5-2)と(5-3)で、これら記述の意味をより明確に示します。ここを明確にすると、西浦論文が非現実的な想定で計算した事情が見えてきます。

(5-2) 西浦論文が採らなかった方法は、計算不可能

(5-1) の論文引用箇所 は、

緊急事態宣言の効果を定量化してモデルに組み込むのではなく、実効再生産数を推定することで、宣言の影響を捉えた。

などと述べています。これは、

緊急事態宣言の効果を定量化してモデルに組み込むこともできたのだが、その方法は選ばずに、実効再生産数を推定して宣言の影響を捉える方法を採った、

とも読めそうですが、そう読むべきではなく、

緊急事態宣言の効果を定量化してモデルに組み込むことはできないので、実効再生産数を推定して宣言の影響を捉える方法を採った、

と読むべきです。

なぜならまず、引用箇所 にある緊急事態宣言の効果の定量化は不可能です。西浦論文は、それが不可能なのを承知の上でこう書いていると思われます。また不可能なのはこれだけではありません。このような定量化の方法(=各種の対策の効果を定量化し、それを積み上げて対策全体の効果を計算する方法)を採るなら、国だけでなく企業や個人が行うあらゆる対策の効果、ならびにその継時的変化を定量化する必要が少なくともあります。これらはさらに不可能です。

(5-3) 西浦論文の「実効再生産数を推定して宣言(など)の影響を捉える方法」は「何も感染症対策を変えなければ」

(5-1) の論文引用箇所 にある、

実効再生産数を推定して宣言(など)の影響を捉える方法

というのは、

推定した実効再生産数になるような社会の状況を仮定する方法

という意味であり、つまり、

推定した実効再生産数になるような宣言などの感染対策を、社会が取ったと仮定する方法

です。つまりこの記述は、西浦論文が計算の前提として

何も感染症対策を変えなければ

を採用したことを示しています。

(5-4) 西浦論文の引用箇所が示す内容

結局、本稿 (5-1) に示した論文の引用箇所 は、

  • 緊急事態宣言(などのあらゆる感染症対策の)効果を定量化してモデルに組み込むことは不可能なので((5-2))、
  • 「実効再生産数を推定して感染症対策の効果を組み込む方法」つまり「何も感染症対策を変えなければ」を採用した((5-3))、

という事情を説明していることになります。

(6)西浦論文を紹介した記事の、不適切な記述

(3) (4) (5) で述べたように、西浦論文は非現実的想定での計算であり、特殊な計算です。 この論文を報じた記事がありました。

2023-11-16 コロナワクチンで死者9割以上減 京都大チームが推計
https://nordot.app/1097790003285230554

この記事は、論文が計算した期間の実際の感染と死者数を述べ、その後に論文の内容を「ワクチンがなければ、それぞれ約6330万人と約36万人に達した恐れがあるとしている。」と説明しています。

この記述は、不適切です。
「恐れがある」という表現では、論文の想定が非現実的だという意味が伝わらないからです。この点でこれは誤った記事と言えますし、少なくとも誤解を招くものです。

「コロナワクチンで死者9割以上減 京都大チームが推計」という記事タイトルも少なくとも誤解を招くものです。「9割以上減」というのは、非現実的想定での計算上の感染被害と、現実の感染被害とを比較しています。しかし既に述べたように、計算上の感染被害は決して実現することのない非現実なのですから、これを現実の感染被害とこのような簡潔な記事の中で比較すべきではありません。
両者が比較されているのを読んだ読者が、西浦論文の計算結果(実現することのない非現実)を、実現するかも知れない状況だと誤解する可能性があります。

私はこの記事は少なくとも不適切だと思います。

(7)西浦氏による、不適切な記述のある記事の紹介

西浦氏は、前項に示した不適切な記述のある記事を以下のように紹介しました。

共同通信からの配信です。ぜひ多くの方に読んでいただきたいニュースです「新型コロナワクチンで死者9割以上減 京都大推計」
https://twitter.com/nishiurah/status/1725101388987019748 (2023-11-16)

共同通信からの配信です。ぜひ多くの方に読んでいただきたいニュースです「新型コロナワクチンで死者9割以上減 京都大推計」 https://sankei.com/article/20231116-3ICUWGMN3FMVZNSKJJGKPB4ARE/ @Sankei_news から

そもそも西浦論文が非現実的な特殊な想定をしていることは、(3) (4) (5) に述べた通りです。
実現することがないと分かっている状況を計算機の中に仮想的に作り、そこでの感染者数や死者数を計算したのです。

西浦氏が「知的体操」と説明したように、この研究は専門家同士で議論する特殊なものであり、関連の報道を多くの人の向けて「ぜひ読んでいただきたい」 などとするものではありません。
しかもこの記事の内容には、(6) に述べたように不適切なところがあるのですから、なおさら多くの人にこの記事を紹介すべきではありません。

その事情を最も良く理解しているはずの西浦氏が、上のように多くの人の向けて記事を紹介したのは、非常に不可解であり、不適切です。

(8)まとめ

本稿で述べたのは以下です。

本稿は、私が twitter.com の以下のスレなどで書いたことをまとめて加筆したものです。ただし(4)(5)は、新たに書いた内容です。

(99)更新履歴

  • 2023-12-27
    公開。
  • 2023-12-29
    全体の体裁と細かな語句を修正。
  • 2024-01-05
    (3-5) を追加。関連する事項を加筆。
  • 2024-01-06、2024-01-08
    構成を変更し、(3-5) を(4)に。(5)を追加。関連する事項を修正、加筆。
  • 2024-01-09
    (3-1) の一節、「例えば「ワクチン接種に出掛ける」との行動についても・・・」に取り消し線を引きました。引いた理由は (3-1) に書きました。


東洋経済の実効再生産数の簡易計算は、一般的な実効再生産数とは異なる数値を計算している


私がネット上でしていることの まとめ は、こちらに。
https://sarkov28.hatenablog.com/entry/2022/03/29/160915



東洋経済の実効再生産数の簡易計算は、一般的な定義の実効再生産数とは異なる数値を計算しています。

(1)本エントリで示すこと

東洋経済サイトが計算している実効再生産数の簡易計算は、一般的な定義の実効再生産数とは異なる数値を計算しています。この簡易計算は、京都大学大学院の西浦博教授の監修によるものです。

以下に、両者が異なる数値であることを示します。

また、両者の関係式を導出します。

また、この関係式が正しいことを数値計算による検算で確認します。

(2)「一般的な定義の実効再生産数」と「東洋経済サイトの簡易計算における実効再生産数」

「一般的な定義の」というのは、ご存じのような、

  • (a1) 何らかの感染対策を実行しているか、あるいは、
  • (a2) 社会集団の一部が感染症への免疫を保持している場合に、
  • (a3) 1人の感染者が、
  • (a4) 感染力を保持する間に、
  • (a5) 平均的に何人の新規感染者を生むのか、

を示す指標としての実効再生産数です。

単に  Rt と書くと区別ができなくなるので、以下ではこれを  sRt と書きます。

東洋経済サイトが計算している実効再生産数の簡易計算を  toyoRt とします。

 toyoRt の定義は、後ほど示します。

(3)以下の議論の前提

以下の議論は、感染拡大初期を前提とします。

シンプルな SIR モデル https://sarkov28.hatenablog.com/entry/2021/01/04/113012 では、感染拡大初期の感染者数  I(t) は、時刻  t の指数関数、式(1) で近似できます。

 I(t) = C e^{kt} ...... (1)

ここで、 C = I(0)、k = \beta S(0) - \gamma、\beta = R_0 / (N T) です。...... (2)

 I(0)、S(0) t=0 における感染者数と感受性者数、 \beta、\gamma は SIR モデルにおける定数、 R_0 は基本再生産数、 N は総人口、 T は平均世代時間です。

感染拡大初期の感染者数  I(t) を指数関数、式(1) で近似する件については、「マルサス法則の周辺について。」とした連ツイ https://twitter.com/sarkov28/status/1405849050524524552 で検討しました。
感染者数の指数関数近似の検討

式(1)(2) の導出もこの連ツイの途中、https://twitter.com/sarkov28/status/1405851265272205317 以下に(数理疫学者の稲葉氏による説明 [1] とともに)示してあります。
感染者数の指数関数近似に関する稲葉氏の説明

式(1)(2) は、私が導いたものではなく、上に示した稲葉氏の説明にある式

 I(t) = I(0) e^{(\beta S(0) - \gamma)t}

を分けて書いたものです。

(4)東洋経済の簡易計算の定義

 t における新規感染者数を、 newI(t) とします。

 C_a(t) を、

 \displaystyle C_a(t) = \sum_{d=0}^6{newI(t-d)} ...... (3)

とすると、東洋経済サイトが示している簡易計算  toyoRt は、定義

toyoRt=(直近7日間の新規陽性者数/その前7日間の新規陽性者数)^(平均世代時間/報告間隔)、(報告間隔= 7、平均世代時間= T

から、

 \displaystyle toyoRt = \left(\frac{C_a(t)}{C_a(t-7)}\right)^{\frac{T}{7}} ...... (4)

と書くことができます。

(5) toyoRt と sRt の関係式の導出

感染拡大初期なので、新規感染者数は感染者数の差で近似できます。

 newI(t) \fallingdotseq I(t) - I(t - 1) = C e^{kt} - C e^{k(t - 1)} = C e^{k(t - 1)}(e^{k} - 1) ...... (5)

この近似は、「感染拡大初期の新規回復者数はゼロと見なすことができる」としたものです。この近似に関する補足を(9)に書きました。

式(3) に式(5) を使うと、

 C_a(t)

 \displaystyle = \sum_{d=0}^6{newI(t-d)}

 \fallingdotseq I(t) - I(t - 7)  = C e^{kt} - C e^{k(t - 7)}  = C e^{kt} (1 - e^{-7k}) ...... (6)

となります。式(4) に式(6) を使うと、

 toyoRt

 \displaystyle = \left(\frac{C_a(t)}{C_a(t-7)}\right)^\frac{T}{7}

 \displaystyle \fallingdotseq \left(\frac{C e^{kt}(1 - e^{-7k})}{C e^{k(t - 7)}(1 - e^{-7k})}\right)^\frac{T}{7}  = (e^{7k})^\frac{T}{7} = e^{kT} ...... (7)

となります。

検討しているような SIR 系では、 \gamma = 1 / T なので、式(2) から  k = \beta S(0) - 1/T \beta = R_0 /(N T) です。したがって、

 \displaystyle kT = \beta S(0) T - 1 = \frac{S(0) R_0}{N} - 1 ...... (8)

です。

さらに流行初期なので  S(0) = N とできるので、式(8) は

 kT = R_0 - 1 ...... (9)

となります。

式(7) に 式(9) を使うと、

 toyoRt = e^{R_0 - 1} ...... (10)

となります。

流行初期の実効再生産数  sRt は、 R_0 に一致するので、 R_0 = sRt とできます。

式(10) で  R_0 = sRt とすれば、 toyoRt sRt の関係式、

 toyoRt = e^{sRt - 1} ...... (11)

が得られます。

なお、式(11) を  sRt について解くと、

 sRt = log_e(toyoRt) + 1 ...... (12)

となります。

式(11) のグラフは以下になります。

toyoRt と sRt の関係
  graph.1a  toyoRt sRt の関係のグラフ(修正版)
感受性者がほとんどを占める状態における、私が正しいと考える実効再生産数  sRt(横軸)と、東洋経済の実効再生産数の簡易計算  toyoRt(縦軸)の関係、ならびに検算(横軸  R0、縦軸  toyoRt)の結果を示している。
青線:私が導出した関係式のグラフ  toyoRt = exp(sRt - 1) sRt について解くと  sRt = log_e(toyoRt) + 1 となる。
赤丸:SIR モデルを使った検算の結果(青線にほぼ一致)
緑線: toyoRt = R0 のグラフ。 toyoRt が実効再生産数の簡易計算として正しければ、検算の赤丸はこの上(付近)に乗る。しかし、赤丸は青線の上に乗っている。
感受性者がほとんどを占める状態にいては、東洋経済の実効再生産数(西浦氏監修の式)は、数値が高めに出る。

(6)まとめ

 toyoRt が、 sRt の簡易計算なら、 toyoRt = sRt か、これに近い関係が必要です。

しかし実際には、上記の導出に用いた式(1) と式(5) の近似が成立する場合には、式(11)

 toyoRt = e^{sRt - 1}

が成立していて、 toyoRt = sRt やそれに近い関係はありません。

東洋経済サイトが計算している実効再生産数の簡易計算  toyoRt は、一般的な定義の実効再生産数  sRt とは異なる数値を計算していることが分かりました。

なお、以下が言えます。

  • (6-1)  sRt = 1 の場合には、式(11) から  toyoRt = 1 になるので、 toyoRt = sRt となります。
  • (6-2)  sRt toyoRt が 1 に近い場合(例えば、 0.5~1.5 程度)には、式(11) から両者の値は近くなります。
  • (6-3)  toyoRt が 1 に近くない場合、例えば、新型コロナの基本再生産数  R_0 = 2.5 付近の  sRt = 2.5 では、 toyoRt = 4.48 となり、両者はとても近い値とは言えません。

(6-1)、(6-2)、(6-3)は、下のグラフ graph.1a で確認できます。

toyoRt と sRt の関係
  (再掲)graph.1a  toyoRt sRt の関係のグラフ(修正版)

幾つかの  sRt における  toyoRt の理論値( = e^{sRt - 1})を示しておきます。

 sRt 1.0 1.2 1.4 1.6 1.8 2.0
 toyoRt 1.00 1.22 1.49 1.82 2.23 2.72
 sRt 2.2 2.4 2.6 2.8 3.0
 toyoRt 3.32 4.06 4.95 6.05 7.39
 sRt 1.7 2.5
 toyoRt 2.01 4.48

[表 1] 幾つかの  sRt における  toyoRt の値( toyoRt = e^{sRt - 1}

1.7 と 2.5 は、西浦氏が 2020-03-02 と 2020-03-19 の専門家会議で示した日本のコロナの基本再生産数の値なので計算してあります。

逆の表も示しておきます。

 toyoRt 1.0 1.2 1.4 1.6 1.8 2.0
 sRt 1.00 1.18 1.34 1.47 1.59 1.69
 toyoRt 2.2 2.4 2.6 2.8 3.0
 sRt 1.79 1.88 1.96 2.03 2.10
 toyoRt 1.7 2.5
 sRt 1.53 1.92

[表 2] 幾つかの  toyoRt における  sRt の値( sRt = log_e(toyoRt) + 1

(7)まとめの補足その1(式(11)が正しいことは数値計算による検算で確認されている)

式(11)

 toyoRt = e^{sRt - 1} ...... (11)

の導出には、幾らかの近似(式(1) と式(5))を用いています。したがってこの近似が適切なのかという疑問があります。

しかし式(11) が正しい事は、実際に計算して確認済です。

理論値(式(11))は、計算値に一致するのです。

このことは、既に示してあるグラフですが、

https://twitter.com/sarkov28/status/1405108027938594822

で述べた「検算」でご確認頂けます。

toyoRt と sRt の関係
  (再掲)graph.1a  toyoRt sRt の関係のグラフ(修正版)

このグラフで「緑の直線」が  sRt です。

「青い曲線」が  toyoRt の理論値、式(11) です。

複数の「赤い〇」は計算値であり、以下の要領で計算したものです。

これらの「赤い〇」が「青い曲線」の上に乗っていること(縦軸の値が一致していること)は、「式(11) が正しい事」を示しています。

「検算」すなわち「赤い〇」の計算要領:
例として、 sRt = 2.0 の場合の「赤い〇」の値を計算します。

  • 実効再生産数を  2.0 として、シンプルな SIR モデル
    https://sarkov28.hatenablog.com/entry/2021/01/04/113012
    で時間経過による感染者数の拡大を計算します。
  • 得られた計算結果に基づいて、東洋経済サイトで示されている  toyoRt の定義(式(4))の通りに、 toyoRt を計算します。
  • 得られた  toyoRt の値が、「 sRt = 2.0 の場合の「赤い〇」の値」です。

同様の計算を、別の  sRt の値に対して行った結果が、グラフにある複数の「赤い〇」です。

(8)まとめの補足その2(東洋経済の実効再生産数簡易計算は、感染状況を過大に見せている)

まとめに書いたように、 toyoRt は高めの値になります。

このため、東洋経済の簡易計算で計算された実効再生産数  toyoRt を見ていると、感染状況を過大に見積もってしまう弊害があったことを特に指摘しておきます。

以下のような状況を考えます。

2020年の春。西浦氏が示すように、新型コロナの基本再生産数が  R_0 = 2.5 だとします。

ある時、何らかの感染対策を行った結果、感染拡大を2割削減できたとします。

この時の実効再生産数は、定義から、あるいは西浦氏自身が説明しているところから、

 Rt = 2.5 \times (2割削減) = 2.5 \times 0.8 = 2.0 ...... (13)

になります。東洋経済の実効再生産数の簡易計算においても、この値(または近い値)が計算されるべきです。

しかし、東洋経済の簡易計算は、これに近い値を計算しません。

式(13) で  2.0 と計算した  Rt は、一般的な定義の((2)の(a1)~(a5) の意味での)実効再生産数なので、これは  sRt です。

この時、東洋経済の簡易計算  toyoRt は、式(11) から、

 toyoRt = e^{sRt - 1} = e^{2.0 - 1} = 2.71

になります。

字義通り「簡易計算」なら、東洋経済の簡易計算は、式(13) に示すように  2.0 かそれに近い値を計算してくるべきです。ところが実際には、 2.71 を出してくるのです。

非常に不適切なことに、ここでは大小関係が逆転しています。

「感染対策により、感染拡大を2割削減できている状況」であるにも関わらず、「東洋経済の簡易計算は、対策なしにおける数値( 2.5)よりも高い数値( 2.71)を示します。

感染拡大初期に東洋経済の簡易計算を見ていると、感染対策の評価を(過小評価する方向に)誤ってしまうのです。(上記の「2割削減」の例では「過小評価」どころか「逆効果的評価」になっています。)(なお現在は、それなりに免疫保持者割合が増加しているので、本エントリが想定している「感染拡大初期」には該当しないと思われます。)

東洋経済はこの不適切さをどのようにお考えなのでしょうか。

(9)式(5)の近似の補足

感染拡大初期なので、新規感染者数は感染者数の差で近似できる、というのが式(5) でした。

 newI(t) \fallingdotseq I(t) - I(t - 1) = C e^{kt} - C e^{k(t - 1)} = C e^{k(t - 1)}(e^{k} - 1) ...... (5)

(9-1) 式(5) の成立のためにどういう仮定を置いたのか

式(5) は、「感染拡大初期の新規回復者数はゼロと見なすことができる」との仮定から導かれるものです。

時刻  t と 時刻  t - 1 の感染者数については、以下の式が成立します。

 I(t) = I(t - 1) + newI(t) - newR(t) ...... (14)

ここで  newR(t) は、時刻  t における新規回復者数です。

「感染拡大初期の新規回復者数はゼロと見なすことができる」のであれば、(今考えている  t は感染拡大初期なので、) newR(t) = 0 とできます。

式(14) で、 newR(t) = 0 として、式変形すると、

 newI(t) = I(t) - I(t - 1) ...... (15)

となります。

「感染拡大初期の新規回復者数はゼロと見なすことができる」との仮定から、式(5) を導きました。

(9-2) 式(5) の成立のための仮定は妥当なのか

次に、「感染拡大初期の新規回復者数はゼロと見なすことができる」との仮定が妥当であることを示します。

  • 感染拡大初期なので、感染者数は比較的小さいです。
  • 回復者は、感染者の一部です。
  • 感染者数が比較的小さいので、感染者数の一部である回復者数はさらに小さいです。
  • したがって、「感染拡大初期の新規回復者数はゼロと見なすことができる」との仮定は妥当です。

どの程度までが「感染拡大初期」と言えるのかについては、別に検討する予定です。

(10)実効再生産数の定式化

(2)に示した「一般的な定義の実効再生産数」は定式化されています。まず、(2)で示した定義を再掲します。

  • (a1) 何らかの感染対策を実行しているか、あるいは、
  • (a2) 社会集団の一部が感染症への免疫を保持している場合に、
  • (a3) 1人の感染者が、
  • (a4) 感染力を保持する間に、
  • (a5) 平均的に何人の新規感染者を生むのか、

これは、参考文献 [3] に説明されている2つの実効再生産数のうちの1つ、瞬間的再生産数(instantaneous reproduction number)に該当しています。

瞬間的再生産数は、[3] では  R(t) と表現され、「1人の感染者あたりが時刻  t に生み出す新規感染者数の平均値」と説明されています。

 \displaystyle R(t) = \frac{i(t)}{\int_{0}^{\infty} i(t - \tau) g(\tau) d\tau} ...... (16)

 i(t) は、時刻  t における新規感染者数。

 g(\tau) は、各感染源にとっての感染後の経過時刻  \tau(=感染齢)における2次感染発生の相対的頻度。

(式(16)と、 i(t) g(\tau) の定義は、[3] の (10-5)式とその周辺から。)

(同様の記述は、参考文献 [2] の (2.113)式とその周辺にもありますが、[2] の説明はやや明確さを欠きます。)

(11)twitter における関連事項

ここに示す各連ツイは、このエントリに関連する事項を述べています。

[0.目次]
https://twitter.com/sarkov28/status/1405140064561029120

[1.概要]
https://twitter.com/sarkov28/status/1402934465261756416

[2.東洋経済の計算式はどこが間違っているのか]
https://twitter.com/sarkov28/status/1405099777688170498

[3. 西浦氏による例示その1(東洋経済の計算式で不適切になる例その1)]
https://twitter.com/sarkov28/status/1437747938713686025

[4. 西浦氏による例示その2(東洋経済の計算式で不適切になる例その2)]
https://twitter.com/sarkov28/status/1440880535002234890

[5. 西浦氏による例示その3(東洋経済の計算式を使っている例)]
https://twitter.com/sarkov28/status/1440972220856078341

[6. (本エントリを紹介する予定の連ツイ)]
(まだですが、たぶんもうすぐ)

[7. (「2つの実効再生産数と、東洋経済の簡易計算」みたいな、一連のまとめっぽいことを書く予定の連ツイ)]
(まだ)

(12)導出した関係式などは、報告間隔や平均世代時間に依存しない

本エントリでは、西浦氏監修の東洋経済の実効再生産数の簡易計算に沿って論じたため、報告間隔を  7 としていました((4)の式(4))。

しかし明らかに、本エントリが導出した関係式である式(11) や式(12) やそれに関連する事項は、他の報告間隔を採用した同様の計算でも成立します。これらは報告間隔に依存しません。

また、本エントリでは、平均世代時間を  T としました((4)の式(4))が、本エントリが導出した関係式である式(11) や式(12) やそれに関連する事項は、平均世代時間に依存しません。

参考文献

  • [1] 稲葉寿 2008 微分方程式と感染症数理疫学
    https://www.researchgate.net/publication/341121176_weifenfangchengshitoganranzhengshuliyixue
  • [2] 稲葉寿編著「感染症の数理モデル 増補版」, 培風館, 2020
    上の(10)で言及の 2章は、西浦氏による。
  • [3] 西浦博編著「感染症疫学のためのデータ入門」, 金芳堂, 2021
    上の(10)で言及の Chapter 10 は、スンモク・ジョン氏と西浦氏の共著。

修正履歴

大きく修正した場合は、ここに書いておきます。

  • 2023-07-19
    公開。公開後に幾らか修正しました。
  • 2023-07-26
    「(9)式(5)の近似の補足」を追加し、関連する修正を行いました。
    (6)の議論の後半を整理し、(6-1)(6-2)(6-3) という形にしました。
    「本エントリで用いている近似が、式(1) と式(5) であること」を、明確にしました。
    この他、細かな表記を修正しました。
  • 2023-07-28
    「(10)実効再生産数の定式化」を追加し、関連する修正を行いました。
    参考文献に [2]、[3] を追加しました。
  • 2023-07-30
     sRtを「一般的な定義の」「通常言われている定義の」などと不統一な表現で説明していました。これを「一般的な定義の」という表現に統一しました。
  • 2024-03-28
    • 初出の graph.1a の下に、図中に書いたのと同様の説明を加筆しました。
    • (6)の、 sRt の値に対応する  toyoRt の値の表を、[表 1] としました。
      [表 2] として  toyoRt の値に対応する  sRt の表を新たに作成しました。
  • 2024-04-06
    • (4)で  toyoRt の定義を twitter への参照で示していましたが、直接書くようにしました。内容は同一です。
    • (12)を追加しました。


基本再生産数 R0 が小さくなると集団免疫閾値が小さくなるのはなぜか(集団免疫閾値の公式の導出)


私がネット上でしていることの まとめ は、こちらに。
https://sarkov28.hatenablog.com/entry/2022/03/29/160915



西浦氏の「42万死亡推計」などの SIR 系のモデルでは、

(a) 小さな基本再生産数  R_0 では、集団免疫閾値が小さな値になります。

本エントリでは、なぜそうなるのかを四則演算だけを使って説明します。
(私が独自に見出した説明をするのではなく、他で示されている説明を繰り返すだけです。詳細は(4)に。)

用語を確認しておきます。

集団免疫閾値とは、「(集団内に免疫保持者が増え、)その時点での感染者数がそれ以上増えなくなった時の、免疫保持者の割合」です。(「累積感染者数がそれ以上増えなくなった時」ではありません。)

基本再生産数  R_0 とは、「感染拡大開始時(=全ての国民が免疫を保持していない時)に、一人の感染者が感染させる人数の平均値」です。

説明の前提などは後ろの方に書いてあります。


目次


(1)集団免疫閾値とは「どの辺」でのことなのか

集団免疫閾値に到達するのは、感染の波における「どの辺」でのことなのでしょうか。
波の終盤で集団免疫閾値に到達するのではありません。
中盤で到達します。

その「終盤」「中盤」という感じを見て頂くために、グラフを示します。

集団免疫閾値というのは、「(集団内に免疫保持者が増え、)その時点での感染者数がそれ以上増えなくなった時」なので、感染者数のピークで到達します。
「終盤」ではなく、割と「中盤」です。

fig.1 「何もしなければ」での集団免疫閾値到達は、波の「中盤」です。

このグラフ付近には、集団免疫閾値に関連する事項が他にもありますので、もう少し説明します。

上の fig.1 は、感染者数だけのグラフでした。これに感受性者数などを加えます。

感受性者とは、「免疫を持っていなくて、感染していない人」です。

ここでは、流行開始時点では、(ほぼ)全人口が感受性者だと考えます。
感染した人は感受性者ではありません。
感染して回復した人は感受性者ではありません。免疫があると考えるので。

fig.1 に感受性者などを加えたのが、fig.2 です。fig.1 での感染者数の赤曲線は、同じ色の「小さな波」で描かれています。
「こんなに小さい波なのか」と思われるかも知れませんが、この赤曲線は「最大の時点では、全人口の2割(0.2)以上が感染者となる」と主張しているのですから、この波は小さくありません。

fig.2 「新規感染者数=新規回復者数」の時が、集団免疫閾値への到達です。

この fig.2 で重要な部分は下の方にあるので、そこを拡大して fig.3 とします。

fig.3 「新規感染者数=新規回復者数」の時が、集団免疫閾値への到達です。

このグラフでは、〇で強調した箇所が縦に3つ並んでいます。
これは、 R_0=2.5 の場合に、

  • (1) 新規感染者数と、新規回復者数が等しくなる時、
  • (2) 感染者数がピークになる時、
  • (3) 感受性者割合が 0.4 になる時(免疫保持者割合が 0.6 になる時)、

が一致していることを示しています。
fig.3 の (3) で示した「縦の黒点線と青曲線との交点」は、「集団免疫閾値が 0.6、つまり 60% であること」を示しています。
西浦氏も 2020-04-07 の twitter 動画 で、「 R_0=2.5 での集団免疫閾値は60%」という旨を説明しています。
「[tex: R_0=2.5] での集団免疫閾値は60%」との西浦氏による説明
2020-04-07 twitter 動画 から。)

次の(2)で、西浦氏が書いている公式「  集団免疫閾値=1-1/R_0」を導出します。

(2)四則演算だけを使った、集団免疫閾値の公式の導出

「基本再生産数と、集団免疫閾値の関係」つまり「 集団免疫閾値=1-1/R_0」の公式を、四則演算だけを使って導出します。

(2-1)実効再生産数と基本再生産数に関する式の説明

まず前段階として、実効再生産数と基本再生産数に関する式を説明します。

基本再生産数は、

基本再生産数  R_0 とは、「感染拡大開始時(=全ての国民が免疫を保持していない時)に、一人の感染者が感染させる人数の平均値」

でした。類似した用語に、実効再生産数があります。

実効再生産数  R_t とは、「一人の感染者が感染させる人数の平均値」

です。

基本再生産数と、実効再生産数の定義から、以下が言えます。

(a) 流行開始時点の、(ほぼ)全人口が感受性者の時の実効再生産数  R_t は、基本再生産数  R_0 に一致します。

また、感染は、感染者と感受性者の接触で起きるので、感受性者が減れば、感染は減ります。
実効再生産数が感染者一人当たりの数値であることを考えると、以下が言えます。

(b) 感受性者が(流行開始時点の大きさから)1% 減ると、実効再生産数  R_t は(流行開始時点の大きさから)1% 減ります。

この(a)(b)を満たす実効再生産数  R_t は、以下の式になります。

(c)  R_t = (1-ε)R_0、( ε は、人口中の免疫保持者の割合)

(c)が(a)(b)を満たすことを確認します。

まず、(c)が(a)を満たすこと。

流行開始時点では、ほぼ全人口  N が感受性だと考えるので、免疫保持者はいません。したがって免疫保持者の割合  ε は、 ε=0 になります。
すると、(c) は、 R_t = (1-0)×R_0=R_0 となり、 R_0 に一致します。
つまり、(c)は(a)を満たしています。

次に、(c)が(b)を満たすこと。

(c)  R_t = (1-ε)R_0 の中の  (1-ε) というのは、「免疫保持者ではない人の割合」なので、「感受性者の割合」です。
「感受性者が 1% 減る」というのは、 (1-ε) 0.01 減るということです。
この場合、(c) の  R_t は 1% 減ります。
つまり、(c)は(b)を満たしています。

以上により、

(c)  R_t = (1-ε)R_0、( ε は、人口中の免疫保持者の割合)

が、(a)(b) を満たしていることが分かりました。
(この式は、私独自のものなどではありません。(4)をご参照下さい。)

(2-2)集団免疫閾値の公式の導出

基本再生産数と集団免疫閾値に関する公式を導出します。

(2-1)での検討の結果、

(c)  R_t = (1-ε)R_0、( ε は、人口中の免疫保持者の割合)

であることが分かっています。

集団免疫閾値というのは、「その時点での感染者数がそれ以上増えなくなった時」の「免疫保持者の割合」でした。
「その時点での感染者数がそれ以上増えなくなった時」の実効再生産数は 1 です。つまり  R_t=1 です。
つまり (c) から、
 R_t = (1-ε)R_0 = 1
の時に、免疫保持者の割合がどうなっているかを考えればいいのです。

 (1-ε)R_0 = 1 を変形していくと、
 1-ε = 1/R_0
 ε = 1-1/R_0
となります。
これは、「その時点での感染者数がそれ以上増えなくなった時」の「免疫保持者の割合」を示しています。

つまり、これが求めたかった公式、

 \displaystyle 集団免疫閾値=1-\frac{1}{R_0}

です。

西浦氏がホワイトボードに書いていた式はこれです。(西浦氏は、免疫保持者の割合として  ε ではなく  p を、実効再生産数として  R_t ではなく  R_e を使っています。)

「[tex: R_0=2.5] での集団免疫閾値は60%」との西浦氏による説明
2020-04-07 twitter 動画 から。)

(3)小さな基本再生産数  R_0 では、集団免疫閾値が小さな値になること

小さな基本再生産数  R_0 では、集団免疫閾値が小さな値になることを示します。

(2)での検討により、公式

 \displaystyle 集団免疫閾値=1-\frac{1}{R_0}

が分かっています。
幾つかの  R_0 における集団免疫閾値を計算してみましょう。

 R_0 集団免疫閾値
=1-1/R_0
1.4 0.286
1.5 0.333
1.7 0.412
2.0 0.500
2.5 0.600
3.0 0.667
3.5 0.714
4.0 0.750

table.1 幾つかの  R_0 における、集団免疫閾値

公式「  集団免疫閾値=1-1/R_0」と table.1 は、「小さな基本再生産数  R_0 では、集団免疫閾値が小さな値になる」ことを示しています。
本エントリで説明しようしていたのはこのことです。

(4)稲葉寿氏による、(2)(3)と同じ説明

上の(2)(3)は、稲葉寿氏と同じ説明をしているだけですので、貼っておきます。
稲葉寿氏による(2)(3)と同じ説明。2020-12-18「感染症数理モデルとCOVID-19」から
2020-12-18「感染症数理モデルとCOVID-19」(https://www.covid19-jma-medical-expert-meeting.jp/topic/3925)から。

(90)詳細な事項

他の値でも良かったのですが、グラフの計算は、

  • 基本再生産数  R_0 R_0=2.5
  • 人口は、約  1.27 億人、

と、「42万死亡推計」と同じものを使いました。また、

  • 「何もしなければ」での感染の波を考えています。
  • 計算に使用したモデルは「シンプルな SIR モデル」です。
    異質性のない個体を想定しています。
    西浦氏が「42万死亡推計」に使用した「3世代型 SIR モデル」ではありません。

「シンプルな SIR モデル」は、SIR モデルが説明される時にしばしば最初に示される、最もシンプルな、良く知られたモデルであり、特殊なものではありません。上に挙げた稲葉氏の論文にも例示されています(https://www.covid19-jma-medical-expert-meeting.jp/topic/3925 の p.06、式(1))。

これらモデルについては、

に少し説明を書きました。

また、単に免疫と書いていますが、これは感染予防免疫を指しています。
現実の感染予防免疫は時間が経つにつれて減じることがありますが、本エントリでは議論を単純にするため、時間が経過しても効果が減じることは考えていません。

 ε を「人口中の免疫保持者の割合」としていますが、文脈を見て頂ければ明らかなように、これは回復者だけでなく感染者を含んでいます。

修正履歴

  • 2023-05-27
    公開
  • 2024-03-23
    table.1 の縦と横を入れ替えました。
     R_0 = 1.4 の値を追加しました。2020-03-02 の専門家会議の資料で挙がっていた値だからです。
    table.1 の集団免疫閾値の項目の有効数字を、2桁から3桁にしました。


免疫や感染症対策が、最終的な感染者数(最終規模)や集団免疫閾値に与える影響


私がネット上でしていることの まとめ は、こちらに。
https://sarkov28.hatenablog.com/entry/2022/03/29/160915



免疫や感染症対策が、最終的な感染者数に与える影響を検討します。
この際、「免疫での実効再生産数の低下」と「感染症対策での実効再生産数の低下」の違いを検討します。
また、感染症対策が、集団免疫閾値に与える影響を検討します。

最終的な感染者数(の比率)は「最終規模 final size」と呼ばれるので、以下ではそう書きます [稲葉 2008]。
つまり以下では、免疫や感染症対策が、最終規模に与える影響について検討します。

このブログエントリの全てのグラフは、google colab スクリプト
https://colab.research.google.com/drive/1XUk08QtkniHMoYI7z_YbfYC3MhDauDv6?usp=sharing
で計算、描画しました。

このブログエントリに関連する事項を、twitter に連ツイ https://twitter.com/sarkov28/status/1646049139594088451 として書きました。ただしこれはこのブログエントリの概要であり、こちらに書いてない事項はないと思います。


目次


(1)議論の前提

以下で行う議論の前提です。

以下で「免疫」というのは、「感染予防に有効なワクチンによる免疫」や「感染による免疫」を指します。
2023-03 時点にける新型コロナのワクチンは、「感染予防に有効なワクチン」には該当しません。また新型コロナにおいては、「感染による免疫」も時間経過により効果が減衰すると言われています。
以下では、議論を簡単にするため、こうした点を単純化して考えます。
すなわち「ワクチン接種などの方法で感染予防の免疫が得られる」「得られた免疫の効果は時間経過で減衰することはない」とします。

また以下では、主に「シンプルな SIR モデル https://sarkov28.hatenablog.com/entry/2021/01/04/113012」で考えます。
私は、まだ詳しくは書いていませんが、日本国民全体のような大規模集団を対象として、こうしたモデルで感染シミュレーションを計算することは机上の空論であり不適切だと考えています。このモデルで計算するのは、議論を簡単にするためです。机上の空論であることを承知の上でなら、ここでの計算からも得られるものはあると考えます。
(西浦氏が「42万死亡推計」で採用した「3世代型 SIR モデル https://sarkov28.hatenablog.com/entry/2021/01/25/170000」を使った計算も同様に机上の空論であり不適切です。このモデルを用いても、以下の議論は概ね、同様に成立します。ただし、モデルがやや複雑になる分だけ、議論がやや複雑になってしまいます。)

(2)記号の意味

以下で用いる記号の意味です。

  •  S : 感受性者数。時刻  t では  S(t)
  •  I : 感染者数。時刻  t では  I(t)
  •  R : 回復者数。時刻  t では  R(t)
  •  N : 全人口。( S + I + R = N です。)
  •  R_0 : 基本再生産数。
  •  Re(t) : 時刻  t における実効再生産数。

実効再生産数は、 R_t R(t) と書くこともありますが、本エントリでは  Re(t) と書きます。時刻  t における回復者数  R の表現  R(t) と区別するためです。

(3)基本再生産数が変化した場合の最終規模の理論値

本題に入ります。
基本再生産数が変化した場合の最終規模の理論値を検討します。

「免疫や感染症対策が、最終的な感染者数に与える影響」について検討する準備として、「シンプルな SIR モデル」で、「何もしなければ」の条件下で基本再生産数が変化した場合の最終規模を検討します。(感染症対策がなく、感染拡大開始時には免疫もない場合です。)
この場合の基本再生産数  R_0 と最終規模  p (全人口中の最終的な感染者数の割合)とは、最終規模方程式(final size equation)と呼ばれる式(1)
 1 - p = e^{- R_0 p} ...... (1)
を満たします。
この式から、 R_0 の値に対する最終規模  p の理論値を計算できます。
([稲葉 2008]。式(1)の導出も説明されています。)


適当な範囲の(小刻みな) R_0 と、その  R_0 から計算される最終規模  p とをプロットすると、以下のようなグラフになります。

fig.1-1 R0 を変化させた時の最終規模の変化(最終規模方程式から得た理論値)
fig.1-1  R_0 を変化させた時の最終規模の変化(最終規模方程式から得た理論値)

 p R_0 の増加につれて、単調に増加します。([稲葉 2009] には、fig.1-1 と同様のグラフがあります。)
幾つかの  R_0 に対する  p の値を示しておきます。

基本再生産数  R_0 1.5 2.0 2.5 3.0 3.5
最終規模  p 0.58 0.80 0.89 0.94 0.97

table.1 幾つかの  R_0 における最終規模の値(シンプルな SIR モデルの理論値)

「3世代型 SIR モデル」などの異質性を導入したモデルには、式(1) で示したような理論値は一般に存在しないと思います。

(4)基本再生産数が変化した場合の最終規模の計算値

基本再生産数が変化した場合の最終規模の計算値です。

前節と同じ条件における最終規模を、SIR 方程式をシミュレーションすることによって計算してみます。
具体的には、

  • ある範囲の(小刻みな) R_0 において、
  • 時刻 0 での初期値から、 S I R をそれぞれの時間 t に対して計算し、
  •  t が十分に経過した時点での  p = R / N(=(感染した人数)/(全人口))を計算して、
  •  R_0 p をプロットする、

ことになります。
fig.1-3 が得られます。(比較のため、再掲の fig.1-1 と共に示します。)

fig.1-1 R0 を変化させた時の最終規模の変化(最終規模方程式から得た理論値)と fig.1-3 R0 を変化させた時の最終規模の変化(SIR 方程式の計算結果)
fig.1-1  R_0 を変化させた時の最終規模の変化(最終規模方程式から得た理論値)と
fig.1-3  R_0 を変化させた時の最終規模の変化(SIR 方程式の計算結果)

(図の番号は連続していません。図の番号は、グラフを作成した google colab スクリプト [sarkov28 2023] における番号を示しているためです。)

fig.1-3 は、fig.1-1 と一致していて、視覚上の差異はありません。

こちらでも、幾つかの  R_0 に対する p の値を示しておきます。

基本再生産数  R_0 1.5 2.0 2.5 3.0 3.5
最終規模  p 0.58 0.80 0.89 0.94 0.97

table.2 幾つかの  R_0 における最終規模の値(シンプルな SIR モデルの計算値)

(計算方法が全く違うのに)前節 table.1 の値に一致しています。(もっと下の桁では不一致もありますが、計算上の誤差と考えます。)


実効再生産数  Re(t) は、流行開始時には  R_0 に一致していますが、感受性者が減少するにつれて低下していきます。

 Re(t) = (S(t) / N) \ R_0 ...... (2)

 Re(t) = 1 となる時がいわゆる集団免疫の状態であり、二次感染者数は増加しなくなります。

 Re(t) = (S(t) / N) \ R_0 = 1 なのですから、

 S(t) / N = 1 / R_0 ...... (3)

であり、感受性者の比率  S(t) / N 1 / R_0 になった時です。
集団免疫に関するよく見かける式は、感受性者の比率ではなく、(接触削減で達成される場合もある)感受性者の減少率で論じたものであり、1 から式(3)の両辺を引いて、

 (接触削減などによる感受性者の減少幅) = (N - S(t)) / N = 1 - 1 / R_0 ...... (4)

を集団免疫に必要な接触削減量として示したものです。
式(4) は、 R_0 = 2.5 であれば、 1 - 1 / 2.5 = 0.6 となります。

(5)感染症対策によって実効再生産数が低下した場合の、最終規模の計算値

感染症対策によって実効再生産数が低下した場合の、最終規模の計算値を検討します。

西浦氏らは「接触削減などの感染症対策は実効再生産数を低下させる」と想定しています。様々な異論があるところですが、ここではこの想定を受け入れて、検討を続けます。

「接触削減などにより、再生産数が相対的に  q 減少する」とした場合の実効再生産数  Re(t) の時刻  t = 0 での値、 Re(0) は、以下になります。

 Re(0) = (1 - q) \ R_0 ...... (5)

一般の時刻  t では、式(2) から

 Re(t) = (1 - q) \ (S(t) / N) \ R_0 ...... (6)

です。(流行期間中で  q が一定であるとの非現実的な想定をしています。)
この  Re(0) で感染拡大が開始した時の最終規模は、前節「基本再生産数が変化した場合の最終規模の計算値」と同じグラフ fig.1-3 になります。同じですのでグラフは略します。

(6)免疫によって実効再生産数が低下した場合の、最終規模の計算値

免疫によって実効再生産数が低下した場合の、最終規模の計算値を検討します。

感染拡大開始時点において、

  • 全人口が感受性者ならば実効再生産数は、基本再生産数になる( Re(0) = R_0)が、
  • 全人口の一部に免疫があり、その分、実効再生産数は低いところから開始する、

という状況を考えます。
(ここでの免疫というのは、「感染を予防する免疫」という意味です。)
(上の記述は、基本再生産数、実効再生産数という言葉の本来の定義からは外れていますが、検討の便宜上、このように書きます。)


前節で示したように、「接触削減などにより、再生産数が相対的に  q 減少する」とした場合の実効再生産数  Re(t) は、時刻  t = 0 と一般の時刻  t では、

 Re(0) = (1 - q) \ R_0 ...... 再掲(5)

 Re(t) = (1 - q) (S(t) / N) \ R_0 ...... 再掲(6)

となります。

ある割合が免疫を持っているために、時刻  t = 0 での感受性人口が  u \ N だとします( u t = 0 での感受性者比率であり、 0.0 \leq u \leq 1.0)。

この場合の実効再生産数  Re(t) は、時刻  t = 0 と一般の時刻  t で、

 Re(0) = (1 - q) \ u \ R_0 ...... (7)

 Re(t) = (1 - q) \ u \ (S(t) / N) \ R_0 ...... (8)

となります。

この  Re(0) での最終規模を、 u を変化させながら計算したのが以下のグラフです。
(「 t = 0 で(ワクチンなどにより)ある割合が免疫を持っている」としたのに、その後の免疫獲得は感染によるものだけを考える(=ワクチンによるものは考えない)という非現実的な想定をしています。)
fig.1-3 と比較したいので  R_0 は同じ値  R_0 = 2.5 としました。
接触削減などの対策の影響がない場合を見たいので  q = 0 としました。

fig.1-5 感受性者比率の初期値 u を変化させた時の最終規模の変化(SIR 方程式の計算結果)
fig.1-5 感受性者比率の初期値  u を変化させた時の最終規模の変化(SIR 方程式の計算結果)

幾つかの  u に対する p の値を示しておきます。

感受性者比率の初期値  u 0.6 0.7 0.8 0.9 1.0
最終規模  p 0.35 0.50 0.64 0.77 0.89

table.3 幾つかの感受性者比率の初期値  u における最終規模の値

 u = 1.0 の時、 p = 0.89 です。
table.1、table.2 では、 R_0 = 2.5 の時、 p = 0.89 でした。
2つの  p の値が一致しているのは偶然ではありません。
table.1 の 0.89 は、「全人口が感受性者だった場合の  R_0 = 2.5 での最終規模」です。
table.3 の 0.89 は、「 R_0 = 2.5 で、全人口のうちの、 u = 1.0 が(=全部が)感受性者の時の最終規模」です。
この2つは同じですから、一致するのが当然です。

しかし気になることがあります。
この fig.1-5 を、上の fig.1-3(fig.1-1 は fig.1-3 と同じ)を比べると、曲がり方が違います。
次節では、この違いについて検討します。

(7)fig.1-3 と fig.1-5 の違いについて

fig.1-3 と fig.1-5 の違いについて検討します。

前節で述べたように、全人口の中のある割合が免疫を持っていて、感染拡大開始時点での感受性人口が  u \ N である場合、実効再生産数  Re(t) の時刻  t = 0 での値、 Re(0) は、以下になります。

 Re(0) = (1 - q) \ u \ R_0 ...... 再掲 (7)

接触削減などによる効果  q と感受性者比率の初期値  u が以下の(条件a)と(条件b)の場合、この  Re(0) は同じ値  0.6 \ R_0 になります。

 q  u  Re(0) = (1 - q) \ u \ R_0 最終規模
(table.2 での)条件a  0.4  1.0  0.6 × 1.0 × R_0 = 0.6 \ R_0 0.58
(table.3 での)条件b  0.0  0.6  1.0 × 0.6 × R_0 = 0.6 \ R_0 0.35

table.4 同じ  Re(0) になる  q u の値の例

前節(6)では、 R_0 = 2.5 で考えましたので、この値を使うと  Re(0) = 0.6 \ R_0 = 1.5 となります。
興味深いのは、(条件a)と(条件b)で最終規模が違うことです。感染拡大開始時の  Re(0) は同じなのに。
(条件a)では 0.58 であり、(条件b)では 0.35 です。この違いが、本節の検討対象である「fig.1-3 と fig.1-5 の違い」になっています。


(条件a) と (条件b) で SIR 方程式を数値計算し、この点を確認しました。

fig.2-1 感受性者比率の時間推移
fig.2-1 感受性者比率の時間推移

fig.2-2 感受性者比率の [tex: t = 0] からの増減の時間推移
fig.2-2 感受性者比率の  t = 0 からの増減の時間推移

fig.2-2 は、fig.2-1 の2つの曲線を下に平行移動したものです。
「fig.2-2 の右端(t = 250 付近)」と「table.4」とは、以下のように対応しています。

fig.2-2(t = 250 付近) table.4
赤曲線の右端、感受性者比率の減少が 0.58 条件a の最終規模が 0.58
青曲線の右端、感受性者比率の減少が 0.35 条件b の最終規模が 0.35

table.5 「fig.2-2(t = 250 付近)」と「table.4」との対応

fig.2-2 右端の(最終的な)感受性者比率の減少と、table.4 の最終規模(最終的な感染者数の比率)は予測されたことですが一致しています。

fig.2-1、fig.2-2 は (条件a) と (条件b) の比較でしたが、他のパラメータでも同様のことが起こります。
fig.1-3 と fig.1-5 の違いは、こうした違いのために発生しています。

この違いを言葉で書きますと、

(条件b) では、(条件a) に比べて、感染拡大開始時( t = 0)で既に感受性者比率が下がっている。このため、その後の感受性者比率の低下による集団免疫的効果が (条件a) よりも強く生じる。このため、最終規模(最終的な感染者数の比率)は小さくなる。

となります。(もっと明確に論証したかったのですが、とりあえずこれでご容赦ください。)
ここで「集団免疫的効果」というのは、いわゆる集団免疫の時点だけを指すのではなく、(それ以外の時点でも)感受性者数が少ないことにより実効再生産数が低下することを指しています。

(8)「3世代型 SIR モデル」における最終規模

「3世代型 SIR モデル」における最終規模について少し検討します。

これまでの検討は、全て「シンプルな SIR モデル」https://sarkov28.hatenablog.com/entry/2021/01/04/113012 によるものでした。
「3世代型 SIR モデル」というのは、西浦氏が「42万死亡推計」で用いたモデル https://sarkov28.hatenablog.com/entry/2021/01/25/170000 です。

「3世代型 SIR モデル」について、fig.1-3、fig.1-5 と同様のグラフを描画しました。
2つのモデルの最終規模のグラフは、table.6 のように対応しています。

グラフの内容\計算モデル シンプルな SIR モデル 3世代型 SIR モデル
最終規模の理論値 fig.1-1 fig.1-2
感染症対策で Rt が変化した時の最終規模 fig.1-3 fig.1-4
免疫で計算開始時の Rt が変化した時の最終規模 fig.1-5 fig.1-6

table.6 2種類の SIR モデルによるグラフの対応

fig.1-1~fig.1-6 計算モデルによる最終規模の違い
fig.1-1~fig.1-6 計算モデルによる最終規模の違い

「fig.1-3 と fig.1-4」「fig.1-5 と fig.1-6」のどちらにおいても、左「シンプルな」よりも右「3世代型」の方がグラフが低くなり、最終規模が小さくなり、つまり最終的な感染者数が少なくなっています。
これは「3世代型」が導入した異質性の効果です。

異質性を導入しないモデル(シンプルな SIR モデル)は、「全ての人々は均一に接触する」と仮定しています。 どんな年齢でも、どこに住んでいても、です。
「3世代型モデル」は、異質性を部分的に導入しています。世代毎の接触の違いです。
(世代数が3つでいいのか、この3つの世代でいいとしても西浦氏が導入した感受性の異質性だけでいいのか、世代以外の異質性は考慮しなくていいのか、などの論点は残っています。西浦氏もこのモデルへの異質性の導入が不十分であることを認めています [西浦 2020]。)

異質性を導入すると、感染者と感受性者の接触に偏りが生じます。すると、偏りがない場合(=均一に接触する場合)と比べて、感染者と感受性者の接触は(ある程度)生じにくくなります。
その結果、感染拡大が(ある程度)抑制されるのです。

table.2 では、「シンプルな SIR モデル」における最終規模の数値計算の値を示しましたが、これと比較して「3世代型 SIR モデル」の値を書いておきます。

使用モデルやグラフ\基本再生産数  R_0 1.5 2.0 2.5 3.0 3.5
シンプルな SIR モデル fig.1-3 0.58 0.80 0.89 0.94 0.97
3世代型 SIR モデル fig.1-4 0.52 0.71 0.80 0.85 0.87

table.7 幾つかの  R_0 における最終規模の値
シンプルな SIR モデルでの値と、3世代型 SIR モデルでの値を比較

「3世代型モデル」への異質性の導入については、こちら https://sarkov28.hatenablog.com/entry/2021/01/25/170000 や、西浦氏による 2020-06 の Newsweek の記事 [西浦 2020] に少し説明があります。
私はこの Newsweek 記事の西浦氏の記述の中には、自身の「42万死亡推計」計算を正当化するためと思われる不適切な記述があると考えていますが、この論点については長くなるので、別の機会にします。おそらく twitter( https://twitter.com/sarkov28 )に書きます。

また、この付近には、ここまでの議論をひっくり返すような重要な論点があるのですが、この論点についても長くなるので、別の機会にします。
おそらく twitter( https://twitter.com/sarkov28 )に書きますが、ブログの別エントリにまとめる予定です。

(9)基本再生産数が変化した場合の集団免疫閾値

基本再生産数が変化した場合の集団免疫閾値を計算します。

「シンプルな SIR モデル」と「3世代型 SIR モデル」で計算し、比較します。

これまで検討してきた最終規模(最終的な感染者数(の比率))と深い関係にあるのが、集団免疫閾値です。
この両者は、どちらも集団免疫の効果によるからです。

fig.3-1~fig.3-4 R0 を変化させた時の集団免疫閾値の変化

fig.3-1  R_0 を変化させた時の集団免疫閾値の変化(シンプルな SIR モデル)
fig.3-2  R_0 を変化させた時の集団免疫閾値の変化(3世代型 SIR モデル)
fig.3-3 fig.3-1 の横軸を拡大したもの
fig.3-4 fig.3-2 の横軸を拡大したもの

fig.3-1 と fig.3-2 は、これまで検討してきた最終規模のグラフ、例えば fig.1-1 と横軸範囲を合わせてあります。
少しこれだと横軸範囲が狭いと思ったので、拡大したのが fig.3-3 と fig.3-4 です。

fig.3-3 より fig.3-4 の方が少しグラフが低くなっています。
幾つかの  R_0 での集団免疫閾値を比較してみましょう。

使用モデルやグラフ\基本再生産数  R_0 1.5 2.0 2.5 3.0 3.5
シンプルな SIR モデル(理論値) 0.33 0.50 0.60 0.67 0.71
シンプルな SIR モデル(計算値)fig.3-3 0.33 0.50 0.60 0.67 0.71
3世代型 SIR モデル(計算値)fig.3-4 0.30 0.45 0.53 0.59 0.64

table.8 幾つかの  R_0 における集団免疫閾値の理論値と計算値

「シンプルな SIR モデル(理論値)」というのは、集団免疫閾値の理論値としてよく見るこの式

 \displaystyle (集団免疫閾値) = 1 - \frac{1}{R_0} ...... (9)

の右辺の  R_0 に、表に示した値を代入して得た数値です。
table.8 の「シンプルな SIR モデル(理論値)」でのこの値は、「シンプルな SIR モデル(計算値)」に全て一致しています。

「3世代型 SIR モデル」などの異質性を導入したモデルには、式(9) で示したような理論値は一般に存在しないと思います。

table.8 で、「シンプルな SIR モデル(計算値)」と「3世代型 SIR モデル(計算値)」を比較すると、「シンプルな SIR モデル(計算値)」の方が全てで少しずつ大きな値になっています。

例えば  R_0=2.5 の時、「シンプルな SIR モデル」で考えると、集団免疫に到達するには、0.60 つまり60%の人が感染(して免疫を獲得)する必要があります。
しかしこれを「3世代型 SIR モデル」で考えると、0.53 つまり 53%の人が感染すれば済むことになります。

この付近には、ここまでの議論をひっくり返すような重要な論点があるのですが、この論点については長くなるので、別の機会にします。
おそらく twitter( https://twitter.com/sarkov28 )に書きますが、ブログの別エントリにまとめる予定です。

(10)参考文献

参考文献です。

(11)修正履歴

修正履歴です。大きく修正した場合には、ここに書きます。

  • 2023-04-11
    公開。
  • 2023-04-12
    あちこちに加筆しました。
  • 2023-05-15
    「(9)基本再生産数が変化した場合の集団免疫閾値」を追加しました。
    「(8)「3世代型 SIR モデル」における最終規模」の最後の「また、この付近には」以後を加筆しました。


コロナ関連の資料


私がネット上でしていることの まとめ は、こちらに。
https://sarkov28.hatenablog.com/entry/2022/03/29/160915



コロナ関連で重要と考える資料と、関連するコメントを書いておきます。著者名は敬称を略すことがあります。

twitter に書いている諸々のうちの2つの事項に関連する資料については、別のエントリに書いています。目的が違うので、このエントリや以下の2つのエントリの内容は、重複するかも知れません。


目次


2008 稲葉寿 微分方程式感染症数理疫学

微分方程式感染症数理疫学, 数理科学 46 (4), 19-25, 2008-04, サイエンス社

この論文を参照するために使っていたリンクが切れていたので、ここに改めて所在を書いておきます。

この論文は、私が以下のツイ他で参照しているものです。


2020-06-11 西浦博「「8割おじさん」の数理モデルとその根拠」(Newsweek 2020-06-09 号)

Newsweek 記事、「「8割おじさん」の数理モデルとその根拠──西浦博・北大教授」です。

これは特に重要な記事です。以下2つは同内容だと思います。

以下、ネット記事に基づいて書きます。

主な内容は、西浦氏による「42万死亡推計」の解説です。
p.03 脚注で、github にある「42万死亡推計」の詳細を示しています。
また、p.05 脚注で、「42万死亡推計」に関係のある重要論文として、異質性に関する2つの論文を示しています。

  • Gomes GM, et al., "Individual Variation in Susceptibility or Exposure to SARS-Cov-2 Lowers the Herd Immunity Threshold," medRxiv, posted on April 27, 2020
  • Britton T, et al., "The disease-induced herd immunity level for Covid-19 is substantially lower than the classical herd immunity level," medRxiv, posted on May 6, 2020

西浦氏はこの記事で

  • (a)「42万死亡推計」が計算時点(2020-03-19)、発表時点(2020-04-15)において正当だったこと
  • (b) 2020-05 後半から 2020-06 上旬の日本の感染状況が、かなり抑えられていること

という、両立が難しそうな2つの事項の説明を試み、一応はこれを行っています。
ただし私は、このうち少なくとも(a)の説明には異論があります。この異論についてはいずれ(おそらく twitter に)書く予定です。


2020-08 西浦博「感染症数理モデル元年に機構と外挿の狭間に立つ」

西浦博「感染症数理モデル元年に機構と外挿の狭間に立つ」です。
数学セミナー 2020-09 の発行は、2020-08 です。

この論考は、以下の3つにあります。
(a)(b)(c)はほぼ同内容ですが、(a)(b)にある冒頭の一節が、(c)にはありません。
この一節は長くありませんが(雑誌紙面で、ページ半分弱)、「(西浦氏が)コードを全開示している」旨の(私にとっては)重要な表明が含まれています。(これ以外の一致を詳細に確認した訳ではありません。)

ムック本の (b) は、(a) の特集を集めたもの。(2023-12-08 時点で、入手可能なのは電子版だけのようです。)
単行本の (c) は、西浦氏らが雑誌などで発表した論文を集めたものであり、その中にこの論考が収められています。

この論考には2つのポイントがあると思います。
一つ目のポイントは、「モデルには、機構と外挿がある」ということです。

「機構」というのは「感染機構を組み込んだモデル」というほどの意味であり、つまり SIR 系のモデルを念頭に置いていると思われます。
(「感染機構を組み込んだ」というのは以下のような意味です。例えば、シンプルな SIR モデル(https://sarkov28.hatenablog.com/entry/2021/01/04/113012)における、新規感染者数を表す項である  βSI は、「新規感染者数が、感染しやすさに比例する定数  β と、感受性者数  S と、感染者数  I とに比例する」と主張しています。これは感染機構、つまり「感染は感受性者と感染者との接触によって起こる」という事情を表現していると言えます。)

西浦氏は「外挿」を「外挿的なモデル(extrapolation model)は誤解を恐れずに対比して書くならば,現象論として観察データを説明するのに「正しいっぽい曲線」を当てはめたものに相当する」と説明しています。例えば中野貴志氏の「K値」の計算はこれに該当すると私は思います。(西浦氏が「K値が」と具体的に書いている訳ではありません。)

「感染機構を組み込んだモデルの方が優れたところがある」という旨で書かれています。
(西浦氏が「外挿」のモデルを全く扱わないという意味ではありません。「(自分も)ロジスティック曲線を使うことがある」旨を説明しています。)

私も、

  • 実測データへの当てはまりの良さが同程度の2つのモデルがあり、
  • 一方が感染機構を組み込んだモデル、
  • 他方が外挿によるモデル、

であるなら、感染機構を組み込んだモデルの方が有望そうだとは思います。

二つ目のポイントは、「厚生労働省HIV 感染者数の推定に用いていたモデルは、外挿によるものだった。これは適切なモデルとは言えなかった」ということです。
HIV 感染者数の推定のモデルについてのその後の事情は、論考に明記がありません。厚労省は、モデル構築において西浦氏に「敗北」した可能性もあります。例えば現在も、西浦氏が関与したモデルを使っているのかも知れません。
もしそうであれば、新型コロナで西浦氏が重用されたこと、重用され続けていることに繋がりそうです。


2023-03-13 瀬名秀明ら「知の統合は可能か: パンデミックに突きつけられた問い」

瀬名秀明ら「知の統合は可能か: パンデミックに突きつけられた問い」です。
https://www.amazon.co.jp/dp/4788718693
この amazon によると、著者は瀬名秀明渡辺政隆、押谷仁、小坂健の各氏であり、発売日は 2023-03-13 です。


2023-03-18 瀬名秀明ら「知の統合は可能か: パンデミックに突きつけられた問い」文献リスト

「知の統合は可能か: パンデミックに突きつけられた問い」の文献リストです。

これは pdf ファイルとしてネット上に示されています。
非常に大量の文献のリストです。
重要と思いますので、別の項目にしました。

2023-03-18 瀬名秀明 新型コロナウイルス感染症に関する書籍リスト Ver.1.0
https://bookpub.jiji.com/files/BookList_2023_0318(1).pdf


2023-09-22 尾身茂「1100日間の葛藤 新型コロナ・パンデミック、専門家たちの記録」

尾身氏の手記、「1100日間の葛藤 新型コロナ・パンデミック、専門家たちの記録」です。

https://www.amazon.co.jp/dp/4296202553/


更新履歴

  • 2023-03-09
    公開。
    西浦博「感染症数理モデル元年に機構と外挿の狭間に立つ」を書きました。
  • 2023-03-23
    瀬名秀明ら「知の統合は可能か: パンデミックに突きつけられた問い」と、「文献リスト」を追加しました。
  • 2023-05-07
    「「8割おじさん」の数理モデルとその根拠」の項を追加しました。
  • 2023-09-24
    尾身茂「1100日間の葛藤 新型コロナ・パンデミック、専門家たちの記録」の項を追加しました。
  • 2023-12-08
    感染症数理モデル元年に機構と外挿の狭間に立つ」の記述を変更しました。
    まず最初の部分、(a)(b)(c)の異同についての記述を変更しました。変更前は、「(a)(b)(c)は実質的に同内容だ」という旨でした。また「外挿」の説明を変更した他、細かな修も加えました。