DAGと因果グラフで学んだように、因果関係を有向非巡回グラフ(DAG)で表現すると、変数間の因果構造を体系的に分析できます。しかし、DAGから因果効果を「計算」するための完全な方法論はまだ紹介していませんでした。
do計算(do-calculus)は、ジュディア・パール(Judea Pearl)が開発した、DAGの構造から介入確率 $P(Y \mid \text{do}(X = x))$ を観測確率に変換するための3つの推論ルールの体系です。do計算は完全(complete)であることが証明されており、DAGから因果効果が識別可能であれば、do計算のルールを適用することで必ず観測量で表現できます。
本記事では、do計算の3つのルール、バックドア基準とフロントドア基準を解説し、Pythonで調整公式の正しさを検証します。
本記事の内容
- do演算子の定義と意味
- do計算の3つのルール
- バックドア基準と調整公式の導出
- フロントドア基準
- Pythonシミュレーションによる調整公式の検証
前提知識
do演算子の定義
介入の数学的表現
因果推論入門で導入した $\text{do}$ 演算子を、DAGの言葉で正式に定義します。
$\text{do}(X = x)$ は、DAGにおいて $X$ への入射辺をすべて切断し、$X$ の値を $x$ に固定する操作を意味します。この操作をDAGのグラフ的な変換として表すと、元のDAG $G$ から $X$ への矢印をすべて除去した修正DAG $G_{\overline{X}}$ が得られます。
$$ \begin{equation} P(Y \mid \text{do}(X = x)) = P_{G_{\overline{X}}}(Y \mid X = x) \end{equation} $$
右辺は修正DAG $G_{\overline{X}}$ のもとでの $Y$ の分布です。
グラフ手術の詳細
do演算子の定義をより具体的に理解するために、「グラフ手術」(graph surgery)の手順を段階的に見ていきましょう。
元のDAG $G$ における各変数 $V_i$ の構造方程式が $V_i = f_i(\text{pa}_i, U_i)$ で与えられているとします。ここで $\text{pa}_i$ は $V_i$ の親変数の集合、$U_i$ は外生的な擾乱項です。
$\text{do}(X = x)$ の操作は、$X$ の構造方程式 $X = f_X(\text{pa}_X, U_X)$ を $X = x$ という定数関数に置き換えることに対応します。グラフ上では、$X$ への入射辺をすべて除去し、$X$ の値を $x$ に固定します。他のすべての変数の構造方程式はそのまま維持されます。
この手術の結果、修正後のモデルでは以下が成り立ちます。
- $X$ の値は常に $x$ である(外部から固定されている)
- $X$ の親変数は $X$ の値に影響を与えない(入射辺が切断されている)
- $X$ の子孫は、$X = x$ という固定値を通じて影響を受ける
- $X$ の非子孫は、$X$ への介入の影響を受けない
この「手術」の比喩は非常に有用です。外科医が体の一部を切除して人工的な部品に置き換えるように、do演算子は因果メカニズムの一部を破壊して外部からの値に置き換えるのです。
構造方程式モデル(SCM)との関係
do演算子の定義は、パールの構造因果モデル(Structural Causal Model, SCM)の枠組みにおいて最も自然に理解されます。SCMは以下の要素から構成されます。
- 外生変数 $\bm{U}$: モデルの外部から決定される変数
- 内生変数 $\bm{V}$: モデル内で決定される変数
- 構造方程式 $\bm{F}$: 各内生変数 $V_i$ を親変数と外生変数の関数として定義する方程式の集合
- 外生変数の確率分布 $P(\bm{U})$
SCMが与えられると、内生変数の同時分布 $P(\bm{V})$ が一意に決まります。そして、do演算子は構造方程式の一部を書き換える操作として定義され、介入後の同時分布 $P(\bm{V} \mid \text{do}(X = x))$ が得られます。
SCMの重要な特徴は、同じSCMからレベル1(観察)、レベル2(介入)、レベル3(反事実)のすべてのクエリに回答できることです。DAGはSCMの「因果構造」を視覚化したものであり、do計算はDAGの情報だけからレベル2のクエリに回答する方法を提供します。
観察と介入の違い(DAGの観点)
条件付き確率 $P(Y \mid X = x)$ と介入確率 $P(Y \mid \text{do}(X = x))$ の違いをDAGで理解しましょう。
$Z \to X \to Y$ かつ $Z \to Y$ というDAGを考えます。
観察 $P(Y \mid X = x)$: $X = x$ であるサンプルを選ぶ。このとき、$X$ の値は $Z$ の影響を受けているため、$X = x$ という条件付けにより $Z$ の分布も変わります。結果として、$Y$ の値は $X$ の因果効果と $Z$ を通じた交絡の両方の影響を受けます。
介入 $P(Y \mid \text{do}(X = x))$: $X$ の値を外部から $x$ に強制設定。$Z \to X$ の矢印が切断されるため、$Z$ の分布は変わりません。$Y$ の値は $X$ の因果効果のみの影響を受けます。
介入と条件付けの根本的な違い
do計算を理解するうえで、介入(intervention)と条件付け(conditioning)の違いを直感的に理解しておくことが決定的に重要です。
条件付けの意味
$P(Y \mid X = x)$ は「$X = x$ であるデータだけを抽出して $Y$ の分布を見る」ことです。この操作はデータのフィルタリングに過ぎず、世界の因果構造を変えるものではありません。$X = x$ という条件のもとでは、$X$ の原因となる変数の分布が変わります(なぜならば、$X$ の値を知ることでその原因について情報が得られるからです)。
介入の意味
$P(Y \mid \text{do}(X = x))$ は「$X$ の値を外部から $x$ に強制設定したときの $Y$ の分布」です。この操作は $X$ を生み出すメカニズムそのものを破壊し、新しい値で置き換えます。$X$ の原因となる変数は $X$ の値から何の情報も得られません(因果の矢印が切断されたため)。
なぜ区別が重要か
この区別が重要な理由は、交絡が存在する場合に $P(Y \mid X = x)$ と $P(Y \mid \text{do}(X = x))$ が異なる値を取るからです。例えば、「バリスタがいる街区は地価が高い」($P(\text{地価高い} \mid \text{バリスタあり})$ は高い) からといって、「バリスタを配置すれば地価が上がる」($P(\text{地価高い} \mid \text{do}(\text{バリスタ配置}))$ は高い)とは言えません。バリスタと地価の両方に影響する交絡因子(街区の富裕度など)が存在するためです。
do計算は、この「介入の効果」を「観測データ」だけから計算するための体系的な方法論を提供するのです。
do計算の3つのルール
ルール1: 不要な観察の挿入/削除
$$ \begin{equation} P(Y \mid \text{do}(X), Z, W) = P(Y \mid \text{do}(X), W) \quad \text{if } (Y \perp\!\!\!\perp Z \mid X, W)_{G_{\overline{X}}} \end{equation} $$
$G_{\overline{X}}$ は $X$ への入射辺をすべて除去したグラフです。このグラフで $Y$ と $Z$ がd分離されているならば、$Z$ の条件付けは不要です。
ルール2: 介入の観察への変換
$$ \begin{equation} P(Y \mid \text{do}(X), \text{do}(Z), W) = P(Y \mid \text{do}(X), Z, W) \quad \text{if } (Y \perp\!\!\!\perp Z \mid X, W)_{G_{\overline{X}, \underline{Z}}} \end{equation} $$
$G_{\overline{X}, \underline{Z}}$ は $X$ への入射辺と $Z$ からの出射辺の両方を除去したグラフです。この条件が成り立つとき、$Z$ への介入を観察に置き換えられます。
ルール3: 介入の削除
$$ \begin{equation} P(Y \mid \text{do}(X), \text{do}(Z), W) = P(Y \mid \text{do}(X), W) \quad \text{if } (Y \perp\!\!\!\perp Z \mid X, W)_{G_{\overline{X}, \overline{Z(S)}}} \end{equation} $$
ここで $Z(S)$ は $G_{\overline{X}}$ で $Z$ の祖先でないノードの集合から $Z$ を除外したものです。この条件が成り立つとき、$Z$ への介入を完全に削除できます。
do計算の完全性
do計算の3つのルールは完全です。つまり、DAGから因果効果が識別可能であれば、これらのルールを繰り返し適用することで $P(Y \mid \text{do}(X))$ を観測確率のみで表現できます。逆に、これらのルールで変換できなければ、そのDAGからは因果効果を識別できません。
この完全性はHuang and Valtorta(2006)、およびShpitser and Pearl(2006)によって証明されました。完全性の証明は、do計算の3つのルールが因果推論のための「十分な」ツールセットであることを保証するものです。つまり、これらのルール以外の追加ルールは不要であり、すべての識別可能な因果効果はこの3つのルールで計算できます。
do計算の適用戦略
3つのルールを実際に適用する際の一般的な戦略を整理しましょう。目標は $P(Y \mid \text{do}(X))$ のdo演算子をすべて除去し、観測確率のみの式に変換することです。
- 確率の法則で展開する: 周辺化($\sum_z$ を導入)や条件付けにより式を展開する
- ルール2で介入を観察に変換: 可能な箇所で $\text{do}(Z)$ を $Z$ に置き換える
- ルール3で不要な介入を削除: $Y$ に影響しない介入を除去する
- ルール1で不要な観察を削除: 余分な条件付けを除去する
この戦略を体系的に適用するのがIDアルゴリズムであり、手動での計算が困難な複雑なDAGに対しても自動的に識別式を導出できます。
バックドア基準と調整公式
バックドア基準の定式化
変数集合 $\bm{Z}$ が処置 $X$ からアウトカム $Y$ への因果効果に対してバックドア基準を満たすとは、以下の2条件が成り立つことです。
- $\bm{Z}$ は $X$ の子孫を含まない
- $\bm{Z}$ は $G_{\overline{X}}$($X$ への入射辺をすべて除去したグラフ)において $X$ と $Y$ を結ぶすべてのパスをd分離する
調整公式の導出
バックドア基準を満たす $\bm{Z}$ が存在するとき、因果効果は次の調整公式(adjustment formula / backdoor formula)で計算されます。
$$ \begin{equation} P(Y \mid \text{do}(X = x)) = \sum_{\bm{z}} P(Y \mid X = x, \bm{Z} = \bm{z})\,P(\bm{Z} = \bm{z}) \end{equation} $$
この式の導出をdo計算のルールを使って示します。
因果マルコフ条件と修正DAG $G_{\overline{X}}$ の同時分布から出発します。
$$ P(Y \mid \text{do}(X = x)) = \sum_{\bm{z}} P(Y \mid X = x, \bm{Z} = \bm{z}, \text{do}(X = x)) \cdot P(\bm{Z} = \bm{z} \mid \text{do}(X = x)) $$
バックドア基準の条件2により、$G_{\overline{X}}$ で $\bm{Z}$ と $Y$ がd分離されないため、$P(Y \mid X, \bm{Z}, \text{do}(X)) = P(Y \mid X, \bm{Z})$ となります(ルール2の適用)。
条件1($\bm{Z}$ は $X$ の子孫でない)により、$P(\bm{Z} \mid \text{do}(X)) = P(\bm{Z})$ となります(介入が $\bm{Z}$ に影響しない)。
これらを組み合わせると調整公式(5)が得られます。
バックドア基準が満たされない例
バックドア基準がすべてのDAGで使えるわけではないことを確認しましょう。以下のDAGを考えます。
$$ U \to X, \quad U \to Y, \quad X \to Y $$
ここで $U$ は未観測です。バックドアパスは $X \leftarrow U \to Y$ であり、$U$ で条件付ければブロックされますが、$U$ は観測できません。他に条件付けられる変数がないため、バックドア基準を満たす変数集合は存在しません。
このようなDAGでは、因果効果はバックドア基準だけでは識別不可能です。しかし、中間変数 $M$($X \to M \to Y$)が観測可能であれば、フロントドア基準が使える場合があります。このように、do計算の強みはバックドア基準を超えた識別を可能にすることにあります。
複数の有効な調整集合
バックドア基準を満たす変数集合は一般に一意ではありません。複数の変数集合がバックドア基準を満たす場合、どれを選ぶかは推定の効率性に影響します。
最小十分調整集合: バックドア基準を満たす最小の変数集合。余分な変数を含まないため、次元の呪いの影響が小さい。
最適調整集合: Henckel, Perovic, and Maathuis(2022)は、推定の漸近分散を最小化する最適な調整集合をDAGから導出する方法を提案しました。一般に、アウトカム $Y$ の予測に寄与する変数を含めることで効率が向上します。
実践的な指針としては、バックドア基準を満たしつつ、できるだけアウトカムの説明力が高い変数を含めることが推奨されます。
do計算のルールの直感的理解
ルール1の直感
ルール1は「余分な観測は無害に削除(または挿入)できる」ことを述べています。修正グラフ $G_{\overline{X}}$ で $Y$ と $Z$ がd分離されているとは、$X$ に介入した世界では $Z$ が $Y$ に追加情報を提供しないことを意味します。したがって $Z$ を条件に加えても外しても結果は変わりません。
例えば、$X \to Y$ かつ $Z$ が孤立($X$, $Y$ のいずれとも無関連)のDAGでは、$G_{\overline{X}}$ でも $Z$ と $Y$ はd分離されるため、$P(Y \mid \text{do}(X), Z) = P(Y \mid \text{do}(X))$ です。
ルール2の直感
ルール2は「介入を観察に置き換えてよい条件」を示します。$G_{\overline{X}, \underline{Z}}$ は $X$ への矢印と $Z$ からの矢印を除去したグラフです。このグラフで $Y$ と $Z$ がd分離されるということは、$Z$ からの因果パスを除いても $Y$ と $Z$ の関連が消えること、すなわち $Z$ と $Y$ の間の関連が $Z$ の因果効果のみに由来していることを意味します。このとき、$Z$ への介入と $Z$ の観察は同じ効果を持ちます。
直感的には、「$Z$ に交絡がない」場合に相当します。$Z$ と $Y$ の関連が純粋に $Z \to Y$ の因果効果によるものであれば、$Z$ を観察しても介入しても同じ情報が得られるのです。
ルール3の直感
ルール3は「介入を完全に削除できる条件」を示します。ある変数 $Z$ への介入が $Y$ に影響を及ぼさない場合、その介入を無視できます。
具体的には、$G_{\overline{X}, \overline{Z(S)}}$ で $Y$ と $Z$ がd分離されるとき、$Z$ への介入は $Y$ の分布に影響しないため、$\text{do}(Z)$ を削除できます。
do計算の適用例: フロントドア基準の導出
フロントドア基準の公式は、do計算の3つのルールを順に適用することで導出されます。DAG $U \to X$, $U \to Y$, $X \to M$, $M \to Y$($U$ は未観測)において、$P(Y \mid \text{do}(X))$ を求めたいとします。
ステップ1: $P(Y \mid \text{do}(X)) = \sum_m P(Y \mid \text{do}(X), M=m) P(M=m \mid \text{do}(X))$(確率の法則)
ステップ2: $P(M \mid \text{do}(X)) = P(M \mid X)$($X$ から $M$ にバックドアパスがないため、ルール2の適用)
ステップ3: $P(Y \mid \text{do}(X), M=m) = P(Y \mid \text{do}(M=m))$(ルール2と3の適用。介入 $\text{do}(X)$ は $M$ が固定されていれば $Y$ に影響しない)
ステップ4: $P(Y \mid \text{do}(M=m)) = \sum_{x’} P(Y \mid M=m, X=x’) P(X=x’)$($M$ から $Y$ へのバックドアパスを $X$ で調整。バックドア基準の適用)
組み合わせ: これらを組み合わせると、フロントドア公式
$$ P(Y \mid \text{do}(X=x)) = \sum_m P(M=m \mid X=x) \sum_{x’} P(Y \mid X=x’, M=m) P(X=x’) $$
が得られます。このように、do計算は複雑な識別問題を体系的に解く強力なツールです。
識別可能性と識別アルゴリズム
識別不可能なDAG
すべてのDAGで因果効果が識別可能であるわけではありません。典型的な識別不可能の例は、処置 $X$ とアウトカム $Y$ の間に未観測の交絡 $U$ があり、他に調整や媒介に使える変数がない場合です。
$$ U \to X, \quad U \to Y, \quad X \to Y $$
このDAGでは、$P(Y \mid \text{do}(X))$ を観測量だけでは表現できません。バックドアパス $X \leftarrow U \to Y$ をブロックするための変数が観測されておらず、フロントドア基準も適用できないためです。
ID アルゴリズム
TianとPearl(2002)は、IDアルゴリズム(Identification Algorithm)を提案しました。これはDAGが与えられたとき、因果効果が識別可能かどうかを判定し、識別可能な場合は観測量による表現を出力する完全なアルゴリズムです。
IDアルゴリズムの基本的なアイデアは、DAGのc-component($U$ の辺で接続された成分)に分解し、再帰的に識別を試みることです。アルゴリズムが失敗(FAIL)を返した場合、その因果効果はDAGの構造のもとでは識別不可能です。
IDアルゴリズムは causaleffect パッケージ(R)や DoWhy ライブラリ(Python)で実装されており、実際のDAGに対して自動的に識別可能性を判定し、識別式を出力できます。
条件付き因果効果の識別
IDアルゴリズムは $P(Y \mid \text{do}(X))$ だけでなく、共変量 $\bm{W}$ で条件付けた条件付き因果効果 $P(Y \mid \text{do}(X), \bm{W})$ の識別にも拡張できます。条件付き因果効果は因果効果の異質性(処置効果が個体の特性によって異なること)の分析に不可欠です。
例えば、「薬を飲むと血圧が下がるか」($P(Y \mid \text{do}(X))$)だけでなく、「年齢層ごとに薬の効果はどう異なるか」($P(Y \mid \text{do}(X), \text{Age})$)を問うことができます。
条件付き因果効果の識別には、IDアルゴリズムの拡張であるcIDアルゴリズム(conditional ID algorithm)が使われます。条件付き変数 $\bm{W}$ が処置 $X$ の子孫である場合、条件付き因果効果の識別は非自明であり、do計算の枠組みが不可欠です。
一般化調整基準
バックドア基準を一般化した一般化調整基準(generalized adjustment criterion, Perkovic et al., 2018)は、バックドア基準よりも広いクラスの変数集合で因果効果を識別できます。バックドア基準では「処置の子孫を含まない」という条件がありましたが、一般化調整基準ではこの条件を緩和しつつ、正しい識別を保証します。
フロントドア基準
バックドア基準が使えない場合
バックドア基準は、すべてのバックドアパスをブロックする変数集合が観測されている場合にのみ適用可能です。もし交絡変数が観測されていなければ、バックドア基準は使えません。
しかし、do計算はバックドア基準よりも強力です。フロントドア基準(front-door criterion)は、未観測の交絡が存在する場合でも因果効果を識別できる場合があることを示す重要な例です。
フロントドア基準の設定
以下のDAGを考えます。$U$ は未観測の交絡変数です。
$$ U \to X, \quad U \to Y, \quad X \to M \to Y $$
$X$ から $Y$ への因果パスは $X \to M \to Y$ です。バックドアパスは $X \leftarrow U \to Y$ であり、$U$ が観測されないためバックドア基準は使えません。
しかし、中間変数 $M$ が観測されていれば、フロントドア基準により因果効果が識別可能です。
フロントドア基準の条件
変数 $M$ がフロントドア基準を満たすための条件を正式に定義しましょう。変数集合 $\bm{M}$ が処置 $X$ からアウトカム $Y$ への因果効果に対してフロントドア基準を満たすとは、以下の3条件が成り立つことです。
- $X$ から $Y$ への因果パスがすべて $\bm{M}$ を経由する($\bm{M}$ は $X$ から $Y$ へのすべての有向パスを遮断する)
- $X$ から $\bm{M}$ へのバックドアパスが存在しない
- $\bm{M}$ から $Y$ へのすべてのバックドアパスが $X$ で条件付けることでブロックされる
条件1は、$X$ が $Y$ に影響を与えるのはすべて $\bm{M}$ を通じてであること($X$ から $Y$ への直接効果がないこと)を要求します。条件2は、$X$ から $\bm{M}$ への因果効果がバックドア基準なしで識別可能であることを保証します。条件3は、$\bm{M}$ から $Y$ への因果効果が $X$ で条件付けることで識別可能であることを保証します。
フロントドア基準の適用例
フロントドア基準の適用例として、喫煙とがんの関係を考えましょう。
喫煙($X$)が肺がん($Y$)に影響するかを調べたいが、遺伝的素因($U$、未観測)が喫煙傾向と肺がんリスクの両方に影響する交絡因子である場合を考えます。中間変数としてタール沈着量($M$)を測定できるとします。
DAGは $U \to X$, $U \to Y$, $X \to M \to Y$ です。このとき、タール沈着量 $M$ はフロントドア基準を満たします。なぜなら、喫煙が肺がんに影響するのはタール沈着を通じてのみ(条件1)、喫煙からタール沈着へのバックドアパスはない(条件2)、タール沈着から肺がんへのバックドアパスは喫煙で条件付けることでブロックされる(条件3)からです。
ただし、フロントドア基準が実際に適用可能な場面は限られています。条件1(すべての因果パスが中間変数を経由する)は特に厳しい仮定であり、処置からアウトカムへの直接効果が存在しないことを要求します。
フロントドア公式
$$ \begin{equation} P(Y \mid \text{do}(X = x)) = \sum_m P(M = m \mid X = x) \sum_{x’} P(Y \mid X = x’, M = m)\,P(X = x’) \end{equation} $$
この公式は2つのステップで因果効果を計算しています。内側の和は $P(Y \mid \text{do}(M = m))$ を計算し、外側の和は $P(M \mid \text{do}(X = x))$ との合成で $P(Y \mid \text{do}(X = x))$ を計算しています。
フロントドア公式の直感的理解
フロントドア公式がなぜ機能するかを、2つの段階に分けて直感的に理解しましょう。
第1段階: $X$ が $M$ に与える因果効果の推定
$P(M = m \mid X = x)$ は単なる条件付き確率ですが、フロントドア基準の条件2により、$X$ から $M$ へのバックドアパスが存在しないため、この条件付き確率は因果効果 $P(M = m \mid \text{do}(X = x))$ に等しくなります。つまり、未観測交絡 $U$ が存在しても、$X$ から $M$ への因果効果は直接推定できるのです。
これは直感的にも理解できます。タール沈着量は喫煙量によって直接決まり、遺伝的素因がタール沈着量に直接影響することはないためです。
第2段階: $M$ が $Y$ に与える因果効果の推定
$M$ から $Y$ への因果効果を推定する際には、バックドアパス $M \leftarrow X \leftarrow U \to Y$ が存在します。しかし、このパスは $X$ で条件付けることでブロックされます(条件3)。したがって、$X$ で条件付けた $P(Y \mid M = m, X = x’)$ を $P(X = x’)$ で加重平均することで、$P(Y \mid \text{do}(M = m))$ が推定できます。
2つの段階の組み合わせ
第1段階で $X \to M$ の因果効果を、第2段階で $M \to Y$ の因果効果を推定し、これらを合成することで、$X \to M \to Y$ という因果連鎖全体の効果を推定しているのです。未観測交絡 $U$ はこの推定のどの段階にも影響しないため、正しい因果効果が得られます。
do計算と因果推論ソフトウェア
計算の自動化
do計算の3つのルールを手動で適用するのは、複雑なDAGに対しては困難です。近年では、DAGを入力として因果効果の識別式を自動的に導出するソフトウェアが整備されています。
PythonのDoWhyライブラリは、DAGの指定から因果効果の識別、推定、反駁検証までの一連の因果推論ワークフローを提供します。内部的にはIDアルゴリズムを実装しており、与えられたDAGに対して因果効果が識別可能かどうかを自動的に判定し、識別可能な場合は推定に進みます。
Rのcausaleffectパッケージは、IDアルゴリズムの直接的な実装を提供し、DAGから識別式をシンボリックに出力します。研究者はこの識別式をもとに具体的な推定を行うことができます。
dagitty(R / Webアプリ)は、DAGの描画、d分離の判定、調整集合の列挙など、DAGベースの因果分析のための包括的なツールを提供します。バックドア基準を満たすすべての最小調整集合を自動的に列挙する機能は、実務的に非常に有用です。
これらのツールを活用することで、研究者はdo計算の数学的な複雑さを気にせずに、DAGの構造的な仮定と推定の結果に集中できます。
Pythonシミュレーション
調整公式の検証
import numpy as np
import matplotlib.pyplot as plt
np.random.seed(42)
n = 50000
# --- DAG: Z -> X, Z -> Y, X -> Y ---
Z = np.random.normal(0, 1, n) # 交絡変数(観測可能)
X = 0.5 + 0.8 * Z + np.random.normal(0, 1, n)
beta_true = 2.0
Y = 1.0 + beta_true * X + 1.5 * Z + np.random.normal(0, 1, n)
# --- 介入分布の計算(シミュレーション)---
x_val = 1.0
# 方法1: 実際の介入(do(X=1)をシミュレーション)
Z_do = np.random.normal(0, 1, n) # Zは介入の影響を受けない
Y_do = 1.0 + beta_true * x_val + 1.5 * Z_do + np.random.normal(0, 1, n)
E_Y_do_true = Y_do.mean()
# 方法2: 調整公式
n_bins = 50
Z_bins = np.linspace(Z.min(), Z.max(), n_bins + 1)
Z_mids = (Z_bins[:-1] + Z_bins[1:]) / 2
Z_widths = Z_bins[1] - Z_bins[0]
E_Y_adjusted = 0.0
for k in range(n_bins):
mask = (Z >= Z_bins[k]) & (Z < Z_bins[k+1])
if mask.sum() > 0:
# P(Y | X=x, Z=z_k)の推定
mask_x = mask & (np.abs(X - x_val) < 0.3)
if mask_x.sum() > 5:
E_Y_given_XZ = Y[mask_x].mean()
else:
continue
# P(Z=z_k)の推定
P_Z = mask.sum() / n
E_Y_adjusted += E_Y_given_XZ * P_Z
# 方法3: ナイーブな条件付き(バイアスあり)
mask_naive = np.abs(X - x_val) < 0.3
E_Y_naive = Y[mask_naive].mean()
print(f"E[Y | do(X={x_val})] (介入シミュレーション): {E_Y_do_true:.4f}")
print(f"E[Y | do(X={x_val})] (調整公式): {E_Y_adjusted:.4f}")
print(f"E[Y | X={x_val}] (ナイーブ条件付き): {E_Y_naive:.4f}")
print(f"\nバイアス (ナイーブ): {E_Y_naive - E_Y_do_true:.4f}")
print(f"バイアス (調整公式): {E_Y_adjusted - E_Y_do_true:.4f}")
# --- 可視化 ---
fig, axes = plt.subplots(1, 2, figsize=(13, 5.5))
# (a) E[Y|do(X=x)]の比較
ax = axes[0]
x_range = np.linspace(-2, 3, 20)
do_vals = []
naive_vals = []
adj_vals = []
for xv in x_range:
# 真の介入
Y_do_x = 1.0 + beta_true * xv + 1.5 * Z_do + np.random.normal(0, 1, n)
do_vals.append(Y_do_x.mean())
# ナイーブ
mask = np.abs(X - xv) < 0.3
if mask.sum() > 10:
naive_vals.append(Y[mask].mean())
else:
naive_vals.append(np.nan)
# 回帰ベースの調整
from sklearn.linear_model import LinearRegression
reg = LinearRegression().fit(np.column_stack([X, Z]), Y)
adj_vals.append(reg.predict([[xv, 0]])[0])
ax.plot(x_range, do_vals, "g-", linewidth=2, label=r"$E[Y|\mathrm{do}(X)]$ (true)")
ax.plot(x_range, naive_vals, "r--", linewidth=2, label=r"$E[Y|X]$ (naive)")
ax.plot(x_range, adj_vals, "b:", linewidth=2, label=r"Adjusted (backdoor)")
ax.set_xlabel("x", fontsize=12)
ax.set_ylabel("E[Y]", fontsize=12)
ax.set_title(r"$E[Y|\mathrm{do}(X=x)]$ vs $E[Y|X=x]$", fontsize=12)
ax.legend(fontsize=10)
ax.grid(True, alpha=0.3)
# (b) バイアスの図示
ax = axes[1]
bias_naive = np.array(naive_vals) - np.array(do_vals)
bias_adj = np.array(adj_vals) - np.array(do_vals)
ax.plot(x_range, bias_naive, "r-", linewidth=2, label="Naive bias")
ax.plot(x_range, bias_adj, "b-", linewidth=2, label="Adjusted bias")
ax.axhline(0, color="black", linewidth=0.5)
ax.set_xlabel("x", fontsize=12)
ax.set_ylabel("Bias", fontsize=12)
ax.set_title("Bias comparison", fontsize=12)
ax.legend(fontsize=10)
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.savefig("do_calculus_verification.png", dpi=150, bbox_inches="tight")
plt.show()
このシミュレーション結果から、do計算と調整公式の正しさが確認できます。
-
ナイーブな条件付き $E[Y \mid X=x]$ は介入分布 $E[Y \mid \text{do}(X=x)]$ からずれている(赤の破線が緑の実線より上)。これは交絡変数 $Z$ の影響が混入しているためです
-
調整公式(バックドア基準)による推定は介入分布に一致しています(青の点線が緑の実線に重なる)。$Z$ で条件付けて周辺化することで、交絡の影響が正しく除去されています
-
バイアスの図(右): ナイーブなバイアスは一定の正の値を取りますが、調整公式のバイアスはほぼゼロです
調整公式の回帰表現
回帰分析とバックドア調整
調整公式は層別化(stratification)や周辺化の形で表現されましたが、実務では回帰分析を使って調整を行うのが一般的です。線形モデルのもとでは、バックドア調整と回帰分析の関係が明瞭に理解できます。
DAG $Z \to X$, $Z \to Y$, $X \to Y$ において、バックドア基準は $Z$ での調整を要求します。線形モデル
$$ \begin{equation} Y = \alpha + \beta X + \gamma Z + \epsilon \end{equation} $$
を最小二乗法で推定すると、$\hat{\beta}$ がバックドア調整された因果効果の推定量になります。これは、回帰分析が $Z$ を「制御」して $X$ の純粋な効果を取り出すためです。
しかし、非線形モデルや効果の修飾(effect modification)が存在する場合は、回帰分析と調整公式が一致しないことがあります。調整公式は
$$ P(Y \mid \text{do}(X = x)) = \sum_z P(Y \mid X = x, Z = z) P(Z = z) $$
であり、$P(Y \mid X, Z)$ を各層ごとに推定してから $P(Z)$ で加重平均する「g-computation」とも呼ばれます。Robins(1986)は、g-computationが非線形モデルでも正しく因果効果を推定できることを示しました。
逆確率重み付け(IPW)
調整公式のもう一つの実装方法が逆確率重み付け(Inverse Probability Weighting, IPW)です。IPWは傾向スコア $e(Z) = P(X = 1 \mid Z)$ を使って、
$$ \begin{equation} E[Y \mid \text{do}(X = 1)] = E\left[\frac{X \cdot Y}{e(Z)}\right] \end{equation} $$
と推定します。各個体を「その処置を受ける確率の逆数」で重み付けすることで、ランダム化されたデータを「擬似的に」再現するのです。IPWの詳細は傾向スコア法で解説しています。
二重ロバスト推定量
g-computation とIPWを組み合わせた二重ロバスト推定量(Doubly Robust estimator, DR)は、アウトカムモデル $E[Y \mid X, Z]$ と傾向スコアモデル $P(X \mid Z)$ のどちらか一方が正しければ、因果効果を一致推定できます。
$$ \begin{equation} \hat{\tau}_{\text{DR}} = \frac{1}{n}\sum_{i=1}^n \left[\hat{\mu}_1(Z_i) – \hat{\mu}_0(Z_i) + \frac{X_i(Y_i – \hat{\mu}_1(Z_i))}{\hat{e}(Z_i)} – \frac{(1-X_i)(Y_i – \hat{\mu}_0(Z_i))}{1-\hat{e}(Z_i)}\right] \end{equation} $$
二重ロバスト推定量は、do計算のバックドア調整の実装として最も信頼性の高い方法の1つであり、因果推論の実務で広く推奨されています。
do計算と反事実
反事実クエリとdo計算の関係
パールの因果推論の階層(Pearl’s causal hierarchy / ladder of causation)では、因果推論を3つのレベルに分類します。
レベル1(連想): $P(Y \mid X = x)$(データの観察から)
レベル2(介入): $P(Y \mid \text{do}(X = x))$(do計算で扱う)
レベル3(反事実): $P(Y_x = y \mid X = x’, Y = y’)$(「もし $X$ を $x$ にしていたら $Y$ はどうなっていたか」)
do計算はレベル2の介入クエリを扱います。レベル3の反事実クエリはdo計算だけでは一般に答えられず、構造方程式モデル(SCM)の完全な仕様が必要です。
例えば、「この患者が薬を飲まなかったら回復していたか」(反事実)と「薬を飲ませたら平均的にどの程度回復するか」(介入)は異なるレベルの問いです。前者は個体レベルの反事実であり、後者は集団レベルの因果効果です。
自然直接効果と間接効果
媒介分析で扱う自然直接効果(NDE)と自然間接効果(NIE)は、レベル3の反事実量を含みます。例えば、NDE は $E[Y_{x, M_{x’}} – Y_{x’, M_{x’}}]$ と定義され、同一個体に対して処置と非処置の両方の反事実世界を参照するため、do計算のみでは表現できません。
これらの反事実量の識別には、SCMの構造的仮定(例えば、交差世界の独立性仮定)が追加で必要になります。
因果推論の梯子の意味
パールの因果推論の梯子(ladder of causation)は、統計学と因果推論の根本的な違いを浮き彫りにします。
レベル1の問い(例えば「この薬を飲んでいる患者はどの程度回復しているか」)は、純粋に観測データの集計で回答でき、通常の統計学の範囲内です。
レベル2の問い(例えば「この薬を飲ませたら患者は回復するか」)は、介入の効果を問うものであり、交絡がある場合は観測データの単純な集計では回答できません。do計算はこのレベルの問いに回答するための方法論です。
レベル3の問い(例えば「この患者が薬を飲んでいなかったら回復していたか」)は、同一個体の反事実的な状態を問うものであり、最も情報量が多い推論ですが、最も強い仮定を必要とします。
重要なのは、データ量をいくら増やしても、より低いレベルの情報からより高いレベルの情報を導くことは一般にはできないという点です。大量の観測データ(レベル1)があっても、因果効果(レベル2)は自動的には分かりません。因果推論には、データに加えて「因果構造の仮定」(DAG)が不可欠なのです。この認識こそが、パールの因果推論フレームワークの核心的な主張です。
まとめ
本記事では、パールのdo計算とバックドア基準・フロントドア基準を解説しました。
- do演算子はDAG上で介入を表現する操作であり、処置変数への入射辺を切断することに対応する
- do計算の3つのルールは介入確率を観測確率に変換するための完全な推論体系
- バックドア基準は最も一般的な因果効果の識別法であり、調整公式を導く
- フロントドア基準は未観測交絡が存在する場合でも因果効果を識別できる場合がある
次のステップとして、以下の記事も参考にしてください。
- 媒介分析 — 直接効果と間接効果の分解
- 選択バイアスとヘックマン補正 — 選択に関するバイアスの制御
- 因果推論の感度分析 — 未観測交絡の影響評価