Pythonで始める因果推論入門

昨今、ビジネスでは、データに基づいた意思決定がますます求められています。しかし、データから得られるのは、あくまで相関関係です。因果関係を正しく理解することで、我々はより確実な予測を行い、効果的な対策を立てることができるようになります。

例えば、新しい薬の開発において、因果推論は薬の効果を正確に評価し、副作用を最小限に抑えるために不可欠です。また、経済政策の立案においても、政策の効果を事前に予測し、より良い社会を実現するために因果推論が活用されています。

今回は、因果推論の基礎からPython実装までを、初学者向けにできるだけ分かりやすく解説していきます。

相関関係 vs. 因果関係

相関関係は、2つの変数間にどのような関係性があるのかを把握するための第一歩です。しかし、相関関係だけでは、その関係性が因果関係によるものなのか、それとも擬似相関によるものなのかを判断することはできません。因果関係を明らかにするためには、より厳密な統計手法や実験設計が必要となります。

相関関係と因果関係の違い

データ分析において、2つの変数間の関係性を調べることは非常に重要です。この際、しばしば「相関関係」と「因果関係」という二つの概念が用いられます。

  • 相関関係: 2つの変数が同時に変化する度合いを表します。つまり、一方の変数が変化すると、もう一方の変数も変化する傾向があることを示します。
  • 因果関係: 1つの変数が変化することで、もう一つの変数が変化するという、原因と結果の関係を表します。

2変数間の相関関係を示す尺度として、一般的に統計学では相関係数と呼ばれるものが使用されます。

相関係数\(\rho_{xy}\)の具体的な定義は次の通りです。

$$
\rho_{xy} = \frac{\sigma_{xy}}{\sigma_x \sigma_y} = \frac{\frac{1}{N} \sum_{i=1}^N (x_i – \bar{x})(y_i – \bar{y})}{\sqrt{\frac{1}{N} \sum_{i=1}^N (x_i – \bar{x})^2}\sqrt{\frac{1}{N} \sum_{i=1}^N (y_i – \bar{y})^2}}
$$

ここで、\(\sigma_x\)と\(\sigma_y\)は、それぞれ変数\(x\)と\(y\)の標準偏差、\(\sigma_{xy}\)は\(x\)と\(y\)の共分散を示しています。共分散\(\sigma_{xy}\)は、\(x\)と\(y\)のそれぞれの偏差の積の平均値を表し、2つの変数の変動の揃い具合を示す指標です。また、\(\bar{x}\)と\(\bar{y}\)は、それぞれ変数\(x\)と\(y\)の平均値、\(N\)はデータの個数を表します。

この相関係数は、\(-1.0 \leq \rho_{xy} \leq 1.0\)の範囲の値をとり、\(0.0\)より大きいとき正の相関があるといい、\(0.0\)より小さいとき負の相関があるといいます。

相関係数の絶対値が\(0.0 \leq |\rho_{xy}| \leq 0.2\)のとき「ほとんど相関関係がない」、\(0.2 \le |\rho_{xy}| \leq 0.4\)のとき「やや相関関係がある」、\(0.4 \le |\rho_{xy}| \leq 0.7\)のとき「かなり相関関係がある」、\(0.7 \le |\rho_{xy}| \leq 1.0\)のとき「強い相関関係がある」などと表現されます。

なぜ相関関係だけでは不十分なのか?

相関関係があるからといって、必ずしも因果関係があるとは限りません。このことを理解するために、擬似相関という概念を導入します。

擬似相関とは、2つの変数間に強い相関関係が見られるにも関わらず、実際には因果関係が存在しない、あるいは因果関係が逆であったり、第3の変数が影響しているような状況を指します。

擬似相関の種類

擬似相関には、主に以下の3つの種類があります。

  • 因果の関係が逆: 2つの変数の因果関係が逆になっている場合です。例えば、気温上昇がアイスクリーム売上の増加をもたらしますが、その逆は成り立ちません。
  • 交絡因子: 2つの変数に同時に影響を与える第3の変数が存在し、それによって誤った因果関係を推測してしまう可能性がある状況です。例えば、気温上昇がプールの利用増加とアイスクリーム売上の増加の両方に影響を与えているため、これらの変数間に直接的な因果関係があるように誤解されることがあります。
  • 合流点バイアス: 2つの変数が共通の原因を持つため、一方の変数の変化がもう一方の変数の変化をもたらすように誤って解釈される可能性がある状況です。例えば、数学と英語の成績が大学の合否に関係している場合、数学の成績が高いことだけが直接的に合格率を上げる要因であると結論づけることはできません。
図1. 疑似相関の種類

因果推論の考え方

因果推論は、ある事象が別の事象に与える因果関係を定量的に評価するための統計学的な手法です。相関関係と異なり、因果関係は「ある事象が変化したときに、別の事象がどのように変化するか」という、より深い洞察を提供します。

反実仮想

因果推論の基礎となる概念の一つが反実仮想(Counterfactual)です。これは、ドナルド・ルービン(1972)によって提唱された考え方で、ある事象に対して「もし、こうだったらどうなっていただろうか」という、現実には起こらなかった状況を仮想的に考えることを指します。

例えば、ある薬を飲んだ場合と飲まなかった場合、どちらの状態がその人にどのような影響を与えたかを考える際、実際に薬を飲んだ状態と飲まなかった状態は同時に観測することはできません。しかし、因果効果を評価するためには、両方の状態を比較する必要があります。このとき、実際に観測されなかった状態を「反実仮想」と呼びます。

より具体的に説明します。サンプルデータ\(i\)について、ある処置を受けたとき(\(z=1\))の結果を\(Y_1^i\)​、受けなかったとき(\(z=0\))の結果を\(Y_0^i\)​とします。この\(Y_1^i\)​と\(Y_0^i\)​は、個体\(i\)にとっての潜在的結果変数と呼ばれます。潜在的結果変数は、現実に存在し得る結果ですが、「\(Y_1^i\)​と\(Y_0^i\)は同時に存在し得ない」という特徴があります。なぜなら、個体\(i\)はどちらか一方の処置しか受けられないためです。

実際に観測できた変数(例えば、薬を飲んだ後の症状の変化など)を結果変数と呼び、潜在的結果変数との対比で理解することができます。

反実仮想の概念は、因果効果を評価する上で非常に重要です。因果効果とは、ある処置を受けた場合と受けなかった場合の潜在的結果変数の差として定義されるためです。

反実仮想の例

  • ある新薬の効き目を評価する場合:新薬を飲んだ場合の症状の変化(\(Y_1\)​)と、プラセボ薬を飲んだ場合の症状の変化(\(Y_0\)​)を比較する。
  • 新しいマーケティング戦略の効果を評価する場合:新しい戦略を採用した場合の売上(\(Y_1\))と、従来の戦略を続けた場合の売上(\(Y_0\)​)を比較する。

反実仮想を扱う手法

  • DiCE (Microsoft):反事実的な質問に対して、機械学習モデルを用いて回答するツール。
  • What-if分析: ある変数を変化させた場合に、他の変数がどのように変化するかをシミュレーションする手法。

因果推論では最終的に何を計算したいのか?

因果推論では、ある処置が結果変数に与える影響を定量的に評価することを目指します。この際、様々な種類の因果効果が考えられます。

因果効果の種類定義説明
ITE: Individual Treatment Effect
(個体レベルの因果効果)
\(ITE = Y_1^i – Y_0^i\)個体\(i\)において、処置を受けた場合の結果\(Y_1^i\)と、受けなかった場合結果\(Y_0^i\)の差。個々の単位に焦点を当てた因果効果です。しかし、すべての個体のITEを正確に観測することは困難です。
ATE: Average Treatment Effect
(集団レベルの因果効果)
\(ATE = \frac{1}{N} \sum_{i=1}^N (Y_1^i – Y_0^i) = \frac{1}{N} \sum_{i=1}^N Y_1^i – \frac{1}{N} \sum_{i=1}^N Y_0^i\)全体の集団における平均的な因果効果です。多くの場合、このATEの推定が目標となります。
ATT: Average Treatment effect on the Treated
(処置群における平均処置効果)
\(ATT = E(Y_1 – Y_0 | Z = 1)\)処置を受けた集団に焦点を当てた因果効果です。例えば、ある薬の効果を評価する際、実際にその薬を飲んだ人々に対する平均的な効果を知りたい場合に用います。
ATU: Average Treatment effect on the Untreated
(対照群における平均処置効果)
\(ATU = E(Y_1 – Y_0 | Z = 0)\)処置を受けていない集団に焦点を当てた因果効果です。例えば、ある政策の効果を評価する際、その政策を受けていない人々にもしその政策を受けていたらどの程度の効果があったのかを知りたい場合に用います。
CATE: Conditional Average Treatment Effect
(条件付き平均処置効果)
\(CATE(x) = E(Y_1 – Y_0 | x)\)特定のサブグループにおける因果効果です。例えば、年齢や性別などの背景変数によって、処置の効果が異なる場合があるため、これらの変数を条件とした因果効果を推定したい場合に用います。

介入とdoオペレーター

因果推論において、介入とは、ある変数の値を人為的に固定し、その結果他の変数がどのように変化するかを調べる操作を指します。これは、反実仮想的な思考実験を形式的に表現するための重要な概念です。

介入の概念は、ジューディア・パール(1995)によって提唱されたdoオペレータを用いて数式的に表現されます。例えば、変数\(Z\)の値を\(1\)に固定する介入は、\(do(Z=1)\)と表記されます。これは、「\(Z\)を\(1\)に固定した場合」という意味を表します。

図2. 介入とdoオペレーター

介入と処置群における平均処置効果(ATT)の違い

介入と、これまで説明した処置群における平均処置効果(ATT)は、一見似ているように思えますが、重要な違いがあります。

  • ATT: 実際に処置を受けた集団(\(Z=1\))における平均的な因果効果を表します。これは、観察データから直接計算できる量です。
  • 介入: ある変数の値を人為的に固定した状況下での因果効果を表します。これは、反実仮想的な概念であり、直接観測することはできません。

なぜ、ATTと介入が異なるのか?

これは、交絡因子の存在が原因です。例えば、ある薬の効果を評価する場合、薬を飲んだ人(\(Z=1\))と飲まなかった人(\(Z=0\))を比較する際、健康状態や年齢など、薬の効果に影響を与える他の変数(交絡因子)が両群で異なっている可能性があります。この場合、ATTは、薬の効果だけでなく、これらの交絡因子の影響も含まれてしまうため、純粋な薬の効果を評価できないことがあります。

一方、介入\(do(Z=1)\)は、これらの交絡因子の影響を排除した状態での薬の効果を表します。つまり、すべての個体に強制的に薬を飲ませた場合、どのような効果が得られるかを表します。

介入と平均処置効果(ATE)の関係

一方、介入と集団全体の平均処置効果(ATE)は、以下のような関係があります。

$$
ATE = E(Y_1) – E(Y_0) = E(Y_1 | do(Z=1)) – E(Y_0 | do(Z=0))
$$

これは、ATEが、介入によって得られる平均的な因果効果と等しいことを意味します。つまり、集団全体で見た場合、介入による効果とATEは一致します。

介入の重要性

  • 因果関係の特定: 介入を用いることで、複数の変数間の因果関係をより明確に特定することができます。
  • 交絡因子の調整: 介入は、交絡因子の影響を排除し、純粋な因果効果を推定するための強力なツールです。
  • 政策評価: ある政策の効果を評価する際、介入の概念は非常に有用です。例えば、ある政策を導入した場合と導入しなかった場合を比較することで、その政策の効果を定量的に評価することができます。

介入操作を表すdoオペレーターは、通常の数学の演算子とは異なるため、一般的な計算式で直接扱うことができません。そこで、この介入操作を数学的に表現するために調整化公式というものが導入されています。この公式は、ある変数に介入した際の他の変数の値の確率(期待値)を計算するためのもので、具体的には、以下の式で表されます。

$$
E(Y=y|do(Z=z)) = \sum_x E(Y=y | Z = z, X = z) E(X=x)
$$

構造方程式モデリングと因果ダイアグラム

構造方程式モデリングと因果ダイアグラムは、複雑な因果関係を分析するための強力なツールです。SEMは、因果関係を数式で表現し、因果ダイアグラムは、因果関係を視覚的に表現します。これら2つの手法を組み合わせることで、より深く因果関係を理解することができます。

構造方程式モデリング

構造方程式モデリング(Structural Equation Modeling, SEM)は、複数の変数間の因果関係を数式を用いて表現し、その関係性を推定する統計モデリング手法です。SEMでは、観測変数だけでなく、直接観測できない潜在変数も扱うことができ、複雑な因果関係をモデル化することができます。

例えば、3つの変数\((x, y, z)\)の関係を以下の式で表すことができます。

$$
x = f_x (y, z, e_x) = 0.0 \cdot y + 0.0 \cdot z + e_x \\
y = f_y (x, z, e_y) = 2.0 \cdot x + 0.0 \cdot z + e_y \\
z = f_z (x, y, e_z) = 1.0 \cdot x + 0.5 \cdot y + e_z
$$

ここで、\((e_x, e_y, e_z)\) は誤差変数を表しています。また、この例のように、変数の1次成分のみで表現されるモデルを線形構造方程式モデルと呼びます。線形構造方程式モデルは、解析的な扱いやすさから、多くの場合で利用されてます。

因果ダイアグラム

因果ダイアグラムは、線形構造方程式モデルにおける変数間の因果関係を視覚的に表現するための図です。各変数をノードで表し、因果関係を矢印(エッジ)で結ぶことで、複雑な因果構造を簡潔に表現することができます。

一般に、因果推論で扱う因果ダイアグラムは有向非循環グラフ(Directed Acyclic Graph, DAG)となります。つまり、循環的な因果関係は存在せず、ある変数から別の変数へ一方通行で因果関係が流れるケースを扱います。また、自己フィードバックが他の因果よりも十分に小さく、無視することができるケースを扱います(間接的な自己フィードバックを含むものも原則として取り扱いません)。

図3. 因果ダイアグラムの例

d分離

d分離は、因果ダイアグラム上で変数間の独立性を確認するための重要な手法です。この手法を用いることで、複雑な因果構造の中から、本当に知りたい因果効果を抽出することができます。

バックドアパスとバックドア基準

因果ダイアグラムでは、変数間の影響は直接的なものだけでなく、間接的な経路を通じて生じることがあります。例えば、「\(z \to y\)」という直接的な因果の他に、「\(z \to a \to y\)」という間接的な影響が存在する場合があります。このような間接的な経路をバックドアパスと呼び、途中にある変数 \(a\) は中間変数と呼ばれます。

さらに、因果推論において注意すべき交絡因子が存在する場合、これを調整することで因果関係の推定が可能になります。この交絡因子を調整する条件をバックドア基準と呼びます。

因果推論では、以下のような操作を行うことが重要です:

  • 不要なバックドアパスが新たに生じないようにする。
  • 既存のバックドアパスを適切に閉じる。

加えて、介入(doオペレータ)を用いて交絡因子の影響を遮断することで、正確な因果関係の推定が可能になります。

図4. バックドアパスとバックドア基準

d分離の手順

因果ダイアグラムにおける正しい因果推論を行うためには、変数間の不要な影響を取り除く必要があります。この目的で行われる操作がd分離(”d”は”directional”の略)です。d分離は、因果ダイアグラムを適切に整理し、因果推論を容易にする方法です。

以下は、d分離を行う際の基本的な手順です:

  1. 中間変数を無視する: バックドアパスを閉じるために、不要な中間変数を調整または無視します。
  2. 注目している変数より下流にある変数を無視する: 因果の影響を受ける側の変数(下流の変数)は、注目している変数の因果推論には直接関係しないため無視します。
  3. 注目している変数より上流の合流点にある変数を無視する: 因果の発生源となる上流の変数のうち、合流点で影響を分離できるものは無視します。
  4. バックドア基準を満たす交絡因子を介入(doオペレータ)で遮断する: 必要な交絡因子を調整することで、交絡の影響を排除します。
図5. d分離の例

因果推論の各種手法

ここからは、Pythonを用いた因果推論の手法を、具体的な実装例と合わせて紹介します。因果推論の実行には、DoWhyEconMLという2つのパッケージを使用します。これらのパッケージに加えて、因果グラフを視覚化する際に便利なPyGraphvizもインストールしましょう。

なお、PyGraphvizのインストールには、事前に以下のソフトウェアをインストールしておく必要があります。

上記をインストール後、Pythonの仮想環境などで以下のコマンドを実行し、必要なパッケージをインストールしてください。

# DoWhy と EconML をインストール
$ pip install dowhy econml

# PyGraphviz のインストール
# 詳細は https://pygraphviz.github.io/documentation/stable/install.html を参照のこと
$ pip install --config-settings="--global-option=build_ext" --config-settings="--global-option=-IC:\Program Files\Graphviz\include" --config-settings="--global-option=-LC:\Program Files\Graphviz\lib" pygraphviz

回帰分析による因果推論

回帰分析は、観測データをもとに因果効果を推定する基本的な手法の1つです。この手法は、英語では “Regression Adjustment” と呼ばれ、平均処置効果(ATE: Average Treatment Effect)を推定する際によく用いられます。

d分離を用いた変数選択

回帰分析を因果推論に適用する際、重要なのは使用する変数の選択です。観測・取得した全データをそのまま使ってしまうと、単なる機械学習にとどまり、因果推論としての意味が薄れてしまいます。ここで役立つのがd分離です。

因果ダイアグラムを基に、d分離の考え方を用いて、介入変数 \(Z\) と目的変数 \(Y\) の因果関係を特定するために必要な説明変数 \(X_i\) を選びます。これにより、交絡因子の影響を取り除きつつ、適切なモデルを構築できます。

回帰モデルの設定

回帰分析を因果推論に用いる場合、以下のような回帰モデルを構築します。

$$
Y = \beta + \alpha_Z Z + \sum_{i=1}^N \alpha_i X_i + \phi_{\text{noise}}
$$

ここで、

  • \(Y\):目的変数(アウトカム)
  • \(Z\):介入変数(処置を受けたかどうか)
  • \(X_i\):説明変数(d分離で選ばれた変数など)
  • \(\alpha_Z\):介入効果(因果効果)
  • \(\alpha_i\):説明変数の回帰係数
  • \(\phi_{\text{noise}}\):ノイズ項

回帰モデルを学習させることで、介入効果 \(\alpha_Z\) を推定できます。この \(\alpha_Z\) が因果推論における平均処置効果(ATE)として解釈されます。

なお、線形回帰モデルの場合、処置群における平均処置効果ATT対照群における平均処置効果ATUは、平均処置効果ATEに一致します。

回帰分析による因果推論のサンプルコードは、ここをクリックして展開して確認できます。

回帰モデルによる因果効果の推定には、DoWhyと呼ばれるパッケージが便利です。まずは、必要なパッケージの読み込みと、因果推論で使用するデータセットを準備します。今回はダミーのデータセットを用いて、効果推定を行います。なお、平均処置効果が10になるようにデータセットを準備します(dowhy.datasets.linear_datasetメソッドの引数betaで因果効果の程度を設定できる)。

import dowhy.datasets
import pandas as pd
from dowhy import CausalModel

pd.set_option("display.unicode.east_asian_width", True)

# 生成されたデータを確認
# カラムの内容: Y=結果(目的変数), v0=介入(処置)、Wi=説明変数 (i=0, 1, ...)
data = dowhy.datasets.linear_dataset(
    beta=10,                # 真の因果効果
    num_common_causes=3,    # 観測された共通原因
    num_samples=3000,       # サンプル数
)

print(data["df"])
            W0        W1        W2     v0          y
0    -0.767973 -0.507762  0.279605  False  -5.809510
1     0.884462 -0.764673 -1.449482  False  -0.367094
2    -1.204831 -0.653101 -0.287883  False  -8.979550
3    -1.353726  0.230209 -2.500752  False  -7.264057
4    -0.215305 -1.454302 -0.072621  False  -7.796488
...        ...       ...       ...    ...        ...
2995  1.019666  1.282668 -0.278459   True  20.572007
2996  0.737367  0.575113  1.046534   True  16.955961
2997  1.586630 -2.970984  1.631861  False  -4.938772
2998 -0.534554 -0.056926 -1.571492  False  -3.975299
2999  2.291931 -0.199890 -0.311141   True  19.773264

[3000 rows x 5 columns]

次に、CausalModelクラスを利用して因果モデルを作成し、作成したデータについて因果ダイアグラムを確認してみましょう。

# 因果モデルを作成
model = CausalModel(
    data=data["df"],
    treatment=data["treatment_name"],   # 介入があったデータの変数名
    outcome=data["outcome_name"],       # 出力変数名
    graph=data["gml_graph"]             # グラフ情報
)

# 因果ダイアグラムを表示
model.view_model()

DoWhyを利用して作成した因果モデルは、identify_effectメソッドを使用することで簡単にd分離を行うことができます。

# d分離
estimand = model.identify_effect(proceed_when_unidentifiable=True)
print(estimand)
Estimand type: EstimandType.NONPARAMETRIC_ATE

### Estimand : 1
Estimand name: backdoor
Estimand expression:
  d
─────(E[y|W1,W0,W2])
d[v₀]
Estimand assumption 1, Unconfoundedness: If U→{v0} and U→y then P(y|v0,W1,W0,W2,U) = P(y|v0,W1,W0,W2)

### Estimand : 2
Estimand name: iv
No such variable(s) found!

### Estimand : 3
Estimand name: frontdoor
No such variable(s) found!

では、d分離したデータに対して回帰モデルによる因果推論をやってみます。次のようにestimate_effectメソッドを使用することで因果推論を実行できますが、引数のmethod_nameで利用するアルゴリズムを選択できます。

estimate = model.estimate_effect(
    identified_estimand=estimand,
    method_name="backdoor.linear_regression",   # アルゴリズムに回帰分析を指定
    target_units="ate"                          # 平均処置効果(ATE)を推定
)
print(f"ATE: {estimate.value}")
ATE: 10.000205514811354

結果を確認すると、推定された因果効果(ATE)は、10.000205514811354となりました。ダミーデータの準備の際に、因果効果を10と設定しているので、うまく評価できていることがわかると思います。

傾向スコアを用いた逆確率重み付け法

傾向スコアを用いた逆確率重み付け法(Inverse Probability Weighting, IPW)は、「処置を受ける傾向」を示す傾向スコア(Propensity Score)を活用して因果効果を推定する手法です。この方法は、処置群と対照群における交絡の影響を調整し、平均処置効果(ATE)を計算します。

傾向スコアとは?

傾向スコアとは、ある個体が介入を受ける確率を指します(例えば、ユーザーに特定の広告が表示される確率)。傾向スコアは、交絡因子を説明変数としてロジスティック回帰を用いて推定されることが一般的です。

例えば、傾向スコアが同一の個体群では、処置群と対照群がランダムに割り振られているとみなすことができます。これにより、交絡因子の影響を適切に除去できます。

逆確率重み付け法の計算方法

傾向スコアを用いて、各サンプルに重みを付けることで、介入群と対照群の期待値を計算します。以下の式で平均処置効果(ATE)を求めます。

$$
\text{ATE} = \frac{1}{N} \sum_{i}^{N} \frac{y_i}{P(Z=1 | X=x_i)}Z_i – \frac{1}{N} \sum_{i}^{N} \frac{y_i}{P(Z=0 | X=x_i)}(1 – Z_i)
$$

ここで、

  • \(P(Z=1 | X=x_i)\):傾向スコア(介入される確率)
  • \(y_i\):個体 \(i\) のアウトカム
  • \(Z_i\):個体 \(i\) が処置を受けたかどうか(\(1\)なら処置群、\(0\)なら対照群)

この重み付けにより、介入群と対照群の差分から平均的な因果効果を推定できます。

逆確率重み付け法のサンプルコードは、ここをクリックして展開して確認できます。

IPWによる因果推論もDoWhyパッケージを使用します。基本的な実装方法は前述した回帰分析による因果推論と同じですので、以下にコード全体をまとめて紹介します。異なるのは、estimate_effectメソッドの引数method_nameに指定するアルゴリズム名のみです。

以下のコードで推定されるATEは、回帰分析による因果推論よりも精度が悪く見えるという結果が出ていますが、これは、回帰モデルに合うようにダミーデータが生成されているためです。そのため、この結果のみからIPWの精度が低いと結論づけることはできません。

# IPWによる因果推論

import dowhy.datasets
import pandas as pd
from dowhy import CausalModel

pd.set_option('display.unicode.east_asian_width', True)

def create_dataset() -> dict:
    # Y:結果(アウトカム)、 v0:介入(処置)、W:共通原因
    data = dowhy.datasets.linear_dataset(
        beta=10,                # 真の因果効果
        num_common_causes=3,    # 観測された共通原因
        num_samples=3000,       # サンプル数
    )

    return data


if __name__ == "__main__":
    # 因果推論用のダミーデータ生成
    data = create_dataset()

    # 生成されたデータを確認
    # カラムの内容: Y=結果(目的変数), v0=介入(処置)、Wi=説明変数 (i=0, 1, ...)
    print(data["df"])

    # 因果モデルを作成
    model = CausalModel(
        data=data["df"],
        treatment=data["treatment_name"],   # 介入があったデータの変数名
        outcome=data["outcome_name"],       # 出力変数名
        graph=data["gml_graph"]             # グラフ情報
    )

    # 因果ダイアグラムを表示
    model.view_model(size=(8, 4))

    # d分離
    estimand = model.identify_effect(proceed_when_unidentifiable=True)
    print(estimand)

    # IPWによる因果効果の推定
    estimate = model.estimate_effect(
        identified_estimand=estimand,
        method_name="backdoor.propensity_score_weighting",  # アルゴリズムにIPWを指定
        target_units="ate"                                  # 平均処置効果(ATE)を推定
    )
    print(f"ATE: {estimate.value}")

Meta Learners

Meta Learnersは、機械学習と因果推論の考え方を組み合わせて、条件付き平均処置効果(CATE: Conditional Average Treatment Effect)を推定する手法の総称です。この手法は、共変量 \(X\) に基づいて処置効果が異なる非線形な関係を捉えることができます。

Meta Learnersの必要性

従来の逆確率重み付け法(IPW)は、目的変数 \(Y\) が各変数の線形和では表せない場合や、共変量 \(X\) によって効果が異なる非線形な効果には対応できません。この課題を解決するためにMeta Learnersが開発されました。

Meta Learnersでは、未観測データ(反実仮想)を機械学習モデルを使って予測し、その結果を因果推論に活用します。主に使用される機械学習モデルには、Random Forestや勾配ブースティングなどのアンサンブルモデルがあります。

Meta Learnersの種類

モデル特徴メリットデメリット
S-Learner1つのモデルで介入変数 \(Z\) を含めた予測を行う動作が速く、シンプルなモデルで実装可能介入効果が小さい場合、過小評価のリスクがある
T-Learner介入群と対照群でそれぞれ独立したモデルを構築介入効果の過小評価を抑制データ量が少ない場合、精度が低下する可能性がある
X-LearnerT-Learnerの推定値を傾向スコアで補正し精度を向上データがインバランスな場合に優れた性能を発揮モデル構築が複雑で計算負荷が高い
Meta Learnersのサンプルコードは、ここをクリックして展開して確認できます。

Meta-Learnersは、EconMLパッケージに実装されています。まずは、必要なパッケージの読み込みとデータセットの準備を行います。

今回、実験用のデータセットにはLaLonde Datasetを使用します。LaLondeデータセットは、1970年代のアメリカで行われた就労支援プログラムであるNational Support Work (NSW)の実験データです。このプログラムは、失業中の若年層に対して就労訓練や職業紹介などの支援を行い、その効果を検証することを目的として作成されてます。詳しくは、公式ページを参照してください。

まずは、パッケージの読み込みとデータセットの準備を行います。

# パッケージ
import pandas as pd
from econml.metalearners import SLearner, TLearner, XLearner
from sklearn.ensemble import RandomForestRegressor
from sklearn.linear_model import LogisticRegression

# データセットのダウンロード先
DATA_URLS = [
    "http://www.nber.org/~rdehejia/data/nswre74_treated.txt",
    "http://www.nber.org/~rdehejia/data/nswre74_control.txt",
    "http://www.nber.org/~rdehejia/data/psid_controls.txt",
    "http://www.nber.org/~rdehejia/data/psid2_controls.txt",
    "http://www.nber.org/~rdehejia/data/psid3_controls.txt",
    "http://www.nber.org/~rdehejia/data/cps_controls.txt",
    "http://www.nber.org/~rdehejia/data/cps2_controls.txt",
    "http://www.nber.org/~rdehejia/data/cps3_controls.txt",
]

# データセットのカラム名
COLUMN_NAMES = [
    "training",     # 処置割り当て指標
    "age",          # 参加者の年齢
    "education",    # 教育年数
    "black",        # 黒人であるかどうかの指標
    "hispanic",     # ヒスパニックであるかどうかの指標
    "married",      # 結婚しているかどうかの指標
    "nodegree",     # 高校卒業資格を持っていないかどうかの指標
    "re74",         # 1974年の実質収入(研究参加前)
    "re75",         # 1975年の実質収入(研究参加前)
    "re78",         # 1978年の実質収入(研究終了後)
]

# データセットのダウンロード
files = [pd.read_csv(url, sep=r"\s+", header=None, names=COLUMN_NAMES) for url in DATA_URLS]
lalonde = pd.concat(files, ignore_index=True)
lalonde = lalonde.sample(frac=1.0, random_state=42)  # Shuffle

print(lalonde)
       training   age  education  black  hispanic  married  nodegree       re74         re75         re78
16827       0.0  26.0       13.0    0.0       0.0      0.0       0.0     58.778     50.12903     31.03226
5412        0.0  27.0       12.0    0.0       0.0      1.0       0.0  16297.180  13429.21000  19562.14000
15399       0.0  26.0       12.0    0.0       0.0      0.0       0.0   5217.527   3174.24200  25564.67000
13077       0.0  38.0       16.0    0.0       0.0      1.0       0.0  23713.010   9178.98400  18814.41000
2189        0.0  55.0        8.0    0.0       0.0      1.0       1.0      0.000      0.00000      0.00000
...         ...   ...        ...    ...       ...      ...       ...        ...          ...          ...
11964       0.0  29.0       18.0    0.0       0.0      1.0       0.0  10329.250   3580.64500    678.27650
21575       0.0  17.0       10.0    1.0       0.0      0.0       1.0      0.000      0.00000   1053.61900
5390        0.0  29.0       12.0    1.0       0.0      1.0       0.0  17308.160   5473.01600  21731.45000
860         0.0  22.0       12.0    0.0       0.0      1.0       0.0  12051.450  16025.17700  13151.76700
15795       0.0  21.0       12.0    1.0       0.0      0.0       0.0  24394.830  16839.77000  22065.41000

[22106 rows x 10 columns]

次に、目的変数、説明変数、介入変数をそれぞれ定義していきます。

# 目的変数
COLUMN_Y = "re78"

# 介入変数
COLUMN_T = "training"

# 説明変数
COLUMN_X = ["age", "education", "black", "hispanic", "married", "nodegree", "re74", "re75"]

では、まずはS-Learnerを使用して条件付き平均処置効果(CATE: Conditional Average Treatment Effect)を推定していきます。

# CATEの推定に使用する機械学習モデル
# 今回はランダムフォレスト使用するが、勾配ブースティングやLightGBMなどのモデルで代替してもよい
overall_model = RandomForestRegressor()

# S-Learnerの訓練
s_learner = SLearner(overall_model=overall_model)
s_learner.fit(
    Y=lalonde[COLUMN_Y],    # 目的変数
    T=lalonde[COLUMN_T],    # 介入変数
    X=lalonde[COLUMN_X],    # 説明変数
)

# 因果効果の推定
effects = s_learner.effect(X=lalonde[COLUMN_X])
print(f"CATE: {effects.mean()}")
CATE: 317.31242275798223

次は、T-Learnerで因果効果を推定してみます。

# CATEの推定に使用する機械学習モデル
# 勾配ブースティングやLightGBMなどのモデルで代替してもよい
model = RandomForestRegressor(n_estimators=100, max_depth=5)

# T-Learnerの訓練
t_learner = TLearner(models=model)
t_learner.fit(
    Y=lalonde[COLUMN_Y],    # 目的変数
    T=lalonde[COLUMN_T],    # 介入変数
    X=lalonde[COLUMN_X],    # 説明変数
)

# 因果効果の推定
effects = t_learner.effect(X=lalonde[COLUMN_X])
print(f"CATE: {effects.mean()}")
CATE: -5682.151436708611

最後に、X-Learnerを試してみましょう。

# 傾向スコアを推定するモデル
# 通常、ロジスティック回帰が使用される
propensity_model = LogisticRegression(max_iter=10_000)

# CATEの推定に使用する機械学習モデル
# 勾配ブースティングやLightGBMなどのモデルで代替してもよい
model = RandomForestRegressor()

# X-Learnerの訓練
x_learner = XLearner(models=model, propensity_model=propensity_model)
x_learner.fit(
    Y=lalonde[COLUMN_Y],    # 目的変数
    T=lalonde[COLUMN_T],    # 介入変数
    X=lalonde[COLUMN_X],    # 説明変数
)

effects = x_learner.effect(X=lalonde[COLUMN_X])
print(f"CATE: {effects.mean()}")
CATE: -3979.270936181339

以上が、Meta-Learnersを使用した因果推論の実装例になります。

まとめ

今回は、因果推論の基礎的な考え方から、代表的な手法までを解説しました。因果推論は、日々進化を続けている分野であり、機械学習との融合や、より複雑な因果構造の推定など、新たな研究が盛んに行われています。

因果推論を学ぶことで、より深い洞察を得て、データからより多くの価値を引き出すことができるようになります。

More Information: