0
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?

因果推論 Day 13/全30回 フロントドア基準、交絡を測れないときの抜け道

0
Posted at

この連載について

因果推論を「本を読んだ」で終わらせず、自分の言葉で説明でき、コードで再現できる状態まで落とす30日連載です。前回のDay 12では、「観察する」と「介入する」の違いをdo演算子で式にして、交絡因子が観測できていれば調整化公式で介入分布を計算できることを確かめました。今日はその大前提を外します。交絡因子がどうしても測れないとき、それでも因果効果を点で当てられる抜け道、フロントドア基準です。真値との答え合わせをしたいので、今日は全編合成データで機構をゼロから組みます。


TL;DR(3行)

  • フロントドア基準=バックドアが塞げないときの迂回路。 未観測の交絡$U$(体質)が喫煙$X$と肺がん$Y$の両方に効く世界では、素朴な差 0.420 は真値 0.239 の倍近く。$U$を測れない以上、バックドア調整は打つ手がありません。
  • 鍵は媒介$M$(タール)を「通り道の計測器」に使うこと。 フロントドア式は2回のバックドア調整の合成 $P(y|do(x))=\sum_m P(m|x)\sum_{x'}P(y|x',m)P(x')$ で、手組みの推定は 0.242。$U$を一度も使わずに真値へ当たります。
  • ただし前提は厳しい。 $X \to Y$の直接経路を足すとフロントドアは0.240のまま真値 0.390 を取り逃し、$U$が$M$にも効くと0.275 vs 真値 0.180 に崩れます。構造を買えるなら点識別、買えないなら感度分析(Day 9)です。

今日の問い

バックドア基準は「交絡因子を測って塞ぐ」だった。では、交絡因子がどうしても測れなかったら?
未観測の交絡をまたいだまま、因果効果を「点で」当てられる構造があるとしたら、それはどんな形か?

答えは「ある」です。今日のゴールは、その構造の条件(フロントドア基準)を言葉にし、フロントドア式を「2回のバックドア調整の合成」として1行ずつ導き、合成データで本当に真値が復元されることと、前提が壊れたときの崩れ方まで、コードで確かめることです。


概念① フロントドア基準:塞げないなら、通り道を測る

舞台は1950年代の喫煙論争です。喫煙が肺がんを増やすという観察結果に対して、「喫煙も肺がんも引き起こす遺伝的な体質$U$があるだけではないか」という反論がありました。体質は測れません。DAGで書くと、バックドアパス $X \leftarrow U \rightarrow Y$ が観測できない変数でできていて、Day 11の意味で塞ぎようがない状況です。

ここでPearlが持ち出すのが機構の知識です。喫煙が肺がんを増やすのは、肺へのタール沈着$M$を介してだと考えます。しかも体質$U$は、タールが溜まるかどうかには直接関係しません。この構造を描いたのが次のDAGです。

このとき$M$は、$X$から$Y$への因果の「唯一の通り道」に立っています。フロントドア基準は、媒介$M$がこの役を務められる条件を3つにまとめたものです。

  1. $M$が$X \to Y$のすべての有向経路を媒介する。 $X$の因果は$M$を通らずに$Y$へ届かない(完全媒介)。
  2. $X \to M$にバックドアパスがない。 $X$から$M$への効果は、そのまま観測できる。
  3. $M \to Y$のすべてのバックドアパスを$X$が塞ぐ。 $M$から$Y$への効果は、$X$で調整すれば観測できる。

直感はこうです。$X \to M$は「きれいな区間」、$M \to Y$は「$X$で掃除できる区間」、そして因果は$M$経由の一本道。ならば2つの区間の効果を別々に測って、つなげばいい。バックドア(裏口)が塞げないなら、正面玄関(フロントドア)から因果が通る様子をそのまま計測する、という発想です。

概念② フロントドア式:2回のバックドア調整の合成

フロントドア式は一見複雑ですが、正体はDay 11とDay 12で組んだバックドア調整を2回使って合成しているだけです。1行ずつ分解します。

第1段。$X \to M$の効果は素通しで読めます。

$$
P(m \mid do(x)) = P(m \mid x)
$$

$X$から$M$へのバックドアパスは $X \leftarrow U \rightarrow Y \leftarrow M$ の1本だけで、$Y$が合流点なので条件付けしなくても閉じています(条件2)。調整セットは空集合。つまり観測された条件付き確率が、そのまま介入分布です。

第2段。$M \to Y$の効果は$X$で塞いで読めます。

$$
P(y \mid do(m)) = \sum_{x'} P(y \mid x', m), P(x')
$$

$M$から$Y$へのバックドアパスは $M \leftarrow X \leftarrow U \rightarrow Y$。$U$は観測できませんが、その手前の$X$で塞げます(条件3)。これはDay 12の調整化公式そのもので、調整セットは${X}$です。

第3段。2つを合成します。

$$
P(y \mid do(x)) = \sum_m P(m \mid do(x)), P(y \mid do(m))
$$

$X$が$Y$に効く経路は$M$を通るものしかありません(条件1)。$do(x)$が$Y$を動かすには、まず$M$を動かし、その動いた$M$が$Y$を動かすしかない。だから「$X$を動かしたときの$M$の反応」と「$M$を動かしたときの$Y$の反応」を掛けて足せば、$X$から$Y$への効果になります。第1段と第2段を代入すると、フロントドア式が出ます。

$$
P(y \mid do(x)) = \sum_m P(m \mid x) \sum_{x'} P(y \mid x', m), P(x')
$$

右辺を見てください。$do$も$U$も残っていません。全部、観測データ$(X, M, Y)$から数えられる量です。因果量 $P(y|do(x))$ がこのように観測分布の式へ書き換えられるとき、その因果量は識別可能であるといいます。ここで大事なのは、識別可能性がデータの性質ではなく「データ+DAGの仮定」の性質だということです。同じ3変数のデータでも、DAGが $X \leftarrow U \rightarrow Y$ だけ($M$なし)なら識別不能、フロントドアDAGなら識別可能。答えを出せるかどうかを決めているのは、データの量ではなく構造の仮定のほうです。

フロントドア基準の図解(左:Uが測れずバックドアX←U→Yが塞げない、中:第1段X→Mはバックドアなしで素通し、右:第2段M→YはXで塞ぐ。合成するとdoもUも消える)


手を動かす①:バックドア調整が本当に「詰み」であることを確認する

まず、上のDAGどおりの世界を合成します。真値は私たちが握ります。$X$が$M$の発生確率を0.8動かし、$M$が$Y$の発生確率を0.3動かすので、$X$の$Y$への真の効果(リスク差)は解析的に $0.3 \times 0.8 = 0.24$ です。答え合わせはdo介入の世界を実際に作って行います。

import numpy as np

rng = np.random.default_rng(13)
N = 1_000_000
u = rng.binomial(1, 0.5, N)                    # 体質(観測できない交絡因子)
x = rng.binomial(1, 0.2 + 0.6 * u)             # 喫煙:体質が強く後押しする
m = rng.binomial(1, 0.1 + 0.8 * x)             # タール:喫煙だけで決まる(Uは効かない)
y = rng.binomial(1, 0.05 + 0.3 * m + 0.3 * u)  # 肺がん:タールと体質で決まる(直接経路なし)

# 真値(オラクル): do(X=1) と do(X=0) の世界を実際に作って比べる
def do_world(x_val, seed=131):
    r = np.random.default_rng(seed)
    u_ = r.binomial(1, 0.5, N)
    m_ = r.binomial(1, 0.1 + 0.8 * x_val, N)
    return r.binomial(1, 0.05 + 0.3 * m_ + 0.3 * u_).mean()

def adjust(x, y, s):
    """sの層の中で処置差を取り、層サイズで加重平均する(層別調整)"""
    return sum((y[(x == 1) & (s == k)].mean() - y[(x == 0) & (s == k)].mean())
               * (s == k).mean() for k in (0, 1))

print("真値 P(Y|do(X=1)) - P(Y|do(X=0)) = %.3f" % (do_world(1) - do_world(0)))
print("素朴な差     = %.3f" % (y[x == 1].mean() - y[x == 0].mean()))
print("Mで層別調整  = %.3f  (媒介変数への悪い調整)" % adjust(x, y, m))
print("Uで層別調整  = %.3f  (オラクル: 現実には不可能)" % adjust(x, y, u))
真値 P(Y|do(X=1)) - P(Y|do(X=0)) = 0.239
素朴な差     = 0.420
Mで層別調整  = 0.179  (媒介変数への悪い調整)
Uで層別調整  = 0.240  (オラクル: 現実には不可能)

素朴な差は 0.420 で、真値0.239の倍近くあります。体質$U$を持つ人は喫煙もしやすく肺がんにもなりやすいので、差が水増しされるわけです。では観測できる変数で調整できるか。手元にあるのは$M$だけですが、$M$で層別すると 0.179。これは因果の通り道そのものを塞いでしまう「悪い調整」(Day 11)で、残った0.179は交絡の成分です。真値に戻せるのは$U$での調整(0.240)だけですが、$U$は観測できません。バックドア調整はこの世界では本当に詰んでいます。


手を動かす②:フロントドア式を手組みして、Uなしで真値を復元する

いよいよ迂回路です。概念②の式をそのまま実装します。使うのは $x, m, y$ の3列だけで、u はコードのどこにも登場しません。

import numpy as np

rng = np.random.default_rng(13)
N = 1_000_000
u = rng.binomial(1, 0.5, N)                    # データ生成は①とまったく同じ(seedも同じ)
x = rng.binomial(1, 0.2 + 0.6 * u)
m = rng.binomial(1, 0.1 + 0.8 * x)
y = rng.binomial(1, 0.05 + 0.3 * m + 0.3 * u)

def front_door(x, m, y):
    """P(y|do(x)) = Σ_m P(m|x) Σ_x' P(y|x',m) P(x') を手組みする。uは一切使わない。"""
    p_x = [(x == 0).mean(), (x == 1).mean()]               # P(x')
    q = []
    for xv in (0, 1):                                      # 介入で選ぶ x
        total = 0.0
        for mv in (0, 1):
            p_m_given_x = (m[x == xv] == mv).mean()        # 第1段: P(m|x)
            p_y_do_m = sum(y[(x == xp) & (m == mv)].mean() * p_x[xp]
                           for xp in (0, 1))               # 第2段: Σ_x' P(y|x',m)P(x')
            total += p_m_given_x * p_y_do_m                # 合成
        q.append(total)
    return q

q0, q1 = front_door(x, m, y)
print("P(Y=1|do(X=1)) = %.3f / P(Y=1|do(X=0)) = %.3f" % (q1, q0))
print("フロントドア推定 = %.3f  (解析的な真値 0.240、誤差 %.3f)" % (q1 - q0, abs(q1 - q0 - 0.24)))
P(Y=1|do(X=1)) = 0.471 / P(Y=1|do(X=0)) = 0.229
フロントドア推定 = 0.242  (解析的な真値 0.240、誤差 0.002)

フロントドア推定は 0.242。観測できない$U$に一度も触れずに、真値0.240を誤差 0.002 で復元しました。コードを見返すと、外側のループ変数 xv が「介入で選ぶ$x$」、第2段の xp が「観測されたXの分布を走る$x'$」で、はっきり別物として扱われています。この2つの区別がフロントドア式の芯です。素朴な差・悪い調整・フロントドアを並べたのが下の図で、真値の赤線に当たるのはフロントドアだけです。

3推定の比較(素朴な差0.420は交絡で上振れ、Mでの調整0.179は因果経路を塞いで交絡だけが残る、フロントドア0.242だけが真値0.239の赤線に当たる)


手を動かす③:前提を壊すと、フロントドアは静かに崩れる

うますぎる話には値段があります。フロントドアの3条件を1つずつ壊して、推定がどうなるかを見ます。壊し方は2通り。$X \to Y$の直接経路を足す(条件1が破れる)と、$U$が$M$にも効くようにする(条件2・3が破れる)です。

import numpy as np

N = 1_000_000

def make(rng, n, x_val=None, direct=0.0, u_to_m=0.0):
    """direct: X→Yの直接経路の強さ/u_to_m: UがMに効く強さ(前提を壊すダイヤル)"""
    u = rng.binomial(1, 0.5, n)
    x = rng.binomial(1, 0.2 + 0.6 * u) if x_val is None else np.full(n, x_val)
    m = rng.binomial(1, 0.1 + (0.8 - u_to_m) * x + u_to_m * u)
    y = rng.binomial(1, 0.05 + 0.3 * m + 0.3 * u + direct * x)
    return x, m, y

def front_door(x, m, y):
    p_x = [(x == 0).mean(), (x == 1).mean()]
    q = [sum((m[x == xv] == mv).mean()
             * sum(y[(x == xp) & (m == mv)].mean() * p_x[xp] for xp in (0, 1))
             for mv in (0, 1)) for xv in (0, 1)]
    return q[1] - q[0]

for name, kw in [("前提どおりの世界   ", {}),
                 ("直接経路X→Yを足す", {"direct": 0.15}),
                 ("UがMにも効く      ", {"u_to_m": 0.2})]:
    x, m, y = make(np.random.default_rng(13), N, **kw)
    truth = (make(np.random.default_rng(131), N, x_val=1, **kw)[2].mean()
             - make(np.random.default_rng(131), N, x_val=0, **kw)[2].mean())
    print("%s: フロントドア = %.3f / 真値 = %.3f" % (name, front_door(x, m, y), truth))
前提どおりの世界   : フロントドア = 0.242 / 真値 = 0.239
直接経路X→Yを足す: フロントドア = 0.240 / 真値 = 0.390
UがMにも効く      : フロントドア = 0.275 / 真値 = 0.180

直接経路を足した世界では、真の効果は0.390に増えたのに、フロントドアは 0.240 のまま。$M$を通らない0.15ぶんの効果を丸ごと取り逃しています。フロントドアが測れるのはあくまで「$M$経由の因果」だけだからです。$U$が$M$にも効く世界では逆に 0.275 と過大になります。第1段の $P(m|x)$ が交絡で水増しされ、それがそのまま推定に乗るためです。怖いのは、どちらの世界でもコードはエラーひとつ吐かず、もっともらしい数字を返してくることです。崩れたことを教えてくれるのは、真値を握ったオラクルだけでした。

前提が壊れたときの崩れ方(左:前提どおりなら0.242で真値0.239に一致、中:直接経路を足すと0.240のまま真値0.390を取り逃す、右:UがMに効くと0.275で真値0.180を過大推定)

ここでDay 9との住み分けを整理しておきます。フロントドアは、DAGという構造の仮定を買う代わりに点識別を得る道です。感度分析(E-value)は、構造が分からないことを認めた上で、結論が隠れ交絡にどこまで耐えるかを報告する道です。そしてフロントドアの前提、とくに完全媒介はかなり贅沢な仮定です。喫煙の例でさえ「タール以外の経路が本当にないのか」は議論の余地があり、教科書の外でこのDAGがそのまま成り立つ場面は正直まれです。前提を買えるなら点で答え、買えないなら耐性を数字で添える。両方を持っていることが報告の幅になります。


つまづき・誤解しやすい点

  • 「Mで調整する」とフロントドアは別物です。 $M$で条件付けるのは因果の通り道を塞ぐ悪い調整で、手を動かす①のとおり交絡だけが残った0.179が出てきます。フロントドアは$M$を調整セットに入れるのではなく、$X \to M$と$M \to Y$という2つの因果効果の計測器として使います。「媒介変数は調整しない、計測する」と覚えると混同しません。
  • 内側の和の$x'$は、介入で選んだ$x$とは別物です。 第2段の $\sum_{x'} P(y|x',m)P(x')$ は$M \to Y$のバックドアを塞ぐための平均で、$x'$は観測されたXの分布全体を走ります。ここに$do$の$x$を代入してしまうのが定番の実装事故です。手を動かす②のコードでループ変数を xvxp に分けているのはそのためです。
  • 前提はデータからは検証しきれません。 $X \to Y$の直接経路や$U \to M$の矢印があるかどうかを、観測分布 $P(x,m,y)$ から判別することは基本的にできません。手を動かす③で見たとおり、壊れた世界でもコードは平然と数字を返します。フロントドアを使えるかは機構の知識で決まるのであって、識別可能性は「データ+仮定」に対して定義される。Day 6の「非交絡性は検証できない」と同じ形の告白が、DAGの世界でも繰り返されるわけです。

GISデータ実務での使い方

  • 施策→行動変容→成果のチェーンが測れているなら、フロントドアの土俵に乗ります。 例えば、漏水リスクの優先度リストを提供する施策$X$、それを受けた職員の調査行動$M$、漏水発見数$Y$。自治体の意欲$U$が測れなくても、行動ログ$M$が取れていて3条件が立つなら、点識別の候補になります。媒介の計測を設計に入れておくことが、後から効いてくる布石です。
  • 完全媒介の仮定は、経路の棚卸しで吟味します。 リストが調査行動以外の経路で成果に効かないか。例えばリストの存在自体が議会や住民の関心を高めて予算が付くなら、それは$M$を通らない直接経路で、手を動かす③の真ん中の崩れ方になります。「$M$を通らない効き方を3つ挙げてみて、全部潰せるか」を仮定の検収項目にします。
  • $U$が$M$に触れないかは、いちばん破れやすい条件です。 意欲の高い自治体は調査行動そのものも多い、というのが典型的な破れで、右端の崩れ方になります。$M$を「施策が指定した行動の実施ログ」に限定して定義するなど、$M$の測り方の設計で条件2・3に寄せる余地はあります。フロントドアは「使える場面を探す」道具というより、3条件をチェックリストにして「何が測れれば識別できるか」を逆算する設計図として使うのが実務的だと感じます。

参考(本棚)

  • Pearl (2009) Causal Inference in Statistics: An Overview. Statistics Surveys 3。第3.3.2節、調整に十分な集合が存在しなくても2段階の調整で識別できる例(式26)。今日の下敷き
  • Pearl, Glymour & Jewell『入門 統計的因果推論』。フロントドア基準の定義と喫煙→タール→肺がんの古典例
  • 黒木学『構造的因果モデルの基礎』。バックドア・フロントドアを含む識別条件の体系的な整理
  • 図はすべて合成データです。真値と未観測交絡$U$を私たちが握ることで、バックドア調整の不可能性とフロントドア式の復元、前提が壊れたときの崩れ方を再現しました。生成スクリプトは元ノート側に添付しています。

次回予告(Day 14)

バックドア、do演算子、フロントドアと、集団への介入 $P(Y|do(X))$ を扱う道具が出そろいました。ただ、介入の言葉で答えられるのはあくまで「みんなにやったらどうなるか」です。次回は最後の階段を登ります。「この患者が投薬を受けていなかったら、回復していたか」。すでに起きたことと反対の世界を個体単位で問う反事実を、SCMのノイズ復元による3ステップ(遡って、書き換えて、予測する)で計算します。関連・介入・反事実という因果のはしごの、いちばん上の段です。

0
0
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
0
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?