拓端tecdat|Stata估算观测数据的风险比
原文链接：http://tecdat.cn/?p=6419
原文出处：拓端数据部落公众号
在分析二元结果时，逻辑回归是分析师对回归建模的默认方法。随机研究中，当然很容易估计比较两个治疗组的风险比。对于观察数据，治疗不是随机分配的，估计治疗效果的风险比有点棘手。
理想情况 - 随机治疗分配
理想情况下，我们首先模拟（在Stata中）一个大型数据集，该数据集可能在随机试验中出现：
gen x = rnormal（）
gen z =（runiform（）<0.5）
gen xb = x + z
gen pr = exp（xb）/（1 + exp（xb））
gen y =（runiform（）<pr）
此代码为10,000个人生成数据集。每个都有一个基线变量x的值，它是从标准N（0,1）分布模拟的。接下来，根据随机研究，我们模拟一个二进制变量z，概率0.5为1，概率0.5为0.然后生成二元结果y，我们从逻辑回归模型生成它，对数几率为1等于x + z。因此，对于x，调整z = 1与z = 0的真实优势比是exp（1）= 2.72。
由于处理是随机分配的，我们可以忽略x并使用带有日志链接的GLM命令估计比较z = 1到z = 0的风险比：
Generalized linear models No. of obs = 10000
Optimization : ML Residual df = 9998Scale parameter = 1
Deviance = 12960.39176 (1/df) Deviance = 1.296298
Pearson = 10000 (1/df) Pearson = 1.0002Variance function: V(u) = u*(1-u) [Bernoulli]
Link function : g(u) = ln(u) [Log]AIC = 1.296439
Log likelihood = -6480.195879 BIC = -79124.59------------------------------------------------------------------------------| OIMy | Risk Ratio Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------z | 1.429341 .0241683 21.13 0.000 1.382748 1.477503_cons | .495886 .0070829 -49.11 0.000 .4821963 .5099643
------------------------------------------------------------------------------风险比估计为1.43，因为数据集很大，95％置信区间非常窄。
估算观测数据的风险比
现在让我们考虑观测数据的情况。为此，我们模拟了一个新的数据集 ：
gen x=rnormal()
gen tr_xb=x
gen tr_pr=exp(tr_xb)/(1+exp(tr_xb))
gen z=(runiform() < tr_pr)
gen xb=x+z
gen pr=exp(xb)/(1+exp(xb))
gen y=(runiform() < pr)如果我们为y运行相同的GLM模型，忽略x，我们得到：
. glm y z, family(binomial) link(log) eformGeneralized linear models No. of obs = 10000
Optimization : ML Residual df = 9998Scale parameter = 1
Deviance = 12014.93919 (1/df) Deviance = 1.201734
Pearson = 10000 (1/df) Pearson = 1.0002Variance function: V(u) = u*(1-u) [Bernoulli]
Link function : g(u) = ln(u) [Log]AIC = 1.201894
Log likelihood = -6007.469594 BIC = -80070.04------------------------------------------------------------------------------| OIMy | Risk Ratio Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------z | 1.904111 .0353876 34.65 0.000 1.836 1.974748_cons | .4101531 .0069811 -52.36 0.000 .396696 .4240667
------------------------------------------------------------------------------
使用对数广义线性模型
最明显的方法是在我们的GLM命令中添加x：
. glm y z x, family(binomial) link(log) eformIteration 0: log likelihood = -9271.4631
Iteration 1: log likelihood = -5833.7585
Iteration 2: log likelihood = -5733.9167 (not concave)
Iteration 3: log likelihood = -5733.9167 (not concave)
...然而，这无法收敛 。
通过逻辑模型估计风险比率
一个相对简单的替代方案是使用逻辑模型来估计调整x的治疗风险比。
Logistic regression Number of obs = 10000LR chi2(2) = 2797.60Prob > chi2 = 0.0000
Log likelihood = -5343.6848 Pseudo R2 = 0.2075------------------------------------------------------------------------------y | Odds Ratio Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------z | 2.947912 .1442639 22.09 0.000 2.678297 3.244668x | 2.693029 .0811728 32.87 0.000 2.538542 2.856918_cons | .9910342 .0322706 -0.28 0.782 .929761 1.056345
------------------------------------------------------------------------------然而，我们可以使用该模型来计算每个人的预测概率，首先假设所有个体都未被治疗（z = 0），然后假设所有个体都被治疗（z = 1）：
replace z=0
predict pr_z0, pr
replace z=1
predict pr_z1, pr此代码首先生成一个新变量zcopy，它保留原始治疗分配变量的副本。然后我们将所有个体z设置为0，并询问y = 1的预测概率。然后我们将所有个体设置为z = 1，并再次计算P（y = 1）。现在估计z = 1与z = 0相比的风险比，我们简单地考虑这两个条件下的边际风险比，即P（y = 1 | z = 1）/ P（y = 1） | Z = 0）：
summ pr_z0 pr_z1Variable | Obs Mean Std. Dev. Min Max
-------------+--------------------------------------------------------pr_z0 | 10000 .4970649 .2060767 .0166056 .9765812pr_z1 | 10000 .7094445 .1765126 .0474181 .9919309di .7094445/.4970649
1.4272673
replace z=zcopy给出了我们估计的风险比，比较z = 1到z = 0，为1.43，与我们第一次模拟数据时估计的风险比相同，其中治疗分配是完全随机的（特别是独立于x）。
置信区间
我们已经找到了风险比的点估计，但我们当然也喜欢置信区间，以指示估计的精确度。
首先，当z设置为0然后设置为1时，我们使用teffects给出我们对二元结果的边际均值的估计（相当于y = 1的概率）：
Iteration 0: EE criterion = 1.898e-21
Iteration 1: EE criterion = 5.519e-34 Treatment-effects estimation Number of obs = 10000
Estimator : regression adjustment
Outcome model : logit
Treatment model: none
------------------------------------------------------------------------------| Robusty | Coef. Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------
POmeans |z |0 | .4957607 .0072813 68.09 0.000 .4814897 .51003181 | .7078432 .0071566 98.91 0.000 .6938164 .7218699
------------------------------------------------------------------------------为了计算风险比和置信区间，我们首先使用teffects ra，coeflegend来查找Stata保存估算值的名称：
Treatment-effects estimation Number of obs = 10000
Estimator : regression adjustment
Outcome model : logit
Treatment model: none
------------------------------------------------------------------------------y | Coef. Legend
-------------+----------------------------------------------------------------
POmeans |z |0 | .4957607 _b[POmeans:0bn.z]1 | .7078432 _b[POmeans:1.z]
------------------------------------------------------------------------------我们现在可以使用nlcom计算风险比及其置信区间。但是，由于这将为我们提供基于Wald的对称置信区间，因此最好找到对数风险比的这个区间，然后将得到的区间反向转换为风险比例：
_nl_1: log(_b[POmeans:1.z] / _b[POmeans:0bn.z])------------------------------------------------------------------------------y | Coef. Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------_nl_1 | .3561291 .0172327 20.67 0.000 .3223536 .3899047
------------------------------------------------------------------------------. di exp(.3561291)
1.4277919. di exp(.3223536)
1.3803728. di exp(.3899047)
1.47684
因此，我们估计风险比为1.43,95％置信区间为1.38至1.48。
非常感谢您阅读本文，有任何问题请在下方留言！
拓端tecdat|Stata估算观测数据的风险比相关推荐
- 拓端tecdat荣获掘金社区入驻新人奖
2021年7月,由掘金发起了"入驻成长礼"颁奖活动.本次活动邀请到知名开发者.服务机构代表等业界人士. 据了解,掘金社区"新入驻创作者礼"主要对已经积累了一定历 ...
- 拓端tecdat荣获2022年度51CTO博主之星
相信技术,传递价值,这是51CTO每一个技术创作者的动力与信念,2022 年度,拓端tecdat 作为新锐的数据分析咨询公司,在51CTO平台上,不断的输出优质的技术文章,分享前沿创新技术,输出最佳生 ...
- 拓端tecdat|bilibili视频流量数据潜望镜
最近我们被客户要求撰写关于bilibili视频流量的研究报告,包括一些图形和统计输出. 最新研究表明,中国有超过7亿人在观看在线视频内容.Bilibili,被称为哔哩哔哩或简称为B站,是中国大陆第二个 ...
- 拓端tecdat|R语言逻辑回归（Logistic回归）模型分类预测病人冠心病风险
最近我们被客户要求撰写关于冠心病风险的研究报告,包括一些图形和统计输出. 相关视频:R语言逻辑回归(Logistic回归)模型分类预测病人冠心病风险 逻辑回归Logistic模型原理和R语言分类预测冠 ...
- 拓端tecdat|R语言用LOESS(局部加权回归)季节趋势分解（STL）进行时间序列异常检测
最近我们被客户要求撰写关于LOESS(局部加权回归)的研究报告,包括一些图形和统计输出. 这篇文章描述了一种对涉及季节性和趋势成分的时间序列的中点进行建模的方法.我们将对一种叫做STL的算法进行研究, ...
- 拓端tecdat|R语言向量误差修正模型 (VECMs)分析长期利率和通胀率影响关系
最近我们被客户要求撰写关于向量误差修正模型的研究报告,包括一些图形和统计输出. 向量自回归模型估计的先决条件之一是被分析的时间序列是平稳的.但是,经济理论认为,经济变量之间在水平上存在着均衡关系,可以 ...
- 拓端tecdat|R语言线性回归和时间序列分析北京房价影响因素可视化案例
最近我们被客户要求撰写关于北京房价影响因素的研究报告,包括一些图形和统计输出. 目的 房价有关的数据可能反映了中国近年来的变化: 人们得到更多的资源(薪水),期望有更好的房子 人口众多 独生子女政策: ...
- R语言广义加性模型GAMs分析温度、臭氧环境数据绘制偏回归图与偏残差图
最近我们被客户要求撰写关于广义加性模型的研究报告,包括一些图形和统计输出. 视频:R语言广义相加模型(GAM)在电力负荷预测中的应用 拓端tecdat:R语言广义相加模型(GAM)在电力负荷预测中的应 ...
- R语言中的广义线性模型（GLM）和广义相加模型（GAM）：多元（平滑）回归分析保险资金投资组合信用风险敞口
最近我们被客户要求撰写关于信用风险敞口的研究报告,包括一些图形和统计输出. 在之前的课堂上,我们已经看到了如何可视化多元回归模型(带有两个连续的解释变量).在此,目标是使用一些协变量(例如,驾驶员的年 ...
- R语言惩罚逻辑回归、线性判别分析LDA、广义加性模型GAM、多元自适应回归样条MARS、KNN、二次判别分析QDA、决策树、随机森林、支持向量机SVM分类优质劣质葡萄酒十折交叉验证和ROC可视化
最近我们被客户要求撰写关于葡萄酒的研究报告,包括一些图形和统计输出. 介绍 数据包含有关葡萄牙"Vinho Verde"葡萄酒的信息.该数据集有1599个观测值和12个变量,分别是 ...
最新文章
- linux 任务计划 权限设置,Linux系统 文件权限＋计划任务＋日志系统
- pytest框架安装（MacOS）
- 【机器学习】数据挖掘算法——关联规则（三），FP-growth算法
- mongodb及其索引的使用例子
- 基于 IdentityServer3 实现 OAuth 2.0 授权服务【密码模式(Resource Owner Password Credentials)】...
- leetcode题解70-爬楼梯
- .NET ：在Visual Studio的不同Tab之间切换
- python下载后安装包在哪里找到_python安装包里idle在哪
- 【2017宁波联考】生成树
- 使用Python对文件进行批量改名
- 如何快速去除图片上的水印？去除图片水印怎么做？
- 数据中心22年基础架构演进史
- Linux替换Docker镜像源
- 搭载3D立体相册网页 加入背景音乐 真香！
- 电子束光刻胶 （正胶SX AR-P 6200，AR-P617PMMA/MA共聚物），负胶AR-N 7700
- android 64位进程,简单科普一下，安卓下的64位和32位
- 谷歌浏览器打开JSP页面依然输出源代码
- android设备变砖,安卓手机变砖了怎么办？
- 为什么黑帽子从不用鼠标，一直在敲键盘？看完长见识了！
- 隐私计算技术｜私有信息检索(PIR)及其应用场景