媒介分析 — 直接効果と間接効果の因果的分解

ある教育プログラムが学生のテスト成績を向上させたとします。しかし、その効果は「どのようにして」もたらされたのでしょうか。プログラムが学生の学習時間を増やし、その増加した学習時間が成績向上につながったのか。それとも、プログラムが学習方法を変えて直接的に学習効率を高めたのか。あるいは、その両方なのか。

因果効果の「大きさ」がわかっても、それだけでは十分ではない場合があります。政策立案者や研究者にとって、効果のメカニズムを理解することは極めて重要です。メカニズムがわかれば、より効果的な介入を設計したり、効果が別の文脈でも再現可能かどうかを予測したりできるからです。

媒介分析(mediation analysis)は、処置効果の「メカニズム」を解明するための手法です。処置 $T$ からアウトカム $Y$ への全体効果(total effect)を、媒介変数(mediator)$M$ を経由する間接効果(indirect effect)と、$M$ を経由しない直接効果(direct effect)に分解します。

媒介分析は以下の幅広い分野で活用されています。

  • 心理学: 認知行動療法の効果が「認知の歪みの修正」を通じて生じるのか、「行動パターンの変化」を通じるのかの検証
  • 疫学: 喫煙の肺がんへの影響が「タールの蓄積」を通じるのか、「免疫機能の低下」を通じるのかの分析
  • 政策評価: 教育プログラムの効果が「学習時間の増加」を介するのか、「学習方法の改善」を介するのかの解明
  • マーケティング: 広告の売上効果が「ブランド認知の向上」を介するのか、「価格感度の低下」を介するのかの分析
  • 社会科学: 人種や性別による賃金格差が「教育水準の違い」を介するのか、「直接的な差別」によるのかの分解

本記事の内容

  • 媒介分析の基本的な設定とDAGによる表現
  • 全体効果、直接効果、間接効果の因果的定義(ポテンシャルアウトカム)
  • Baron-Kenny法の古典的アプローチと限界
  • 因果媒介分析: 自然直接効果(NDE)と自然間接効果(NIE)
  • 識別条件と逐次的無視可能性仮定
  • $M$-$Y$ 交絡の問題と処置誘導交絡
  • 非線形モデルと交互作用がある場合の拡張
  • Pythonシミュレーションによる効果分解の実装

前提知識

この記事を読む前に、以下の知識があると理解が深まります。

媒介分析の設定

問題の構造

媒介分析は、処置 $T$ がアウトカム $Y$ に影響を与える経路(pathway)を解明する枠組みです。核心的なアイデアは、$T$ から $Y$ への効果を「$M$ を経由する部分」と「$M$ を経由しない部分」に分解することです。

具体的な例で考えましょう。新薬($T$)が血圧($Y$)を下げるメカニズムを調べるとします。新薬は体内の酵素の活性($M$)を変化させ、その酵素活性の変化が血圧に影響するかもしれません。同時に、新薬が酵素活性を経由せずに直接的に血管を拡張させる効果もあるかもしれません。

  • 全体効果: 新薬が血圧に与える総合的な影響
  • 間接効果: 新薬 → 酵素活性 → 血圧($M$ を経由する経路)
  • 直接効果: 新薬 → 血圧($M$ を経由しない経路)

DAGによる表現

媒介分析の基本的なDAG(有向非巡回グラフ)は次のように表されます。

$$ T \rightarrow M \rightarrow Y, \quad T \rightarrow Y $$

ここで: – $T$: 処置(treatment) – $M$: 媒介変数(mediator)— 処置とアウトカムの間に位置する変数 – $Y$: アウトカム(outcome)

$T \rightarrow M \rightarrow Y$ が間接パス(mediation path)、$T \rightarrow Y$ が直接パス(direct path)です。

実際の状況では、共変量 $\bm{X}$ が存在し、DAGはより複雑になります。

$$ \bm{X} \rightarrow T, \quad \bm{X} \rightarrow M, \quad \bm{X} \rightarrow Y, \quad T \rightarrow M \rightarrow Y, \quad T \rightarrow Y $$

共変量 $\bm{X}$ は $T$, $M$, $Y$ のいずれにも影響しうる交絡因子です。

直感的な理解: なぜ分解が重要か

効果の分解が重要である理由を、もう少し深く考えてみましょう。

ある運動プログラム($T$)が心臓病リスク($Y$)を低減する場合、その効果が「体重減少」($M_1$)を介するのか、「ストレス軽減」($M_2$)を介するのか、「直接的な心肺機能向上」を介するのかによって、次のアクションが変わります。

  • 体重減少が主な経路であれば、食事療法でも同様の効果が期待できる
  • ストレス軽減が主な経路であれば、ヨガや瞑想でも代替可能かもしれない
  • 直接的な心肺機能向上が主な経路であれば、運動自体が不可欠

このように、効果のメカニズムを理解することで、より効率的な介入の設計やターゲットの特定が可能になります。

では、これらの概念を数学的に厳密に定義しましょう。

効果の因果的定義 — ポテンシャルアウトカムの枠組み

ネスト化ポテンシャルアウトカム

媒介分析をポテンシャルアウトカムの枠組みで定式化するために、ネスト化(入れ子型)ポテンシャルアウトカム(nested potential outcomes)を導入します。

各個体 $i$ に対して、以下を定義します。

  • $M_i(t)$: 処置を $t$ に設定した場合の媒介変数のポテンシャルアウトカム
  • $Y_i(t, m)$: 処置を $t$ に設定し、媒介変数を $m$ に設定した場合のアウトカムのポテンシャルアウトカム
  • $Y_i(t) = Y_i(t, M_i(t))$: 処置を $t$ に設定した場合の(自然な)アウトカム

ここで重要なのは、$Y_i(t, m)$ では処置 $t$ と媒介変数 $m$ を独立に設定できると仮定している点です。実際には処置が媒介変数に影響するため、この「独立に設定する」操作は反事実的です。

具体例で理解しましょう。ある患者に対して: – $Y_i(1, M_i(1))$: 薬を投与し、媒介変数はその結果に自然に従う → 実際の処置下のアウトカム – $Y_i(0, M_i(0))$: 薬を投与せず、媒介変数はその結果に自然に従う → 実際の非処置下のアウトカム – $Y_i(1, M_i(0))$: 薬を投与するが、媒介変数はもし薬を投与しなかった場合の水準に固定 → 反事実的な量

3番目の量 $Y_i(1, M_i(0))$ が媒介分析の核心です。これは「処置は受けるが、媒介変数だけは処置を受けなかった場合の値になる」という、現実には観察不可能な状況を表しています。

全体効果(TE)

全体効果は、通常の平均処置効果と同じです。

$$ \begin{equation} \text{TE} = E[Y_i(1) – Y_i(0)] = E[Y_i(1, M_i(1)) – Y_i(0, M_i(0))] \end{equation} $$

処置を受けた場合と受けなかった場合のアウトカムの差の期待値です。

自然直接効果(NDE)

自然直接効果(Natural Direct Effect, NDE)は、媒介変数を処置なしの水準に固定したまま、処置を変えた場合の効果です。

$$ \begin{equation} \text{NDE} = E[Y_i(1, M_i(0)) – Y_i(0, M_i(0))] \end{equation} $$

NDEは「$M$ を経由しない直接パスの効果」を捉えています。媒介変数を $M_i(0)$(処置なしの場合の値)に固定しているので、処置が $M$ を通じて $Y$ に影響する経路は遮断され、直接パスの効果だけが残ります。

直感的には、「薬の酵素活性への影響を人工的にブロックした上で、薬を投与した場合の効果」に相当します。

自然間接効果(NIE)

自然間接効果(Natural Indirect Effect, NIE)は、処置を受けた状態で、媒介変数だけを処置なしの水準から処置ありの水準に変えた場合の効果です。

$$ \begin{equation} \text{NIE} = E[Y_i(1, M_i(1)) – Y_i(1, M_i(0))] \end{equation} $$

NIEは「$M$ を経由する間接パスの効果」を捉えています。処置は $t = 1$ に固定されているので、アウトカムの変化は純粋に $M$ の変化によるものです。

直感的には、「薬を投与した状態で、酵素活性だけが変化した場合の効果」に相当します。

分解の成立

全体効果はNDEとNIEの和に分解されます。

$$ \begin{align} \text{TE} &= E[Y_i(1, M_i(1)) – Y_i(0, M_i(0))] \\ &= E[Y_i(1, M_i(1)) – Y_i(1, M_i(0))] + E[Y_i(1, M_i(0)) – Y_i(0, M_i(0))] \\ &= \text{NIE} + \text{NDE} \end{align} $$

2行目では $Y_i(1, M_i(0))$ を足して引いています。この項は反事実的な量ですが、数学的な分解には問題ありません。

別の分解: NDE(1)とNIE(0)

上記の分解では $M$ を $M(0)$ に固定したNDE(NDE(0)と書くことがある)と $T$ を1に固定したNIE(NIE(1)と書くことがある)を使いましたが、逆の組み合わせも可能です。

$$ \text{NDE}(1) = E[Y_i(1, M_i(1)) – Y_i(0, M_i(1))] $$ $$ \text{NIE}(0) = E[Y_i(0, M_i(1)) – Y_i(0, M_i(0))] $$

この場合も $\text{TE} = \text{NDE}(1) + \text{NIE}(0)$ が成り立ちます。一般に $\text{NDE}(0) \neq \text{NDE}(1)$ かつ $\text{NIE}(0) \neq \text{NIE}(1)$ ですが、$T$ と $M$ の交互作用がない場合は一致します。

ここまでで因果的な定義を確立しました。次に、これらの量を観測データからどのように推定するかという識別の問題に進みます。

Baron-Kenny法 — 古典的アプローチ

方法の概要

因果媒介分析の理論が確立される以前から、心理学を中心に広く使われてきたのがBaron-Kenny法(Baron and Kenny, 1986)です。この方法は3つの回帰式を用いて間接効果を推定します。

ステップ1: 全体効果の推定

$$ Y = c_0 + cT + \epsilon_1 $$

$c$ は全体効果の推定量です。

ステップ2: $T$ から $M$ への効果の推定

$$ M = a_0 + aT + \epsilon_2 $$

$a$ は処置が媒介変数に与える効果です。

ステップ3: $M$ を制御した上での $T$ の効果の推定

$$ Y = c_0′ + c’T + bM + \epsilon_3 $$

$c’$ は「$M$ をコントロールした後の直接効果」、$b$ は「$T$ をコントロールした後の $M$ から $Y$ への効果」です。

Baron-Kenny法では: – 直接効果 = $c’$ – 間接効果 = $a \times b$($a$ と $b$ の積) – 分解: $c = c’ + ab$(全体効果 = 直接効果 + 間接効果)

Sobelの検定

間接効果 $ab$ の統計的有意性を検定する方法として、Sobelの検定が使われます。間接効果の標準誤差を次のように近似します。

$$ \begin{equation} \text{SE}_{ab} \approx \sqrt{a^2 \text{SE}_b^2 + b^2 \text{SE}_a^2} \end{equation} $$

$z$ 統計量 $z = ab / \text{SE}_{ab}$ を標準正規分布と比較します。ただし、$ab$ の分布は正規分布から大きく逸脱する場合があるため、ブートストラップ法による信頼区間の構成が推奨されます。

Baron-Kenny法の限界

Baron-Kenny法は直感的でわかりやすい反面、因果推論の観点から重要な限界があります。

限界1: 全ての関係が線形であることを仮定

3つの回帰式がすべて線形モデルであるため、$T$ と $M$ の間の非線形な関係や、$M$ と $Y$ の間の非線形な関係を捉えられません。また、$T$ と $M$ の交互作用($T$ によって $M$ の効果が変わる場合)も想定していません。

限界2: $M$-$Y$ の交絡が無視されている

ステップ3で $M$ を回帰に含めていますが、$M$ と $Y$ の間に交絡変数が存在する場合、$b$ の推定にバイアスが生じます。$M$ は処置後の変数(post-treatment variable)であるため、$M$ をコントロールすることは通常の回帰調整よりも複雑な問題を引き起こします。

限界3: 処置誘導交絡を扱えない

処置 $T$ が $M$-$Y$ 間の交絡変数 $L$ に影響する場合($T \rightarrow L \rightarrow Y$ かつ $L \rightarrow M$)、Baron-Kenny法では NDE と NIE を正しく識別できません。

これらの限界を克服するために発展したのが、次に説明する因果媒介分析の枠組みです。

因果媒介分析の識別条件

逐次的無視可能性仮定

Imai, Keele, and Tingley(2010)が定式化した逐次的無視可能性(sequential ignorability)仮定は、以下の2つの条件からなります。

仮定1(処置の無視可能性):

$$ \begin{equation} \{Y_i(t’, m), M_i(t)\} \perp\!\!\!\perp T_i \mid \bm{X}_i = \bm{x} \end{equation} $$

共変量 $\bm{X}$ で条件付けたとき、処置の割り当てがポテンシャルアウトカムおよびポテンシャル媒介変数と独立である。RCTでは自動的に満たされます。

仮定2(媒介変数の無視可能性):

$$ \begin{equation} Y_i(t’, m) \perp\!\!\!\perp M_i(t) \mid T_i = t, \bm{X}_i = \bm{x} \end{equation} $$

処置と共変量で条件付けたとき、媒介変数の値がポテンシャルアウトカムと独立である。

仮定2は仮定1よりもはるかに強い仮定です。仮定1はRCTで保証できますが、仮定2は処置がランダムに割り当てられていても自動的には保証されません。なぜなら、$M$ は処置後の変数であり、$M$ と $Y$ の関係に交絡がないことを要求しているからです。

処置誘導交絡の問題

仮定2が特に問題となるのは、処置誘導交絡(treatment-induced confounding)がある場合です。

DAGで表すと:

$$ T \rightarrow L \rightarrow M, \quad T \rightarrow L \rightarrow Y, \quad T \rightarrow M \rightarrow Y, \quad T \rightarrow Y $$

ここで $L$ は処置 $T$ の影響を受ける変数で、かつ $M$ と $Y$ の両方に影響します。例えば、教育プログラム($T$)がモチベーション($L$)を変え、モチベーションが学習時間($M$)と成績($Y$)の両方に影響する場合です。

$L$ で条件付けると $T$ から $M$ への因果パスの一部がブロックされ、$L$ で条件付けないと $M$-$Y$ 間の交絡が残ります。どちらにしても正しい推定ができないというジレンマが生じます。

処置誘導交絡は、NDEとNIEの非パラメトリック識別を不可能にする深刻な問題です。この場合の対処法としては、パラメトリック仮定を置いた感度分析や、フロントドア基準の活用などが研究されています。

識別のもとでの推定式

逐次的無視可能性仮定と正値性条件(全ての共変量の組み合わせで処置群・対照群の両方が存在し、媒介変数の条件付き分布にオーバーラップがある)のもとで、NDEとNIEは次のように識別されます。

$$ \begin{equation} \text{NDE} = \int \int E[Y | T=1, M=m, \bm{X}=\bm{x}] \cdot f(m | T=0, \bm{X}=\bm{x}) \, dm \, dF(\bm{x}) – E[Y | T=0] \end{equation} $$

$$ \begin{equation} \text{NIE} = \int \int E[Y | T=1, M=m, \bm{X}=\bm{x}] \cdot [f(m | T=1, \bm{X}=\bm{x}) – f(m | T=0, \bm{X}=\bm{x})] \, dm \, dF(\bm{x}) \end{equation} $$

NDEの式の直感は次の通りです。処置群のアウトカムモデル $E[Y | T=1, M=m, \bm{X}]$ に対して、媒介変数を対照群の分布 $f(m | T=0, \bm{X})$ で周辺化しています。これにより、「処置は受けるが、$M$ は処置を受けなかった場合の分布に従う」という反事実的状況を再現しています。

線形モデルの場合の簡略化

全ての関係が線形であり、$T$ と $M$ の交互作用がないと仮定すると、推定は大幅に簡略化されます。

$$ M = \alpha_0 + \alpha_1 T + \bm{\alpha}_X^T \bm{X} + \epsilon_M $$ $$ Y = \gamma_0 + \gamma_1 T + \gamma_2 M + \bm{\gamma}_X^T \bm{X} + \epsilon_Y $$

このモデルのもとで: – $\text{NDE} = \gamma_1$($Y$ の回帰における $T$ の係数) – $\text{NIE} = \alpha_1 \gamma_2$($T \to M$ の効果 $\times$ $M \to Y$ の効果)

これはBaron-Kenny法と同じ結果です。つまり、線形モデルで交互作用がなく、逐次的無視可能性が成り立つ場合、Baron-Kenny法は因果的に正しい推定を与えます。

しかし、交互作用がある場合($T$ と $M$ の相互作用が $Y$ に影響する場合)は、追加的な項が必要になります。

交互作用がある場合の拡張

$T \times M$ 交互作用モデル

実際の多くの状況では、処置と媒介変数の間に交互作用が存在します。例えば、薬の直接効果の大きさが酵素活性のレベルによって変わるような場合です。

交互作用を含むモデルは次のようになります。

$$ \begin{equation} Y = \gamma_0 + \gamma_1 T + \gamma_2 M + \gamma_3 T \times M + \bm{\gamma}_X^T \bm{X} + \epsilon_Y \end{equation} $$

このとき、NDEとNIEは次のように表されます。

$$ \begin{align} \text{NDE}(0) &= \gamma_1 + \gamma_3 E[M(0)] \\ \text{NIE}(1) &= (\gamma_2 + \gamma_3) \alpha_1 \\ \text{NDE}(1) &= \gamma_1 + \gamma_3 E[M(1)] \\ \text{NIE}(0) &= \gamma_2 \alpha_1 \end{align} $$

交互作用がある場合、NDEとNIEは $M$ を固定する水準($M(0)$ か $M(1)$ か)によって異なる値を取ります。

割合媒介(Proportion Mediated)

媒介分析でよく報告される指標に割合媒介(proportion mediated)があります。

$$ \begin{equation} \text{PM} = \frac{\text{NIE}}{\text{TE}} = \frac{\text{NIE}}{\text{NDE} + \text{NIE}} \end{equation} $$

$\text{PM} = 0.6$ であれば、全体効果の60%が媒介変数を経由する間接経路で説明されることを意味します。ただし、PMの信頼区間は広くなりがちであり、解釈には注意が必要です。

ここまでの理論を、Pythonシミュレーションで確認しましょう。

Pythonシミュレーション

シミュレーション1: 線形モデルでの媒介分析(Baron-Kenny法)

まず、全ての関係が線形で交絡がない理想的な状況で、Baron-Kenny法が正しく機能することを確認します。

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

np.random.seed(42)
n = 5000

# --- データ生成(線形、交絡なし)---
T = np.random.binomial(1, 0.5, n)  # 処置(ランダム化)
a_true = 1.5  # T -> M の効果
b_true = 2.0  # M -> Y の効果
c_true = 1.0  # T -> Y の直接効果

M = a_true * T + np.random.normal(0, 1, n)
Y = c_true * T + b_true * M + np.random.normal(0, 1, n)

# 真の効果
true_TE = c_true + a_true * b_true  # = 1.0 + 1.5*2.0 = 4.0
true_NDE = c_true                    # = 1.0
true_NIE = a_true * b_true           # = 3.0

# --- Baron-Kenny法 ---
# ステップ1: Y ~ T(全体効果)
reg1 = LinearRegression().fit(T.reshape(-1, 1), Y)
c_hat = reg1.coef_[0]

# ステップ2: M ~ T(T -> M)
reg2 = LinearRegression().fit(T.reshape(-1, 1), M)
a_hat = reg2.coef_[0]

# ステップ3: Y ~ T + M(直接効果 + M -> Y)
reg3 = LinearRegression().fit(np.column_stack([T, M]), Y)
c_prime_hat = reg3.coef_[0]
b_hat = reg3.coef_[1]

indirect_hat = a_hat * b_hat
direct_hat = c_prime_hat

print("=== 真の効果 ===")
print(f"全体効果 (TE):     {true_TE:.4f}")
print(f"自然直接効果 (NDE): {true_NDE:.4f}")
print(f"自然間接効果 (NIE): {true_NIE:.4f}")
print(f"\n=== Baron-Kenny法の推定 ===")
print(f"全体効果 (c):      {c_hat:.4f}")
print(f"直接効果 (c'):     {direct_hat:.4f}")
print(f"間接効果 (a*b):    {indirect_hat:.4f}")
print(f"分解の検証: c' + a*b = {direct_hat + indirect_hat:.4f}")
print(f"割合媒介: {indirect_hat / c_hat:.4f}")

推定結果から、線形モデルかつ交絡なしの設定では、Baron-Kenny法が真の効果を正確に分解できることが確認できます。全体効果4.0は直接効果1.0と間接効果3.0に正しく分解されています。割合媒介は約0.75であり、効果の75%が媒介変数を経由する間接経路で説明されています。

シミュレーション2: $M$-$Y$ 交絡がある場合のバイアス

次に、$M$ と $Y$ の間に未観測の交絡変数が存在する場合、Baron-Kenny法にバイアスが生じることを確認します。

import numpy as np
from sklearn.linear_model import LinearRegression

np.random.seed(0)
n = 10000

T = np.random.binomial(1, 0.5, n)  # ランダム化処置
U = np.random.normal(0, 1, n)       # M-Y間の未観測交絡

a_true = 1.0
b_true = 1.5
c_true = 2.0
delta_UM = 1.0  # U -> M の効果
delta_UY = 1.5  # U -> Y の効果

M = a_true * T + delta_UM * U + np.random.normal(0, 0.5, n)
Y = c_true * T + b_true * M + delta_UY * U + np.random.normal(0, 0.5, n)

true_NDE = c_true           # = 2.0
true_NIE = a_true * b_true  # = 1.5
true_TE = true_NDE + true_NIE  # = 3.5

# Baron-Kenny法(Uを無視)
reg_M = LinearRegression().fit(T.reshape(-1, 1), M)
a_hat = reg_M.coef_[0]

reg_Y = LinearRegression().fit(np.column_stack([T, M]), Y)
c_prime_hat = reg_Y.coef_[0]
b_hat = reg_Y.coef_[1]

indirect_biased = a_hat * b_hat
direct_biased = c_prime_hat

# 正しい推定(Uを含める)
reg_Y_correct = LinearRegression().fit(np.column_stack([T, M, U]), Y)
c_prime_correct = reg_Y_correct.coef_[0]
b_correct = reg_Y_correct.coef_[1]

indirect_correct = a_hat * b_correct
direct_correct = c_prime_correct

print("=== 真の効果 ===")
print(f"TE = {true_TE:.4f}, NDE = {true_NDE:.4f}, NIE = {true_NIE:.4f}")
print(f"\n=== Baron-Kenny法(Uを無視、バイアスあり)===")
print(f"直接効果: {direct_biased:.4f} (真値: {true_NDE:.4f})")
print(f"間接効果: {indirect_biased:.4f} (真値: {true_NIE:.4f})")
print(f"b (M->Y): {b_hat:.4f} (真値: {b_true:.4f})")
print(f"\n=== Uを制御した推定(バイアスなし)===")
print(f"直接効果: {direct_correct:.4f} (真値: {true_NDE:.4f})")
print(f"間接効果: {indirect_correct:.4f} (真値: {true_NIE:.4f})")
print(f"b (M->Y): {b_correct:.4f} (真値: {b_true:.4f})")

$M$-$Y$ 交絡がある場合の結果から、重要な観察ができます。

  1. $b$ の推定にバイアスが生じる: $U$ を無視すると、$M$ と $Y$ の関係に $U$ の効果が混入するため、$b$ の推定値が真の値から乖離します。$U$ が $M$ と $Y$ の両方に正の影響を与えているため、$b$ は過大推定されます。
  2. 間接効果が過大推定される: $b$ の過大推定は間接効果 $ab$ の過大推定を引き起こします。
  3. 直接効果が過小推定される: 全体効果は正しく推定されるため($T$ はランダム化されているため)、間接効果の過大推定は直接効果の過小推定として現れます。
  4. $U$ を制御すれば正しい推定が得られる: 交絡変数 $U$ を回帰に含めることで、バイアスが消えます。

この結果は、処置のランダム化だけでは媒介分析のバイアスを解消できないという重要な教訓を示しています。$T$ をランダム化しても、$M$-$Y$ 間の交絡は残るのです。

シミュレーション3: ブートストラップによる間接効果の信頼区間

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

np.random.seed(42)
n = 200  # やや小さなサンプルサイズ
n_bootstrap = 5000

T = np.random.binomial(1, 0.5, n)
a_true, b_true, c_true = 0.8, 1.2, 1.5

M = a_true * T + np.random.normal(0, 1, n)
Y = c_true * T + b_true * M + np.random.normal(0, 1, n)

# ブートストラップ
boot_indirect = []
boot_direct = []
boot_total = []

for _ in range(n_bootstrap):
    idx = np.random.choice(n, n, replace=True)
    T_b, M_b, Y_b = T[idx], M[idx], Y[idx]

    reg_m = LinearRegression().fit(T_b.reshape(-1, 1), M_b)
    a_b = reg_m.coef_[0]

    reg_y = LinearRegression().fit(np.column_stack([T_b, M_b]), Y_b)
    c_b = reg_y.coef_[0]
    b_b = reg_y.coef_[1]

    boot_indirect.append(a_b * b_b)
    boot_direct.append(c_b)
    boot_total.append(a_b * b_b + c_b)

boot_indirect = np.array(boot_indirect)
boot_direct = np.array(boot_direct)
boot_total = np.array(boot_total)

# 信頼区間(パーセンタイル法)
ci_indirect = np.percentile(boot_indirect, [2.5, 97.5])
ci_direct = np.percentile(boot_direct, [2.5, 97.5])
ci_total = np.percentile(boot_total, [2.5, 97.5])

print("=== ブートストラップ信頼区間 (95%) ===")
print(f"間接効果: {np.mean(boot_indirect):.4f} [{ci_indirect[0]:.4f}, {ci_indirect[1]:.4f}]"
      f" (真値: {a_true * b_true:.4f})")
print(f"直接効果: {np.mean(boot_direct):.4f} [{ci_direct[0]:.4f}, {ci_direct[1]:.4f}]"
      f" (真値: {c_true:.4f})")
print(f"全体効果: {np.mean(boot_total):.4f} [{ci_total[0]:.4f}, {ci_total[1]:.4f}]"
      f" (真値: {c_true + a_true * b_true:.4f})")

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

for ax, data, label, true_val, color in zip(
    axes,
    [boot_indirect, boot_direct, boot_total],
    ['Indirect Effect (NIE)', 'Direct Effect (NDE)', 'Total Effect (TE)'],
    [a_true * b_true, c_true, c_true + a_true * b_true],
    ['steelblue', 'darkorange', 'forestgreen']
):
    ax.hist(data, bins=60, density=True, alpha=0.6, color=color, edgecolor='white')
    ax.axvline(true_val, color='red', linestyle='--', linewidth=2,
               label=f'True = {true_val:.2f}')
    ci = np.percentile(data, [2.5, 97.5])
    ax.axvline(ci[0], color='black', linestyle=':', linewidth=1.5)
    ax.axvline(ci[1], color='black', linestyle=':', linewidth=1.5,
               label=f'95% CI: [{ci[0]:.2f}, {ci[1]:.2f}]')
    ax.set_xlabel('Effect size', fontsize=11)
    ax.set_ylabel('Density', fontsize=11)
    ax.set_title(label, fontsize=12)
    ax.legend(fontsize=9)
    ax.grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig('mediation_bootstrap.png', dpi=150, bbox_inches='tight')
plt.show()

3つのパネルは、間接効果、直接効果、全体効果のブートストラップ分布を示しています。赤い破線が真の値、黒い点線が95%信頼区間の境界です。

この結果から以下のことがわかります。

  1. 間接効果の分布はやや右に歪んでいる: $ab$ は2つの推定量の積であるため、正規分布からの逸脱が見られます。これがSobel検定(正規近似)よりもブートストラップ法が推奨される理由です。
  2. 信頼区間は真の値を含んでいる: $n = 200$ の比較的小さなサンプルサイズでも、95%信頼区間は真の値をカバーしています。
  3. 間接効果の信頼区間は直接効果の信頼区間よりも幅が広い: 間接効果は2つのパラメータの積であるため、不確実性がより大きくなります。

シミュレーション4: 交互作用がある場合

import numpy as np
from sklearn.linear_model import LinearRegression

np.random.seed(123)
n = 10000

T = np.random.binomial(1, 0.5, n)

# T x M 交互作用あり
a_true = 1.0
gamma1 = 2.0  # T -> Y 直接
gamma2 = 1.5  # M -> Y
gamma3 = 0.8  # T * M -> Y 交互作用

M = a_true * T + np.random.normal(0, 1, n)
Y = gamma1 * T + gamma2 * M + gamma3 * T * M + np.random.normal(0, 1, n)

# 真のNDE, NIE
E_M0 = 0.0   # E[M(0)] = a_true * 0 = 0
E_M1 = a_true  # E[M(1)] = a_true * 1 = 1

true_NDE_0 = gamma1 + gamma3 * E_M0  # 2.0 + 0.8*0 = 2.0
true_NDE_1 = gamma1 + gamma3 * E_M1  # 2.0 + 0.8*1 = 2.8
true_NIE_0 = gamma2 * a_true          # 1.5 * 1.0 = 1.5
true_NIE_1 = (gamma2 + gamma3) * a_true  # (1.5+0.8) * 1.0 = 2.3
true_TE = gamma1 + (gamma2 + gamma3) * a_true  # 2.0 + 2.3 = 4.3

# Baron-Kenny法(交互作用を無視)
reg_M = LinearRegression().fit(T.reshape(-1, 1), M)
a_hat = reg_M.coef_[0]

reg_Y_no_int = LinearRegression().fit(np.column_stack([T, M]), Y)
c_prime_no_int = reg_Y_no_int.coef_[0]
b_no_int = reg_Y_no_int.coef_[1]

# 交互作用を含む推定
TM = T * M
reg_Y_int = LinearRegression().fit(np.column_stack([T, M, TM]), Y)
gamma1_hat = reg_Y_int.coef_[0]
gamma2_hat = reg_Y_int.coef_[1]
gamma3_hat = reg_Y_int.coef_[2]

print("=== 真の効果 ===")
print(f"TE = {true_TE:.4f}")
print(f"NDE(0) = {true_NDE_0:.4f}, NDE(1) = {true_NDE_1:.4f}")
print(f"NIE(0) = {true_NIE_0:.4f}, NIE(1) = {true_NIE_1:.4f}")
print(f"\n=== Baron-Kenny法(交互作用無視)===")
print(f"直接効果: {c_prime_no_int:.4f}")
print(f"間接効果: {a_hat * b_no_int:.4f}")
print(f"\n=== 交互作用を含む推定 ===")
print(f"gamma1 (T): {gamma1_hat:.4f} (真値: {gamma1})")
print(f"gamma2 (M): {gamma2_hat:.4f} (真値: {gamma2})")
print(f"gamma3 (TxM): {gamma3_hat:.4f} (真値: {gamma3})")
print(f"NDE(0) = gamma1 + gamma3 * E[M(0)] = {gamma1_hat + gamma3_hat * 0:.4f}")
print(f"NIE(1) = (gamma2 + gamma3) * a = {(gamma2_hat + gamma3_hat) * a_hat:.4f}")

交互作用がある場合の結果から重要な点が確認できます。

  1. Baron-Kenny法(交互作用を無視)は、直接効果と間接効果の推定にバイアスを生じる: 交互作用項を無視すると、$M$ の係数 $b$ に $T$ の影響が混入します。
  2. 交互作用を含むモデルでは、NDEとNIEが固定するレベルによって異なる値を取る: $\text{NDE}(0) = 2.0$ と $\text{NDE}(1) = 2.8$ は異なり、$\text{NIE}(0) = 1.5$ と $\text{NIE}(1) = 2.3$ も異なります。
  3. 交互作用は「処置が媒介変数の効果を修飾する」ことを意味する: $\gamma_3 > 0$ は、処置を受けている人($T=1$)では媒介変数の効果がより大きくなることを示しています。

まとめ

本記事では、媒介分析の理論を因果推論の枠組みから解説しました。

  • 媒介分析は全体効果を、媒介変数 $M$ を経由する間接効果(NIE)と経由しない直接効果(NDE)に分解する手法
  • ネスト化ポテンシャルアウトカム $Y(t, m)$ を用いてNDEとNIEを因果的に定義する
  • Baron-Kenny法は直感的だが、線形性、交絡なし、交互作用なしの仮定に依存する
  • 逐次的無視可能性がNDEとNIEの識別に必要。特に $M$-$Y$ 間の交絡の制御が重要
  • 処置誘導交絡がある場合、非パラメトリック識別が困難になる
  • 処置のランダム化だけでは不十分: $T$ をランダム化しても $M$-$Y$ 交絡は残りうる
  • 交互作用がある場合: NDEとNIEは固定するレベルに依存し、Baron-Kenny法ではバイアスが生じる

媒介分析は因果推論の中でも特に仮定が強く、結果の解釈には慎重さが求められます。識別条件(特に逐次的無視可能性)が成り立つかどうかを十分に検討し、感度分析を併用することが推奨されます。

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