PSM倾向匹配详细步骤和程序

1.安装psmatch2统计包。
2.数据准备
3.数据分析及命令解读
4.结果解读
- 4.1模型拟合结果，此处无太多实际意义。
- 4.2试验组可匹配的观测概览，按照命令中设定的匹配规则，试验组有8例患者未能匹配到合适对照。
- 4.3结果解读的重点应该是对stata新生成的中间变量的解读。
- 4.4均衡性检验结果
- 4.5匹配结果的图示化
5.Stata命令汇总
最后留下两个问题给大家思考：
题目

试验设计中，匹配的目的在于确保干预效应估计是建立在可比个体之间的不同结果的基础上。最简单的匹配方式是将干预组和对照组中协变量值相同的两个个体进行配对分析。但是，如果协变量并不是某一个变量，而是一组变量时，这种简单的匹配方式也就不再适用，而是采用倾向得分匹配方式进行匹配。倾向性匹配得分（PSM）分析，主流统计学软件SAS、Stata、SPSS（22.0以上版本）、R语言均可实现。但SAS难度较高，不推荐；SPSS虽然操作简便，但是仅能实现1:1匹配，如无特殊需求可以尝试。笔者重点推荐使用Stata或者R语言完成PSM分析。下面笔者将以实例演示的形式讲解Stata软件在倾向性匹配得分中的应用。

1.安装psmatch2统计包。

命令如下：

.ssc install psmatch2

需要在联网状态下键入上述命令，然后软件自动搜索对应的程序包进行安装，成功安装后会有以下提示：

checking psmatch2 consistency and verifying not already installed…
installing into .\ado\plus… installation complete.（出现此提示表示安装完成）

为了验证是否成功安装以及查看psmatch2命令的帮助菜单，可在命令窗口键入

.help psmatch2

如果能顺利弹出帮助文件，表示安装成功，可正常使用。

2.数据准备

数据如下图所示，共有10个变量，614个观测，试验组185例，对照组429例。treat变量即为分组变量，“1”=试验组，“0”=对照组。age, educ, black, hispan, married, nodegree, re74, re75为协变量, re78为结局变量。事实上，倾向性匹配得分分析是要建立一个以分组变量（treat）为因变量，各个协变量（age, educ, black, hispan, married, nodegree, re74, re75）为自变量的回归方程。而结局变量（re78）在PSM过程中几乎不参与建模。

图1. 数据整理

3.数据分析及命令解读

命令窗口键入如下命令：

.gen tmp = runiform()
.sort tmp (以上两步对所有观测值进行随机排序)
.psmatch2 treat age educ black hispan married nodegree re74 re75, out(re78) logit neighbor(1) common caliper(.05) ties
.pstest, both
.psgraph

命令解读：
以下是帮助菜单中psmatch2语法格式，

psmatch2 depvar [indepvars] [if exp] [in range] [, outcome(varlist)
pscore(varname) neighbor(integer) radius caliper(real)
mahalanobis(varlist) ai(integer) population altvariance kernel llr
kerneltype(type) bwidth(real) spline nknots(integer) common trim(real)
noreplacement descending odds index logit ties quietly w(matrix) ate]

简单说就是：psmatch2 因变量协变量，[选择项]。重点解读命令语句中选择项的含义。本例中选择“nearest neighbor matching within caliper”匹配方法。out(re78)指明结局变量。logit指定使用logit模型进行拟合，默认的是probit模型。neighbor(1)指定按照1:1进行匹配，如果要按照1:3进行匹配，则设定为neighbor(3)，本例中因对照组样本量有限，仅适合1:1进行匹配。common强制排除试验组中倾向值大于对照组最大倾向值或低于对照组最小倾向值。caliper(.05)试验组与匹配对照所允许的最大距离为0.05。ties强制当试验组观测有不止一个最优匹配时同时记录。
pstest, both做匹配后均衡性检验，理论上说此处只能对连续变量做均衡性检验，对分类变量的均衡性检验应该重新整理数据后运用χ2检验或者秩和检验。但此处对于分类变量也有一定的参考价值。
psgraph对匹配的结果进行图示。

4.结果解读

4.1模型拟合结果，此处无太多实际意义。

图2. 回归结果

4.2试验组可匹配的观测概览，按照命令中设定的匹配规则，试验组有8例患者未能匹配到合适对照。

图3. 匹配情况概览

4.3结果解读的重点应该是对stata新生成的中间变量的解读。

打开数据编辑窗口，会发现软件自动生成了几个新变量：其中_pscore是每个观测值对应的倾向值；_id是自动生成的每一个观测对象唯一的ID（事实上这列变量即是对_pscore排序）；_treated表示某个对象是否试验组；_n1表示的是他被匹配到的对照对象的_id（如果是1:3匹配，还会生成_n2, _n3）；_pdif表示一组匹配了的观察对象他们概率值的差。为了观察方便可以按照id变量进行排序，排序后结果如下图所示：

图4. 匹配后的数据
匹配后数据整理进行统计分析即可。

4.4均衡性检验结果

图5. 均衡性检验结果
由均衡性检验结果可知，（1）各变量匹配后在试验组和对照组间是均衡的。（2）只有educ这个变量匹配前后试验组较对照组p值无变化，匹配前该变量试验组和对照组就无差别，匹配后不太可能出现差异，因此在建模的时候也可以考虑把educ这个变量排除，事实证明排除这个变量后匹配结果更为理想，读者可自行尝试。需要再次强调的是，此处理论上说只能对连续变量做均衡性检验，对分类变量的均衡性检验应该重新整理数据后运用χ2检验或者秩和检验等方法。

4.5匹配结果的图示化

图6. 匹配结果

5.Stata命令汇总

.ssc install psmatch2 #安装程序包
.use "F:\lalonde.dta" #调用F盘存储数据
.gen tmp = runiform()
.sort tmp #对所有观测随机排序
.psmatch2 treat age educ black hispan married nodegree re74 re75, out(re78) logit neighbor(1) common caliper(.05) ties #PSM分析
.pstest, both #均衡性检验
.psgraph #图示匹配结果

最后留下两个问题给大家思考：

问题1：倾向性匹配得分分析这么牛，是不是可以替代设计良好的随机对照试验？
问题2：PSM既然可以很好的均衡基线特征，统计分析时继续运用多元回归是否还有必要？

题目

凡事有利必有弊！PSM肯定不是完美无缺的（废话，如果PSM完美无缺，那还要随机对照试验干什么？）。大家觉得PSM这种方法有何缺陷呢？当然，也可以顺便谈谈PSM的优点。
查看源图
计量经济圈资深圈友：Inno·静，贡献倾向得分匹配分析代码
第一种：
教程: 倾向匹配分析深度（Propsensity matching analysis）
安装系统包:ssc install psmatch2

统计缺失值:misstable sum smoke2_new qa1age employ2014 cfps2010edu_best qe1_best urban qq301 qg307siops feduc meduc

建模:probit smoke2_new qa1age employ2014 cfps2010edu_best qe1_best urban qq301 qg307siops feduc meduc

计算得分:predict pscore, p

匹配:

psmatch2 smoke2_new, pscore(pscore) noreplacement
attnd med_per fasmoke3_why $xlist, pscore(myscore) comsup boot reps($breps) dots

第二种：

psmatch2  smoke2_new qa1age employ2014 cfps2010edu_best qe1_best urban qq301 feduc meduc, out(income) logit neighbor(1) common caliper(.05) tiespstest, bothpsgraph

common强制排除试验组中倾向值大于对照组最大倾向值或低于对照组最小倾向值。caliper(.05)试验组与匹配对照所允许的最大距离为0.05。ties强制当试验组观测有不止一个最优匹配时同时记录。

pstest, both做匹配后均衡性检验，理论上说此处只能对连续变量做均衡性检验，对分类变量的均衡性检验应该重新整理数据后运用χ2检验或者秩和检验。但此处对于分类变量也有一定的参考价值。

psgraph对匹配的结果进行图示。

第三种：

global ylist logmed_per5global treat fasmoke1_whyglobal xlist falcohol3_why fasmokenum2_why smoke_longth logindincglobal xlist falcohol3_why logindincpscore $treat $xlist, pscore(mypscore) blockid(myblock) detail,if household_hsy==1psgraph, treated($treat)pscore(mypscore)gen logitpscore = log(mypscore/(1-mypscore))sum logitpscorequi psmatch2 $treat, outcome($ylist) pscore(mypscore) caliper(.01138) neighbor(1)pstest $xlist, treated($treat) both graph, if household_hsy==1drop mypscore myblock