上一期帖子我们讲了propensity score matching(PS)的注意事项,有朋友说最好能讲讲怎么做PS,其实关于怎么做PS已经有不少公众号,博主都做过比较详尽的介绍,今天袁老师就以“A Step-by-Step Guide to Propensity Score Matchingin R”这篇文献为蓝本,手把手教大家怎么用R做PS。
What!用R做,我根本不会编程怎么破。不用担心,虽然R是个编程统计软件,但得益于大量免费共享的函数包,很多新的方法都可以很容易的上手,而PS的代码更是只有寥寥几句而已。
这篇文献里用到的是MatchIt包,所以在开始做PS前你需要安装好R软件和MatchIt包,R软件从https://cran.r-project.org/下载免费使用,安装完成后运行:
Install.packages(“MatchIt”)和library(MatchIt)即完成MatchIt包的安装和加载。
接下来就是准备数据了,对于现有的数据你可以通过以下几种方式载入到R软件:
1)直接从电脑指定的文件夹载入数据,如有个位于C盘r文件夹的newyork的csv文件,可以通过运行:mydata
2)通过运行getwd()找到R软件的工作目录,然后把数据拷贝到这个工作目录文件夹里面去。袁老师安装的R软件工作目录是C盘的我的文档,那就把csv文件拷贝到“我的文档”文件夹即可。
getwd()
[1]"C:/Users/Yoga/Documents"
袁老师产生了个随机数据集,包含3个变量age、gender和group,病例组样本量为250,age从均数55,标准差为8的正态分布随机产生,gender从概率分布为0.7和0.3的二项分布随机产生;对照组1000,age从均数50,标准差为8的正态分布随机产生,gender从概率分布为0.5和0.5的二项分布随机产生。group取代表对照组,取1代表病例组,注意在matchit函数中,分组变量group只能取和1,且为对照组,对照组的样本量需要大于等于病例组。两组age分布情况如下,左边是病例组,右边是对照组
对照组gender为482/518,病例组gender为71/179。两组age和gender组间差异性比较结果如下,可以看到都是有统计学差异的。
t.test(age~group,data=mydata)$p.value
[1] 1.14814e-17
chisq.test(table(mydata$gender,mydata$group))$p.value
[1] 2.597042e-08
接下来我们就以最近邻匹配的方式进行1:1匹配,
match.it
summary(match.it)
可以看到匹配后两组age均数相近,gender的比例相近,两组各250例样本。也可以通过图示的方式来看看两组匹配前后倾向性评分的分布情况来直观的了解两组均衡性。
plot(match.it, type = 'jitter', interactive = FALSE)
plot(match.it, type = 'hist')
输出匹配后的数据集,并进行组间差异性比较,结果age和gender都无统计学差异。
df.match
t.test(age~group,data=df.match)$p.value
[1] 0.7140559
chisq.test(table(df.match$gender,df.match$group))$p.value
[1] 1
对于匹配前后均衡性评估除了做差异性检验外,目前比较多文献推荐采用标准化差值standard difference(也有人称之为effect size)来反映,以标准化差值
此外,还记得上期PS帖子提到的注意事项的最后一点吗,别忘了报告下未匹配样本的特征。
领取专属 10元无门槛券
私享最新 技术干货