読了: Downes& Carlin (2020) Mr.P (マルチレベル回帰・層化)ってやっぱり良いの? はい、良いみたいです

 代表性のない標本調査に基づいて母集団特性をうまいこと推測する、という話題のひとつに、Mr.P ことMRP(マルチレベル層化・回帰)というのがある。選挙予測で注目を集めた手法で(たしかYouGovがジョンソン勝利を当てたときに使ったんですよね)、私は面白がってGelman先生たちのarXivの論文のメモをとったりしてたんだけど、国内でもあれよあれよという間に有名になり、周囲から「MRPって知ってます?」と訊かれ、学会発表やそのへんの調査結果リリースでも見かけるようになり、なんと日本の伝統的大手メディアの方まで採用を検討しはじめ、ついには「MRPってできますか」と問い合わせが来たりするようになった。さすがにここまで来ると、個人的な熱は冷めてしまう。

 このたびちょっと思い立って検索してみたら、もはや心理学などの隣接諸分野で啓蒙論文が出るレベルにまで普及したようだ。2020年以降の論文は件数が少なく、とっくにハイプ・カーブの峠を越えたという感じである。
 そのなかから試しに一本読んでみた。疫学分野での手法検証研究である。

Downes, M., Carlin, J.B. (2020) Multilevel Regression and Poststratification Versus Sample Weighting for Estimating Populatoin Quantities in Large Population Health Studies. Americal Journal of Epidemiology, 189(7), 717-725.

(イントロダクション)
 非代表標本からの母集団パラメータ推定のモデルベース・アプローチとしてMRP(マルチレベル回帰・層化)が用いられている。[以下をreferしている: Gelman & Little (1997 Surv.Methodol.), Park, Gelma, & Bafumi (2004 Political Anal.), Wang, Rothschild, Goel, et al.(2015 Int.J.Forecast.), Kieding & Louis (2016 JRSS A).]
[…中略…]

 フォーマルに言うとMRPとはこういうのである。
 単一のアウトカム\(Y\)について、まず、事後層別セル\(j\)における平均(二値変数なら平均のロジット)\(\mu_j\)の線形予測子を指定する。$$ g(\mu_j) = g(E(Y_{j[i]}) = \beta_0 + \mathbf{X}^\top_j \beta + \sum_l a^k_{l[j]}$$ \(l[j]\)はセル番号\(j\)から適切なカテゴリ\(l\)へのマップ。すべてのランダム効果はIID正規 \(a^k_l \sim N(0, \sigma^2_k)\)とする。
 このモデルで\(\hat{\mu}_j\)を得て、重み\(N_j\)をつけて平均する。\(N_j\)はそのセルに対応する母集団サイズとする。$$ \hat{\mu}^{PS} = \frac{\sum_j N_j \hat{\mu}_j}{\sum_j N_j} $$ 同様に、下位母集団\(s\)のパラメータも推定できる。セル番号の下位集合を\(J_s\)として、上の式のサメーションを\(j \in J_s\)について回せばよろしい。

 複雑なhealth studyで母集団の記述統計量を推定したいとき、標本と母集団のずれをどう調整するか、というのは昔からの問題だった。伝統的には標本ウェイティングが用いられていた。ウェイトは標本選択確率の逆数と、無回答の標本ベース調整を組み合わせて求められた。事後層別やレイキングのようなカリブレーションがなされることもあった。たとえば一般化レイキングは、2次元以上の頻度表における母集団の周辺カウントを補足情報として用い、バイアスを取り除き推定量の効率性を高めようとする手法であった。

 我々はすでに、伝統的な標本ウェイティングとMRPの比較について報告している。

  • Downes, et al.(2018 Am.J.Epidemiol.): Ten to Men: The Australian Longitudinal Study on Male Health というのを使った事例研究。ベースライン・ウェイブ[通常の調査回ってことですかね]から35%の非代表標本を抜き出し、3つのアウトカム変数について推定した。MRPはウェイティングに比べて一致性も精度も高かった。
  • Downes & Calin (2020 Biom.J.) MRPとウェイティングの正確性と精度をシミュレーション実験で調べた。[…中略するけど、MRPは良かったんだったってさ]

 しかし、前者では真値が分かんないのでバイアスを評価できなかったし、後者では現実にあるような複雑なデータ生成メカニズムを作るのが難しかった。
 本研究では、もう一度Ten to Men調査を分析するんだけど、真値の代理として2011 Australian Census or Population and Housingを使う。

方法
アウトカム指標
 Ten to MenとCensusの両方で測っている変数のうち以下に注目する: 就労状態(元は6カテゴリだが2値に落とす)、週当たり労働時間。
 [考察のところに、確かにこの研究で注目した指標はhealth関連の指標ではないけれど、だからなに? なにか問題あんの? というくだりがあった。ちょっと笑っちゃった]

共変量
 以下を用いる: リモートネス(大都市, inner regional, outer regionalの3カテゴリ)、州(8)、年代(8)、先住民族か(2)、英語の流暢さ(4)、結婚状態(5)、社会経済的指標(SEIFA)(10)。州xリモートネスと年代xリモートネスの交互作用もいれる。

モデル
 就労状態については、\(\mu_j = Pr(Y_{j[i]} = 1)\)について、$$ \mu_j = \mathrm{logit}^{-1} (\beta_0 + \beta_{先住民族} + a^{リモート}_{l[i]} + a^{州}_{m[i]} + a^{年代}_{o[i]} + a^{英語}_{r[i]} + a^{結婚}_{s[i]} + a^{SEIFA}_{t[i]} + a^{州リモート}_{m[i], l[i]} + a^{年代リモート}_{o[i], l[i]})$$
\(a\)は全部ランダム効果で、iid正規とする。たとえば\(a^{リモート}_l \sim N(0, \sigma^2_{リモート})\)。
 [気持ち悪い書き方だなあ。これ、いまのところは個票のモデルなんですよね? そもそも右辺に\(j\)がないですよね? 左辺は\(Pr(Y_i = 1)\)のほうがよくないですか? この下のパラグラフで、セルが全ての共変量の直積として定義されることが明らかにされるので、ようやくこれがセルレベルのモデルでもある、つまり左辺を\(\mu_j = Pr(Y_{j[i]} = 1)\)とも書けることがわかるわけだけど…]

 労働時間についても同様。[こっちは連続量だから線型モデルにするんだろうけど]

事後層別
 事後層別セルは3x8x8x2x4x5x10 = 74800となる。構造欠損があって実際には60800。母集団サイズはCensusからとれる。

ベイジアン・モデリング
 RStanで推定した。パラメータの事前分布はすべての弱情報コーシー。無事収束しましたよ。

推定値
 unweighted/weightedの推定値と95%信頼区間を、surveyパッケージのsvymeanとsvycipropで得た。
 [あっそうか、Stanで完結するわけじゃないのね。MRPではモデリングとウェイティングが別だから。なるほど]

 Ten to Menにはベースライン・ウェイトがある(地域ごと)。さらに、年代、先住民族、英語、SEIFA、州xリモートを周辺にとったレイキングも行った。[ややこしい…]

州レベル共変量
ついでに、州レベル共変量をいれたモデルもつくった。中等教育完了率と就労者にしめるパートタイム率。結局MRPの州レベル推定値には効かなかったんだけど、以下は入れたモデルについて報告する。

結果
 unewighted, weighted(simple), weighted(raking), MRPによる推定値をCensusと比べる。

就労状態
 国レベルではMRP, raking, simple, unweightedの順に正確。州レベルでも大きな州ではだいたいその順。小さい州では…[略]
 MRPの州レベル推定値は他の方法に比べて変動が小さい。信頼区間も狭い。

労働時間
 国レベルでは、どの方法もoverestimateしているんだけど、MRP、raking, simple, unweightedの順に正確。州レベルでもだいたいそんな感じ…[略]

考察
 […中略…]
 MRPの精度は、以前の研究ほどには他と比べて高くならなかった。モデルが複雑だからだろう。小さい州ではやっぱり利得が大きい。
 州レベル共変量は効かなかったけど、これは今回そうだっただけで効くこともあるだろう。
 モデルをもっと複雑にしてもいいけど、母集団レベル頻度データを手に入れるのが難しいだろう。
 今回はRStanを使ったけど、lmer()やglmer()でもよい。もっともサンプルサイズが小さいと違いが出てくるかも。
 Ten to Menのウェイトを使う方法はもっと工夫すればもう少し性能が上がったかもしれない。でもMRPの精度には追いつかない。
 云々。
——–
 というわけで、やっぱしMRPがいいね、という話でした。なんだよ、つまんないな。(すいません)

 私がポイントを掴み損ねているのかもしれないけれど、この論文、理論的示唆があるのかどうかはよくわからない。センサスを”真”と見立てたとき、アウトカムのデータを使わないデザインベースの推定と、アウトカムのデータをみる(それどころか共変量の”真”の分布まで利用する)モデルベースの推定を比べたら、後者のほうが”真”と近かった、ということですよね? そりゃまあそうなんじゃないでしょうか。
 おそらく、公衆衛生の研究における実務的な寄与があるんだろうな。実際のところ、いま大規模な非代表標本があるとして、MRP以外の補正を考えるとしたら、たしかにその有力候補はデザインベースの補正だろうなと思う。その意味では、ベンチマークとしてデザインベースの方法を使うのは必ずしも藁人形というわけでもないのだろう。