相关文章推荐
绅士的匕首  ·  windows: pip install ...·  昨天    · 
奔跑的水桶  ·  group by ...·  1 年前    · 
活泼的芒果  ·  救救孩子吧，到现在还搞不懂@Componen ...·  1 年前    · 
朝气蓬勃的面包  ·  Helix editor 笔记 ...·  1 年前    · 
鼻子大的人字拖  ·  博士申请 | ...·  1 年前    · 
绅士的水煮肉  ·  从大熊猫数据栏中删除时区 - 问答 - ...·  2 年前    · 
Code  ›  因果推断相关python库学习 （一）
python
https://juejin.cn/post/7178461156489560124
活泼的小蝌蚪
2 年前
import numpy as np from matplotlib import pyplot as plt from sklearn.linear_model import LinearRegression from sklearn.model_selection import train_test_split import statsmodels.api as sm from xgboost import XGBRegressor import warnings from causalml.inference.meta import LRSRegressor from causalml.inference.meta import XGBTRegressor, MLPTRegressor from causalml.inference.meta import BaseXRegressor, BaseRRegressor, BaseSRegressor, BaseTRegressor from causalml. match import NearestNeighborMatch, MatchOptimizer, create_table_one from causalml.propensity import ElasticNetPropensityModel from causalml.dataset import * from causalml.metrics import * # Generate synthetic data using mode 1 y, X, treatment, tau, b, e = synthetic_data(mode= 1 , n= 10000 , p= 8 , sigma= 1.0 )

s-learner、t-learner、x-learner、r-learner 预估 

# Ready-to-use S-Learner using LinearRegression
learner_s = LRSRegressor()
# Calling the Base Learner class and feeding in XGB
learner_t = BaseTRegressor(learner=XGBRegressor())
# Calling the Base Learner class and feeding in XGB
learner_x = BaseXRegressor(XGBRegressor())
# Calling the Base Learner class and feeding in XGB
learner_r = BaseRRegressor(learner=XGBRegressor())
# 估计ATE(以t-learner为例，其他均相同) 
ate_t = learner_t.estimate_ate(X, treatment, y)
# 估计ITE(以t-learner为例，其他均相同)
# fit_predict先后执行fit -> predict两个过程
cate_t = learner_t.fit_predict(X=X, treatment=treatment, y=y)

特征重要性和指标评估


# 特征重要性排序
slearner.plot_importance(X=X, 
                         tau=cate_t, 
                         normalize=True, 
                         method='auto', 
                         features=feature_names)

# 增益曲线和计算auuc
df_preds = pd.DataFrame([cate_t.ravel(),
                         treatment.ravel(), 
                         y.ravel()], 
                        index=["T", "w", "y"]).T
plot_gain(df_preds, normalize=True)
print(auuc_score(df_preds))

二，psmatching


PSM 倾向性得分匹配法的Python代码实操


# 1. 安装
!pip install psmatching
# 2. 导入
import psmatching.match as psm
import pytest
import pandas as pd
import numpy as np
from psmatching.utilities import *
import statsmodels.api as sm
# 3. 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)
# 4. 倾向得分计算
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
# 5. 随机排序干预组，减少原始排序的影响
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。
    n1, n2, g1, g2 = n2, n1, g2, g1
m_order = list(np.random.permutation(groups[groups==1].index))
# 6. 根据倾向性得分，进行匹配
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匹配完成!")
# 7. 将匹配完成的结果合并起来
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)
matched_data = get_matched_data(matches, df)
# 8. 对匹配结果进行变量检验
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变量全部通过匹配")

未完待续
  
 
   
相关推荐
 
  •   
        Andy_Qin
      
     
       简单聊聊因果推断
     
    简单聊聊因果推断这个事。 Thoughtworks 曾与一家国际知名药企合作，解决患者用药依从性低(从药物治疗的角度，药物依从性是指患 者对药物治疗方案的执行程度)的问题。 如果仅着眼于患者本身，其不
     
    •  4910
      •  
     
       因果推断在哈啰的实践探索
     
    因果推断在智能营销场景有着广泛的应用。在哈啰酒店营销业务中，自研Tree-based因果推断模型通过修改分裂准则，使得模型目标和业务目标的一致性更高，在智能补贴模块中取得了较好的效果。




    

     
    •  1182
      •  
     
       因果推断学习笔记.5：结构因果模型 (SCM)
     
    上一篇的综述论文读的比较懵。为此，从这一篇开始记录Judea Pearl的Causal Inference in Statistics: A Primer的阅读笔记，将经典的因果推断理论基础打好。
     
    •  632
      •  
     
       《因果推断：如果》第一章学习笔记
     
    最近好忙，需要看的东西好多，真的有点烦，在同时需要看很多书的时候也许换着看换着写是一种好的提升看书效率的方法，所以也开始更Causal Inferences: What if这本书
     
    •  730
      •  
     
       Python实战|利用Dowhy框架实现因果推断实战（二）
     
    相关性不等于因果性，数据分析师必会的因果推断 这篇文章是基于上篇的基础继续分析，前面数据预处理的过程可以通过阅读全文查看源代码，这篇文章主要聚焦于因果推断的实现。
     
    •  1254
      •  
     
       因果推断学习笔记.7：前门、后门准则
     
    为了解决其他变量/协变量干扰因果效应的问题，老爷子提出了后门准则 & 前门准则。 后门准则 假设在变量图中，想要探究 $X\to Y$的因果效应 但是图中总有一些指向$X$的、看得见（图中有显示）却摸
     
    •  528
      •  
     
       Python日志库logging总结-可能是目前为止将logging库总结的最好的一篇文章
     
    在部署项目时，不可能直接将所有的信息都输出到控制台中，我们可以将这些信息记录到日志文件中，这样不仅方便我们查看程序运行时的情况，也可以在项目出现故障时根据运行时产生的日志快速定位问题出现的位置。 Python 标准库 logging 用作记录日志，默认分为六种日志级别（括号为级…
     
    •  31.8w
      •  
     
       因果推断学习笔记.6：什么是“干预”
     
    因果推断中的干预概念。 上一篇主要讲因果模型的形式化表示方法（结构因果模型）。从这一篇开始，主要围绕如何在现实场景中推断因果的问题展开叙述。
     
    •  479
      •  
     
       美团内部讲座 | 清华大学崔鹏：因果推断技术最新的发展趋势
     
    因果推理模型的出现，给机器学习打开了一个新的思路。美团技术团队特别邀请到清华大学计算机学院长聘副教授崔鹏老师，请他为美团技术团队的同学分享了因果推断技术最新的发展趋势，以及现阶段取得的一些成果。
     
    •  1318
      •  
     
       【1 月最新】前端 100 问：能搞懂 80% 的请把简历给我
     
    半年时间，几千人参与，精选大厂前端面试高频 100 题，这就是「壹题」。 在 2019 年 1 月 21 日这天，「壹题」项目正式开始，在这之后每个工作日都会出一道高频面试题，主要涵盖阿里、腾讯、头条、百度、网易等大公司和常见题型。得益于大家热情参与，现在每道题都有很多答案，提…
     
    •  65.0w
      •  
     
       [译] 使用因果分析优化 Go HTTP/2 服务器
     
    如果你一直都有关注本博客，那么你应该看过这篇介绍因果分析的论文。这种分析方式旨在建立性能消耗周期与性能优化之间的联系。我已经在 Go 语言中实践了这种分析方式。我觉得是时候在一个真正的软件中 —— Go 标准库的 HTTP/2 实现中去实践一下了。 HTTP/2 是我们熟悉并且…
     
    •  5481
      •  
     
       ChatGPT保姆级教程，一分钟学会使用ChatGPT！
     
    最近ChatGPT大火！微软退出首款ChatGPT搜索引擎，阿里等国内巨头也纷纷爆出自家产品，一夜之间，全球最大的科技公司仿佛都回到了自己年轻时的样子！ 然而，ChatGPT这么火，这么好玩的东西，国
     
    •  70.1w
      •  
     
       30 道 Vue 面试题，内含详细讲解（涵盖入门到精通，自测 Vue 掌握程度）
     
    本文以前端面试官的角度出发，对 Vue 框架中一些重要的特性、框架的原理以问题的形式进行整理汇总，意在帮助作者及读者自测下 Vue 掌握的程度。本文章节结构以从易到难进行组织，建议读者按章节顺序进行阅读，当然大佬级别的请随意。希望读者读完本文，有一定的启发思考，也能对自己的 V…
     
    •  65.0w
      •  
     
       中高级前端大厂面试秘籍，为你保驾护航金三银四，直通大厂(上)
     
    当下，正面临着近几年来的最严重的互联网寒冬，听得最多的一句话便是：相见于江湖~🤣。缩减HC、裁员不绝于耳，大家都是人心惶惶，年前如此，年后想必肯定又是一场更为惨烈的江湖厮杀。但博主始终相信，寒冬之中，人才更是尤为珍贵。只要有过硬的操作和装备，在逆风局下，同样也能来一波收割翻盘…
     
    •  66.6w
      •  
     
       因果推理的春天 - 数据挖掘中的Confounding, Collidar, Mediation Bias
     
    如何预估干预对用户的影响？如何从观测数据中推断因果？统计推断和因果推理有何区别？这一章让我们通过The Book of Why的经典案例来了解因果推断，跟随图灵奖得主一起推开因果世界的大门
     
    •  1289
      •   
     因果推理的春天 - 数据挖掘中的Confounding, Collidar, Mediation Bias
      
     
    •  
    友情链接：
     
     
      
    私信
     
    推荐文章
    绅士的匕首  ·  windows: pip install 报错SSLError Can‘t connect to HTTPS URL because the SSL module is not available_p
    昨天
    奔跑的水桶  ·  group by 多个字段（多个条件去分具体多少个组等价于条件分支）_group by分组 多少组-CSDN博客
    1 年前
    活泼的芒果  ·  救救孩子吧，到现在还搞不懂@ComponentScan注解的使用 - 掘金
    1 年前
    朝气蓬勃的面包  ·  Helix editor 笔记 (持续更新) - 知乎
    1 年前
    鼻子大的人字拖  ·  博士申请 | 香港大学倪卓娴教授招收机器学习和智慧物联网方向全奖博士生_PaperWeekly的博客-CSDN博客
    1 年前
    绅士的水煮肉  ·  从大熊猫数据栏中删除时区 - 问答 - 腾讯云开发者社区-腾讯云
    2 年前
    今天看啥   ·   Py中国   ·   codingpro   ·   小百科   ·   link之家   ·   卧龙AI搜索
    删除内容请联系邮箱 2879853325@qq.com
    Code - 代码工具平台
    © 2024 ~ 沪ICP备11025650号