2022年3月14日月曜日

まん防に効果があったのか、Synthetic Control Methodで検証してみた

このエントリーをはてなブックマークに追加
Pocket

経済学者の北村周平氏がまん延防止等重点措置(通称、まん防)に対して、差分の差分法(DID)と呼ばれる統計的因果推論の方法を用いて分析を行い、効果なしと言う結果を、「専門家のみならず、一般の方々にもご理解していただけるように解説」を試みた*1のだが、分析および解釈がよろしくないと統計学と経済学の2つのPh.Dを持つ統計学者に批判されている。

擁護者も出てきて論争しているのだが、「批判するだけではなくて、やってみろ」と言うものもあった。代わりに分析してみよう。

1. 分析方法

分析方法としては、2022年1月からのサンプルに絞って、合成対照法(SCM)を用いて、長めの期間でまん防の効果を評価する。また、被説明変数は累積陽性者数の昨週比*2とする。

疫学者と共著にしたほうが…と話が出ていたが、上記の方針は疫学方面の議論をある程度は参照して決めた。

  1. 感染者数がだんだんと増えてまん防になるだけなら簡単なのだが、まん防→非常事態宣言→まん防のようなケースがあると非常事態宣言に対するまん防の効果が混じることになるし、まん防→解除→まん防で解除期間が短い場合も前回の残効果との比較になってしまうが、2022年1月からののサンプルはこれらの影響がない。また、ウイルス株が異なること(特にオミクロン株)などの影響も、分析期間を絞ることである程度は排除できる。
  2. SIRモデルなどの感染症数理モデルでは、感染症の予防効果は実効再生産数に影響を与えると考えているが、実効再生産数は新規感染者数に線形の影響を与えるのではなく、感染者数に応じて新規感染者数に影響する。新規感染者数(やそれを都道府県人口で割った率)は、感染症数理モデルにあわない。また、明らかに新規陽性数には曜日バイアスがあり、不要なノイズを生じさせる。昨週比をとることで、よりモデルに近いノイズが小さい被説明変数を作ることができる。
  3. 感染から発症、そして検査陽性までは平均12日ぐらいかかると言われている。これは2020年後半に言われていた話ではあるが、今でも病院にいくまで悪化しないと検査を受けない人も多いはずだ。『COVID-19 診療の手引』と言う医療関係者向けの冊子などを読むと、その場合は2週間以上、感染から感染確認まで時間がかかることになる。まん防指定から2週間ぐらいまでは効果が確認できない可能性は高い。

DIDではなくSCMを用いるのは、ここ4年間ぐらい他人にSCMを薦めながら自分で使ったことがなかったので…ではなくて、もう少しマトモな理由がある。DIDは、まん防を実施した都道府県(介入群)の平均と、実施しなかった県(対照群)の平均の、実施前のトレンドが平行であることが前提になるが、この前提を満たすようにするのが難しく苦しくなる。SCMの場合は、介入群のそれぞれに対して対照群のウェイト付き平均を計算するので、平行トレンドの要求が緩くなるハズ。

2. データセット

典型的なSCMは介入実施前までのサブサンプルで、内部的に2種類の行列を交互に収束するまで推定する。被説明変数以外は、片方の推定ではpredictor(定訳不明)と呼ばれる非時系列のデータセット*3、もう片方の推定では対照群の被説明変数と同じ時系列変数を用いる。

predictorとして、令和2年「人口統計」からの人口、14歳未満人口率、65歳以上人口率、「第六十三回日本統計年鑑」からの人口密度、全域に占める人口集中地区の人口割合、デジタル庁からのワクチン2回目接種率を用いた。被説明変数は、内閣官房新型コロナウイルス等感染症対策推進室から得たデータで作った陽性数の昨週比で、介入群は東京都、対照群は愛媛県、岩手県、宮城県、滋賀県、秋田県、鳥取県、徳島県、奈良県、富山県、福井県。まん防の実施期間(東京都は2022/1/21から)については、面倒だったのでWikipediaを参照している('-' )\(--;)BAKI

3. 結果

データセットをつくってしまったら、RのSynthパッケージにぽちっとなとお任せしたら結果が出てくる。一応、Abadie (2021)を読んだりしていたのだが、ほぼ役立たなかった。

SCMの結果は以下。実線が実際の陽性数の昨週比で、破線が反実仮想の陽性数の昨週比。実線が破線を下回った分だけ感染抑制効果があったことを示唆する。

0日がまん防実施開始日になるのだが、-20日から0日までのまん防実施前の当てはまりはよく、東京都の反実仮想は上手くいっているように思える。19日ぐらいに破線の反実仮想よりも実線の観測値の方が小さくなっていくことが分かる。すぐに効果は出ないが、3週間ぐらいでまん防の効果が見えてくると分析からは言える。なお、処置群を埼玉、千葉、神奈川にしても同様の傾向であった。TJO氏がCausalImpactで分析した結果も、だいたい似ている。

ただし、分析の頑強性の評価(placebo test)はよくない。推定モデルによる反実仮想がうまく行っていれば、対照群の観測値と反実仮想値の差は0付近に集まり、まん防に効果があるにして処置群は外れ値になるはずだが、微妙なところである。

SCMは介入群の個体と似ている個体が対照群にないと上手く機能しない。対照群の動きで、介入群をシミュレーションするからだ。しかし、介入群は都会であり、対照群は田舎になるので、無理があるのかも知れない。まん防がメディアを通じて全国にもたらす自粛ムードは、介入群と対照群の差異にならないので検出できないし、まん防実施都道府県が過半であり、まん防実施都道府県の隣接県を対照群から除外すると分析ができなくなるため、ある個体への処置は他の個体へ影響を与えないと言うSUTVAの仮定はかなり怪しく、効果量に過小推定バイアスが入っていることも影響しているであろう。

総じて評すれば、効果がありそうな気がするが、サンプルに限界があるのでなんとも言えないと言う話になる。無念。

4. まとめ

統計解析に慣れている人、統計分析の結果を知りたい人には、うまく行かない情報も有益だと思うが、そうではない人々にこういう不明瞭な結果を見せると混乱しそうである。さらに、まん防の効果を見つけることができませんでしたなどとまとめると、まん防に効果が無いと早とちりしそうだ。データセットの限界からかも知れない点は、強調しておく方が望ましい。既に謎のプレプリント・ペーパーにある憶測が伝播して、混乱している人々もいる。なお、Synthetic DIDなどの別手法を用いたら、もう少し話がはっきりするかも知れない。余力があれば試して見たい。

*1以下のツイート

*2新規陽性者数の昨週比ではなくて、累積陽性数の昨週比でよい理由は次の考えからである。SIRモデルの場合だと、t時点の累積感染者数をJₜを以下のように書ける(関連記事:直近一週間のデータで基本再生産数R₀を推定する方法)。

感染者数が多くなるとJ₀-C₁は相対的に小さくなるので、Jₜ/Jₜ₋₁をexp(βS₀-γ)で近似できる。βS₀-γは基本再生産数に直結する数字であり、これはt時点の感染者数によらない予防効果の指標になる。

追記(2022/03/16 07:54):やはり新規感染者数の昨週比を使うべきではないかと言う指摘を受けた。累積値だとt期のβではなく、0期からt期までのデータを使ったβの推定量を見ることになるわけで、各期のβの変化が緩やかにしか観測されないし、累積開始時点の影響を受ける。新規感染者数であれば、J₀-C₁のようなtに対するオーダーを理由に無視しないといけない項が残らない。理屈の上ではもっともなのだが、ノイズが大きくなるためか反実仮想の予測力が良く無かった。また、このバイアスは全ての都道府県にのるので、反実仮想自体には影響しないはずだ。

*3時系列データは期間全体で平均をとるなどして非時系列化して使う。

0 コメント:

コメントを投稿