この連載について
因果推論を「本を読んだ」で終わらせず、自分の言葉で説明でき、コードで再現できる状態まで落とす30日連載です。前回のDay 8では、LaLondeの職業訓練データでマッチングを組み、観察データが実験ベンチマークにどこまで迫れるかを試しました。マッチングも傾向スコアも、土台は「測った共変量で条件付ければ割付はランダムに見える」という非交絡性の仮定です。そしてDay 6で見たとおり、この仮定はデータから検証できません。今日は、その検証できない仮定にそれでも数字で向き合う道具、感度分析とE-valueを扱います。
TL;DR（3行）
- 感度分析＝検証できない仮定を「覆すのに必要な交絡の強さ」に翻訳する道具。 仮定を全部捨てると区間が広すぎて何も言えず（Manski境界）、信じ切ると「仮定が正しければ」の但し書きが宙に浮きます。その間をつなぐダイヤルが感度分析です。
- E-value は式1本 $E = RR + \sqrt{RR(RR-1)}$。 nhefs実データの禁煙→体重変化 +3.46kg（近似リスク比 1.49）のE-valueは 2.35。測った共変量をすべて調整してなお、処置ともアウトカムともリスク比2.35で結びつく隠れ交絡がない限り、この効果はゼロになりません。
- オラクル実験で体感: 強さ2.5の交絡を仕込むと観測1.49は調整後 1.00 に消え、強さ1.5では 1.35 が残ります。E-value 2.35はちょうどその境界線です。
今日の問い
非交絡性は検証できない。それはDay 6で身をもって確かめた。隠れた交絡が「あるか、ないか」をデータは教えてくれない。
ならば問いを変えよう。いま見えている効果を丸ごと消し去るには、どれだけ強い隠れ交絡が必要か。それなら計算できるのではないか？
答えは「計算できる」です。今日のゴールは、検証不能な仮定を「覆すのに必要な交絡の強さ」という1つの数字に翻訳する発想を言語化し、その代表であるE-valueを、式1本の手組みから実データの適用まで組むことです。仕上げに、E-valueが示す強さの交絡を合成データに本当に仕込み、効果が本当に消えることをオラクル実験で確かめます。
概念① 感度分析：検証できない仮定を、覆すのに必要な強さに翻訳する
Day 7のIPWもDay 8のマッチングも、推定値には「非交絡性が成り立つなら」という但し書きが付いています。この但し書きへの向き合い方は3通りあります。信じて進むか、捨てるか、破れを少しずつ許して結論の壊れ方を測るかです。
仮定を全部捨てるとどうなるかを、Imbens & Rubin 第22章はManskiの**境界（bounds）**で見せます。宝くじ当選が就労に与える効果の分析で、非交絡性を仮定すれば点推定は $-0.134$。ところが仮定を捨てると、二値アウトカムの場合、データがどれだけあっても区間は必ず幅1になり、必ずゼロをまたぎます。宝くじデータでは $[-0.56,\ 0.44]$ でした。仮定ゼロの正直さは「ほぼ何も言えない」とセットで返ってくるわけです。
そこで中間の道が感度分析です。Rosenbaum-Rubin流の定式化では「未観測の共変量 $U$ まで条件付ければ非交絡になる」と考え、$U$ が処置・アウトカムと持つ関連の強さを、推定するのではなくダイヤルとして外から振ります。強さゼロなら非交絡性そのもの、無限大ならManski境界に一致します。感度分析はこの両極をつなぐ中間地帯を埋める道具です。Rosenbaumの $\Gamma$ はその代表で、「同じ $X$ を持つ2人でも、処置を受けるオッズは最大 $\Gamma$ 倍まで違いうる」と許した上で、結論（p値）が生き残る $\Gamma$ の限界を報告します。宝くじの例では、観測共変量の中で最強のもの相当（$\Gamma$ に翻訳して1.12）の隠れ交絡を許すと、p値の上限は0.99まで跳ね、負の効果の証拠は消えます。
概念② E-value：関連を消し去るのに必要な交絡の最小の強さ
Rosenbaum-Rubin流は柔軟ですが、パラメータが複数あって報告が重くなります。VanderWeele & Ding (2017) はこれを1つの数字に圧縮しました。E-valueとは、「測定済みの共変量を調整した上でなお残る処置とアウトカムの関連（リスク比 $RR$）を丸ごと説明し去るために、未測定の交絡因子が処置との関連とアウトカムとの関連の両方に持つ必要のある、リスク比スケールの最小の強さ」です。式は1本です。
$$
\text{E-value} = RR + \sqrt{RR \times (RR - 1)}
$$
背後にはDing & VanderWeele (2016) のbounding factorがあります。隠れ交絡 $U$ が処置と持つ関連を $RR_{EU}$、アウトカムと持つ関連を $RR_{UD}$ とすると、この交絡が作り出せる見かけの関連には上限があります。
$$
B = \frac{RR_{EU} \times RR_{UD}}{RR_{EU} + RR_{UD} - 1}
$$
観測された $RR$ のうち、交絡では $B$ までしか説明できません。つまり真の効果は少なくとも $RR / B$ 残ります。2つの関連を同じ値にそろえて $B = RR$ を解くと、上のE-valueの式が出てきます（2次方程式の解です）。片方が弱ければもう片方はより強くなければならず、そのトレードオフの全体像は手を動かす②の図で確かめます。
単位はリスク比です。回帰係数をそのまま式に入れることはできません。今日のnhefsのような連続アウトカムでは、効果をアウトカムの標準偏差で割った標準化効果量 $d$ に直してから、近似変換 $RR \approx \exp(0.91 \times d)$ を使います。0.91という係数は、$d$ を対数オッズ比に直すChinn (2000) の近似（$\ln OR \approx 1.81 d$）と、まれではないアウトカムでの $RR \approx \sqrt{OR}$ の近似を合わせたものです（$1.81 / 2 \approx 0.91$）。信頼区間には、区間の端のうち1に近い側へ同じ式を適用します。区間が1をまたぐならE-valueは1、つまり交絡を持ち出すまでもなく結論はゼロと区別できません。
手を動かす①：E-value を手組みする
式は1本なので、実装は数行で終わります。まず相場観をつかみます。
import numpy as np
def evalue(rr):
"""点推定のE-value。rr<1（保護的効果）は逆数にしてから使う。"""
rr = max(rr, 1 / rr)
return rr + np.sqrt(rr * (rr - 1))
for rr in [1.2, 1.5, 2.0, 9.0]:
print("RR = %.1f -> E-value = %.2f" % (rr, evalue(rr)))
RR = 1.2 -> E-value = 1.69
RR = 1.5 -> E-value = 2.37
RR = 2.0 -> E-value = 3.41
RR = 9.0 -> E-value = 17.49
読み方の目安ができました。観測されたリスク比が1.5なら、それを丸ごと消すには処置・アウトカムの両方とリスク比 2.37 で結びつく隠れ交絡が要ります。最後の行は喫煙と肺がんです。リスク比およそ9という強い関連に対し、E-valueは 17.49。1950年代、Cornfieldらは「これを説明し去る隠れ交絡は喫煙と桁違いに強く結びつく何かでなければならず、現実的に考えにくい」と論じて、タバコ論争の流れを変えました。E-valueはこの古典的な議論を、どの研究でも使える1行の計算にしたものといえます。
手を動かす②：nhefsの禁煙効果は E-value 2.35
実データに適用します。データはDay 7と同じ nhefs（NHEFS: 米国の全国健康栄養調査の追跡研究）で、1971〜75年時点で喫煙していた成人を1982年まで追った縦断データです。体重変化（kg）が記録された1,566人について、禁煙した人としなかった人を比べます。年齢や喫煙量などの交絡を回帰で調整し、標準化してリスク比に変換し、E-valueまで一気に計算します。
import numpy as np
import statsmodels.formula.api as smf
from causaldata import nhefs
df = nhefs.load_pandas().data.dropna(subset=["wt82_71"]).copy()
# 禁煙(qsmk)→体重変化(wt82_71)。Hernán & Robins と同じ共変量で回帰調整
f = ("wt82_71 ~ 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)")
m = smf.ols(f, data=df).fit()
b, se = m.params["qsmk"], m.bse["qsmk"]
# 連続アウトカムを標準化し、近似リスク比へ変換（VanderWeele & Ding 2017）
d, se_d = b / df["wt82_71"].std(), se / df["wt82_71"].std()
rr = np.exp(0.91 * d)
rr_lo = np.exp(0.91 * d - 1.78 * se_d) # 1.78 = 0.91 × 1.96
def evalue(rr):
rr = max(rr, 1 / rr)
return rr + np.sqrt(rr * (rr - 1))
naive = df.loc[df.qsmk == 1, "wt82_71"].mean() - df.loc[df.qsmk == 0, "wt82_71"].mean()
print("素朴な差 = %.2f kg" % naive)
print("調整後効果 = %.2f kg (95%%CI [%.2f, %.2f])" % (b, b - 1.96 * se, b + 1.96 * se))
print("標準化 d = %.2f -> 近似リスク比 RR = %.2f (CI下限 %.2f)" % (d, rr, rr_lo))
print("E-value: 点推定 %.2f ／ 信頼区間下限 %.2f" % (evalue(rr), evalue(rr_lo)))
素朴な差 = 2.54 kg
調整後効果 = 3.46 kg (95%CI [2.60, 4.32])
標準化 d = 0.44 -> 近似リスク比 RR = 1.49 (CI下限 1.35)
E-value: 点推定 2.35 ／ 信頼区間下限 2.04
まず素朴な差2.54kgが調整後に 3.46kg へ増えている点が面白いところです。禁煙した人は年齢が高めで、年齢が高いほど体重は増えにくいため、この例の交絡は効果を「隠す」方向に働いていました。調整後の3.46kgは標準化すると $d = 0.44$、近似リスク比で 1.49 です。
E-valueの読み方はこうです。年齢・性別・喫煙量など測った共変量をすべて調整した上でなお、禁煙とも体重変化とも、それぞれリスク比2.35以上で結びつく未測定の交絡因子が存在しない限り、この3.46kgをゼロにはできません。信頼区間の下限まで消すだけなら、必要な強さは 2.04 に下がります。点推定と信頼区間下限の両方を報告するのが作法です。
この2.35が「強い」のか「すぐ超えられる」のかは、測れている共変量をものさしにして判断します。Imbens & Rubin 第22章（Imbens 2003に遡る作法）は、観測済み共変量が処置・アウトカムと持つ関連の強さを計算し、隠れ交絡がそれを超えると考えるべき理由があるかを問います。もし手元の最強の共変量でもリスク比2に届かないなら、2.35の隠れ交絡は「観測されたどの要因よりも強い何か」を意味します。
下の図が、概念②で述べたトレードオフの全体像です。境界曲線より上の強さの組み合わせなら1.49を丸ごと説明でき、曲線の下では必ず効果が残ります。片方の関連が4あるなら、もう片方は1.8で足りますが、両方を同じ強さにするなら2.35が最小です。それが対角線との交点、すなわちE-valueです。
手を動かす③：E-valueが示す強さの交絡を、本当に仕込んでみる
E-valueが2.35だと言われても、「強さ2.35の交絡が効果を消す」を実感するのは難しいものです。そこでオラクル実験です。E-valueの土俵であるリスク比の世界で二値アウトカムの合成データを作り、観測される関連がnhefsと同じ1.49になるように、真の効果と交絡の配合だけを変えた2つの世界を用意します。世界Aは真の効果ゼロで、E-valueを超える強さ2.5の交絡 $U$ がすべてを作る世界。世界Bは強さ1.5の弱い交絡しかない世界です。私たちはオラクルなので、本来観測できない $U$ で調整して答え合わせができます。
import numpy as np
def mh_rr(e, y, u):
"""層別（Mantel-Haenszel）リスク比。Uの層の中で比べてから集計する。"""
num = den = 0.0
for k in (0, 1):
m = u == k
n1, n0 = (e[m] == 1).sum(), (e[m] == 0).sum()
if n1 == 0 or n0 == 0:
continue
num += y[m & (e == 1)].sum() * n0 / m.sum()
den += y[m & (e == 0)].sum() * n1 / m.sum()
return num / den
rng = np.random.default_rng(9)
N, p0, target = 1_000_000, 0.10, 1.49 # 観測される関連を nhefs の近似RR 1.49 にそろえる
# 世界A: 真の効果はゼロ。強さ2.5（E-value 2.35 超え）の隠れ交絡Uだけで1.49を作る
s = 2.5
pi1 = s * (target - 1) / ((s - 1) * (s - target)) # 観測RRが1.49になるUの配置を逆算
e = rng.binomial(1, 0.5, N)
u = rng.binomial(1, np.where(e == 1, pi1, pi1 / s)) # Uの保有率の比 RR_EU = 2.5
y = rng.binomial(1, p0 * s**u) # RR_UD = 2.5。eは式に無い＝効果ゼロ
print("世界Aの配置: 処置群の%.0f%%、対照群の%.0f%%がUを持つ" % (100 * pi1, 100 * pi1 / s))
print("世界A（真の効果1.00・交絡2.5）: 素朴RR = %.2f -> U調整後 = %.2f"
% (y[e == 1].mean() / y[e == 0].mean(), mh_rr(e, y, u)))
# 世界B: 弱い交絡（強さ1.5 < E-value）。観測1.49には本物の効果が必要になる
s2, pi1b = 1.5, 0.8
fab = (1 + (s2 - 1) * pi1b) / (1 + (s2 - 1) * pi1b / s2) # 交絡が作れる見かけ = 1.11
rr_true = target / fab # 足りない分は本物の効果 = 1.35
e2 = rng.binomial(1, 0.5, N)
u2 = rng.binomial(1, np.where(e2 == 1, pi1b, pi1b / s2))
y2 = rng.binomial(1, p0 * s2**u2 * rr_true**e2)
print("世界B（真の効果%.2f・交絡1.5）: 素朴RR = %.2f -> U調整後 = %.2f"
% (rr_true, y2[e2 == 1].mean() / y2[e2 == 0].mean(), mh_rr(e2, y2, u2)))
世界Aの配置: 処置群の81%、対照群の32%がUを持つ
世界A（真の効果1.00・交絡2.5）: 素朴RR = 1.49 -> U調整後 = 1.00
世界B（真の効果1.35・交絡1.5）: 素朴RR = 1.49 -> U調整後 = 1.35
世界Aでは、処置群の81%・対照群の32%が $U$ を持つという配置で、真の効果がゼロなのに観測リスク比1.49が出現します。$U$ で層別調整すると 1.00、効果は跡形もなく消えました。E-valueを超える強さの交絡は、観測された関連を本当に丸ごと作れるのです。ちなみに強さをちょうどE-valueの2.35にすると、1.49を作るには「処置群は全員 $U$ を持つ」という極端な配置まで必要になります。E-valueが「最も不利な配置を許したときの最小の強さ」だという意味がここに表れています。
世界Bは対照的です。強さ1.5の交絡が作れる見かけの関連は、bounding factorの上限でも $B = 1.5^2 / (2 \times 1.5 - 1) = 1.13$ 止まりです。観測された1.49と整合するには本物の効果が少なくとも1.32必要で、実際この世界では真の効果 1.35 を仕込まないと1.49になりませんでした。$U$ 調整後も1.35がそのまま残ります。E-value未満の交絡は、配置をどう選んでも結論を消せません。下の図で、この2つの世界が理論曲線の上に着地していることを確認できます。
つまづき・誤解しやすい点
- 「E-valueが大きい＝因果の証明」ではありません。 E-valueが測るのは未測定交絡への耐性だけで、測定誤差・選択バイアス・モデルの誤設定については何も言いません。しかも最悪の配置を想定した計算なので、E-valueと同じ強さの交絡が実在しても、配置次第では結論は覆りません。保証があるのは逆方向、「E-value未満の強さでは決して覆せない」の側です。
- 単位はリスク比です。回帰係数やオッズ比をそのまま式に入れてはいけません。 連続アウトカムは標準化して $RR \approx \exp(0.91 \times d)$、まれではないアウトカムのオッズ比は $\sqrt{OR}$ で近似してから使います。変換はどれも近似なので、E-valueの小数第2位を細かく議論するのは精度の過信です。
- 感度分析は仮定の身代わりではなく、報告の作法です。 E-valueを書いても非交絡性の仮定が不要になるわけではなく、「仮定がどれだけ破れたら結論が崩れるか」を読み手に開示しているだけです。仮定そのものを図として宣言し、何を調整すべきかを設計するのは、Phase 2のDAG（Day 10から）の仕事です。
実務での使いどころ（空間データ／天地人の文脈）
- 衛星で測れない交絡を、数字で語れるようにする。 施策を導入した自治体と未導入の自治体を比べるとき、人口・財政・管路延長・老朽化率は統計や衛星データで調整できても、首長の意欲や現場の体制は測れません。「E-valueは2.3。意欲や体制の違いが、導入の有無とも成果とも、この強さで結びついていない限り効果はゼロにならない」と言えれば、議論は「交絡があるかもしれない」の水掛け論から「その強さは現実的か」へ進みます。
- 効果検証レポートに1行添える運用にする。 点推定と信頼区間の下に「E-value: 点推定2.35／CI下限2.04」を定型で置くだけです。コストは数行のコードで、読み手に仮定への耐性を毎回開示できます。結論が覆りやすい分析ほど、この1行が正直さの証明になります。
- 測れている共変量をものさしにする。 隠れ交絡の強さが現実的かは、観測済み共変量の中で最強のもの（例えば管路の布設年代や自治体規模）が処置・アウトカムと持つ関連と比べて判断します。手を動かす②で触れたImbens & Rubin 第22章の作法で、E-valueに文脈を与える一手間です。
参考（本棚）
- Imbens & Rubin『Causal Inference for Statistics, Social, and Biomedical Sciences』第22章 感度分析と境界。Manski境界とRosenbaum-Rubin感度分析、観測共変量をものさしにする作法。今日の下敷き
- VanderWeele & Ding (2017) Sensitivity Analysis in Observational Research: Introducing the E-Value. Annals of Internal Medicine 167(4)。E-valueの原典。連続アウトカムの変換表もこの論文
- Ding & VanderWeele (2016) Sensitivity Analysis Without Assumptions. Epidemiology 27(3)。bounding factorの導出。E-valueの数学的土台
- Hernán & Robins『Causal Inference: What If』。nhefsデータの出典。禁煙と体重変化の定番教材
- データ: NHEFS（NHANES I Epidemiologic Followup Study の教材版）。
causaldataパッケージ（MITライセンス）経由で取得（2026-08-15）。図と数値はすべて本文コードとシード固定で再現できます
次回予告（Day 10）
Phase 1はこれで一区切りです。潜在的結果から始まり、RCT、非交絡性、傾向スコア、マッチング、そして今日の感度分析まで、「仮定を信じて推定し、仮定への耐性を測って報告する」一連の流れがそろいました。ただここまで、肝心の「どの変数を調整すべきか」は素朴な直感に頼ってきました。次回からのPhase 2では、因果の仮定そのものを1枚の図に描く言語、**DAG（有向非巡回グラフ）**に入ります。Day 1の疑似相関や合流点の罠が、実は3つの基本パターンとして再登場します。調整してよい変数といけない変数を、勘ではなく規則で見分ける準備です。