「アイスクリームの売上が増えると、水難事故も増える」――この2つの変数の間には明確な正の相関があります。では、アイスクリームを販売禁止にすれば水難事故は減るのでしょうか。もちろん、そんなことはありません。夏になると気温が上がり、アイスクリームの消費量も水遊びの機会も増える、というだけのことです。ここには相関はあるけれど因果はありません。

この「相関と因果の区別」は、日常生活だけでなく、科学・政策・ビジネスの意思決定において極めて重要です。新薬が本当に効くのか（臨床試験）、広告費を増やすと売上が伸びるのか（マーケティング）、教育プログラムが学力を向上させるのか（教育政策）——これらの問いはすべて、単なる相関ではなく因果的な効果を知りたい場面です。

因果推論（causal inference）は、観測データや実験データから因果関係を正しく推定するための統計学・計量経済学・疫学にまたがる方法論です。この分野を学ぶと、以下のような応用が開けます。

• 医学・疫学: ランダム化比較試験（RCT）の設計と解析、観察研究からの因果効果の推定
• 経済学・政策評価: 最低賃金の引き上げが雇用に与える効果、教育支援策の効果測定
• 機械学習: 予測モデルから因果モデルへの発展、公平性（fairness）の確保
• ビジネス: A/Bテストの設計、マーケティング施策の因果効果の測定

本記事では、因果推論の最も基本的な問いである「相関と因果はどう違うのか」を数理的に掘り下げます。交絡（confounding）の概念を導入し、シンプソンのパラドックスという直感に反する現象を通じて、なぜ単純な相関分析では因果関係を捉えられないのかを明らかにします。

本記事の内容

• 相関と因果の直感的な違い
• 交絡変数が相関を「偽の因果」に見せかけるメカニズム
• シンプソンのパラドックスの数理的構造
• 因果推論の基本的な枠組み（介入と観察の区別）
• Pythonシミュレーションによる交絡の可視化と因果効果の推定

前提知識

この記事を読む前に、以下の記事を読んでおくと理解が深まります。

• 相関係数 — ピアソンの相関係数の定義と性質
• 条件付き確率 — 条件付き確率とベイズの定理の基礎
• 単回帰分析 — 回帰分析の基本的な枠組み

相関とは何か — 統計的関連性の定量化

相関の直感的な意味

まず「相関がある」とは何を意味するのかを確認しましょう。日常語では「AとBには相関がある」と言えば「AとBは何か関係がありそうだ」くらいの意味ですが、統計学では明確な定義があります。

2つの確率変数 $X$ と $Y$ の間の相関とは、一方の値を知ることで他方の値について何らかの情報が得られるという統計的な依存関係のことです。最も標準的にはピアソンの相関係数で定量化されます。

$$ \begin{equation} r_{XY} = \frac{\text{Cov}(X, Y)}{\sqrt{\text{Var}(X)\,\text{Var}(Y)}} = \frac{E[(X – \mu_X)(Y – \mu_Y)]}{\sigma_X \sigma_Y} \end{equation} $$

ここで $\text{Cov}(X, Y)$ は $X$ と $Y$ の共分散、$\mu_X, \mu_Y$ はそれぞれの期待値、$\sigma_X, \sigma_Y$ はそれぞれの標準偏差です。$r_{XY}$ は $-1$ から $1$ の値を取り、$|r_{XY}|$ が大きいほど線形な関連が強いことを示します。

しかし、ここで極めて重要な点があります。相関は2つの変数が一緒に変動するという事実を記述するだけであり、一方が他方を引き起こしているかどうかについては何も言っていません。アイスクリームの売上 $X$ と水難事故 $Y$ の相関係数がたとえ $r_{XY} = 0.85$ であっても、この数値からは「アイスクリームが水難事故を引き起こす」とも「水難事故がアイスクリーム消費を増やす」とも言えないのです。

相関が生じるメカニズム

では、2つの変数の間に相関が生じるのはどのような場合でしょうか。大きく分けて、3つのメカニズムがあります。

メカニズム1: 直接的な因果関係

$X$ が $Y$ を直接引き起こす場合、$X$ と $Y$ の間には当然相関が生じます。例えば、喫煙量 $X$ と肺がんリスク $Y$ の関係がこれにあたります。$X$ を増やすと $Y$ も増える、という因果的な関係が相関として現れています。

メカニズム2: 逆因果（reverse causation）

$Y$ が $X$ を引き起こしている場合も、$X$ と $Y$ の間には相関が生じます。例えば、「警察官の数が多い地域ほど犯罪率が高い」という相関が観測されたとき、警察官が犯罪を引き起こしているのではなく、犯罪率が高い地域に多くの警察官が配置されている、つまり因果の矢印が逆向きなのです。

メカニズム3: 交絡（confounding）

3つ目のメカニズムが最も厄介で、かつ最も頻繁に遭遇するものです。第3の変数 $Z$（交絡変数、confounder）が $X$ と $Y$ の両方に影響を与えている場合、$X$ と $Y$ の間に $X \to Y$ の因果関係がなくても相関が生じます。

アイスクリームと水難事故の例では、気温 $Z$ が交絡変数です。気温が高いとアイスクリームの売上 $X$ も水難事故 $Y$ も増えるため、$X$ と $Y$ の間に見かけの相関が発生しています。

$$ Z \to X, \quad Z \to Y \quad \Longrightarrow \quad X \text{ と } Y \text{ に相関が生じる} $$

この3つのメカニズムのうち、メカニズム1だけが「$X$ を変えれば $Y$ が変わる」という因果的な意味を持ちます。メカニズム2と3では、$X$ を操作しても $Y$ は変わりません。因果推論の中心的な課題は、観測された相関がどのメカニズムで生じているのかを識別することにあります。

では、交絡変数がどのようにして「偽の因果」を生み出すのか、数理的にもう少し詳しく見ていきましょう。

交絡変数の数理 — なぜ見かけの相関が生まれるのか

交絡の数学的な構造

交絡変数の効果を数学的に理解するために、以下のような線形モデルを考えます。処置変数 $X$、アウトカム $Y$、交絡変数 $Z$ の3変数があり、真のデータ生成メカニズムが次のように与えられているとします。

$$ \begin{equation} X = \alpha_Z Z + \epsilon_X \end{equation} $$

$$ \begin{equation} Y = \beta X + \gamma Z + \epsilon_Y \end{equation} $$

ここで $\epsilon_X$ と $\epsilon_Y$ は互いに独立な誤差項（期待値ゼロ）、$\alpha_Z$ は $Z$ が $X$ に与える影響、$\beta$ は $X$ が $Y$ に与える真の因果効果、$\gamma$ は $Z$ が $Y$ に直接与える影響です。

この設定で、$Z$ を無視して $X$ と $Y$ の単回帰を行うとどうなるでしょうか。$Y$ を $X$ に回帰したときの回帰係数 $\hat{\beta}_{\text{naive}}$ は、次のようになります。

$\hat{\beta}_{\text{naive}}$ の期待値を計算するために、まず $\text{Cov}(X, Y)$ を求めます。(2)式を(3)式に代入すると次のようになります。

$$ Y = \beta(\alpha_Z Z + \epsilon_X) + \gamma Z + \epsilon_Y = (\beta \alpha_Z + \gamma)Z + \beta \epsilon_X + \epsilon_Y $$

したがって、$X$ と $Y$ の共分散を計算すると次のようになります。

$$ \text{Cov}(X, Y) = \text{Cov}(\alpha_Z Z + \epsilon_X, (\beta\alpha_Z + \gamma)Z + \beta\epsilon_X + \epsilon_Y) $$

$Z$, $\epsilon_X$, $\epsilon_Y$ が互いに独立であることを使って展開すると次のようになります。

$$ \text{Cov}(X, Y) = \alpha_Z(\beta\alpha_Z + \gamma)\text{Var}(Z) + \beta\text{Var}(\epsilon_X) $$

一方、$X$ の分散は次のとおりです。

$$ \text{Var}(X) = \alpha_Z^2\text{Var}(Z) + \text{Var}(\epsilon_X) $$

$Y$ を $X$ に回帰した回帰係数の確率極限は次のようになります。

$$ \begin{equation} \hat{\beta}_{\text{naive}} \xrightarrow{p} \frac{\text{Cov}(X, Y)}{\text{Var}(X)} = \beta + \frac{\alpha_Z \gamma \text{Var}(Z)}{\alpha_Z^2\text{Var}(Z) + \text{Var}(\epsilon_X)} \end{equation} $$

右辺の第2項が交絡バイアス（confounding bias）です。この項は $\alpha_Z \neq 0$（$Z$ が $X$ に影響する）かつ $\gamma \neq 0$（$Z$ が $Y$ に影響する）のときゼロではなくなります。つまり、交絡変数 $Z$ を無視して $X$ と $Y$ の関係を推定すると、真の因果効果 $\beta$ からずれた値が得られてしまうのです。

交絡バイアスの直感的な解釈

式(4)の結果を直感的に解釈してみましょう。交絡バイアスの符号は $\alpha_Z \gamma$ の符号で決まります。

• $\alpha_Z > 0$ かつ $\gamma > 0$（または $\alpha_Z < 0$ かつ $\gamma < 0$）の場合: 交絡バイアスは正。真の効果 $\beta$ よりも過大な推定値が得られます。つまり、$X$ が $Y$ に与える効果を「水増し」してしまいます

• $\alpha_Z > 0$ かつ $\gamma < 0$（または $\alpha_Z < 0$ かつ $\gamma > 0$）の場合: 交絡バイアスは負。真の効果 $\beta$ よりも過小な推定値が得られます。場合によっては、効果の符号すら逆転してしまうこともあります

極端なケースとして、真の因果効果がゼロ（$\beta = 0$）であっても、交絡バイアスにより $\hat{\beta}_{\text{naive}} \neq 0$ となることがあります。これがまさに「相関はあるが因果はない」状況の数学的な表現です。

交絡の制御

交絡バイアスを除去する最も直接的な方法は、交絡変数 $Z$ を制御する（control for）ことです。具体的には、$Y$ を $X$ と $Z$ の両方に回帰する重回帰分析を行います。

$$ Y = \beta X + \gamma Z + \epsilon_Y $$

この重回帰モデルにおける $X$ の回帰係数 $\hat{\beta}$ は、$Z$ の影響を除いた後の $X$ と $Y$ の関係、すなわち真の因果効果 $\beta$ に一致します（モデルが正しく指定されている場合）。

ただし、この方法が機能するためには、すべての交絡変数を観測できていて、かつモデルに含めていることが必要です。もし重要な交絡変数が観測されていなければ、重回帰分析を行っても交絡バイアスは残ります。この「未観測の交絡」の問題こそが、因果推論を単なる回帰分析よりもはるかに困難にしている核心的な課題です。

交絡の問題を見てきましたが、交絡がいかに直感に反する現象を引き起こすかを示す有名な例として、シンプソンのパラドックスがあります。次にこれを詳しく見てみましょう。

シンプソンのパラドックス — 交絡が生む直感への裏切り

パラドックスの概要

シンプソンのパラドックス（Simpson’s paradox）は、全体で成り立つ相関の方向が、データをグループに分けると逆転するという現象です。これは交絡変数の影響が引き起こす最も劇的な例であり、「データを見るだけでは因果関係を判断できない」ということを強く示唆しています。

具体的な例で考えましょう。ある大学の入学審査で、男女別の合格率を調べたとします。

全体の合格率:

全体で見ると、女性の合格率（70%）が男性の合格率（60%）を上回っています。一見すると、女性が有利に見えます。

ところが、学部別に分けてみると状況が逆転します。

学部Aの合格率:

学部Bの合格率:

学部Aでは男性の合格率（80%）が女性（70%）を上回り、学部Bでも男性の合格率（53.3%）は女性（70%）を下回っています。つまり、どの学部を見ても男性の合格率が低い（あるいは同等）にもかかわらず、全体では男性の合格率の方が低い——というパラドキシカルな状況ではないのですが、ポイントは「全体の数字」と「層別の数字」が矛盾するように見える点です。

パラドックスの数学的な説明

上の例の数値を注意深く見ると、からくりが分かります。男性は合格率の低い学部Bに多く受験しており（800人中600人）、女性は合格率の高い学部Aに多く受験しています（400人中300人）。つまり、「受験する学部の選択」が交絡変数として作用しているのです。

より一般的に、シンプソンのパラドックスを数学的に定式化してみましょう。処置 $T \in \{0, 1\}$、アウトカム $Y \in \{0, 1\}$、層別変数 $S \in \{1, 2, \ldots, K\}$ があるとします。全体での処置群のアウトカム率と各層でのアウトカム率の関係は、次の重み付き平均で表されます。

$$ \begin{equation} P(Y=1 \mid T=1) = \sum_{k=1}^{K} P(Y=1 \mid T=1, S=k)\,P(S=k \mid T=1) \end{equation} $$

同様に対照群については次のようになります。

$$ \begin{equation} P(Y=1 \mid T=0) = \sum_{k=1}^{K} P(Y=1 \mid T=0, S=k)\,P(S=k \mid T=0) \end{equation} $$

ここで重要なのは、重みが処置群と対照群で異なることです。処置群では $P(S=k \mid T=1)$ で重み付けされ、対照群では $P(S=k \mid T=0)$ で重み付けされます。この重みの違いにより、各層では $P(Y=1 \mid T=1, S=k) < P(Y=1 \mid T=0, S=k)$ であっても、全体では $P(Y=1 \mid T=1) > P(Y=1 \mid T=0)$ となり得るのです。

直感的に言えば、層ごとの効果の大きさと、各層に占める割合の偏りの両方が結果に影響するため、層ごとの結果を単純に合計しても全体の結果とは一致しないのです。

因果推論の観点からの解釈

シンプソンのパラドックスは、「全体のデータを見るべきか、層別のデータを見るべきか」という問いを突きつけます。この問いに対する答えは、因果構造に依存します。

もし層別変数 $S$ が交絡変数（処置 $T$ とアウトカム $Y$ の両方の共通原因）であれば、$S$ で層別した結果を見るべきです。層別することで交絡の影響を除去できるからです。

逆に、$S$ が処置 $T$ の結果として生じる中間変数（mediator）であれば、$S$ で層別すべきではありません。層別すると因果効果を歪めてしまうからです。

どちらのケースに該当するかを判断するには、変数間の因果構造を理解する必要があります。これは純粋にデータだけからは判断できず、因果に関する知識や仮定が不可欠です。因果推論が単なるデータ分析を超えた「因果的な思考」を要求する理由はここにあります。

シンプソンのパラドックスは、交絡の問題がいかに根深いかを示しています。では、因果推論はこの問題にどのようにアプローチするのでしょうか。次に、因果推論の基本的な枠組みを見ていきます。

介入と観察 — 因果推論の基本的な枠組み

「見ること」と「行うこと」の違い

因果推論の最も基本的な区別は、観察（observing）と介入（intervening）の区別です。これは統計学の創始者たちからパール（Judea Pearl）に至るまで、因果推論の理論全体を貫く根本的な概念です。

• 観察: データを受動的に収集する。$X = x$ であるサンプルを見て、$Y$ の値を調べる
• 介入: 能動的に $X$ の値を $x$ に設定して、$Y$ の変化を見る

観察で得られるのは条件付き確率 $P(Y \mid X = x)$ です。一方、介入で得られるのは介入確率 $P(Y \mid \text{do}(X = x))$ です。パールの提唱した $\text{do}$ 記法（do-notation）は、この2つを明確に区別するために導入されました。

なぜこの区別が重要なのでしょうか。具体例で考えてみましょう。

病院のデータから「手術を受けた患者の方が回復率が低い」という相関が観測されたとします。

$$ P(\text{回復} \mid \text{手術あり}) < P(\text{回復} \mid \text{手術なし}) $$

この相関から「手術は回復を妨げる」と結論できるでしょうか。もちろんできません。手術を受ける患者は重症であることが多く、「重症度」が交絡変数として作用しています。重症度を $Z$ とすると、次のような因果構造があります。

$$ Z \to X \text{（重症 → 手術を選択）}, \quad Z \to Y \text{（重症 → 回復率が低い）} $$

観察データから得られる $P(Y \mid X)$ には交絡変数 $Z$ の影響が混入しています。介入確率 $P(Y \mid \text{do}(X))$ は、$X$ の値を外部から強制的に設定した場合の $Y$ の分布であり、$Z \to X$ の矢印を切断した状況に対応します。

介入の数学的な定義

介入 $\text{do}(X = x)$ を数学的に定義するには、構造因果モデル（structural causal model, SCM）の枠組みが便利です。以下の連立方程式で表されるモデルを考えます。

$$ Z = f_Z(\epsilon_Z) $$ $$ X = f_X(Z, \epsilon_X) $$ $$ Y = f_Y(X, Z, \epsilon_Y) $$

ここで $\epsilon_Z, \epsilon_X, \epsilon_Y$ は外生変数（exogenous variables）です。

観察とは、このシステムをそのまま動かして、$X = x$ となったサンプルを選び出すことです。

介入 $\text{do}(X = x)$ とは、2番目の方程式 $X = f_X(Z, \epsilon_X)$ を$X = x$ に強制的に置き換えることです。

$$ Z = f_Z(\epsilon_Z) $$ $$ X = x \quad \text{（強制的に設定）} $$ $$ Y = f_Y(x, Z, \epsilon_Y) $$

介入により $Z \to X$ の因果関係が断ち切られ、$X$ はもはや $Z$ に依存しません。しかし $Z \to Y$ の関係はそのまま残るため、$Z$ の影響は $Y$ に直接作用します。介入確率 $P(Y \mid \text{do}(X = x))$ は、この修正されたモデルのもとでの $Y$ の分布です。

因果効果の定義

因果効果（causal effect）は、介入確率を用いて定義されます。最も基本的なものは平均因果効果（Average Causal Effect, ACE）または平均処置効果（Average Treatment Effect, ATE）と呼ばれ、次のように定義されます。

$$ \begin{equation} \text{ATE} = E[Y \mid \text{do}(X = 1)] – E[Y \mid \text{do}(X = 0)] \end{equation} $$

ここで $X = 1$ は「処置あり」、$X = 0$ は「処置なし」を表します。ATEは、集団全体に処置を施した場合と施さなかった場合のアウトカムの期待値の差であり、因果効果を集団レベルで要約した量です。

重要な点は、$\text{ATE} \neq E[Y \mid X = 1] – E[Y \mid X = 0]$ であることです。右辺は単なる条件付き期待値の差であり、交絡が存在する場合には因果効果を正しく反映しません。因果推論の課題は、観測データ（条件付き確率）から介入確率を推定することにあります。

この目標を達成するための最も直接的な方法がランダム化比較試験ですが、観察データからも一定の条件のもとで因果効果を推定できることがこの分野の面白さです。次に、Pythonシミュレーションを通じて、交絡の問題と因果効果の推定を実際に体験してみましょう。

Pythonで体験する交絡と因果効果の推定

シミュレーション1: 交絡がない場合

まず、交絡変数がない場合を確認します。処置 $X$ がアウトカム $Y$ に直接効果を持つが、交絡変数が存在しないシンプルなケースです。

import numpy as np
import matplotlib.pyplot as plt

np.random.seed(42)
n = 2000

# データ生成: 交絡なし
# X -> Y の直接的な因果関係のみ
X = np.random.binomial(1, 0.5, n)  # 処置をランダムに割り当て
epsilon_Y = np.random.normal(0, 1, n)
beta_true = 2.0  # 真の因果効果
Y = beta_true * X + epsilon_Y

# 単純な条件付き期待値の差で推定
Y_treated = Y[X == 1]
Y_control = Y[X == 0]
naive_estimate = Y_treated.mean() - Y_control.mean()

print(f"真の因果効果 (ATE): {beta_true:.2f}")
print(f"単純推定値 (条件付き期待値の差): {naive_estimate:.4f}")
print(f"バイアス: {naive_estimate - beta_true:.4f}")

# 可視化
fig, axes = plt.subplots(1, 2, figsize=(12, 5))

# 箱ひげ図
ax = axes[0]
bp = ax.boxplot([Y_control, Y_treated], labels=["Control (X=0)", "Treated (X=1)"],
                patch_artist=True)
bp["boxes"][0].set_facecolor("lightblue")
bp["boxes"][1].set_facecolor("lightsalmon")
ax.set_ylabel("Outcome Y", fontsize=12)
ax.set_title("No confounding: simple comparison works", fontsize=12)
ax.grid(True, alpha=0.3)

# ヒストグラム
ax = axes[1]
ax.hist(Y_control, bins=30, alpha=0.5, color="blue", label="Control", density=True)
ax.hist(Y_treated, bins=30, alpha=0.5, color="red", label="Treated", density=True)
ax.axvline(Y_control.mean(), color="blue", linestyle="--", linewidth=2,
           label=f"Control mean: {Y_control.mean():.2f}")
ax.axvline(Y_treated.mean(), color="red", linestyle="--", linewidth=2,
           label=f"Treated mean: {Y_treated.mean():.2f}")
ax.set_xlabel("Outcome Y", fontsize=12)
ax.set_ylabel("Density", fontsize=12)
ax.set_title("Distribution of outcomes", fontsize=12)
ax.legend(fontsize=9)
ax.grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig("causal_no_confounding.png", dpi=150, bbox_inches="tight")
plt.show()

交絡がない場合、処置群と対照群の条件付き期待値の差がそのまま真の因果効果を正しく推定できています。推定値は真の値 $\beta = 2.0$ の近くにあり、バイアスはほぼゼロです。これはランダム化比較試験（RCT）に対応する状況で、処置の割り当てがアウトカムに影響する他の要因と独立であるため、交絡が生じないのです。

次に、交絡変数が存在する場合に何が起こるかを見てみましょう。

シミュレーション2: 交絡がある場合

同じ因果効果 $\beta = 2.0$ の設定に交絡変数 $Z$ を追加します。$Z$ が処置の割り当てにもアウトカムにも影響する状況です。

import numpy as np
import matplotlib.pyplot as plt
from sklearn.linear_model import LinearRegression

np.random.seed(42)
n = 2000

# データ生成: 交絡あり
Z = np.random.normal(0, 1, n)  # 交絡変数
prob_X = 1 / (1 + np.exp(-1.5 * Z))  # Zが大きいほどX=1になりやすい
X = np.random.binomial(1, prob_X, n)

beta_true = 2.0   # 真の因果効果
gamma_true = 3.0   # Zの直接効果
epsilon_Y = np.random.normal(0, 1, n)
Y = beta_true * X + gamma_true * Z + epsilon_Y

# 1. ナイーブな推定（Zを無視）
Y_treated = Y[X == 1]
Y_control = Y[X == 0]
naive_estimate = Y_treated.mean() - Y_control.mean()

# 2. 重回帰による推定（Zを制御）
reg = LinearRegression()
reg.fit(np.column_stack([X, Z]), Y)
adjusted_estimate = reg.coef_[0]

print(f"真の因果効果 (ATE): {beta_true:.2f}")
print(f"ナイーブ推定値 (Zを無視): {naive_estimate:.4f}")
print(f"調整済み推定値 (Zを制御): {adjusted_estimate:.4f}")
print(f"ナイーブ推定のバイアス: {naive_estimate - beta_true:.4f}")
print(f"調整済み推定のバイアス: {adjusted_estimate - beta_true:.4f}")

このコードの出力を確認すると、ナイーブな推定値は真の因果効果 $\beta = 2.0$ から大きくずれており、交絡バイアスが発生していることがわかります。一方、重回帰分析で交絡変数 $Z$ を制御すると、推定値は真の因果効果にほぼ一致します。これは式(4)で理論的に予測した結果と一致しています。

シミュレーション3: 交絡バイアスの可視化

交絡バイアスの構造をより直感的に理解するために、散布図で可視化してみましょう。

import numpy as np
import matplotlib.pyplot as plt

np.random.seed(42)
n = 2000

# データ生成（連続的な処置変数の場合）
Z = np.random.normal(0, 1, n)
epsilon_X = np.random.normal(0, 0.5, n)
X = 1.5 * Z + epsilon_X  # ZがXに影響

beta_true = 1.0
gamma_true = 2.0
epsilon_Y = np.random.normal(0, 0.5, n)
Y = beta_true * X + gamma_true * Z + epsilon_Y

fig, axes = plt.subplots(1, 3, figsize=(16, 5))

# (a) XとYの単純な散布図
ax = axes[0]
scatter = ax.scatter(X, Y, c=Z, cmap="coolwarm", alpha=0.4, s=15)
z_sorted = np.sort(X)
# ナイーブな回帰直線
from numpy.polynomial.polynomial import polyfit
b_naive, a_naive = polyfit(X, Y, 1)
ax.plot(z_sorted, a_naive * z_sorted + b_naive, "k-", linewidth=2,
        label=f"Naive slope: {a_naive:.2f}")
ax.set_xlabel("X (Treatment)", fontsize=12)
ax.set_ylabel("Y (Outcome)", fontsize=12)
ax.set_title("(a) Naive: X vs Y\n(confounding bias)", fontsize=11)
ax.legend(fontsize=10)
ax.grid(True, alpha=0.3)
plt.colorbar(scatter, ax=ax, label="Z (confounder)")

# (b) Zの影響を除いた後のXとYの関係
from sklearn.linear_model import LinearRegression
reg_X = LinearRegression().fit(Z.reshape(-1, 1), X)
reg_Y = LinearRegression().fit(Z.reshape(-1, 1), Y)
X_residual = X - reg_X.predict(Z.reshape(-1, 1))
Y_residual = Y - reg_Y.predict(Z.reshape(-1, 1))

ax = axes[1]
ax.scatter(X_residual, Y_residual, alpha=0.4, s=15, color="steelblue")
b_adj, a_adj = polyfit(X_residual, Y_residual, 1)
x_sorted = np.sort(X_residual)
ax.plot(x_sorted, a_adj * x_sorted + b_adj, "k-", linewidth=2,
        label=f"Adjusted slope: {a_adj:.2f}")
ax.set_xlabel("X residual (Z removed)", fontsize=12)
ax.set_ylabel("Y residual (Z removed)", fontsize=12)
ax.set_title("(b) After removing Z\n(causal effect recovered)", fontsize=11)
ax.legend(fontsize=10)
ax.grid(True, alpha=0.3)

# (c) 推定値の比較
ax = axes[2]
labels = ["True effect", "Naive estimate", "Adjusted estimate"]
values = [beta_true, a_naive, a_adj]
colors = ["green", "red", "blue"]
bars = ax.bar(labels, values, color=colors, alpha=0.7, edgecolor="black")
ax.axhline(beta_true, color="green", linestyle="--", linewidth=1.5)
ax.set_ylabel("Estimated causal effect", fontsize=12)
ax.set_title("(c) Comparison of estimates", fontsize=11)
for bar, val in zip(bars, values):
    ax.text(bar.get_x() + bar.get_width() / 2, bar.get_height() + 0.05,
            f"{val:.2f}", ha="center", fontsize=11)
ax.grid(True, alpha=0.3, axis="y")

plt.tight_layout()
plt.savefig("causal_confounding_visualization.png", dpi=150, bbox_inches="tight")
plt.show()

print(f"真の因果効果: {beta_true:.2f}")
print(f"ナイーブ推定値: {a_naive:.4f}")
print(f"調整済み推定値: {a_adj:.4f}")

この3パネルの図から、交絡バイアスの構造が視覚的に明確になります。

1. 図(a): $X$ と $Y$ の散布図を単純にプロットすると、色（交絡変数 $Z$ の値）でグラデーションがかかっています。$Z$ が大きい（赤い点）ほど $X$ も $Y$ も大きくなるため、$X$ と $Y$ の見かけの関係（黒い直線）は真の因果効果よりも急な傾きを示しています。ナイーブな推定値が真の値 $\beta = 1.0$ よりも大きいのは、交絡バイアスが正であるためです

2. 図(b): $X$ と $Y$ のそれぞれから $Z$ の影響を回帰で除去した残差同士をプロットすると、$Z$ の影響が取り除かれ、純粋な $X \to Y$ の関係が現れます。この調整済みの傾きは真の因果効果 $\beta = 1.0$ にほぼ一致しています。これは重回帰分析でフリッシュ・ウォーの定理（Frisch-Waugh-Lovell定理）として知られる結果の視覚的な確認です

3. 図(c): 棒グラフで比較すると、ナイーブな推定値が真の効果から大きくずれている一方、交絡変数を制御した推定値が真の効果を正確に回復していることが一目瞭然です

シミュレーション4: シンプソンのパラドックスの再現

最後に、シンプソンのパラドックスをシミュレーションで再現し、因果推論の重要性を確認します。

import numpy as np
import matplotlib.pyplot as plt

np.random.seed(42)

# 2つのグループ（層）を持つデータを生成
# 層 S: 1（軽症）, 2（重症）
n = 2000

# グループの割り当て
S = np.random.choice([1, 2], size=n, p=[0.5, 0.5])

# 処置の割り当て（重症ほど処置を受けやすい = 交絡）
prob_T = np.where(S == 1, 0.3, 0.7)  # 軽症30%、重症70%が処置
T = np.random.binomial(1, prob_T, n)

# アウトカム生成（処置の真の効果は正: +0.2）
# ただし重症ほどアウトカムが悪い（-0.4の効果）
base_rate = np.where(S == 1, 0.7, 0.3)  # 軽症は基本回復率70%、重症は30%
treatment_effect = 0.2  # 真の処置効果（正）
prob_Y = np.clip(base_rate + treatment_effect * T, 0, 1)
Y = np.random.binomial(1, prob_Y, n)

# 全体の処置効果（ナイーブ）
overall_treated = Y[T == 1].mean()
overall_control = Y[T == 0].mean()

# 層別の処置効果
mild_treated = Y[(T == 1) & (S == 1)].mean()
mild_control = Y[(T == 0) & (S == 1)].mean()
severe_treated = Y[(T == 1) & (S == 2)].mean()
severe_control = Y[(T == 0) & (S == 2)].mean()

print("=== シンプソンのパラドックス ===")
print(f"\n全体:")
print(f"  処置群の回復率: {overall_treated:.4f}")
print(f"  対照群の回復率: {overall_control:.4f}")
print(f"  ナイーブな効果: {overall_treated - overall_control:.4f}")
print(f"\n軽症群 (S=1):")
print(f"  処置群の回復率: {mild_treated:.4f}")
print(f"  対照群の回復率: {mild_control:.4f}")
print(f"  層別効果: {mild_treated - mild_control:.4f}")
print(f"\n重症群 (S=2):")
print(f"  処置群の回復率: {severe_treated:.4f}")
print(f"  対照群の回復率: {severe_control:.4f}")
print(f"  層別効果: {severe_treated - severe_control:.4f}")
print(f"\n真の処置効果: {treatment_effect:.2f}")

# 可視化
fig, axes = plt.subplots(1, 3, figsize=(15, 5))

# (a) 全体
ax = axes[0]
bars = ax.bar(["Control", "Treated"],
              [overall_control, overall_treated],
              color=["lightblue", "lightsalmon"], edgecolor="black")
ax.set_ylabel("Recovery rate", fontsize=12)
ax.set_title(f"(a) Overall\nNaive effect: {overall_treated - overall_control:.3f}",
             fontsize=11)
ax.set_ylim(0, 1)
for bar in bars:
    ax.text(bar.get_x() + bar.get_width()/2, bar.get_height() + 0.02,
            f"{bar.get_height():.3f}", ha="center", fontsize=11)
ax.grid(True, alpha=0.3, axis="y")

# (b) 軽症群
ax = axes[1]
bars = ax.bar(["Control", "Treated"],
              [mild_control, mild_treated],
              color=["lightblue", "lightsalmon"], edgecolor="black")
ax.set_ylabel("Recovery rate", fontsize=12)
ax.set_title(f"(b) Mild (S=1)\nEffect: {mild_treated - mild_control:.3f}",
             fontsize=11)
ax.set_ylim(0, 1)
for bar in bars:
    ax.text(bar.get_x() + bar.get_width()/2, bar.get_height() + 0.02,
            f"{bar.get_height():.3f}", ha="center", fontsize=11)
ax.grid(True, alpha=0.3, axis="y")

# (c) 重症群
ax = axes[2]
bars = ax.bar(["Control", "Treated"],
              [severe_control, severe_treated],
              color=["lightblue", "lightsalmon"], edgecolor="black")
ax.set_ylabel("Recovery rate", fontsize=12)
ax.set_title(f"(c) Severe (S=2)\nEffect: {severe_treated - severe_control:.3f}",
             fontsize=11)
ax.set_ylim(0, 1)
for bar in bars:
    ax.text(bar.get_x() + bar.get_width()/2, bar.get_height() + 0.02,
            f"{bar.get_height():.3f}", ha="center", fontsize=11)
ax.grid(True, alpha=0.3, axis="y")

plt.suptitle("Simpson's Paradox: confounding reverses the apparent effect",
             fontsize=13, y=1.02)
plt.tight_layout()
plt.savefig("causal_simpson_paradox.png", dpi=150, bbox_inches="tight")
plt.show()

このシミュレーション結果から、シンプソンのパラドックスの本質が確認できます。

1. 全体の結果（図a）: ナイーブに見ると、処置群の回復率が対照群よりも低いように見える可能性があります。これは処置が有害であるかのような印象を与えますが、それは誤りです

2. 層別の結果（図b, c）: 軽症群・重症群のそれぞれで見ると、処置群の回復率は対照群よりも高くなっています。これが真の因果効果です。両方の層で処置の効果は正（回復を促進する）であり、真の処置効果 $0.2$ にほぼ一致しています

3. パラドックスの原因: 重症患者は処置を受ける確率が高い（70%）一方、基本的な回復率が低い（30%）です。そのため、処置群には重症患者が多く含まれ、回復率が全体として引き下げられています。この「重症度」という交絡変数が全体の結果を歪めているのです

このパラドックスは、データの集約の仕方によって結論が正反対になりうることを示しています。正しい結論を導くためには、因果構造の理解が不可欠なのです。

因果推論の3つのアプローチ — 概観

ここまでで「相関と因果の違い」と「交絡の問題」を理解しました。では、因果推論はこの問題にどのようにアプローチするのでしょうか。大きく分けて3つの枠組みがあります。

1. ポテンシャルアウトカム（潜在結果）の枠組み

ドナルド・ルービン（Donald Rubin）によって体系化されたこの枠組みでは、各個体について「処置を受けた場合のアウトカム」$Y(1)$ と「処置を受けなかった場合のアウトカム」$Y(0)$ の2つの潜在的なアウトカム（potential outcomes）を考えます。個体レベルの因果効果は $Y(1) – Y(0)$ と定義されます。

しかし、各個体は処置を受けるか受けないかのどちらか一方しか経験できないため、$Y(1)$ と $Y(0)$ を同時に観測することは不可能です。これを因果推論の根本問題（fundamental problem of causal inference）と呼びます。ポテンシャルアウトカムの枠組みでは、この問題を統計的な識別条件（ignorability, unconfoundedness）のもとで集団レベルの因果効果（ATE）として推定します。

2. 構造因果モデル（SCM）とDAG

ジュディア・パール（Judea Pearl）が発展させた枠組みで、変数間の因果関係を有向非巡回グラフ（Directed Acyclic Graph, DAG）で表現します。DAG上での $d$-分離（d-separation）という概念を用いて、どの変数を条件付けすれば因果効果を識別できるかを判定します。

パールの $\text{do}$-計算（do-calculus）は、DAGの構造から介入確率 $P(Y \mid \text{do}(X))$ を観測確率から計算可能かどうかを判定し、可能であれば具体的な計算式を導く体系的な方法を提供します。

3. 自然実験と準実験デザイン

ランダム化が不可能な状況で、自然に生じた「実験のような状況」を利用して因果効果を推定する方法です。差の差法（Difference-in-Differences, DID）、操作変数法（Instrumental Variables, IV）、回帰不連続デザイン（Regression Discontinuity Design, RDD）などが代表的な手法であり、主に計量経済学の分野で発展してきました。

これら3つのアプローチは対立するものではなく、異なる視点から因果推論の問題に取り組む相補的な枠組みです。ポテンシャルアウトカムは「何を推定するか」を明確にし、DAGは「いつ推定可能か」を判定し、準実験デザインは「どのように推定するか」の具体的な戦略を提供します。

まとめ

本記事では、因果推論の最も基本的な問い——「相関と因果はどう違うのか」——を数理的に掘り下げました。

• 相関は因果を含意しない（correlation does not imply causation）。相関が生じるメカニズムには、直接因果、逆因果、交絡の3種類がある
• 交絡変数は処置とアウトカムの両方に影響する第3の変数であり、交絡を無視すると因果効果の推定にバイアスが生じる。交絡バイアスの大きさと方向は数学的に計算できる
• シンプソンのパラドックスは交絡が引き起こす劇的な現象で、全体の相関の方向と層別の相関の方向が逆転しうる。正しい解釈には因果構造の理解が不可欠
• 観察と介入の区別は因果推論の根幹。パールの $\text{do}$ 記法は条件付き確率 $P(Y \mid X)$ と介入確率 $P(Y \mid \text{do}(X))$ を明確に区別する
• Pythonシミュレーションにより、交絡バイアスの可視化、交絡の制御による因果効果の回復、シンプソンのパラドックスの再現を確認した

因果推論は「なぜ」を問う科学の根源的な営みに対して、数学的に厳密な枠組みを与える分野です。本記事は因果推論シリーズの出発点であり、この先の記事で各手法を詳しく解説していきます。

次のステップとして、以下の記事も参考にしてください。

• ポテンシャルアウトカムの枠組み — ルービン因果モデルの詳細
• ランダム化比較試験（RCT） — 因果効果を推定する「ゴールドスタンダード」
• 交絡バイアスとその制御法 — 交絡への対処法の体系的な解説
• DAGと因果グラフ — 因果関係をグラフで表現する方法