标签：得分 匹配 PSM python matches 干预 df 实操 groups
前言
大家好，我是顾先生，PSM倾向性得分匹配法的Python代码实操终于来啦！
对于PSM原理不太熟悉的同学可以看看前一篇文章：PSM倾向得分匹配法【上篇：理论篇】
目前网上PSM实操的相关文章都是R语言、SPSS和STATA实现的，少数Python版本代码不全，可读性有限（有些甚至要钱。。。）
所以我想出一版可读性强、能迅速复用的Python版本PSM，让各位同学既能看懂又能快速上手。
本次Python代码实操主要参考了psmatching的源码，并做了一定的修改，地址在文末参考资料。
这次我把每段代码都做了保姆级的注释，相信每位同学都能理解到位，当然肯定有注释不对的地方，也欢迎后台私信我。
本文的代码和数据集可关注我的公众号“顾先生聊数据”，后台发送“psm”后领取~
数据集介绍
为了更好地进行演示，我现编了一个数据集，该数据集以电商场景为基础，判断给用户发送优惠券PUSH是否会影响到用户。
数据集的类别主要由事实层标签（年龄、性别等）和行为层标签（最近一次购买diff、之前使用优惠券情况等）。
重申一下！数据集是我randbetween现编的，不用太较真具体内容。
完整代码
下面的代码我做了尽可能详细的注释，复用时需要修改的地方我也做了标注，如有不合理的地方欢迎后台私聊我哦~
安装psmatching包。
!pip install psmatching
import psmatching.match as psm
import pytest
import pandas as pd
import numpy as np
from psmatching.utilities import *
import statsmodels.api as sm
path及model设置。
#地址
path = "./data/psm/psm_gxslsj_data.csv"
#model由干预项和其他类别标签组成，形式为"干预项~类别特征+列别特征。。。"
model = "PUSH ~ AGE + SEX + VIP_LEVEL + LASTDAY_BUY_DIFF + PREFER_TYPE + LOGTIME_PREFER + USE_COUPON_BEFORE + ACTIVE_LEVEL"
#想要几个匹配项，如k=3，那一个push=1的用户就会匹配三个push=0的近似用户
k = "3"
m = psm.PSMatch(path, model, k)
获得倾向性匹配得分。
df = pd.read_csv(path)
#将用户ID作为数据的新索引
df = df.set_index("ID")
print("\n计算倾向性匹配得分 ...", end = " ")
#利用逻辑回归框架计算倾向得分，即广义线性估计 + 二项式Binomial
glm_binom = sm.formula.glm(formula = model, data = df, family = sm.families.Binomial())
#拟合拟合给定family的广义线性模型
#https://www.w3cschool.cn/doc_statsmodels/statsmodels-generated-statsmodels-genmod-generalized_linear_model-glm-fit.html?lang=en
result = glm_binom.fit()
# 输出回归分析的摘要
# print(result.summary)
propensity_scores = result.fittedvalues
print("\n计算完成!")
#将倾向性匹配得分写入data
df["PROPENSITY"] = propensity_scores
计算倾向性匹配得分 …
计算完成!
groups是干预项，propensity是倾向性匹配得分，这里要分开干预与非干预，且确保n1<n2。
注意：将PUSH替换成自己的干预项。
groups = df.PUSH
propensity = df.PROPENSITY
#把干预项替换成True和False
groups = groups == groups.unique()[1]
n = len(groups)
#计算True和False的数量
n1 = groups[groups==1].sum()
n2 = n-n1
g1, g2 = propensity[groups==1], propensity[groups==0]
#确保n2>n1，,少的匹配多的，否则交换下
if n1 > n2:
n1, n2, g1, g2 = n2, n1, g2, g1
随机排序干预组，减少原始排序的影响
m_order = list(np.random.permutation(groups[groups==1].index))
根据倾向评分差异将干预组与对照组进行匹配
注意：caliper = None可以替换成自己想要的精度
matches = {}
k = int(k)
print("\n给每个干预组匹配 [" + str(k) + "] 个对照组 ... ", end = " ")
for m in m_order:
# 计算所有倾向得分差异,这里用了最粗暴的绝对值
# 将propensity[groups==1]分别拿出来，每一个都与所有的propensity[groups==0]相减
dist = abs(g1[m]-g2)
array = np.array(dist)
#如果无放回地匹配，最后会出现要选取3个匹配对象，但是只有一个候选对照组的错误，故进行判断
if k < len(array):
# 在array里面选择K个最小的数字，并转换成列表
k_smallest = np.partition(array, k)[:k].tolist()
# 用卡尺做判断
caliper = None
if caliper:
caliper = float(caliper)
# 判断k_smallest是否在定义的卡尺范围
keep_diffs = [i for i in k_smallest if i <= caliper]
keep_ids = np.array(dist[dist.isin(keep_diffs)].index)
else:
# 如果不用标尺判断，那就直接上k_smallest了
keep_ids = np.array(dist[dist.isin(k_smallest)].index)
# 如果keep_ids比要匹配的数量多，那随机选择下，如要少，通过补NA配平数量
if len(keep_ids) > k:
matches[m] = list(np.random.choice(keep_ids, k, replace=False))
elif len(keep_ids) < k:
while len(matches[m]) <= k:
matches[m].append("NA")
else:
matches[m] = keep_ids.tolist()
# 判断 replace 是否放回
replace = False
if not replace:
g2 = g2.drop(matches[m])
print("\n匹配完成!")
给每个干预组匹配 [3] 个对照组 …
匹配完成!
这里提一下抽取规则：
抽取规则分无放回匹配和有放回匹配两种。如对照组个体不多，可选择有放回匹配，但重复选择对照组的样本进行匹配，会降低最终匹配样本的样本量，估计精度下降。大家可依据自身样本情况修改代码。
再提一下对应关系：
对应关系分一对一和一对多两种。一对一匹配样本少，估计方差较大，且每个匹配都是最近的，故偏差小； 一对多匹配样本多，估计精度高，但由于干预组个体匹配多个，故相似度降低，偏差会增加。本文使用的是一对多，大家可依据自身样本情况修改代码。
将匹配完成的结果合并起来
matches = pd.DataFrame.from_dict(matches, orient="index")
matches = matches.reset_index()
column_names = {}
column_names["index"] = "干预组"
for i in range(k):
column_names[i] = str("匹配对照组_" + str(i+1))
matches = matches.rename(columns = column_names)
从全部data中提取干预组和匹配对照的数据。
这里直接调用get_matched_data，注意输入的matches是匹配结果，df是全部数据。
matched_data = get_matched_data(matches, df)
将结果数据写入到文档
注意：可以自己更改地址
print("将倾向性匹配得分写入到文档 ...", end = " ")
save_file = path.split(".")[0] + "_倾向性匹配得分.csv"
df.to_csv(save_file, index = False)
print("完成!")
print("将匹配结果写入到文档 ...", end = " ")
save_file = path.split(".")[0] + "_匹配结果.csv"
matches.to_csv(save_file, index = False)
print("完成!")
将倾向性匹配得分写入到文档 … 完成!
将匹配结果写入到文档 … 完成!
对匹配结果进行变量检验
#注意：将PUSH替换成自己的干预项
variables = df.columns.tolist()[0:-2]
results = {}
print("开始评估匹配 ...")
#注意：将PUSH替换成自己的干预项
for var in variables:
# 制作用于卡方检验的频率计数交叉表
crosstable = pd.crosstab(df['PUSH'],df[var])
if len(df[var].unique().tolist()) <= 2:
# 计算 2x2 表的卡方统计量、df 和 p 值
p_val = calc_chi2_2x2(crosstable)[1]
else:
# 计算 2x2 表的卡方统计量、df 和 p 值
p_val = calc_chi2_2xC(crosstable)[1]
results[var] = p_val
print("\t" + var + '(' + str(p_val) + ')', end = "")
if p_val < 0.05:
print(": 未通过")
else:
print(": 通过")
if True in [i < 0.05 for i in results.values()]:
print("\n变量未全部通过匹配")
else:
print("\n变量全部通过匹配")
开始评估匹配 … AGE(0.9904): 通过
SEX(0.6688): 通过
VIP_LEVEL(0.0089): 未通过
LASTDAY_BUY_DIFF(0.5134): 通过
PREFER_TYPE(0.7107): 通过
LOGTIME_PREFER(0.2442): 通过
USE_COUPON_BEFORE(0.2961): 通过
ACTIVE_LEVEL(0.7934): 通过
变量未全部通过匹配
这里提一句，为啥P值<0.05未通过，P值>0.05反而通过呢？
这是因为这里的"通过"与假设检验里面的"通过"不太一样。
我们做PSM，是为了在对照组中找到与干预组类似的个体。
那单个个体中的各个变量，应该和干预项"PUSH"是相互独立的。
所以上面的检验中，VIP_LEVEL(0.0089)显示未通过。
后记
ok，PSM的python实操就介绍到这里了，相信看到这里的同学肯定已经学会了如何使用PSM了吧~
本文的代码和数据集可关注我的公众号“顾先生聊数据”，后台发送“psm”后领取~
如有不清楚的地方，也欢迎后台私信我一起讨论。
参考资料：
⼩洛学习群——因果推断（第⼀期）
陈强 高级计量经济学及stata应用（第二版）
https://github.com/rickydangc/psmatching
https://mattzheng.blog.csdn.net/
https://www.jianshu.com/p/34dd19ebe475
标签：得分,匹配,PSM,python,matches,干预,df,实操,groups 来源： https://blog.csdn.net/weixin_42847656/article/details/122355421
本站声明： 1. iCode9 技术分享网（下文简称本站）提供的所有内容，仅供技术学习、探讨和分享； 2. 关于本站的所有留言、评论、转载及引用，纯属内容发起人的个人观点，与本站观点和立场无关； 3. 关于本站的所有言论和文字，纯属内容发起人的个人观点，与本站观点和立场无关； 4. 本站文章均是网友提供，不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属；如您发现该文章侵犯了您的权益，可联系我们第一时间进行删除； 5. 本站为非盈利性的个人网站，所有内容不会用来进行牟利，也不会利用任何形式的广告来间接获益，纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。