因果推断实战：基于Propensity得分预测实验组与对照组的效果差异

1次阅读

共计 2270 个字符，预计需要花费 6 分钟才能阅读完成。

在因果推断领域，Propensity 得分匹配（PSM）通过模拟随机实验的条件，解决了观察性研究中的选择偏差问题。然而，匹配后直接比较组间差异仍存在两个关键缺陷：

条件独立假设的局限性 ：即使经过 PSM 调整，未观测变量仍可能导致 $Y(1), Y(0) \perp!!!\perp T | X$ 不成立
模型误设风险 ：传统机器学习模型会混淆协变量 $X$ 与处理变量 $T$ 的关系，导致效应估计有偏

方法	核心思想	计算复杂度	适用场景
T-Learner	分别训练实验组 / 对照组子模型	$O(2n)$	处理效应异质性明显时
X-Learner	引入反事实预测校正	$O(3n)$	样本量不平衡时
R-Learner	基于残差拟合的因果森林	$O(n^2)$	高维协变量场景

from sklearn.linear_model import LogisticRegression
from sklearn.model_selection import train_test_split

# 类型注解增强可读性
def calculate_propensity_score(X: pd.DataFrame, 
                             treatment: pd.Series) -> np.ndarray:
    """
    计算倾向性得分
    :param X: 协变量矩阵 (n_samples, n_features)
    :param treatment: 处理变量 (1= 实验组, 0= 对照组)
    :return: propensity scores (n_samples,)
    """ps_model = LogisticRegression(penalty='l2', C=1.0)
    ps_model.fit(X, treatment)
    return ps_model.predict_proba(X)[:, 1]

from sklearn.ensemble import GradientBoostingRegressor

def t_learner_estimate(X_train: pd.DataFrame,
                      y_train: pd.Series,
                      treat_train: pd.Series) -> tuple:
    """
    训练双模型架构
    :return: (模型 1, 模型 0, 平均处理效应 ATE)
    """
    # 实验组模型
    model1 = GradientBoostingRegressor(n_estimators=100).fit(X_train[treat_train==1], y_train[treat_train==1])

    # 对照组模型
    model0 = GradientBoostingRegressor(n_estimators=100).fit(X_train[treat_train==0], y_train[treat_train==0])

    # 计算 ATE
    ate = np.mean(model1.predict(X_train) - model0.predict(X_train))
    return model1, model0, ate

import seaborn as sns

def plot_overlap(ps_score: np.ndarray, treatment: pd.Series):
    """可视化倾向得分分布"""
    plt.figure(figsize=(10,6))
    sns.kdeplot(ps_score[treatment==1], label='Treated', shade=True)
    sns.kdeplot(ps_score[treatment==0], label='Control', shade=True)
    plt.xlabel('Propensity Score')
    plt.legend()
    plt.show()

预筛选机制 ：
使用 Double LASSO 选择与结果和处理都相关的变量
计算协变量 $X_j$ 的重要性得分：
$$I_j = |\beta_{y,j}| \cdot |\beta_{t,j}|$$
正则化调整 ：
在 Propensity 模型中加入 ElasticNet 惩罚项：
$$\min_\beta \left[\frac{1}{n}\sum_{i=1}^n \mathcal{L}(t_i,x_i^T\beta) + \lambda(\alpha|\beta|_1 + \frac{1-\alpha}{2}|\beta|_2^2)\right]$$

import shap

def plot_cate_heterogeneity(model, X_eval):
    """用 SHAP 展示处理效应异质性"""
    explainer = shap.TreeExplainer(model)
    shap_values = explainer.shap_values(X_eval)
    shap.summary_plot(shap_values, X_eval)