まだ生態学に本格導入されていない統計的因果推論手法の紹介：傾向スコア、回帰分断デザイン、操作変数法を中心に

以下のセミナーでの発表内容のアーカイブです。
---
第 21 回　つくばE3セミナー　
「それってあなたの感想ですよね?」にならないための因果推論フェス in 生態学
日時：2022年2月22日(火) 14:00-17:30　オンライン開催（Zoom）

演題
まだ生態学に本格導入されていない統計的因果推論手法の紹介：傾向スコア、回帰分断デザイン、操作変数法を中心に

要旨
他分野で広く用いられてきている統計的因果推論手法の中には、未だ生態学分野にはあまり導入されていないものが幾つかある。それらの手法を無理矢理に生態学に導入する必要はないが、例えば傾向スコア法などは生態学においても有用な局面は多いと思われる。また、従来の生態学でも（暗黙のうちに）傾向スコア法に質的に類似したアプローチが採られている場合があり、そうしたアプローチの可能性とその限界を反省的に捉える上でも、傾向スコア法を理解する意義はあるだろう。本講演では、生態学分野に本格導入されていない手法のうちの代表的なものとして、傾向スコア法、回帰分断デザイン、操作変数法について簡単な解説を行う。

関連論文
1. Ramsey et al. (2018) Using propensity scores for causal inference in ecology: Options, considerations, and a case study. Methods Ecol Evol. 2019;10:320–331. DOI: 10.1111/2041-210X.13111
→Ramseyらによる生態学分野における傾向スコア法の解説

Transcript

1. まだ⽣態学に本格導⼊されていない 統計的因果推論⼿法の紹介： @国⽴環境研究所社会システム領域 傾向スコア、回帰分断デザイン、操作変数法を中⼼に 林岳彦
2. 今⽇は諸解析法のほんの概略をお伝えできれば 詳細や実装等については 最近に良本が出たので そちらをご参照おすすめ ⾼橋将宜著『統計的因果推論の理論と実装』
3. （簡単に）傾向スコア法の考え⽅の説明
4. 「相関」と「因果」がズレるとき 適⽤なし 適⽤あり 各調査地点における獣害施策Tと個体数Yの関係 獣害の 被害数 獣害対策T ？ 例：柵の設置とか 対策してる⽅が 被害が多い
5. 「相関」と「因果」がズレるとき 作物A 作物B メイン作物 獣害対策T 被害数Y 適⽤なし 適⽤あり 各調査地点における獣害施策Tと個体数Yの関係 獣害の 被害数 獣害対策T ⽐較対象間で背景要因が揃っていない（“バランシングしていない”） 例：柵の設置とか
6. シンプルな解決策：層別化により”背景”を揃える “背景が揃った”層へと 層別化して解析できれば 交絡の影響なく 因果効果が推定可能 適⽤なし 適⽤あり 獣害の 被害数 獣害対策T 適⽤なし 適⽤あり 獣害の 被害数 獣害対策T 適⽤なし 適⽤あり 獣害の 被害数 獣害対策T ？ ! ! 作物Aのみ 作物Bのみ ⽐較したい(=異なる処置を受けている)集団間 の”背景条件を揃える”ことが、統計的因果推論 の基本的なアプローチ
7. 何が困るか：”背景”が単純ならいいけれど･･･ 2値の背景要因がm個あるとき、2 個の層を考える必要がある m 調整すべき要因が 多数あると、各層が 細切れになりすぎて 層別化の際にサンプルサ イズが追いつかない 施⽤あり 種Aの 個体数Y 農薬T 施⽤なし 施⽤あり 種Aの 個体数Y 農薬T 施⽤なし 施⽤あり 種Aの 個体数Y 農薬T 施⽤なし ･･･ C1=0, C2=0, C3=0, C4=0, C5=0, C6=0 C1=1, C2=0, C3=0, C4=0, C5=0, C6=0 C1=1, C2=1, C3=1, C4=1, C5=1, C6=1 また、そもそも 背景要因が連続量の 場合にはそのままでは 層別化できない 全64層 ･･･ 適⽤なし 適⽤あり 獣害の 被害数 獣害対策T ？ いわゆる”次元の呪い” 例えば2値の背景要因が6個あると ”背景が揃った”層の 組み合わせは2 =64個 6
8. 何が困るか：回帰モデルでの調整も可能ではあるが･･･ 回帰モデルでの調整はモデルの妥当性に強く依存する 本当に全ての共変量に対して線形モデルが妥当なの？ モデルが妥当でない（実態と合っていない）と 因果効果の推定量にバイアスがかかることが知られている 共変量が多くなるほどモデル⾃体の妥当性の吟味や正当化は⼀般に難しい ＆（実データがない）外挿的な推定となる空間が多⼤となり「モデル頼みの推定」となる 𝑌 = 𝛽! 𝑇 + 𝛽" 𝐶" + 𝛽# 𝐶# +𝛽$ 𝐶$ +𝛽% 𝐶% +𝛽& 𝐶& + 𝛽' 𝐶' +…
9. 傾向スコア法の考え⽅：”背景”を⼀次元に縮約する C1 処置T 結果Y C2 Cm … C1 処置T 結果Y C2 Cm … 傾向スコアe これら全部いちいち個別に 調整するの⾯倒アンド困難 「処置Tの予測」の 観点から エイヤっと ⼀次元に縮約 しちゃおう！ ＊傾向スコア法では基本的に 処置Tは⼆値（少なくとも離散） 変数であることが前提 ここで傾向スコアeiは、共変量C1,…,Cmからモデルをもとに計算 される処置Tの予測スコア（個体i がT=1となる確率） 傾向スコアeの計算には 多くの場合ロジスティック 回帰モデルが⽤いられる ⽣態学分野だと連続量の原因を扱うケース も多く、ここがネックの⼀つとなりがち logit(𝑇) = 𝛼 + 𝛽!𝐶! + ⋯ + 𝛽"𝐶" ⽣態学者はみんな⼤好きGLM ⽬的変数は処置Tであることに注意 logit(𝑒#) = 0 𝛼 + 1 𝛽!𝐶!,# + ⋯ + 1 𝛽"𝐶",# Propensity Score 「T=1となる確率」の予測値を計算
10. 傾向スコア法の何が嬉しいのか 似たような傾向スコアeの値を持つ個体同⼠は 諸背景要因が似ている（背景要因が揃っている）と期待できる 傾向スコアを利⽤した マッチング 傾向スコアを利⽤した 層別解析 傾向スコアを利⽤した 逆確率重み付け(IPW) 多様な背景要因をまとめてエイヤっと揃えた上で 処置Tの因果効果を推定できる ＊実際に揃っているかはpost-hoc的に確認できる←なにげにこれが傾向スコア法の⼀番嬉しい点かも ＊ 他にも回帰分析の説明変数として ⽤いる⽅法や⼆重ロバスト法 （IPW+回帰）などもあり
11. 傾向スコア法の適⽤に必要な理論的＆実際的要件（の⼀部） 交絡の調整に必要な（=バックドア基準を満たす）変数が全て観測されている＆ 傾向スコアの推定の際に利⽤されている 傾向スコアの推定の際に、処置Tの下流にある変数を⽤いてはいけない 傾向スコアの推定が概ねうまくいっている＆調整により共変量のバランシングが 実際に改善している（←後で説明します） ⽐較したい群間で、共変量の分布の重なり（common support）が⼀定以上ある 観測データが豊富にあることと、背景の因果構造に関する⼀定以上のドメイン知識が既に あることが上記2つの条件が満たされる（かどうかが判断可能となる）ためには必要 実際に解析をしてみて、あんまりうまくいかないというケースも多い あんまりこういう条件を真⾯⽬に検討していないケースも多い 競争的資⾦の研究計画とかに「傾向スコア法を使う」と書きにくい理由の⼀つ この後紹介するケーススタディもちょっとそういうところある ⽐較したい集団間の隔たりが⼤きすぎると統計的な調整では追いつかない
12. ケーススタディ (Ramsey et al. 2018 in Methods. Ecol. Evol.) ＊論⽂のURLはE3セミナーの⾴の要旨にあります https://doi.org/10.1111/2041-210X.13111
13. ケーススタディ (Ramsey et al. 2018 in Methods. Ecol. Evol.) ニュージーランドにおける⽣物多様性保全のためのBiodiversity Monitoring and Reporting Systemsデータ（公共保護地での8kmグリッドデータ） フクロギツネのブラウジングによる在来樹⽊の減少が問題視されている ニュージーランドの公共保護地における、移⼊種であるフクロギツネ への管理施策（毒⼊りベイト剤の空中散布）の因果効果の推定 各保護対象樹⽊種に対する管理施策導⼊グリットと⾮導⼊グリットの数 樹⽊種 施策ありサイト 施策なしサイト Melicytus ramiflorus 32 145 Metrosideros umbellata 38 81 Raukaua simplex 45 106 Weinmannia racemosa 101 221
14. ケーススタディ (Ramsey et al. 2018 in Methods. Ecol. Evol.) 結果変数は「各保護樹⽊の林冠被覆度」、処置変数は「ベイト散布の有無」 結果変数と処置変数の両者に影響を与えうる共変量（背景要因）である 13変数を⽤いて傾向スコアを推定 ベイト散布T 被覆度Y Total basal area of stems Total number of stems Total number of possum-preferred food plant species Elevation of plot above sea level Slope Plot distance from improved pasture Rainfall Potential evapotranspiration Soil acid soluble phosphorus Calcium Psize MAS Temp ＊⼀部、ベイト散布の下流となりうる変数が含まれており 変数選択が不適切だと思うがひとまずスルーします これらの背景要因を 処置群間で揃えたい 傾向スコアe 傾向スコアで ⼀括調整 フクロギツネ⾃体の個体数等をアウトカムにしているわけではない
15. ケーススタディ (Ramsey et al. 2018 in Methods. Ecol. Evol.) ＊SUTVAが満たされていなそうな気もするので推定値の解釈は正直やや微妙なところはある ＊Ramseyら提供のRコードを利⽤して作図 傾向スコアを⽤いた調整なし／ありでの因果効果の推定結果 調整 なし IPW IPW+ 回帰 調整 なし IPW IPW+ 回帰 調整 なし IPW IPW+ 回帰 調整すると有意！ 傾向スコアを⽤いた交絡の調整法 あまり変わらず 効果量がよりNULL https://commons.wikimedia.org/wiki/File:Rat a-tree.jpg#filelinks https://commons.wikimedia.org/wiki/File: Botany_of_Antarctica-PL012-0029.jpg https://commons.wikimedia.org/wiki/File:KamahiFoliage.jpg 被度への 因果効果 の推定値
16. ケーススタディ (Ramsey et al. 2018 in Methods. Ecol. Evol.) 背景要因が調整によりどれくらい揃ったかを確認する（ラブプロット） 処置あり群-処置なし群間でのMean Standardised Difference 群間での偏り⼤ 群間での偏り⼤ 群間で 良バランス PS調整前 PS調整後 ＊Ramseyら提供のRコードを利⽤して作図 傾向スコアによる調整後 に全体的にバランスが 改善している https://commons .wikimedia.org/wi ki/File:Rata- tree.jpg#filelinks
17. 考察：⽣態学的には⽣息適地分析的な⽂脈に少し似てるかも 例題：「ある農薬Aの使⽤の有無T」が「ある⽔⽥依存種Bの⽣息Y」に影響を 与えているかを広域モニタリングデータから調べたいと思ったらどうする？ 「農薬Aの使⽤の無いサイト」の中には、「市街地などそもそも種Bの⽣息に適さないサイト」が 多く含まれているため、農薬の使⽤有り-無し群間の単純な⽐較では背景因⼦が揃わない 「⽣息確率の予測スコアが同じサイト」でマッチング／層別化して、農薬Xの使⽤の有無Tが 異なるサイトにおける⽣息Yを⽐較する ⼀つの考え⽅として、多変量データに基づく⽣息適地分析からの という⽅法がありうるかもしれない こういう使い⽅は傾向スコア法に概念的にけっこう近いと思う （ただしここでのモデリング対象はアウトカムであることが⼤きな違い） 尚、統計的因果推論の理論からは、こういう使い⽅をする場合には、予測スコアの構築 においては農薬Aの下流の変数を含んではいけないことが含意されることとなる
18. ⽣態学×傾向スコア法のまとめと雑感 「処置効果が⼆値（離散）であり、処置変数の因果効果の推定のみが⽬的である」 ときには⽣態学でも本来は第⼀選択的な⼿法となる場合は多いだろう 「傾向スコア法を知っていたら、⽇の⽬を⾒たデータになってたのに･･･」みたいな可能性 も考えると、ツール／概念としては知っておいたほう絶対よい！ 違う⾔い⽅をすると、上記の⽬的の場合には、回帰モデルの適⽤は不利な（より強い仮定 を必要とする）場合が多い ⽣態学分野では⼆値変数の因果効果推定が⽬的となるケースは相対的には多くはないかも 実際にやろうとすると必要となる条件が合わない場合も多い ⽣態学だと（既存データでは）⼀連の共変量データまで揃ってることはあまりない いざやろうとすると、群間での共変量の重なりがない等で悲しい思いをすることも多い いずれにしろ、もしハマればつよつよの⼿法であるのは確かである！（TRY!）
19. (ごくごく簡単に)回帰分断デザインの説明 細かい実装等を説明する前の段階の、そもそも本解析の前提がハマる状況があるのかというところに ⽣態学的には⼤きなハードルがあるため、細かい話は抜きにして基本的なアイデアのみを説明します
20. 回帰分断デザインの基本的なイメージ 連続量のどこかで施策が切り替わるときの「切り替わり際」に着⽬する 仮想例：クーラー導⼊の労働⽣産性への効果 「午前10時の外気温25℃以上の⽇だけクーラー使⽤可」というルールが あるとする 午前10時の外気温 労 働 ⽣ 産 性 25℃ 15℃ 30℃ 切り替わり近傍でのクーラーの有 無以外の条件はほぼ”揃っている” と想定可能 もちろん諸々の前提の もとで成り⽴つ解析です ＊ この差 ＝クーラーの因果効果
21. ⽣態学で「切り替わり」を利⽤できるってことある？ ⾃治体等の境界で施策が切り替わるケースとか？ 県Bのみで保全施策X が実施されている 県A 県B 県境をまたいだ近傍の ⽐較により施策効果が わかるかも 県境 「トリガー⽔位」を基に故障の⽔管理施策が切り替わるケースとか？
22. (ごくごく簡単に)操作変数法の説明 細かい実装等を説明する前の段階の、そもそも本解析の前提がハマる状況があるのかというところに ⽣態学的には⼤きなハードルがあるため、細かい話は抜きにして基本的なアイデアのみを説明します
23. 未観測交絡要因 操作変数法の基本的なイメージ 操作変数IV 処置T 結果Y C2 Cm … C1 交絡要因となる共変量群とは別のルートで処置Tの有無に影響を与える 変数があると、それを利⽤して因果効果を推定できる場合がある 操作変数（IV）の条件 操作変数が以上の条件を満たすとき T→Yの因果効果 ＝IVとYの共分散／ IVとTの共分散 Instrumental Variable (1) 処置Tと相関がある （尚、⾼い相関でないと推定値は不安定になる） (2) 処置Tを通してのみ結果Yに影響を与える （唯⼀経路条件） (3) IVとYの間にバックドアパスが開いていない （IVは外⽣的である／誤差項と相関がない） 𝜌!"→$ 𝜌$→% 𝜌!"→% =𝜌!"→$×𝜌$→% 𝜌$→% =𝜌!"→% ∕ 𝜌!"→$ パス係数の積から
24. ⽣態学でそんな都合のよい”操作変数” ってある？ 環境団体の抗議でベイト剤による施策が中⽌になるとか？ 台⾵のような外的要因で「柵」が倒れるとか？ 台⾵IV 柵の設置T イノシシ⾷害Y C2 Cm … C1 環境団体に よる抗議IV ベイト剤T 被覆度Y C2 Cm … C1 施策Tが外部要因IVの関数となっていることが肝要 未観測交絡要因 未観測交絡要因
25. 「⽣態学×本発表での因果推論⼿法」のまとめと雑感 今回紹介した⼿法が⽣態学で未だあまり使われていない理由とその対策 単純に知られていないので使われていない テーマ的に使いどころがハマらない 質・量ともにデータが⾜りない 査読者がその意義を理解できない （頑張っても苦労が報われない） 頑張って説明すれば最終的には多分なんとかなる （もちろん新規性で⾼評価が付く場合もある） ⾃分たちでデータを取る／構築するしかないかも 今⽇知ったので使えるよ！ ハマりうるテーマをスルーしないようアンテナを張ろう よく分からなかったら国環研チームにいつでも相談してね！ （我々もそんなに詳しくないけど）みんなで考えよう！