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

西浦氏の 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 に関する部分を削除しました。