この連載について
因果推論を「本を読んだ」で終わらせず、自分の言葉で説明でき、コードで再現できる状態まで落とす30日連載です。前回のDay 13では、交絡因子がどうしても観測できなくても、処置からアウトカムへの影響を媒介変数のチェーンで丸ごと捉えられれば因果効果を点で識別できるという抜け道、フロントドア基準を確かめました。ここまで扱ってきた do は、集団への前向きの問いです。今日はPearlのはしごの最上段、観測済みの個体を後ろ向きに問う反事実に登ります。図と数値はすべて合成データで、SCMを私たちが握った状態で答え合わせをします。
TL;DR（3行）
- 因果のはしごは3段。関連（見る）・介入（する）・反事実（振り返る）。 do は「これから全員に介入したら」という集団の前向きの問い、反事実は「観測済みのこの個体が違う選択をしていたら」という個体の後ろ向きの問い。証拠で個体のノイズを縛るかどうかが分かれ目。
- 反事実は3ステップで解けます。アブダクション（観測からノイズ $U$ を復元）→ 介入（方程式を書き換え）→ 予測（同じ $U$ で再計算）。 血圧76で観測された患者の「服薬しなかった世界」は 86。集団の介入平均 70 では答えられない、個体の量です。
- ただし反事実はSCMの関数形に依存します。 観測分布も介入分布も同一の2つのモデルで、「処置のおかげで回復した割合」PNは 0.623 / 0.754 に割れます。識別には、介入より強い仮定（今日の例では単調性）が要ります。
今日の問い
「この患者は投薬を受けて、回復した。もし受けていなかったら、回復していただろうか」。
Day 12 の do は「これから全員に投薬したら」という問いに答える道具だった。いま問われているのは、すでに結果を見てしまった1人の、起きなかった世界のほうだ。
データではなく、この個体の「もしあの時」を計算するとはどういうことか？
今日のゴールは2つです。まず「do と反事実はどう違うのか」に自分の言葉で答えられるようになること。次に、反事実を解く3ステップ（アブダクション→介入→予測）をnumpyで素直に実装し、介入では答えられない問いに実際に答えてみることです。
概念① 因果のはしご：関連・介入・反事実の3段
Pearlは因果的な問いを3段のはしごに整理しました。第1段は関連、$P(Y \mid X)$ です。服薬した人の血圧を「見る」段で、観測データだけで登れます。第2段は介入、$P(Y \mid do(X))$ です。全員に服薬「させたら」どうなるかで、Day 11・12で見たとおりDAGと調整（あるいは実験）が要ります。そして第3段が今日の反事実です。「服薬して血圧76だったこの患者は、服薬しなければいくつだったか」。起きたことを見た上で、起きなかった世界を「振り返る」段です。
do と反事実の違いは、この記事でいちばん大事な1点なので先に言い切ります。do は前向きで集団の問い、反事実は後ろ向きで個体の問いです。$P(Y \mid do(X=0))$ は「これから母集団全員を処置なしにしたら」を問うており、式のどこにも個体の観測は現れません。一方、反事実 $P(Y_{X=0} \mid X=1, Y=76)$ は、添字に $X=0$、条件に $X=1$ が同居しています。同じ変数の2つの値が1つの式に共存するのが第3段の印で、この形は do 記法では書けません。実際に起きたこと（証拠）でその個体の素性を縛ってから、別の世界を計算し直す。この「証拠で縛る」工程が、介入には無く反事実にはある唯一の違いです。
概念② SCMと3ステップ：証拠でノイズを縛ってから、別の世界を計算する
第3段の入場料が**SCM（構造的因果モデル）**です。SCMは、各変数を「親の関数＋外生ノイズ」で書いた連立方程式です。今日の例では次の形を使います。
$$
Z = U_Z, \qquad X = f_X(Z, U_X), \qquad Y = 70 + 8Z - 10X + U_Y
$$
$Z$ は重症度、$X$ は服薬、$Y$ は収縮期血圧です。DAGは「矢印があるか」しか言いませんが、SCMは矢印の中身（関数形）まで指定します。そして外生ノイズ $U$ が、方程式に現れない個体差、つまりその人の「個性」を全部背負います。関数形とノイズの分布まで指定すれば、関連も介入も反事実も、すべて計算できるようになります。
反事実の計算は、Pearlが定式化した3ステップに従います。
- アブダクション（abduction）：観測された事実と矛盾しないように、その個体の外生ノイズ $U$ の値（または事後分布）を復元します。
- 介入（action）：問いに合わせて方程式を書き換えます。$X$ の式を捨てて $X = 0$ に差し替えます。
- 予測（prediction）：復元した $U$ のまま、書き換えた方程式で計算し直します。
do は、このうちステップ2と3だけを実行したものです。$U$ は母集団の分布のまま、方程式だけを書き換えます。反事実はそこにステップ1が加わります。一言でまとめると、縛らずに介入するのが do、証拠で縛ってから介入するのが反事実です。
手を動かす①：1個体の「もしあの時」を3ステップで解く
重症度 $z=1.5$ の患者が服薬し（$x=1$）、血圧 $y=76$ が観測されたとします。この患者が服薬していなかった世界を、3ステップそのままのコードで解きます。
# 構造方程式（この関数形を信じることが、反事実計算の前提になる）
def f_y(x, z, u_y):
return 70 + 8 * z - 10 * x + u_y # 収縮期血圧: 処置Xで-10、重症度Zで+8
# ある患者の観測: 重症度 z=1.5 で処置を受け (x=1)、血圧は y=76 だった
z_obs, x_obs, y_obs = 1.5, 1, 76.0
# ステップ1 アブダクション: 観測と矛盾しないように、この個体のノイズ U_Y を復元する
u_y = y_obs - f_y(x_obs, z_obs, 0.0)
# ステップ2 介入: 方程式の X を 0 に書き換える（do(X=0)）
x_cf = 0
# ステップ3 予測: 復元した U_Y のまま、書き換えた方程式で再計算する
y_cf = f_y(x_cf, z_obs, u_y)
print("モデルの予測 f(1, z) = %.1f" % f_y(x_obs, z_obs, 0.0))
print("復元されたノイズ U_Y = %+.1f" % u_y)
print("観測された血圧 = %.1f (処置あり)" % y_obs)
print("反事実の血圧 = %.1f (もし処置を受けていなかったら)" % y_cf)
モデルの予測 f(1, z) = 72.0
復元されたノイズ U_Y = +4.0
観測された血圧 = 76.0 (処置あり)
反事実の血圧 = 86.0 (もし処置を受けていなかったら)
モデルが $z=1.5$ の服薬者に予測する血圧は72です。ところが観測は76でした。この差 +4 がアブダクションの成果で、測っていない体質や生活習慣の束が、この個体を平均より4だけ高くしていると分かります。次に方程式の $X$ を0に書き換え、復元した+4を持たせたまま再計算すると、反事実の血圧は 86 です。集団の介入平均 $E[Y \mid do(X=0)] = 70$ とも、重症度1.5の層の介入平均82とも違います。証拠が「この個体は平均より4高い」と教えてくれたので、その4を反事実の世界へ持ち越せたのです。
ひとつ正直に付け加えます。この線形・定数効果のSCMでは、反事実は結局「観測値＋10」という単純な形になります。それが言えるのは、線形・加法ノイズという関数形を信じたからです。関数形への依存がどれほど本質的かは、手を動かす③で数値で確かめます。
手を動かす②：介入分布と反事実分布は別物
次は集団レベルです。Day 6と同じ交絡構造、つまり重症な人ほど服薬しやすい観察データをSCMで作ります。その上で「処置なしの世界」を2通りの方法で問います。1つは全員の方程式を $X=0$ に書き換える介入分布 $P(Y \mid do(X=0))$。もう1つは、実際に処置を受けた人たちだけを取り出し、1人ずつアブダクションしてから $X=0$ で再計算する反事実分布 $P(Y_{X=0} \mid X=1)$ です。
import numpy as np
rng = np.random.default_rng(14)
N = 200_000
# SCM: 重症度Zが処置XとアウトカムYの両方を動かす（Day 6 と同じ交絡構造）
u_z, u_x, u_y = rng.normal(0, 1, N), rng.normal(0, 1, N), rng.normal(0, 3, N)
z = u_z
x = (1.2 * z + u_x > 0).astype(int) # 重症な人ほど処置を受けやすい
y = 70 + 8 * z - 10 * x + u_y # 観測されるアウトカム
# 介入分布 P(Y|do(X=0)): 全員の方程式の X を 0 に書き換える（ノイズはそのまま）
y_do0 = 70 + 8 * z + u_y
# 反事実分布 P(Y_{X=0}|X=1): 処置を受けた人だけを取り出し、
# アブダクションで U_Y を復元してから（zは観測済み）、X=0 で再計算する
t = x == 1
u_rec = y[t] - (70 + 8 * z[t] - 10) # ステップ1: このSCMでは厳密に戻る
y_cf0 = 70 + 8 * z[t] + u_rec # ステップ2+3: 書き換えて再計算
print("介入分布 E[Y|do(X=0)] = %.1f" % y_do0.mean())
print("反事実分布 E[Y_{X=0}|X=1] = %.1f" % y_cf0.mean())
print("処置群の重症度 E[Z|X=1] = %.2f" % z[t].mean())
print("復元誤差の最大値 = %.1e" % np.abs(u_rec - u_y[t]).max())
介入分布 E[Y|do(X=0)] = 70.0
反事実分布 E[Y_{X=0}|X=1] = 74.9
処置群の重症度 E[Z|X=1] = 0.61
復元誤差の最大値 = 7.1e-15
同じSCM、同じ「処置なしの世界」なのに、平均は 70.0 と 74.9 で差が 4.9 あります。理由は証拠です。「処置を受けた」という事実が、その人はもともと重症（$E[Z \mid X=1] = 0.61$）だと教えてくれるので、処置群の反事実分布は高い側へ寄ります。「処置を受けた人たちにとっての、処置がなかった世界」という量（Day 8で見たATTの反事実側、Pearlの言うETTの部品）は、介入分布をいくら精密に推定しても単独では出てきません。なお復元誤差が $10^{-15}$ 台であるとおり、連続アウトカム＋加法ノイズのこのSCMではアブダクションが1人ずつ厳密に決まっています。
手を動かす③：介入では答えられない問い、必要性確率PN
最後に、冒頭の患者の問いに戻ります。二値の回復データで「処置を受けて回復した人のうち、処置のおかげだった（処置がなければ回復しなかった）割合」を計算します。Pearlが**必要性確率（PN: probability of necessity）**と呼ぶ、$P(Y_{X=0}=0 \mid X=1, Y=1)$ です。薬害訴訟や帰責の議論で法廷が問うのは、まさにこの量です。今回はRCTを仕込みます。無作為化してもなお、この問いには第3段の道具が要ることを見るためです。
import numpy as np
rng = np.random.default_rng(14)
N = 200_000
u = rng.uniform(0, 1, N) # 各人の体質（外生ノイズ）。回復のしやすさを決める
x = rng.binomial(1, 0.5, N) # RCT: 処置はコイン投げで割り付ける
# モデルA（単調）: 処置なしなら u<0.3 の30%が自然回復、処置ありなら u<0.8 の80%が回復
y0_A = (u < 0.3).astype(int)
y1_A = (u < 0.8).astype(int)
y_A = np.where(x == 1, y1_A, y0_A)
# 処置を受けて回復した人のうち、処置のおかげ（＝処置なしなら回復しなかった）の割合
m_A = (x == 1) & (y_A == 1)
pn_A = (y0_A[m_A] == 0).mean()
# モデルB（非単調）: 自然回復するのは u∈[0.6, 0.9) の30%。
# 回復率はAと同じだが、u∈[0.8, 0.9) の人には処置がむしろ害になる
y0_B = ((0.6 <= u) & (u < 0.9)).astype(int)
y1_B = y1_A
y_B = np.where(x == 1, y1_B, y0_B)
m_B = (x == 1) & (y_B == 1)
pn_B = (y0_B[m_B] == 0).mean()
for name, y_, y0_ in [("A", y_A, y0_A), ("B", y_B, y0_B)]:
print("モデル%s: 回復率 処置群 %.3f / 対照群 %.3f"
% (name, y_[x == 1].mean(), y_[x == 0].mean()))
print("PN(モデルA・単調) = %.3f (理論値 (0.8-0.3)/0.8 = 0.625)" % pn_A)
print("PN(モデルB・非単調) = %.3f (同じ観測・介入分布なのに答えが変わる)" % pn_B)
モデルA: 回復率 処置群 0.802 / 対照群 0.296
モデルB: 回復率 処置群 0.802 / 対照群 0.303
PN(モデルA・単調) = 0.623 (理論値 (0.8-0.3)/0.8 = 0.625)
PN(モデルB・非単調) = 0.754 (同じ観測・介入分布なのに答えが変わる)
モデルAは単調、つまり処置が害になる人がいない世界です。処置を受けて回復した人のうち 62.3% が処置のおかげで、残りは放っておいても回復した人たちでした。この値は超過リスク比 $(0.8-0.3)/0.8 = 0.625$ と一致します。RCTと単調性がそろえば、PNは介入分布から識別できるのです（Pearlの定理の特別な場合です）。
問題はモデルBです。回復率は処置群0.802、対照群0.303と、モデルAと観測上区別がつきません。介入分布も同じです。ところがPNは 75.4% に跳ねます。Bでは体質 $u$ が $[0.8, 0.9)$ の1割の人にとって処置が害で、「誰が自然回復するのか」の割り当てが違うからです。RCTを何度繰り返してもAとBは区別できないのに、反事実の答えは変わる。反事実はSCMの関数形に依存し、その識別には介入より強い仮定（ここでは単調性）を買う必要がある。これが第3段の力と代償です。
もうひとつ、この例はアブダクションの別の顔も見せています。証拠 $(X=1, Y=1)$ が絞れるのは「$u < 0.8$」という範囲までで、①②のように点では戻りません。だから個体の反事実は点でなく確率、すなわちPNという形になります。
つまづき・誤解しやすい点
- 「do に個体の情報を条件付けすれば反事実になる」わけではありません。 $P(Y \mid do(X=0), Z=1.5)$ は「重症度1.5の人たち全員に介入したら」という集団の量で、答えは82です。手を動かす①の反事実86は、観測 $y=76$ という証拠でノイズを縛ってから介入した量です。絞り込んでから書き換えるのではなく、縛ってから書き換える。証拠にアウトカム自身を使える点が決定的に違います。
- 反事実はDAGだけでは決まりません。 手を動かす③では、同じDAG・同じ観測分布・同じ介入分布の2つのモデルでPNが0.623と0.754に割れました。第2段まではDAGとデータで登れますが、第3段は関数形という追加の仮定を買います。だから反事実の報告には「このSCMを信じるなら」という但し書きが必須です。仮定の破れにどこまで耐えるかを測る発想は、Day 9の感度分析と同じ精神です。
- アブダクションはいつも一意に戻るとは限りません。 連続アウトカム＋加法ノイズなら誤差 $10^{-15}$ で厳密に復元できましたが、二値のように証拠が粗いと $U$ は範囲までしか絞れず、反事実は点でなく確率になります。ノイズの事後分布を丸ごと扱うこの見方は、Day 29のベイズ因果推論で再登場します。
GISデータ実務での使い方
- 「この管路を昨年更新していたら、今回の漏水は起きなかったか」は第3段の問いです。 事故調査や再発防止の議論は、気づかないうちに反事実の言語で行われています。管路更新の平均効果（do）をどれだけ精密に推定しても、目の前の1件が防げたかには答えられません。会議の問いがはしごのどの段かを最初に判定するだけで、使うべき道具と答えられる限界がはっきりします。
- 反事実を語るには、劣化の構造方程式が要ります。 材質・布設年代・土壌・水圧から漏水リスクへの関数形を仮定して初めて「防げたか」が計算できます。物理モデルや劣化モデルの蓄積がある水道分野は、実はSCMと相性の良い領域です。ただし手を動かす③のとおり、モデルが違えば「防げたか」の答えも変わります。事故報告に反事実の主張を書くなら、前提にしたモデルを必ず添えます。
- 前向きは do、後ろ向きは反事実、と使い分けます。 来年度どの管路を更新すべきかは集団への前向きの問いなので、Day 12までの道具で足ります。起きた事故の帰責・検証は個体への後ろ向きの問いで、今日の道具の出番です。この2つは1つの会議で平気で混ざるので、はしごの段で切り分ける癖をつけておくと議論が濁りません。
参考（本棚）
- Pearl (2009) Causal Inference in Statistics: An Overview. Statistics Surveys 3。第5節 Counterfactuals at Work。必要性確率PN・ETTと反事実の識別条件。今日の下敷き
- Pearl, Glymour & Jewell『入門 統計的因果推論』第4章。反事実の3ステップ（アブダクション・介入・予測）の教科書的な導入
- Robert Osazuwa Ness『因果AI』第6章。SCMでの反事実実装。pyroの実装例の考え方を、今日はnumpyに置き換えました
- 図はすべて合成データです。反事実は定義上どのデータにも観測されないため、SCMを私たちが握り、全個体のノイズと両方の潜在的結果をオラクルとして持つことでしか答え合わせできません。生成スクリプトは元ノート側に添付しています
次回予告（Day 15）
Phase 2はこれで登り切りました。DAGで仮定を描き、d分離とバックドアで調整集合を選び、doで介入を定義し、フロントドアで抜け道を探し、今日はしごの最上段に立ちました。Rubin流の傾向スコアとマッチング、Pearl流のDAGとdo。流派の道具は出そろいましたが、実務ではこれらがバラバラに使われがちです。次回Day 15は、DoWhyでモデル化→識別→推定→反証を1本のワークフローに組みます。Day 7と同じnhefsデータを、今度は端から端まで1本のパイプラインで通します。