因果推断与实验评估：从潜在结果到双重稳健估计

模型能告诉你「看了广告的用户购买率高」，但不能告诉你「是广告让他们买的」。前者是预测，后者是因果。在定价、推荐、风控这些「改变一个动作会影响结果」的场景里，因果推断才是决策的正确语言。

相关不等于因果：问题动机

先看一个经典的反例。

辛普森悖论

分组看，用药组都更好；合计看，用药组反而更差。这就是辛普森悖论：分组结论与合计结论相反，根源在于重症患者更多地被分到用药组（用药组中重症占 50%，对照组只占 33%），重症本身拉低了用药组的整体存活率。

混淆变量的作用

上例中「病情严重程度」同时影响「是否用药」与「存活率」，这类变量叫混淆变量（Confounder）。它制造了虚假的关联，让朴素的组间比较失效。

因果推断的核心任务，就是在存在混淆的情况下，识别并估计真实的因果效应。

因果之梯

Pearl 提出的三层因果之梯：

机器学习模型基本停在第一层。要上到第二层，必须引入因果假设。

潜在结果框架与因果效应

潜在结果（Potential Outcomes） 框架，又称 Rubin 因果模型，是量化因果效应的标准语言。

潜在结果与个体处理效应

对每个个体 i，定义两个潜在结果：

• Y_i(1)：接受处理（如看到广告）时的结果。
• Y_i(0)：未接受处理时的结果。

个体处理效应（ITE） 是 τ_i = Y_i(1) - Y_i(0)。

根本难题是：我们只能观察到其中一个——因果推断的基本问题。观察到的是 Y_i = T_i·Y_i(1) + (1-T_i)·Y_i(0)，另一个永远是反事实。

ATE、ATT 与 CATE

既然个体效应不可得，就退而求其次估计群体平均：

ATE = E[Y(1) - Y(0)]                    # 全体平均处理效应
ATT = E[Y(1) - Y(0) | T = 1]            # 处理组的平均效应
CATE = E[Y(1) - Y(0) | X = x]           # 条件平均效应（异质效应）

三者回答不同问题：ATE 用于「要不要全面推广」，ATT 用于「已投放的效果如何」，CATE 用于「该给谁投放」。

三大识别假设

估计因果效应必须声明假设：

1. 一致性（Consistency）：观测结果等于对应处理下的潜在结果。
2. 可交换性（Exchangeability）：给定协变量 X 后，处理分配与潜在结果独立，{Y(0),Y(1)} ⊥ T | X。这等价于「无未观测混淆」。
3. 重叠性（Positivity）：任何 X 下都有正概率接受或不接受处理，0 < P(T=1|X) < 1。

第三条最容易被忽视：如果某个年龄段从来没有用户被投放广告，那这个年龄段的效果无法识别——再好的模型也救不了。

混淆、碰撞与有向无环图

识别假设不能靠拍脑袋，需要工具来判断「该控制哪些变量」。

三类变量的区别

碰撞偏差（Collider Bias） 尤其反直觉：控制一个碰撞变量，会人为地制造出原本不存在的关联。

DAG 与后门准则

有向无环图（DAG） 把变量间假设画成箭头图，让识别问题可视化。

混淆（需控制 Z）：
    Z
   / \
  v   v
  T → Y

碰撞（不可控制 C）：
  T → C ← Y

中介（不可控制 M）：
  T → M → Y

后门准则（Backdoor Criterion）：若变量集合 Z 阻断了 T 到 Y 的所有「后门路径」（以箭头指向 T 开头的路径），则 P(Y|do(T)) = Σ_z P(Y|T,z)P(z)。这就是「控制 Z 就能识别因果效应」的严格表述。

用 DAG 指导特征选择

实践中，先画 DAG 再决定控制哪些协变量，比「把所有特征都塞进去」安全得多——因为控制碰撞变量和中介变量都会引入偏差，特征越多反而越错。工具上可以用 pgmpy 的 is_dconnected 检验给定集合下的 d-分离，从而判断某个协变量集合是否满足后门准则。

随机对照试验与 A/B 测试

随机化是解决混淆最彻底的办法：随机分配处理，让混淆变量在处理组与对照组间期望上均衡。

随机化为何有效

随机化切断了混淆变量到处理分配的箭头：Z 与 T 独立，后门路径自动阻断，此时简单的组间均值差 τ̂ = E[Y|T=1] - E[Y|T=0] 就是 ATE 的无偏估计。

A/B 测试的统计基础

线上 A/B 测试是随机试验的工业实现。核心要素：

• 假设检验：原假设「两组无差异」，用 t 检验或 z 检验计算 p 值。
• 样本量：由效应量、显著性水平 α、统计功效 1-β 决定。
• 置信区间：比 p 值更有信息量，直接给出效应范围。

import numpy as np
from scipy import stats

def ab_test(control, treatment, alpha=0.05):
    """两组均值对比，返回效应量与置信区间"""
    n_c, n_t = len(control), len(treatment)
    m_c, m_t = np.mean(control), np.mean(treatment)
    v_c, v_t = np.var(control, ddof=1), np.var(treatment, ddof=1)
    se = np.sqrt(v_c / n_c + v_t / n_t)
    diff = m_t - m_c
    z = diff / se
    p = 2 * (1 - stats.norm.cdf(abs(z)))
    lo, hi = diff - 1.96 * se, diff + 1.96 * se
    return {"diff": diff, "p_value": p,
            "ci": (lo, hi), "significant": p < alpha}

# 样本量估算：检测 1% 的相对提升，基线转化率 10%
def sample_size(baseline, mde_rel, alpha=0.05, power=0.8):
    p1 = baseline
    p2 = baseline * (1 + mde_rel)
    p_bar = (p1 + p2) / 2
    z_a = stats.norm.ppf(1 - alpha / 2)
    z_b = stats.norm.ppf(power)
    n = ((z_a + z_b) ** 2 * 2 * p_bar * (1 - p_bar)) / (p2 - p1) ** 2
    return int(np.ceil(n))

常见陷阱

• 偷看数据（Peeking）：边跑边看，一旦显著就停，会让假阳性率远超 α。对策是固定样本量或用序贯检验。
• 多重比较：同时测 20 个指标，至少有一个偶然显著。用 Bonferroni 或 FDR 校正。
• 样本比例失衡（SRM）：分组比例偏离预期，说明分流逻辑有 bug，此时结果不可信。上线前必查 SRM。
• 网络效应：社交、双边市场里处理组会溢出影响对照组，破坏 SUTVA 假设，需用聚类随机化或切换实验。

倾向得分与匹配

观察性数据没有随机化，只能用调整的方法模拟随机化。倾向得分是最经典的工具。

倾向得分的定义

倾向得分（Propensity Score） 是个体接受处理的概率：

e(X) = P(T = 1 | X)

Rosenbaum 与 Rubin 证明：给定倾向得分后，协变量分布与处理分配独立。也就是说，把多维协变量压成一维得分，就能完成调整——这被称为倾向得分的「降维性质」。

匹配与加权

倾向得分匹配（PSM）：为每个处理组个体找一个得分相近的对照组个体配对。

逆概率加权（IPW）：用倾向得分的倒数加权，构造一个「伪总体」：

ATE_IPW = (1/n) Σ [ T_i·Y_i / e(X_i) - (1-T_i)·Y_i / (1-e(X_i)) ]

直觉：得分低却接受了处理的人（e 小，1/e 大）代表了那些「不太可能被处理却仍被处理」的群体，给他们更大权重，让样本重新平衡。

import numpy as np
from sklearn.linear_model import LogisticRegression

def ipw_ate(X, T, Y):
    """逆概率加权估计 ATE"""
    ps_model = LogisticRegression(max_iter=1000).fit(X, T)
    e = np.clip(ps_model.predict_proba(X)[:, 1], 1e-6, 1 - 1e-6)
    # 稳定化权重：乘以处理概率的边际
    w1 = T / e
    w0 = (1 - T) / (1 - e)
    ate = np.mean(w1 * Y) / np.mean(w1) - np.mean(w0 * Y) / np.mean(w0)
    return ate, e

重叠性与极端权重

e(X) 接近 0 或 1 时，权重会爆炸，方差失控。诊断与对策：

• 画倾向得分分布：两组分布重叠区越窄，估计越不可靠。
• 截断权重：把权重上限设为 10 或 20。
• 修剪样本：丢弃得分超出 [0.05, 0.95] 的样本。

def check_overlap(e, T, lo=0.05, hi=0.95):
    in_range = (e >= lo) & (e <= hi)
    print(f"重叠区内样本占比：{in_range.mean():.3f}")
    print(f"处理组得分均值 {e[T==1].mean():.3f}，对照组 {e[T==0].mean():.3f}")

双重稳健估计与双机器学习

IPW 依赖倾向得分模型正确，回归调整依赖结果模型正确。双重稳健（Doubly Robust） 方法只要两者之一正确，估计就一致。

AIPW 估计量

增广逆概率加权（AIPW）：

τ̂_AIPW = (1/n) Σ [ μ_1(X_i) - μ_0(X_i)
                    + T_i(Y_i - μ_1(X_i))/e(X_i)
                    - (1-T_i)(Y_i - μ_0(X_i))/(1-e(X_i)) ]

其中 μ_1, μ_0 是两个结果模型（处理组与对照组各自的条件期望）。前半部分是回归调整，后半部分用 IPW 修正回归的残差。

双重稳健的含义：e(X) 对而 μ 错 → 一致；μ 对而 e(X) 错 → 一致；两者都对 → 更高效。

from sklearn.ensemble import GradientBoostingRegressor, GradientBoostingClassifier
import numpy as np

def aipw_ate(X, T, Y):
    # 结果模型
    m1 = GradientBoostingRegressor().fit(X[T == 1], Y[T == 1])
    m0 = GradientBoostingRegressor().fit(X[T == 0], Y[T == 0])
    mu1 = m1.predict(X)
    mu0 = m0.predict(X)
    # 倾向得分模型
    ps = GradientBoostingClassifier().fit(X, T)
    e = np.clip(ps.predict_proba(X)[:, 1], 0.01, 0.99)

    tau = (mu1 - mu0
           + T * (Y - mu1) / e
           - (1 - T) * (Y - mu0) / (1 - e))
    # 影响函数给出标准误
    se = tau.std(ddof=1) / np.sqrt(len(tau))
    return tau.mean(), se

双机器学习

双机器学习（Double Machine Learning, DML） 用交叉拟合（cross-fitting）消除过拟合带来的偏差：

1. 把数据分成 K 折。
2. 对每折，用其余 K-1 折训练 μ 与 e 模型。
3. 用留出折的预测做残差回归，得到效应估计。

交叉拟合的关键作用是防止用同一份数据既训练模型又估计效应导致的偏差。DML 是当前观察性因果估计的主流方法之一。

from sklearn.model_selection import KFold

def dml_ate(X, T, Y, n_splits=5):
    kf = KFold(n_splits=n_splits, shuffle=True, random_state=42)
    res_y = np.zeros(len(Y))
    res_t = np.zeros(len(T))
    for tr, te in kf.split(X):
        m_y = GradientBoostingRegressor().fit(X[tr], Y[tr])
        m_t = GradientBoostingClassifier().fit(X[tr], T[tr])
        res_y[te] = Y[te] - m_y.predict(X[te])
        res_t[te] = T[te] - m_t.predict_proba(X[te])[:, 1]
    # 残差对残差回归，斜率即效应
    tau = np.sum(res_t * res_y) / np.sum(res_t * res_t)
    return tau

工具变量与断点回归

有些混淆变量无法观测，只能用准实验设计绕开。

工具变量

工具变量（Instrumental Variable, IV） 是与处理相关、但只通过处理影响结果的变量。

两个条件：

• 相关性：Z 与 T 强相关。
• 外生性（排他性）：Z 只能通过 T 影响 Y，且与未观测混淆无关。

经典例子：用「距离最近大学的远近」作为「受教育年限」的工具变量。距离影响上学决策，但不直接影响收入（除通过教育）。

两阶段最小二乘（2SLS）：

第一阶段：T = α + βZ + ε        → 得到 T 的预测值 T̂
第二阶段：Y = γ + τ T̂ + δ       → τ 即因果效应

# 用 linearmodels 做 2SLS
from linearmodels.iv import IV2SLS
import pandas as pd

df = pd.DataFrame({"y": Y, "t": T, "z": Z, "x": X[:, 0]})
res = IV2SLS.from_formula("y ~ 1 + x + [t ~ z]", data=df).fit()
print(res.summary)

弱工具变量是最大风险：第一阶段 F 统计量低于 10 时，IV 估计会有严重偏差。

断点回归

断点回归（Regression Discontinuity, RDD） 利用一个阈值规则：阈值两侧的个体几乎相同，唯一差别是是否接受处理。

典型场景：奖学金按分数线发放，分数刚好在线上和线下的学生几乎无差别，比较他们的后续表现即可估计奖学金的因果效应。

Y = α + τ·1(X ≥ c) + f(X - c) + ε

f(X-c) 是断点两侧的平滑趋势，τ 是断点处的跳跃——即因果效应。

import statsmodels.formula.api as smf
import numpy as np

def rdd_estimate(X, Y, cutoff):
    df = {"y": Y, "x": X - cutoff,
          "treated": (X >= cutoff).astype(int)}
    import pandas as pd
    d = pd.DataFrame(df)
    # 局部线性回归，两侧各自拟合
    d["x_pos"] = np.where(d["treated"] == 1, d["x"], 0)
    d["x_neg"] = np.where(d["treated"] == 0, d["x"], 0)
    m = smf.ols("y ~ treated + x_pos + x_neg", data=d).fit()
    return m.params["treated"], m.bse["treated"]

关键假设是个体无法精确操纵阈值变量。若学生能花钱改分数，断点就失效了——可用 McCrary 密度检验排查。

双重差分与合成控制

面板数据提供了另一个识别来源：时间维度上的变化。

双重差分

双重差分（Difference-in-Differences, DID） 用「处理组前后差」减去「对照组前后差」：

τ̂_DID = (Ȳ_treat,post - Ȳ_treat,pre) - (Ȳ_ctrl,post - Ȳ_ctrl,pre)

对照组承担了「如果没有处理，处理组会怎么变化」的角色，从而扣掉时间趋势。

import numpy as np

def did_estimate(y, treat, post):
    """y: 结果  treat: 是否处理组  post: 是否处理后时期"""
    return ((y[(treat == 1) & (post == 1)].mean()
             - y[(treat == 1) & (post == 0)].mean())
            - (y[(treat == 0) & (post == 1)].mean()
               - y[(treat == 0) & (post == 0)].mean()))

核心假设是平行趋势：处理前两组的时间趋势一致。检验方法是画处理前的趋势图，或做事件研究（event study）看处理前的系数是否显著异于零。

合成控制

当对照组难以找到（如「某城市实施了政策」），合成控制法（Synthetic Control） 用多个未处理单元的加权组合，构造一个反事实的「合成对照组」：

目标：min_w || X_treat - Σ w_j X_j ||²
约束：Σ w_j = 1, w_j ≥ 0

权重 w 让合成单元在处理前的路径尽可能贴近处理单元。处理后两者的差距即效应。

方法选型

因果机器学习与异质效应

传统因果方法估计的是平均效应。业务真正想要的是「对谁有效」——异质处理效应（CATE）。

元学习器

元学习器（Meta-Learners） 用现成的机器学习模型拼装 CATE 估计：

• S-Learner：把 T 当作一个普通特征，训练单一模型，CATE 为同一 X 下 T=1 与 T=0 的预测差。简单，但 T 容易被淹没。
• T-Learner：对处理组与对照组各训一个模型，CATE 为两者预测差。样本少时易过拟合。
• X-Learner：用对照组模型估计处理组的反事实，反之亦然，再对两个效应估计做加权。适合处理组与对照组样本量极不平衡的场景。
• R-Learner：用残差化方法，与 DML 同源，稳健性好。

def t_learner(X, T, Y):
    m1 = GradientBoostingRegressor().fit(X[T == 1], Y[T == 1])
    m0 = GradientBoostingRegressor().fit(X[T == 0], Y[T == 0])
    def cate(X_new):
        return m1.predict(X_new) - m0.predict(X_new)
    return cate

因果森林

因果森林（Causal Forest） 是随机森林的因果版本：分裂准则不是「预测更准」，而是让左右子节点的效应估计差异更大。它天然给出 CATE 的置信区间，是异质效应估计的标杆方法。

# econml：微软的因果机器学习库
from econml.dml import CausalForestDML

est = CausalForestDML(
    model_y=GradientBoostingRegressor(),
    model_t=GradientBoostingClassifier(),
    n_estimators=500,
    discrete_treatment=True,
)
est.fit(Y, T, X=X, W=W)              # X 是效应修饰变量，W 是纯混淆变量
cate = est.effect(X_test)
lo, hi = est.effect_interval(X_test, alpha=0.05)

注意 X 与 W 的区分：X 是想要看效应如何随其变化的变量（如用户分层特征），W 是只用来去混淆、不关心其效应异质性的变量。

增益模型与 Uplift

推荐与营销场景把 CATE 直接叫 Uplift，建模目标是「增量转化」而非「转化」：

Uplift 模型的价值在于识别后两类人群，省下预算并避免负效应。评估不能用 AUC，要用 Qini 曲线 或 AUUC。

def qini_curve(uplift_score, T, Y):
    order = np.argsort(-uplift_score)          # 按 uplift 降序
    T, Y = T[order], Y[order]
    xs, gains, cum_t, cum_c = [], [], 0, 0
    for i in range(len(Y)):
        cum_t += Y[i] if T[i] == 1 else 0
        cum_c += Y[i] if T[i] == 0 else 0
        n_t = T[:i + 1].sum()
        n_c = (i + 1) - n_t
        if n_c > 0:
            gains.append(cum_t - cum_c * n_t / n_c)
            xs.append((i + 1) / len(Y))
    return np.array(xs), np.array(gains)

因果推断在推荐与风控的应用

因果方法在工业界有明确的高价值场景。

推荐系统的偏差

推荐数据天然是有偏的：模型只推荐了它认为好的物品，用户只反馈了被推荐的物品。三类偏差：

• 位置偏差：排在前面的物品点击率高，与物品本身质量无关。
• 曝光偏差：未被曝光的物品完全没有反馈。
• 选择偏差：用户自选择进入某些页面。

逆倾向得分（IPS） 是去偏的通用手段：用曝光概率的倒数给样本加权，把有偏的观测分布校正回均匀分布。

def ips_estimate(rewards, propensities):
    """rewards: 观测奖励  propensities: 记录下来的曝光概率"""
    w = 1.0 / np.clip(propensities, 1e-6, 1.0)
    return np.sum(w * rewards) / np.sum(w)

更稳健的变体是 SNIPS（自归一化 IPS）和 Doubly Robust 版本，能显著降低方差。

风控与反欺诈

风控场景的因果问题：

• 策略评估：新的拦截规则上线后，坏账率下降多少？不能只看拦截量，因为「被拦截的用户本来可能也不会违约」。
• 中介分析：某特征通过什么路径影响违约。
• 反事实推断：若没拦截这个用户，会发生什么。

风控里随机试验代价高（拦截了就可能损失好客户），因此观察性因果方法（PSM、DML） 使用更广。

从预测到决策的转变

把因果纳入决策的标准流程：

1. 明确决策变量：是要「预测谁会转化」还是「决定给谁发券」。
2. 画出 DAG：识别混淆、中介、碰撞，确定要控制的变量。
3. 选择识别策略：能随机化就随机化，不能则用 PSM/DML/IV/DID。
4. 估计 CATE：用因果森林或 Uplift 模型找到响应人群。
5. 决策与评估：按增量收益而非绝对指标分配预算。

生产踩坑清单

• 把所有特征都塞进模型：控制中介或碰撞变量会引入偏差。先画 DAG 再选变量。
• A/B 测试中途偷看：假阳性率失控。用序贯检验或固定样本量。
• 忽略 SRM：分流比例异常时结果不可信，上线前必须检查。
• 倾向得分极端值不处理：权重爆炸导致方差失控。做截断与重叠区修剪。
• 单模型做双重稳健：只用 IPW 或只用回归，放弃了一致性保障。用 AIPW 或 DML。
• DID 不做平行趋势检验：处理前趋势不一致时估计完全无效。
• 弱工具变量：第一阶段 F 值太低，IV 估计偏差比朴素回归还大。
• Uplift 用 AUC 评估：AUC 衡量的是预测排序，不是增量排序。必须用 Qini/AUUC。
• 忽略 SUTVA 溢出：社交与双边市场里处理会传染，需聚类随机化。
• 把观测关联当作因果写进报告：这是最常见也最昂贵的错误。

总结

因果推断的工程主线是「先声明假设，再选择识别策略，最后才估计效应」：潜在结果框架给出 ATE/ATT/CATE 的定义与三大识别假设，DAG 与后门准则告诉你该控制哪些变量，随机试验与 A/B 测试是黄金标准，观察性数据则靠倾向得分、双重稳健与双机器学习逼近随机化的效果，无法观测混淆时用工具变量、断点回归、双重差分与合成控制等准实验设计绕开。业务侧，因果森林与 Uplift 模型把平均效应细化到人群级别，让预算投向真正产生增量的人。一句话：预测模型告诉你世界会怎样，因果模型告诉你改变什么能让世界变得更好——后者才是决策系统真正需要的能力。