因果推論 Day 15／全30回　DoWhyの4ステップ、仮定を書き推定し反証する

Posted at 2026-08-27

この連載について

因果推論を「本を読んだ」で終わらせず、自分の言葉で説明でき、コードで再現できる状態まで落とす30日連載です。直前のDay 14では、SCMのノイズを個体ごとに復元して「もしあの時」を解く、反事実の3ステップを実装しました。今日はPhase 2（DAGと識別）の締めくくりです。Day 10から積んだDAG・バックドアと、Day 7までに組んだ傾向スコアが同じ土台に乗ることを、Day 7と同じ実データ nhefs（禁煙と体重変化）のend-to-endで確かめます。道具は DoWhy です。

TL;DR（3行）

DoWhy のワークフローは model → identify → estimate → refute の4ステップ。 本体は推定ではなく仮定の明示化です。DAGを書くと、nhefsでは9共変量のバックドア調整集合と「未観測交絡があれば崩れる」という仮定文が自動で出てきます。
同じ識別結果の上で推定量を差し替えられます。 素朴な差2.54 kgに対し、傾向スコア重み付けで 3.34 kg、回帰で 3.38 kg。Day 7 の2乗項入り傾向スコアを差し込むと 3.44 kg で手組みIPWと小数点まで一致します。
refute は「仮定の証明」ではなく「壊れ方の検査」。 処置列をシャッフルすると効果は 3.38 → 0.00（p=0.96）に消え、ダミー共通原因を足しても 3.38 のまま動きません。合格しても非交絡性（Day 6）は保証されず、その定量化は感度分析（Day 9）の仕事です。

今日の問い

DAG（Day 10）、バックドア基準（Day 11）、do演算子（Day 12）、そして傾向スコア（Day 7）。道具をバラバラに学んできた。
これらを実務の1本のワークフローに組むと、どういう形になるのか？

答えとして今日は DoWhy というライブラリの4ステップを nhefs で通します。ただし主役はライブラリの使い方ではありません。「因果推論の成果物は推定値ではなく、仮定と、その仮定の下での推定と、仮定の検査のセットである」という考え方を、コードの形で言語化するのが今日のゴールです。

概念①　4ステップ：分析の本体は仮定の宣言

DoWhy（Sharma & Kiciman 2020、現在は PyWhy プロジェクト）は、因果分析を4つのステップに分けて強制するライブラリです。

順番が大事です。①と②はデータの数値を1つも使いません。DAGはデータから出てこない仮定の宣言であり（Day 10）、識別はその仮定からの純粋に論理的な帰結です（Day 11、12）。推定量の選択はその後、③で初めて登場します。つまり因果推論の骨格は「仮定 → 識別 → 推定」の一方通行で、機械学習の「モデルを試して精度で選ぶ」ループとは向きが逆です。

私がこのライブラリを気に入っているのは、賢いからではなく仮定を書かされるからです。sklearn は交絡について何も聞かずに fit() してくれます。DoWhy はDAGを渡さないと始まらず、識別結果には「もし未観測の $U$ が処置と結果の両方に矢印を持つなら、この式は成り立たない」という仮定文が印字されます。分析の一番弱い部分が、コードとログに残る形で表に出ます。

概念②　refute：証明ではなく、壊れ方の検査

4ステップ目の refute（反証テスト）は誤解されやすいので、先に位置づけをはっきりさせます。たとえばプラセボ処置テストは、処置列をランダムにシャッフルして同じパイプラインを流し、「効果がゼロに消えるか」を見ます。ランダム共通原因テストは、無関係な乱数の共変量をDAGに足し、「推定値が動かないか」を見ます。

これらが検出できるのは、パイプラインの実装ミスや、たまたまの紐付きに乗った脆い推定です。検出できないのは未観測交絡です。 シャッフルしても乱数を足しても、測っていない変数は測っていないままだからです。Day 6 で見たとおり非交絡性はデータから検証できず、それはDoWhyを通しても1ミリも変わりません。だからrefuteの合格は「因果が正しい」の証明ではなく、健康診断の「今日の検査項目では異常なし」です。仮定がどれだけ崩れると結論が覆るかを定量化する仕事は、E-value（Day 9）が担当します。役割が違う2つの検査だと理解するのが正確です。

なお DoWhy 自身にも add_unobserved_common_cause という感度分析系の refuter があり、仮想の未観測交絡を強さを変えて注入し、推定値がどこまで動くかを見られます。Day 9 で手組みした発想のライブラリ版です。今日は「壊れ方の検査」に役割を絞った3つを使います。

手を動かす①：model と identify、データを見る前に決まること

データは Day 7 と同じ nhefs（Hernán & Robins の教材データ。1971年時点の喫煙者を11年追跡し、禁煙 qsmk が体重変化 wt82_71 を何kg動かすかを問う）です。Day 7 で使った9共変量が、それぞれ禁煙のしやすさと体重変化の両方に効く、という仮定をDAGに書きます。

import warnings

import numpy as np
import pandas as pd
from causaldata import nhefs
from dowhy import CausalModel

warnings.filterwarnings("ignore")   # dowhy内部のsklearn/pandas由来の警告を抑制して出力を読みやすく保つ

d = nhefs.load_pandas().data
d = d.dropna(subset=["wt82_71"]).copy()   # Day 7 と同じ：アウトカム欠測の63人を除外

covs = ["sex", "race", "age", "education", "smokeintensity",
        "smokeyrs", "exercise", "active", "wt71"]

# DAG を GML 文字列で書く：9共変量それぞれが処置(qsmk)と結果(wt82_71)の両方の原因、
# qsmk -> wt82_71 が推定したい因果の矢印
nodes = " ".join('node [id "%s" label "%s"]' % (n, n) for n in ["qsmk", "wt82_71"] + covs)
edges = " ".join('edge [source "%s" target "qsmk"] edge [source "%s" target "wt82_71"]' % (c, c)
                 for c in covs)
gml = 'graph [directed 1 %s edge [source "qsmk" target "wt82_71"] %s]' % (nodes, edges)

model = CausalModel(data=d, treatment="qsmk", outcome="wt82_71", graph=gml)
estimand = model.identify_effect()
print("バックドア調整集合:", sorted(estimand.get_backdoor_variables()))

バックドア調整集合: ['active', 'age', 'education', 'exercise', 'race', 'sex', 'smokeintensity', 'smokeyrs', 'wt71']

identify_effect() がDAGを解析し、バックドア基準（Day 11）を満たす調整集合を自動で列挙しました。今日のDAGは素直な形なので9共変量がそのまま出ますが、中間変数や合流点を含む複雑なDAGなら「調整してはいけない変数」を除いた集合が返ります。print(estimand) で全文を出すと、調整化公式（Day 12）と一緒に Unconfoundedness の仮定文、つまり「もし $U \to \text{qsmk}$ かつ $U \to \text{wt82_71}$ なる未観測 $U$ があれば、この等式は成り立たない」という告白まで印字されます。ここまで、体重のデータは1行も集計していません。

DAGの渡し方はGML文字列のほかに networkx の DiGraph も受け付けます。ノードが増えて文字列がつらくなったら、networkx で組み立てて渡す方が保守しやすいです（今日の9共変量なら上の数行で足ります）。

手を動かす②：estimate、同じ識別結果の上で推定量を並べる

識別が済んだので、推定量を選びます。同じ estimand に対して、傾向スコア重み付け（Day 7 のIPW）と回帰調整の2本を流します。sklearn の LogisticRegression は既定でL2正則化がかかるので、Day 7 の教訓どおり素の最尤推定に切り替えて渡します。

import statsmodels.formula.api as smf
from sklearn.linear_model import LogisticRegression

est_ipw = model.estimate_effect(
    estimand, method_name="backdoor.propensity_score_weighting",
    method_params={"weighting_scheme": "ips_weight",
                   "propensity_score_model": LogisticRegression(penalty=None, max_iter=1000)})
est_reg = model.estimate_effect(estimand, method_name="backdoor.linear_regression")
print("IPW（DoWhy既定の共変量仕様）= %.2f kg" % est_ipw.value)
print("回帰調整                    = %.2f kg" % est_reg.value)

# 答え合わせ：Day 7 と同じ2乗項入りスペックの傾向スコアを propensity_score 列で差し込む
d2 = d.copy()
d2["propensity_score"] = smf.logit(
    "qsmk ~ sex + race + age + I(age**2) + C(education)"
    " + smokeintensity + I(smokeintensity**2) + smokeyrs + I(smokeyrs**2)"
    " + C(exercise) + C(active) + wt71 + I(wt71**2)", data=d2).fit(disp=0).predict(d2)
model7 = CausalModel(data=d2, treatment="qsmk", outcome="wt82_71", graph=gml)
est_day7 = model7.estimate_effect(
    model7.identify_effect(), method_name="backdoor.propensity_score_weighting",
    method_params={"weighting_scheme": "ips_weight"})
print("IPW（Day 7 のスコア差し込み）= %.2f kg" % est_day7.value)

IPW（DoWhy既定の共変量仕様）= 3.34 kg
回帰調整                    = 3.38 kg
IPW（Day 7 のスコア差し込み）= 3.44 kg

素朴な差2.54 kg（Day 7）に対し、IPWが 3.34 kg、回帰が 3.38 kg。教科書の答え約3.4 kgのまわりにそろいます。3本目が今日の答え合わせです。DoWhyは、データに propensity_score 列があると自前の推定をスキップしてその列を使います。そこで Day 7 と同じ2乗項入りロジスティック回帰のスコアを statsmodels で作って差し込むと、3.44 kg。Day 7 で手組みした自己正規化IPWと小数点まで一致しました。Phase 1 の手組みとPhase 2 のライブラリは、同じ土台の上の同じ計算です。ちなみに 3.34 と 3.44 の差は、DoWhy既定の傾向スコアが9共変量を線形のまま使う（2乗項なし）ことから来ます。識別が同じでも、推定段には仕様の自由度が残るわけです。

手を動かす③：refute、自分の結果を壊しにいく

最後に反証テストです。回帰の推定に対して3種類を流します。プラセボ処置（処置列をシャッフル。効果は0に消えるべき）、ランダム共通原因（ダミー交絡を足す。効果は動かないべき）、データサブセット（無作為の80%で再推定。効果は動かないべき）です。

common = dict(num_simulations=100, random_seed=15)   # シード固定。実行が重ければ回数を絞ってよい（本文の数値は100回）

ref_placebo = model.refute_estimate(
    estimand, est_reg, method_name="placebo_treatment_refuter",
    placebo_type="permute", **common)
ref_rcc = model.refute_estimate(
    estimand, est_reg, method_name="random_common_cause", **common)
ref_subset = model.refute_estimate(
    estimand, est_reg, method_name="data_subset_refuter",
    subset_fraction=0.8, **common)

for name, r in [("プラセボ処置（処置列をシャッフル）", ref_placebo),
                ("ランダム共通原因（ダミー交絡を追加）", ref_rcc),
                ("データサブセット（80%で再推定）", ref_subset)]:
    print("%s: 元の効果 %.2f → 反証後 %.2f（p=%.2f）"
          % (name, r.estimated_effect, r.new_effect, r.refutation_result["p_value"]))

プラセボ処置（処置列をシャッフル）: 元の効果 3.38 → 反証後 0.00（p=0.96）
ランダム共通原因（ダミー交絡を追加）: 元の効果 3.38 → 反証後 3.38（p=0.90）
データサブセット（80%で再推定）: 元の効果 3.38 → 反証後 3.38（p=0.94）

処置をシャッフルすると効果は 0.00 に消え、ダミー交絡を足しても80%に間引いても 3.38 のまま動きません。p値は「期待される値（プラセボなら0、他は元の推定値）が反証後の分布から普通に出る確率」で、大きいほど無事です。3つとも期待どおりの壊れ方をした、つまり今日のパイプラインは処置とアウトカムの偶然の紐付きや実装ミスには乗っていない、と読めます。繰り返しますが、これは非交絡性の証明ではありません。禁煙の裏にある「健康への意識の変化」のような未観測交絡が仮にあっても、この3テストは全部合格します。それを見張るのはDAGを書いた人間と感度分析（Day 9）です。

実行時間の感覚も書いておきます。3テストは各100回の再推定で、回帰推定量なら手元のマシンで合計10秒ほどでした。推定量が重い場合は num_simulations を絞っても結論の傾向は変わりません（シードだけ固定して、絞ったことを明記します）。

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

identify が通っても、DAGが正しいことにはなりません。 identify_effect() が検査するのは「このDAGを信じるなら識別できるか」という論理だけで、DAGそのものの正しさは検査していません。喫煙年数から体重への矢印を1本忘れたDAGでも、識別は平然と通ります。DAGは仮定の宣言（Day 10）であり、その品質を保証するのはドメイン知識とレビューです。
同じ識別結果でも、推定値は推定量の仕様で動きます。 今日は 3.34（DoWhy既定の線形スペック）と 3.44（2乗項入り）で 0.1 kg の幅が出ました。どちらも正しい手続きです。だから報告するときは「識別の仮定」と「推定量の仕様」を分けて書き、可能なら複数推定量を並べます。1本の数字だけを出すと、この幅が見えなくなります。
dowhy は手元のデータフレームに列を書き込みます。 私の環境（dowhy 0.14）では、傾向スコア系の推定を実行すると propensity_score などの列が元の d に追加され、しかもその推定に refute をかけると古い傾向スコア列が再利用されます。プラセボ処置なのに効果が0に消えず 0.70 前後が残る、という不可解な結果はこれが原因でした。今日refuteを回帰推定に対して流したのはこのためです。傾向スコア系をrefuteするなら、先に列を落とすか、データのコピーで新しいモデルを作ってください。

GISデータ実務での使い方

分析レポートの定型を「DAG＋識別結果＋推定値＋反証テスト」にします。 衛星データ導入の施策評価なら、人口・財政力・管路老朽化率から導入と成果への矢印を張ったDAGを冒頭に置きます。仮定が図になった瞬間、統計の専門家でない自治体側や営業側も「この矢印はおかしい」と指差しでレビューできます。推定値だけのレポートでは、この会話が起きません。
refute はデータ分析版のテストコードです。 ノートブックで1回流して終わりにせず、プラセボとランダム共通原因をパイプラインに組み込めば、データの更新や前処理の変更で推定が壊れたことに気づけます。実装ミスで「効果があるように見えていた」事故は、モデルの数式ではなくこういう仕組みが拾います。
推定量は1本に絞らず並べて出します。 Fig 1 のように素朴な差・IPW・回帰を1枚に並べると、「調整で2.54が3.4に動く（交絡の大きさ）」と「調整済み同士は0.1しか違わない（仕様への頑健さ）」が同時に伝わります。この形はそのまま意思決定者向けの1枚図になります。
分析を出す前の点検項目としてのワークフローは、Day 27 のチェックリストで15の問いに展開します。

参考（本棚）

Robert Osazuwa Ness『Causal AI』第11章。DoWhyの4ステップワークフローの解説と評価。今日の下敷き
Sharma & Kiciman (2020) DoWhy: An End-to-End Library for Causal Inference。「仮定を明示させる」設計思想の原典
DoWhy 公式ドキュメント（PyWhy）。推定量と反証テストの一覧、GML/networkxでのDAG指定
Hernán & Robins『Causal Inference: What If』第12章。nhefsの本家。約3.4 kgの答え合わせは今日も一致
データ: nhefs（ causaldata パッケージ経由、MITライセンス、取得日2026-08-15）。Day 7 と同一の抽出・欠測処理

次回予告（Day 16）

Phase 2 はDAGで仮定を描き、識別し、ワークフローに束ねるところまで来ました。次回からPhase 3、因果推論のためのデザインと計量経済の道具に入ります。初回は意外かもしれませんが、実務で一番よく使う道具、昔ながらの線形回帰です。深層学習の時代に、なぜ因果推論の実務は回帰から始まるのか。回帰係数が層別比較の加重平均になっていること、そしてFWL定理の「残差×残差」の幾何を、手を動かして確かめます。

You get articles that match your needs
You can efficiently read back useful information
You can use dark theme

What you can do with signing up

因果推論 Day 15／全30回 DoWhyの4ステップ、仮定を書き推定し反証する