异质性处理效应：为什么平均效应可能误导决策

问题背景

想象一下，一家制药公司开发了一种新的降压药。在临床试验中，该药物显示出显著的平均降压效果。但如果深入研究，可能会发现：

• 对年轻患者效果显著
• 对老年患者效果微弱甚至无效
• 对特定基因型的患者有副作用

如果仅依据平均效应做出"该药物有效"的结论并推广至所有患者，将导致资源浪费和潜在的健康风险。这就是研究异质性处理效应的现实意义。

I. 理论基础与核心概念

1.1 平均处理效应的局限性

**平均处理效应(ATE)**定义为：

$$ATE = E[Y(1) - Y(0)]$$

其中$Y(1)$和$Y(0)$分别表示接受处理和不接受处理的潜在结果。

ATE的局限性主要体现在：

1.2 异质性处理效应的定义

**异质性处理效应(HTE)**指处理效应在不同个体或子群体之间的系统性差异。数学上表示为：

$$\tau_i = Y_i(1) - Y_i(0)$$

$$\tau(x) = E[Y(1) - Y(0) | X=x]$$

其中$X$为协变量，$\tau(x)$为条件平均处理效应(Conditional Average Treatment Effect, CATE)。

1.3 异质性的来源

处理效应的异质性可能来源于多种机制：

II. 实例分析：教育干预项目的异质性效应

2.1 案例背景

考虑一个教育科技公司开展的在线辅导项目，旨在提高学生的数学成绩。公司进行了随机对照试验(RCT)，收集了1000名学生的数据，包括：

• 预处理特征：基础数学水平、学习动机、家庭支持、学校资源等
• 处理变量：是否接受在线辅导(随机分配)
• 结果变量：期末数学考试成绩

2.2 平均效应分析

首先，我们计算项目的平均处理效应：

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
from sklearn.ensemble import RandomForestRegressor
from sklearn.model_selection import train_test_split
from sklearn.metrics import mean_squared_error
import warnings
warnings.filterwarnings('ignore')

# 设置随机种子确保结果可重现
np.random.seed(42)

# 生成模拟数据
n = 1000  # 样本量

# 生成学生特征
base_math = np.random.normal(60, 10, n)  # 基础数学水平
motivation = np.random.normal(50, 15, n)  # 学习动机
family_support = np.random.binomial(1, 0.4, n)  # 家庭支持(二分)
school_resources = np.random.normal(50, 12, n)  # 学校资源

# 随机分配处理组和对照组
treatment = np.random.binomial(1, 0.5, n)

# 生成异质性处理效应
# 处理效应取决于基础数学水平和学习动机
treatment_effect = (0.5 * (base_math - 60) + 
                   0.3 * (motivation - 50) + 
                   2 * family_support +
                   np.random.normal(0, 2, n))

# 生成潜在结果
Y0 = (0.8 * base_math + 
      0.5 * motivation + 
      3 * family_support +
      0.4 * school_resources +
      np.random.normal(0, 5, n))

Y1 = Y0 + treatment_effect

# 观察结果
Y_obs = Y0 * (1 - treatment) + Y1 * treatment

# 创建数据集
data = pd.DataFrame({
    'base_math': base_math,
    'motivation': motivation,
    'family_support': family_support,
    'school_resources': school_resources,
    'treatment': treatment,
    'Y_obs': Y_obs,
    'treatment_effect': treatment_effect  # 在实际中不可观测
})

print("数据生成完成！")
print(f"样本量: {len(data)}")
print(f"处理组: {treatment.sum()}, 对照组: {len(data) - treatment.sum()}")

现在计算平均处理效应：

# 计算平均处理效应(ATE)
ate = data[data['treatment'] == 1]['Y_obs'].mean() - data[data['treatment'] == 0]['Y_obs'].mean()
true_ate = data['treatment_effect'].mean()

print("=== 平均处理效应分析 ===")
print(f"真实ATE: {true_ate:.2f}")
print(f"估计ATE: {ate:.2f}")
print(f"ATE 95% 置信区间: [{ate - 1.96 * data['Y_obs'].std()/np.sqrt(n):.2f}, {ate + 1.96 * data['Y_obs'].std()/np.sqrt(n):.2f}]")

# 分组统计
treatment_group = data[data['treatment'] == 1]
control_group = data[data['treatment'] == 0]

print(f"\n处理组平均成绩: {treatment_group['Y_obs'].mean():.2f}")
print(f"对照组平均成绩: {control_group['Y_obs'].mean():.2f}")

2.3 异质性探索分析

虽然平均效应显示项目整体有效，但我们需要深入探究效应如何随学生特征变化：

# 按基础数学水平分组分析
data['math_level'] = pd.cut(data['base_math'], 
                           bins=[0, 50, 70, 100], 
                           labels=['低基础', '中基础', '高基础'])

# 按学习动机分组
data['motivation_level'] = pd.cut(data['motivation'], 
                                 bins=[0, 40, 60, 100], 
                                 labels=['低动机', '中动机', '高动机'])

# 计算各子组的处理效应
subgroup_effects = []

for math_lev in data['math_level'].unique():
    for motiv_lev in data['motivation_level'].unique():
        subgroup = data[(data['math_level'] == math_lev) & 
                       (data['motivation_level'] == motiv_lev)]
        if len(subgroup) > 0:
            effect = (subgroup[subgroup['treatment'] == 1]['Y_obs'].mean() - 
                     subgroup[subgroup['treatment'] == 0]['Y_obs'].mean())
            true_effect = subgroup['treatment_effect'].mean()
            subgroup_effects.append({
                '数学基础': math_lev,
                '学习动机': motiv_lev,
                '样本量': len(subgroup),
                '估计处理效应': effect,
                '真实处理效应': true_effect
            })

subgroup_df = pd.DataFrame(subgroup_effects)
print("\n=== 各子组处理效应分析 ===")
print(subgroup_df.round(2))

2.4 可视化异质性效应

# 绘制处理效应的异质性
fig, axes = plt.subplots(2, 2, figsize=(15, 12))

# 1. 按数学基础分组的处理效应
math_effects = data.groupby('math_level').apply(
    lambda x: x[x['treatment'] == 1]['Y_obs'].mean() - x[x['treatment'] == 0]['Y_obs'].mean()
).reset_index()
math_effects.columns = ['数学基础', '处理效应']

axes[0, 0].bar(math_effects['数学基础'], math_effects['处理效应'], color=['#ff9999', '#66b3ff', '#99ff99'])
axes[0, 0].set_title('按数学基础分组的处理效应', fontsize=14)
axes[0, 0].set_ylabel('处理效应')

# 2. 按学习动机分组的处理效应
motiv_effects = data.groupby('motivation_level').apply(
    lambda x: x[x['treatment'] == 1]['Y_obs'].mean() - x[x['treatment'] == 0]['Y_obs'].mean()
).reset_index()
motiv_effects.columns = ['学习动机', '处理效应']

axes[0, 1].bar(motiv_effects['学习动机'], motiv_effects['处理效应'], color=['#ff9999', '#66b3ff', '#99ff99'])
axes[0, 1].set_title('按学习动机分组的处理效应', fontsize=14)
axes[0, 1].set_ylabel('处理效应')

# 3. 家庭支持的影响
family_effects = data.groupby('family_support').apply(
    lambda x: x[x['treatment'] == 1]['Y_obs'].mean() - x[x['treatment'] == 0]['Y_obs'].mean()
).reset_index()
family_effects.columns = ['家庭支持', '处理效应']
family_effects['家庭支持'] = family_effects['家庭支持'].map({0: '无支持', 1: '有支持'})

axes[1, 0].bar(family_effects['家庭支持'], family_effects['处理效应'], color=['#ff9999', '#66b3ff'])
axes[1, 0].set_title('按家庭支持分组的处理效应', fontsize=14)
axes[1, 0].set_ylabel('处理效应')

# 4. 处理效应分布
axes[1, 1].hist(data['treatment_effect'], bins=30, alpha=0.7, color='skyblue', edgecolor='black')
axes[1, 1].axvline(x=true_ate, color='red', linestyle='--', linewidth=2, label=f'平均效应: {true_ate:.2f}')
axes[1, 1].set_title('个体处理效应分布', fontsize=14)
axes[1, 1].set_xlabel('个体处理效应')
axes[1, 1].set_ylabel('频数')
axes[1, 1].legend()

plt.tight_layout()
plt.show()

III. 异质性处理效应的估计方法

3.1 传统方法与机器学习方法对比

估计异质性处理效应的方法可以分为传统计量经济学方法和现代机器学习方法：

3.2 因果森林实现

下面我们使用因果森林来估计异质性处理效应：

from sklearn.ensemble import RandomForestRegressor
from sklearn.model_selection import cross_val_score

# 准备特征矩阵和结果变量
features = ['base_math', 'motivation', 'family_support', 'school_resources']
X = data[features]
W = data['treatment']  # 处理变量
Y = data['Y_obs']      # 观察结果

# 分割训练集和测试集
X_train, X_test, W_train, W_test, Y_train, Y_test = train_test_split(
    X, W, Y, test_size=0.3, random_state=42
)

print("=== 数据准备完成 ===")
print(f"训练集大小: {len(X_train)}")
print(f"测试集大小: {len(X_test)}")

# 实现简单的因果森林（基于R-learner思想）
class CausalForest:
    def __init__(self, n_estimators=100, max_depth=10, min_samples_leaf=10):
        self.n_estimators = n_estimators
        self.max_depth = max_depth
        self.min_samples_leaf = min_samples_leaf
        self.propensity_forest = None
        self.outcome_forest = None
        self.tau_forest = None
        
    def fit(self, X, W, Y):
        # 第一步：估计倾向得分
        self.propensity_forest = RandomForestRegressor(
            n_estimators=self.n_estimators,
            max_depth=self.max_depth,
            min_samples_leaf=self.min_samples_leaf,
            random_state=42
        )
        self.propensity_forest.fit(X, W)
        propensity_scores = self.propensity_forest.predict(X)
        
        # 第二步：估计结果模型
        self.outcome_forest = RandomForestRegressor(
            n_estimators=self.n_estimators,
            max_depth=self.max_depth,
            min_samples_leaf=self.min_samples_leaf,
            random_state=42
        )
        self.outcome_forest.fit(X, Y)
        outcome_predictions = self.outcome_forest.predict(X)
        
        # 第三步：计算残差并估计处理效应
        # Y_residual = Y - E[Y|X]
        Y_residual = Y - outcome_predictions
        
        # 构建处理效应估计的特征：添加倾向得分作为特征
        X_tau = X.copy()
        X_tau['propensity'] = propensity_scores
        X_tau['W'] = W
        
        # 训练处理效应森林
        self.tau_forest = RandomForestRegressor(
            n_estimators=self.n_estimators,
            max_depth=self.max_depth,
            min_samples_leaf=self.min_samples_leaf,
            random_state=42
        )
        self.tau_forest.fit(X_tau, Y_residual)
        
        return self
    
    def predict_ite(self, X):
        # 预测倾向得分
        propensity_scores = self.propensity_forest.predict(X)
        
        # 构建处理效应预测的特征
        X_tau = X.copy()
        X_tau['propensity'] = propensity_scores
        
        # 预测W=1时的处理效应
        X_tau_1 = X_tau.copy()
        X_tau_1['W'] = 1
        tau_1 = self.tau_forest.predict(X_tau_1)
        
        # 预测W=0时的处理效应
        X_tau_0 = X_tau.copy()
        X_tau_0['W'] = 0
        tau_0 = self.tau_forest.predict(X_tau_0)
        
        # 个体处理效应 = tau_1 - tau_0
        ite = tau_1 - tau_0
        
        return ite

# 训练因果森林模型
print("开始训练因果森林模型...")
cf_model = CausalForest(n_estimators=200, max_depth=15, min_samples_leaf=5)
cf_model.fit(X_train, W_train, Y_train)

# 预测个体处理效应
ite_pred = cf_model.predict_ite(X_test)

# 添加到测试集数据中
test_data = X_test.copy()
test_data['W'] = W_test.values
test_data['Y_obs'] = Y_test.values
test_data['ite_pred'] = ite_pred

# 由于我们有真实数据，可以计算预测精度
if 'treatment_effect' in data.columns:
    # 获取测试集的真实处理效应
    test_indices = X_test.index
    test_data['ite_true'] = data.loc[test_indices, 'treatment_effect'].values
    
    # 计算预测精度
    mse_ite = mean_squared_error(test_data['ite_true'], test_data['ite_pred'])
    correlation = np.corrcoef(test_data['ite_true'], test_data['ite_pred'])[0, 1]
    
    print(f"\n=== 模型评估结果 ===")
    print(f"个体处理效应预测MSE: {mse_ite:.4f}")
    print(f"预测值与真实值相关性: {correlation:.4f}")

3.3 模型评估与验证

# 模型性能可视化
fig, axes = plt.subplots(1, 2, figsize=(15, 6))

# 1. 预测vs真实处理效应散点图
if 'ite_true' in test_data.columns:
    axes[0].scatter(test_data['ite_true'], test_data['ite_pred'], alpha=0.6)
    axes[0].plot([test_data['ite_true'].min(), test_data['ite_true'].max()], 
                [test_data['ite_true'].min(), test_data['ite_true'].max()], 
                'r--', linewidth=2)
    axes[0].set_xlabel('真实处理效应')
    axes[0].set_ylabel('预测处理效应')
    axes[0].set_title('预测vs真实处理效应')
    
    # 添加性能指标文本
    axes[0].text(0.05, 0.95, f'相关性: {correlation:.3f}', 
                transform=axes[0].transAxes, fontsize=12,
                verticalalignment='top', bbox=dict(boxstyle='round', facecolor='white', alpha=0.8))

# 2. 特征重要性分析
feature_importance = cf_model.tau_forest.feature_importances_
feature_names = features + ['propensity', 'W']
importance_df = pd.DataFrame({
    'feature': feature_names,
    'importance': feature_importance
}).sort_values('importance', ascending=False)

axes[1].barh(importance_df['feature'], importance_df['importance'])
axes[1].set_xlabel('特征重要性')
axes[1].set_title('因果森林特征重要性')
axes[1].grid(axis='x', alpha=0.3)

plt.tight_layout()
plt.show()

# 按预测处理效应分组的实际效果验证
test_data['pred_effect_group'] = pd.qcut(test_data['ite_pred'], 5, 
                                        labels=['最低20%', '较低20%', '中间20%', '较高20%', '最高20%'])

# 计算各组的平均预测效应和实际效应
validation_results = []
for group in test_data['pred_effect_group'].unique():
    group_data = test_data[test_data['pred_effect_group'] == group]
    
    # 计算该组的平均预测处理效应
    pred_effect = group_data['ite_pred'].mean()
    
    # 计算该组的实际处理效应（在RCT中可通过分组比较得到）
    actual_effect = (group_data[group_data['W'] == 1]['Y_obs'].mean() - 
                    group_data[group_data['W'] == 0]['Y_obs'].mean())
    
    validation_results.append({
        '效应分组': group,
        '样本量': len(group_data),
        '平均预测效应': pred_effect,
        '实际观察效应': actual_effect
    })

validation_df = pd.DataFrame(validation_results)
print("\n=== 模型验证结果 ===")
print(validation_df.round(3))

IV. 决策优化与应用策略

4.1 基于异质性效应的资源分配

利用预测的个体处理效应，我们可以优化资源分配策略：

# 基于处理效应的资源分配优化
def optimize_allocation(data, budget_constraint=0.3):
    """
    基于预测处理效应优化处理分配
    
    参数:
    data: 包含预测处理效应的数据集
    budget_constraint: 处理比例约束 (0-1)
    
    返回:
    优化分配方案
    """
    # 按预测处理效应降序排列
    sorted_data = data.sort_values('ite_pred', ascending=False)
    
    # 根据预算约束选择处理对象
    n_treat = int(len(sorted_data) * budget_constraint)
    treatment_decision = [1] * n_treat + [0] * (len(sorted_data) - n_treat)
    sorted_data['optimal_treatment'] = treatment_decision
    
    # 计算不同分配策略的效果
    # 1. 随机分配
    random_treatment = np.random.choice([1, 0], size=len(data), 
                                      p=[budget_constraint, 1-budget_constraint])
    random_effect = data.loc[random_treatment == 1, 'ite_true'].mean() if 'ite_true' in data.columns else None
    
    # 2. 基于预测的最优分配
    optimal_effect = sorted_data[sorted_data['optimal_treatment'] == 1]['ite_true'].mean() if 'ite_true' in data.columns else None
    
    # 3. 实际随机分配（RCT）
    actual_effect = data[data['W'] == 1]['ite_true'].mean() if 'ite_true' in data.columns else None
    
    return {
        'random_allocation_effect': random_effect,
        'optimal_allocation_effect': optimal_effect,
        'actual_allocation_effect': actual_effect,
        'allocation_data': sorted_data
    }

# 应用优化算法
optimization_results = optimize_allocation(test_data, budget_constraint=0.4)

print("=== 资源分配优化结果 ===")
print(f"随机分配期望效果: {optimization_results['random_allocation_effect']:.3f}")
print(f"最优分配期望效果: {optimization_results['optimal_allocation_effect']:.3f}")
print(f"实际分配平均效果: {optimization_results['actual_allocation_effect']:.3f}")

improvement = (optimization_results['optimal_allocation_effect'] - 
               optimization_results['random_allocation_effect']) / optimization_results['random_allocation_effect'] * 100
print(f"优化分配提升: {improvement:.1f}%")

# 可视化分配策略比较
strategies_comparison = pd.DataFrame({
    '策略': ['随机分配', '实际RCT', '优化分配'],
    '期望效果': [
        optimization_results['random_allocation_effect'],
        optimization_results['actual_allocation_effect'], 
        optimization_results['optimal_allocation_effect']
    ]
})

plt.figure(figsize=(10, 6))
bars = plt.bar(strategies_comparison['策略'], strategies_comparison['期望效果'], 
               color=['lightblue', 'lightgreen', 'coral'])
plt.ylabel('期望处理效应')
plt.title('不同分配策略的效果比较')
plt.grid(axis='y', alpha=0.3)

# 在柱状图上添加数值标签
for bar, value in zip(bars, strategies_comparison['期望效果']):
    plt.text(bar.get_x() + bar.get_width()/2, bar.get_height() + 0.1, 
             f'{value:.2f}', ha='center', va='bottom')

plt.tight_layout()
plt.show()

4.2 个性化决策规则

基于异质性处理效应，我们可以制定个性化的决策规则：

# 制定个性化决策规则
def personalization_rules(data):
    """
    基于特征制定个性化决策规则
    """
    rules_data = data.copy()
    
    # 规则1: 基于数学基础和学习动机
    rules_data['rule1_treatment'] = ((rules_data['base_math'] < 65) & 
                                   (rules_data['motivation'] > 45)).astype(int)
    
    # 规则2: 基于预测处理效应的阈值
    effect_threshold = rules_data['ite_pred'].quantile(0.6)  # 选择前40%
    rules_data['rule2_treatment'] = (rules_data['ite_pred'] > effect_threshold).astype(int)
    
    # 规则3: 基于随机森林特征重要性的复合规则
    rules_data['rule3_treatment'] = ((rules_data['base_math'] < 70) & 
                                   (rules_data['family_support'] == 1) & 
                                   (rules_data['motivation'] > 40)).astype(int)
    
    # 评估各规则效果
    rules_evaluation = []
    
    for i in range(1, 4):
        rule_col = f'rule{i}_treatment'
        rule_effect = rules_data[rules_data[rule_col] == 1]['ite_true'].mean() if 'ite_true' in rules_data.columns else None
        coverage = rules_data[rule_col].mean()
        
        rules_evaluation.append({
            '规则': f'规则{i}',
            '覆盖率': coverage,
            '期望效果': rule_effect,
            '处理人数': rules_data[rule_col].sum()
        })
    
    return pd.DataFrame(rules_evaluation), rules_data

# 应用个性化规则
rules_eval_df, rules_data = personalization_rules(test_data)

print("\n=== 个性化决策规则评估 ===")
print(rules_eval_df.round(3))

# 可视化规则比较
fig, axes = plt.subplots(1, 2, figsize=(15, 6))

# 规则效果对比
axes[0].bar(rules_eval_df['规则'], rules_eval_df['期望效果'], color=['lightblue', 'lightgreen', 'coral'])
axes[0].set_ylabel('期望处理效应')
axes[0].set_title('不同决策规则的期望效果')
axes[0].grid(axis='y', alpha=0.3)

# 规则覆盖率对比
axes[1].bar(rules_eval_df['规则'], rules_eval_df['覆盖率'], color=['lightblue', 'lightgreen', 'coral'])
axes[1].set_ylabel('覆盖率')
axes[1].set_title('不同决策规则的覆盖率')
axes[1].grid(axis='y', alpha=0.3)

plt.tight_layout()
plt.show()

4.3 成本效益分析与敏感性测试

# 成本效益分析
def cost_benefit_analysis(data, treatment_cost=100, benefit_per_unit=10):
    """
    进行成本效益分析
    """
    analysis_data = data.copy()
    
    # 不同分配策略的成本效益
    strategies = ['random', 'optimal', 'rule1', 'rule2', 'rule3']
    results = []
    
    for strategy in strategies:
        if strategy == 'random':
            treatment_decision = np.random.choice([1, 0], size=len(analysis_data), p=[0.4, 0.6])
        elif strategy == 'optimal':
            # 基于预测效应选择前40%
            sorted_data = analysis_data.sort_values('ite_pred', ascending=False)
            treatment_decision = [1] * int(0.4 * len(sorted_data)) + [0] * (len(sorted_data) - int(0.4 * len(sorted_data)))
            treatment_decision = pd.Series(treatment_decision, index=sorted_data.index).loc[analysis_data.index]
        else:
            treatment_decision = analysis_data[f'{strategy}_treatment']
        
        # 计算总成本
        total_cost = treatment_cost * treatment_decision.sum()
        
        # 计算总收益（使用真实处理效应）
        total_benefit = benefit_per_unit * analysis_data.loc[treatment_decision == 1, 'ite_true'].sum() if 'ite_true' in analysis_data.columns else 0
        
        # 计算净收益
        net_benefit = total_benefit - total_cost
        
        # 计算投资回报率
        roi = (net_benefit / total_cost) * 100 if total_cost > 0 else 0
        
        results.append({
            '策略': strategy,
            '总成本': total_cost,
            '总收益': total_benefit,
            '净收益': net_benefit,
            '投资回报率(%)': roi,
            '处理人数': treatment_decision.sum()
        })
    
    return pd.DataFrame(results)

# 执行成本效益分析
cba_results = cost_benefit_analysis(rules_data)

print("\n=== 成本效益分析 ===")
print(cba_results.round(2))

# 敏感性分析：不同成本参数的影响
def sensitivity_analysis(data, cost_range=[50, 100, 150, 200]):
    """
    对处理成本进行敏感性分析
    """
    sensitivity_results = []
    
    for cost in cost_range:
        cba = cost_benefit_analysis(data, treatment_cost=cost)
        cba['处理成本'] = cost
        sensitivity_results.append(cba)
    
    sensitivity_df = pd.concat(sensitivity_results, ignore_index=True)
    
    # 可视化敏感性分析
    plt.figure(figsize=(12, 8))
    
    for strategy in sensitivity_df['策略'].unique():
        strategy_data = sensitivity_df[sensitivity_df['策略'] == strategy]
        plt.plot(strategy_data['处理成本'], strategy_data['投资回报率(%)'], 
                marker='o', label=strategy, linewidth=2)
    
    plt.xlabel('处理成本')
    plt.ylabel('投资回报率 (%)')
    plt.title('处理成本敏感性分析')
    plt.legend()
    plt.grid(True, alpha=0.3)
    plt.show()
    
    return sensitivity_df

# 执行敏感性分析
sensitivity_df = sensitivity_analysis(rules_data)

V. 实践建议与注意事项

5.1 实施异质性分析的最佳实践

基于我们的分析和实验结果，我们总结出以下实施异质性处理效应分析的最佳实践：

5.2 常见陷阱与解决方案

在实践中，异质性处理效应分析可能遇到多种挑战：

# 常见问题诊断函数
def diagnose_hte_problems(data, model, features):
    """
    诊断异质性处理效应分析的常见问题
    """
    diagnostics = {}
    
    # 1. 检查特征平衡
    treatment_group = data[data['treatment'] == 1]
    control_group = data[data['treatment'] == 0]
    
    balance_stats = []
    for feature in features:
        t_mean = treatment_group[feature].mean()
        c_mean = control_group[feature].mean()
        t_std = treatment_group[feature].std()
        c_std = control_group[feature].std()
        
        # 标准化均值差异
        smd = abs(t_mean - c_mean) / np.sqrt((t_std**2 + c_std**2) / 2)
        
        balance_stats.append({
            '特征': feature,
            '处理组均值': t_mean,
            '对照组均值': c_mean,
            '标准化差异': smd
        })
    
    balance_df = pd.DataFrame(balance_stats)
    diagnostics['特征平衡'] = balance_df
    
    # 2. 检查预测效应的分布
    diagnostics['效应分布'] = {
        '均值': data['ite_pred'].mean(),
        '标准差': data['ite_pred'].std(),
        '偏度': data['ite_pred'].skew(),
        '峰度': data['ite_pred'].kurtosis()
    }
    
    # 3. 检查极端预测值
    q1 = data['ite_pred'].quantile(0.01)
    q99 = data['ite_pred'].quantile(0.99)
    extreme_low = data[data['ite_pred'] < q1]
    extreme_high = data[data['ite_pred'] > q99]
    
    diagnostics['极端预测'] = {
        '极端低值数量': len(extreme_low),
        '极端高值数量': len(extreme_high),
        '极端值比例': (len(extreme_low) + len(extreme_high)) / len(data)
    }
    
    return diagnostics

# 执行诊断
if 'ite_pred' in data.columns:
    hte_diagnostics = diagnose_hte_problems(data, cf_model, features)
    
    print("=== 异质性分析诊断结果 ===")
    print("1. 特征平衡检查:")
    print(hte_diagnostics['特征平衡'].round(3))
    
    print("\n2. 效应分布检查:")
    for key, value in hte_diagnostics['效应分布'].items():
        print(f"   {key}: {value:.3f}")
    
    print("\n3. 极端预测检查:")
    for key, value in hte_diagnostics['极端预测'].items():
        print(f"   {key}: {value}")

5.3 伦理考量与公平性

在应用异质性处理效应分析时，必须考虑伦理和公平性问题：

# 公平性评估
def fairness_assessment(data, sensitive_features):
    """
    评估异质性处理效应在不同敏感群体间的公平性
    """
    fairness_results = {}
    
    for feature in sensitive_features:
        # 检查敏感特征的处理效应差异
        unique_values = data[feature].unique()
        
        group_effects = []
        for value in unique_values:
            group_data = data[data[feature] == value]
            group_effect = group_data['ite_pred'].mean()
            group_size = len(group_data)
            
            group_effects.append({
                '群体': f"{feature}={value}",
                '平均预测效应': group_effect,
                '样本量': group_size
            })
        
        group_df = pd.DataFrame(group_effects)
        
        # 计算群体间最大差异
        max_effect = group_df['平均预测效应'].max()
        min_effect = group_df['平均预测效应'].min()
        effect_range = max_effect - min_effect
        
        fairness_results[feature] = {
            '群体效果': group_df,
            '最大差异': effect_range,
            '是否存在显著差异': effect_range > group_df['平均预测效应'].std()
        }
    
    return fairness_results

# 假设家庭支持是敏感特征
sensitive_features = ['family_support']
fairness_results = fairness_assessment(data, sensitive_features)

print("\n=== 公平性评估 ===")
for feature, results in fairness_results.items():
    print(f"敏感特征: {feature}")
    print(results['群体效果'])
    print(f"最大群体间差异: {results['最大差异']:.3f}")
    print(f"是否存在潜在公平性问题: {results['是否存在显著差异']}")
    print()