倾向得分匹配模型中控制组重复匹配实现方法问询

Hey there! 针对你提到的想要让处理组每个观测值都能匹配到最接近的控制组样本，同时允许控制组样本被重复使用的需求，这在倾向得分匹配（PSM）里其实是个很常见的场景，对应的方法叫做有放回匹配（matching with replacement）。下面我从概念、具体实现步骤到代码示例给你详细拆解：

简单来说，这种匹配方式允许控制组样本被多次选中（相当于“放回”重新参与匹配池），每个处理组样本都能独立找到倾向得分最接近的控制组样本，不用受“控制组只能用一次”的限制。这种方法的优势是能最大化单个处理组样本的匹配精度，尤其适合控制组样本量不足，或者你更看重每个处理组样本匹配质量的场景。

1. 先计算倾向得分 #

首先得给每个样本计算倾向得分——也就是该样本被分配到处理组的概率，通常用逻辑回归拟合，也可以用随机森林、XGBoost等模型提升精度。这里用Python的sklearn举个基础例子：

from sklearn.linear_model import LogisticRegression
import pandas as pd

# 假设你的数据集data里，treat是处理标识（1=处理组，0=控制组），剩下的是协变量
X = data.drop('treat', axis=1)
y = data['treat']

# 拟合逻辑回归模型计算倾向得分
lr_model = LogisticRegression(max_iter=1000)
lr_model.fit(X, y)
data['propensity_score'] = lr_model.predict_proba(X)[:, 1]  # 取预测为1的概率

2. 执行有放回的最邻近匹配 #

这里有两种方式：手动实现逻辑，或者用现成的因果推断库简化操作。

方式一：手动实现（灵活可控）

手动写逻辑的好处是能完全自定义匹配规则，比如加卡尺限制、重复次数限制等：

# 拆分处理组和控制组
treated_group = data[data['treat'] == 1].reset_index(drop=True)
control_group = data[data['treat'] == 0].reset_index(drop=True)

matched_results = []
# 遍历每个处理组样本，找最接近的控制组样本（允许重复）
for idx, treated_row in treated_group.iterrows():
    # 计算当前处理组样本与所有控制组样本的倾向得分差
    control_group['score_diff'] = abs(control_group['propensity_score'] - treated_row['propensity_score'])
    # 找到得分差最小的控制组样本
    best_match = control_group.loc[control_group['score_diff'].idxmin()]
    # 记录配对信息
    matched_results.append({
        'treated_index': idx,
        'control_index': best_match.name,
        'treated_ps': treated_row['propensity_score'],
        'control_ps': best_match['propensity_score'],
        'score_diff': best_match['score_diff']
    })

# 转成DataFrame方便后续分析
matched_df = pd.DataFrame(matched_results)

如果想限制控制组样本的重复次数，比如最多被匹配3次，可以在每次匹配后给控制组样本加个计数，超过次数就跳过：

# 初始化控制组匹配次数
control_group['match_count'] = 0

matched_results = []
for idx, treated_row in treated_group.iterrows():
    # 只选匹配次数<3的控制组样本
    available_control = control_group[control_group['match_count'] < 3]
    if available_control.empty:
        print(f"处理组样本{idx}无可用控制组匹配")
        continue
    available_control['score_diff'] = abs(available_control['propensity_score'] - treated_row['propensity_score'])
    best_match = available_control.loc[available_control['score_diff'].idxmin()]
    # 更新匹配次数
    control_group.loc[best_match.name, 'match_count'] += 1
    matched_results.append({
        'treated_index': idx,
        'control_index': best_match.name,
        'treated_ps': treated_row['propensity_score'],
        'control_ps': best_match['propensity_score'],
        'score_diff': best_match['score_diff']
    })

方式二：用现成库简化操作

如果不想手动写逻辑，可以用causalml这类专门做因果推断的库，它们直接支持有放回匹配：

from causalml.matching import NearestNeighborMatch

# 初始化匹配器，replacement=True就是允许有放回
matcher = NearestNeighborMatch(replacement=True, caliper=0.05)  # caliper是卡尺，限制得分差不超过0.05
# 传入倾向得分和处理标识，得到匹配后的控制组索引
matched_control_indices = matcher.match(data['propensity_score'], data['treat'])

# 把匹配后的控制组数据提取出来
matched_control = data.loc[matched_control_indices[data['treat']==1].values]

3. 匹配后的平衡性校验 #

匹配完成后一定要检查协变量的平衡性，确保处理组和匹配后的控制组协变量分布没有显著差异，常用的指标是标准化均值差（SMD），一般SMD<0.1就算平衡良好：

def calculate_smd(df, treat_col, covariate):
    treated_mean = df[df[treat_col]==1][covariate].mean()
    control_mean = df[df[treat_col]==0][covariate].mean()
    treated_std = df[df[treat_col]==1][covariate].std()
    control_std = df[df[treat_col]==0][covariate].std()
    pooled_std = ((treated_std**2 + control_std**2)/2)**0.5
    return abs(treated_mean - control_mean)/pooled_std

# 遍历所有协变量计算SMD
for covar in X.columns:
    smd_value = calculate_smd(data, 'treat', covar)
    print(f"协变量{covar}的SMD值：{smd_value:.3f}")

• 避免过度重复：如果某个控制组样本被大量处理组样本匹配，会导致估计结果的方差被低估，建议设置最大重复次数或者卡尺来限制。
• 样本量考量：如果控制组样本量远大于处理组，有放回匹配的效果会很好；但如果控制组样本极少，即使有放回，匹配质量也可能不佳，这时候可以考虑结合倾向得分加权等方法。
• 模型选择：如果协变量和处理分配的关系是非线性的，用随机森林、XGBoost等模型计算倾向得分会比逻辑回归更准确。

内容的提问来源于stack exchange，提问作者d142c681