読了: Sutton, Almquest, & Wakefield (2025) COVIDワクチンの小地域接種率をMRPと空間MRPで推定してみたけどからっきしダメでした

Sutton, A., Almquist, Z., Wakefield, J. (2025) Evaluating Multilevel Regression and Poststratification with Spatial Priors with a Big Data Behavioural Survey. arXiv:2503.05915.

先月arXivにあがったプレプリント。普通のMRPと空間MRPを実データで比べるという、関心にジャストフィットな話だったので、試しに読んでみた。JRSS Series Aに投稿中だと書いてあるが、ほんとだろうか。単にlatexのテンプレートがそうなっているだけなのかも。

1. イントロダクション
 古典的MRPでは未知の地域レベル推定値を部分プーリングで予測するけれど、地域間の空間関係は使わない。Gao et al.(2021)は「空間MRP」というのを提案していて、一次空間隣接行列によるBesag-York-Millieモデル(BYM2)を使っている。シミュレーションではバイアスを減らすことが示されているが実データでどうかはわからない。空間ベイジアン法はMCMCが必要でとても複雑だったが最近はINLAという手が出てきた。
 本論文では、古典的MRPと空間MRPを比べます。後者ではBYM2でデモグラカテゴリの内部と地域間の平滑化を行う[???]。INLAを使う。調査はCOVID-19 Trends and Impact Survey (CTIS)。

2. 方法
 関心あるアウトカムは、2021/06/30時点でのカリフォルニアの郡別の初回COVID-19ワクチン接種数。年代、学歴で層別する。性別にモデルを作る。

2.1 データ
[メモ省略]

2.2 モデリング・アプローチ: 空間MRP
 年代(3水準)を\(j\), 学歴(6水準)を\(k\), 群(58水準)を\(i\)と書く。人口を\(n_{ijk}\)とし、\(n_i = \sum_j \sum_k n_{ijk}, N = \sum_i \sum_j \sum_k n_{ijk}\)とする。接種者を\(y_{ijk}\)とする。
 モデルは次の4個。

  1. 古典的MRP、その1。$$ y_{ijk} | \theta_{ijk} \sim \mathrm{BetaBinomial}(n_{ijk}, \theta_{ijk}, d)$$ $$ logit(\theta_{ijk}) = \alpha + \gamma_j + \beta_k + \epsilon_i$$ $$ \epsilon_i \sim_{iid} \mathrm{Normal}(0, \sigma^2_\epsilon) $$ \(d\)はoverdispersionパラメータである。年代の効果は固定、学歴の効果は固定、郡の効果はiidである。
  2. その2。デモグラ層間の平滑化を入れる。\(\beta_k\)を一次のランダム効果としてランダム・ウォーク・モデルを考えよう。$$ y_{ijk} | \theta_{ijk} \sim \mathrm{BetaBinomial}(n_{ijk}, \theta_{ijk}, d)$$ $$ logit(\theta_{ijk}) = \alpha + \gamma_j + \beta_k + \epsilon_i$$ $$ [\beta_1, \ldots, \beta_K] \sim RW1(\sigma^2_\beta) $$ $$ \epsilon_i \sim_{iid} \mathrm{Normal}(0, \sigma^2_\epsilon) $$ \(\alpha\)を入れているからランダム・ウォークモデルには合計ゼロの制約を入れる。超事前分布はpenalized complexity事前分布とする。年代の効果は固定、学歴の効果はRW1, 郡の効果はiid。
  3. 空間的モデル、ひとつめ。$$ y_{ijk} | \theta_{ijk} \sim \mathrm{BetaBinomial}(n_{ijk}, \theta_{ijk}, d)$$ $$ logit(\theta_{ijk}) = \alpha + \lambda_j + \beta_k + \gamma_i$$ $$ [\beta_1, \ldots, \beta_K] \sim RW1(\sigma^2_\beta) $$ $$ [\gamma_1, \ldots, \gamma_I] \sim BYM2 (\sigma^2_\gamma, \phi) $$ 年代の効果は固定、学歴の効果はBYM2、郡の効果は暗黙的にiidである。[しれっと記号が変わっている点に注意。何考えてんだ。BYM2についてもなんの説明もない。不親切にもほどがあるぞ]
  4. ふたつ目。学歴の群ごとに平滑化する。$$ y_{ijk} | \theta_{ijk} \sim \mathrm{BetaBinomial}(n_{ijk}, \theta_{ijk}, d)$$ $$ logit(\theta_{ijk}) = \alpha + \lambda_j + \gamma_{ik}$$ $$ [\gamma_{1k}, \ldots, \gamma_{Ik}] \sim BYM2(\sigma^2_\gamma, \phi_k)$$ 年代の効果は固定、学歴の効果はRW1, 郡の効果はBYM2である。

 それぞれのモデルについて、\(\theta_{ijk}\)をドローして、3x6x58行の人口データを持ってきてかけて、群ごとの推定値を出す。これを1000回繰り返す。

3. 結果
 それぞれのモデルの結果を平均LCPO(logarithmic conditional predictive ordinate)で評価する。[恥ずかしながら初めて聞く名前だが、検索したところ自分のブログ記事がヒットした。要するに交差妥当化で求めた周辺尤度のことみたいだ]
 男女いずれも、モデル1が勝ってしまった。
 実はCDCによる正解的なデータを持っているので、突き合わせると、どのモデルの群レベル推定値も同じくらいに過小であった。

4. 考察
 MRPは実装という観点からはエキサイティングだが、目標母集団のアウトカムをモデルが正しく表現しているときにのみ有効であって、その検証はユーザ任せである。
 今回の問題についていうと、ワクチン接種へのためらいは歴史的に周辺化されたエスニシティで強く、また2021年初頭のCOVID-19ワクチン接種はカリフォルニアでは学歴と関係が強かった。
 モデルをどうやったら改善できただろうか? Kastellec et al.(2015 J.Politics)いわく、事後層別表を補足変数で拡張すべし。補足変数として地域レベルの生態学的・人口諸特性を入れ、調査標本と目標母集団の違いを減らす。今回の場合でいえばエスニシティ。
 Kennedy & Gelman (2021 Psych.Methods)曰く、汝、汝の母集団を知れ。標本を超えて一般化できるような結果を生むモデルを知れ。今回の場合でいえばCTISの標本選択プロセス。CTISはFacebook上での調査なので、アクティブ・ユーザがオーバーサンプリングされていたはずである。Facebookのような会社を説得して、ベンチマークにできるような指標を提供してもらうのが望ましい。
 […後略…]

5. 結論
 [パス]
————–
 要するに、正解的なデータがあるような状況で普通のMRPと空間MRPを比較してみたら、そもそもどっちもダメだった。という論文である。いさぎよいはなしだ。
 まじで当てに行く気なら、郡レベル共変量を入れたり、学歴の郡内での傾きをランダム効果にして郡間で変動させたりすればよかったんじゃない? と思うのだが、ま、余計なお世話であろう。ほんとに必要な説明変数が手に入んないときって、どんなにジタバタしてもダメなんでしょうね。