Stata绘图：回归系数可视化-论文更出彩

本文转载自公众号连享会( ID：lianxh_cn )

作者： 万莉 (北京航空航天大学)

E-Mail:  [e-mail]

Source: 世界银行 stata-visual-library^[1]，项目主页^[2] 提供了大量 Stata 绘图可视化案例和代码。

注：本推文相关数据，do file 及资料获取方式如下：

• 百度云盘：https://pan.baidu.com/s/1XkPfCMA6rECZc_M-1xwzKg?pwd=k417 提取码: k417

长按扫码，获取资源

• 项目主页：https://github.com/worldbank/stata-visual-library

目录

• 1.简介
• 2.范例一：比较不同模型的边际效应
• 3.范例二：图表可视化回归系数
• 4.范例三：比较不同组别的边际效应
• 5.范例四：展示回归的预测结果/拟合值
• 6.总结

1. 简介

本文参考 GitHub 项目 stata-visual-library^[3] 提供的数据及代码，节选四个 Stata 范例来说明如何“用图展示回归结果”。

在学习本文案例时，不用太纠结弄懂原始数据中变量的定义（原作者并没提供足够的数据说明），只需借鉴案例中用到的数据处理思路和绘图命令。每节的基本结构为：展示完整的 Stata 范例（代码+部分注释）-- 展示输出的图 -- 重点强调并解释部分语法。

若你还不太熟悉 Stata 绘图命令，可先参考 一文看尽 Stata 绘图 这篇推文，了解基本语法结构。

2. 范例一：比较不同模型的边际效应

该小节的代码和数据来自 stata-visual-library/Library/reg-models^[4]。

* Figure: Comparison of marginal effects from linear and logistic specifications

* 1. 将图像的基本设置保存在全局暂元 graph_opts 里
  global graph_opts title(, justification(left) color(black) span pos(11)) ///
    graphregion(color(white)) ylab(,angle(0) nogrid) ///
    xtit(,placement(left) justification(left)) ///
    yscale(noline) xscale(noline) legend(region(lc(none) fc(none)))

* 2. 加载数据
  use "data\reg-models.dta" , clear

* 3. 创建画图需要的数据
  cap mat drop theResults 
  local theLabels ""
  local x = 15.5

  qui foreach var of varlist ??_correct refer med_any med_class_any_6 med_class_any_16 {
  // ??表示匹配任意两个字符
    local theLabel : var label `var'
    local theLabels `"`theLabels' `x' "`theLabel'""'
    local x = `x' - 2
    
    reg `var' facility_private i.case_code //线性回归
    
 //可输入return list查看 Stata 命令的返回值
    mat a = r(table)  //用矩阵保存线性回归命令的返回值
       
    local b = a[1,1]  //提取估计系数
    local ll = a[5,1] //提取左置信区间
    local ul = a[6,1] //提取右置信区间
    mat a = [`b',`ll',`ul',1] //横向合并矩阵
    mat rownames a = "`var'"  //以变量名为矩阵的行名
        
    logit `var' facility_private i.case_code //logit回归
    margins , dydx(facility_private)
        
    mat b = r(table)  //用矩阵保存logit回归命令的返回值
         
    local b = b[1,1]  //提取估计系数
    local ll = b[5,1] //提取左置信区间
    local ul = b[6,1] //提取右置信区间
    mat b = [`b',`ll',`ul',2] //横向合并矩阵
    mat rownames b = "`var'" //以变量名为矩阵的行名
        
    mat theResults = nullmat(theResults) \ a \ b  
 // \表示纵向合并矩阵
  }
        
  mat colnames theResults = "b" "ll" "ul" "type" 
  //修改矩阵的列名
  matlist theResults //查看矩阵元素

  clear //清除原始数据
  svmat theResults , names(col) 
  //将矩阵转化为变量
    
  gen n = _n
  replace n = 17-n //逆序排列

* 4. 绘图
  tw ///
   (rcap ll ul n if type == 1 , hor lc(navy)) ///
   (scatter n b if type == 1 , mc(black)) ///
   (rcap ll ul n if type == 2 , hor lc(maroon)) ///
   (scatter n b if type == 2 , mc(black)) ///
   , ///
   ${graph_opts} ///
   ylab(`theLabels') ///
   ytit(" ") ///
   xlab(-1 "-100p.p." -.5 `""-50p.p." "{&larr} Favors Public""' ///
   0 "No Effect" .5 `""+50p.p." "Favors Private {&rarr}""' 1 "+100p.p.") ///
   xline(0 , lc(gs12) lp(dash)) ///
   legend(order(2 "Marginal Effect" 1 "Linear Model" 3 "Logistic Model") r(1))
   
/*
由于局部暂元(local macro)运行一次即自动删除，需要连续运行上述代码
qui 表示悄悄处理命令，不显示结果
/// 表示换行，避免一行内容过多
tw 为 twoway 的缩写
hor 是horizontal的缩写，表示水平绘制；默认为竖直绘制
你可以把 $graph_opts 删掉看看差异
*/

语法要点：

• 利用 twoway rcap 绘制置信区间；
• twoway rcap 用于以|——|图示数据范围，左右两端分别为最小、最大值。

3. 范例二：图表可视化回归系数

该小节的代码和数据来自 stata-visual-library/Library/reg-chartable^[5]。基于案例一，该案例增加了表格形式。

* Figure: Impact of provider qualifications on main standardized patient outcomes

* 1. 运行adofiles
// 该adofile仅适用于本案例的数据，不适合推广到任何数据（需要自己修改）
// 运行该adofiles还需利用findit安装firthlogit estadd xml_tab 这三个命令。
* firthlogit, in package firthlogit from http://fmwww.bc.edu/RePEc/bocode/f
* estadd, in package st0085_2 from http://www.stata-journal.com/software/sj14-2
* xml_tab, in package dm0037 from http://www.stata-journal.com/software/sj8-3
findit firthlogit
findit estadd
findit xml_tab

qui do "ado\chartablesheet.ado"
version 13

* 2. 加载数据
use "data\reg-chartable.dta", clear

* 3. 绘图
* Sample Restriction Logic
// 选取部分样本
global cxr "& sp_case <3"
global sputum "& sp_case <3"
global sp_drugs_tb "& sp_case >2"
        
* Output
chartable correct cxr sputum s5_referral ///
          sp_drugs_tb sp_drugs_antibio sp_drugs_quin ///
          using "results.xlsx" ///
          , ///
          c(xi: firthlogit) ///
          rhs(q_mbbs pro_age pro_male patients_waiting_in i.sp_id) or p globalif

语法要点：

• 本案例使用了作者自己编写的外部命令 chartable ；
• 本案例不具有普适性，仅适用于本案例中的数据，若要推广则需手动修改作者的 adofile。

4. 范例三：比较不同组别的边际效应

该小节的代码和数据来自 stata-visual-library/Library/reg-het^[6]。

*FIGURE: Marginal effects by heterogeneity gorups

* 1. 加载系统自带数据
sysuse auto, clear

* 2. 预处理
ssc install coefplot, replace //安装命令

local label1 Foreign
local label0 Domestic
 
*在合并图像中，使其共用一个y轴，设置y轴
local yopts1 ylab(, labcolor(white)) /// 
             ytitle("") ///
             yscale(noline)
     
gr drop _all

*3. 利用循环语句绘制单个图
* Loop over heterogeneity variables
forvalues foreign = 0/1 {

* 运行线性回归得到边际效应
* -----------------------------------
  reg price mpg headroom if foreign == `foreign'   

* 绘制不同组别的边际效应
* -----------------------------
  coefplot ///
    , ///
    drop(_cons) /// 不显示常数项
    vertical ///
    yline(0, lpattern(dash)) /// 判断显著性
    ciopt(color(gs10 gs10) ///
    recast(rcap)) ///
    levels(90 95) /// 展示90%和95%置信区间
    ylab(, noticks glcolor(gs15)) ///
    xlab(, noticks labsize(small)) ///
    xtitle("") ///
    ytitle("Marginal effect on price") ///
    title(`label`foreign'', box bexpand bcolor(gs15)) /// 
    color(black) ///
    name(f`foreign') ///
    graphregion(color(white)) ///
   `yopts`foreign'' 
}
// title(...,box ...)添加标题，并添加边框线
 
* 4. 合并图形
gr combine  f0 f1, ///
   ycommon ///
   graphregion(color(white))
// gr是graph的缩写

语法要点：

• 可视化回归系数的外部命令为 coefplot；

• coefplot 的具体使用教程可参考 Stata可视化：让他看懂我的结果！。

5. 范例四：展示回归的预测结果/拟合值

该小节的代码和数据来自 stata-visual-library/Library/reg-predicted^[7]。

* 1. 加载系统自带数据
sysuse auto, clear

* 2. 创建画图需要的数据
reg price ibn.foreign //线性回归
// ibn.表示不设置基准组
  
test 1.foreign == 0.foreign
local pval : di %5.3f `r(p)'
local pval = trim("`pval'")
//去 trim 掉两端的空格
//用暂元保存p-value，用于画图
  
* 计算不同种类的边际效应
margins foreign

* 3. 绘图
marginsplot, ///
   recast(bar) /// 
   plotopts(barwidth(0.5) bargap(10)) ///
   ciopts(recast(rcap) color(gs10)) ///
   ylabel(0(2000)8000) ///
   title("Adjusted prediction of car price" "with 95% confidence interval") ///
   xlabel(, noticks) ///
   xtitle(Model origin) ///
   graphregion(color(white)) ///
   note(Note: P-value of F-test for coefficient equality across categories is `pval'.)
// recast(bar) 是将散点变成条形图形式。

语法要点：

其一，ibn.foreign, 1.foreign, 0.foreign 的用法涉及因子变量(Factor Variable)。有关因子变量的使用请参考连享会推文 Stata：因子变量全攻略 。

图示边际效应的命令为 marginsplot，具体使用教程可参考以下文章：

• 专题：回归分析^[8]
• 专题：交乘项-调节^[9]

6. 总结

本推文较为详细地介绍了 4 个 Stata 案例，展示了如何将估计系数结果简洁、直观的呈现出来。

注：本推文相关数据，do file 及资料获取方式如下：

• 百度云盘：https://pan.baidu.com/s/1XkPfCMA6rECZc_M-1xwzKg?pwd=k417 提取码: k417

长按扫码，获取资源

• 项目主页：https://github.com/worldbank/stata-visual-library

友情链接

[1]

stata-visual-library:  https://github.com/worldbank/stata-visual-library

[2]

项目主页:  https://worldbank.github.io/stata-visual-library/

[3]

stata-visual-library:  https://github.com/worldbank/stata-visual-library

[4]

stata-visual-library/Library/reg-models: https://github.com/worldbank/stata-visual-library/tree/master/Library/do

[5]

stata-visual-library/Library/reg-chartable: https://github.com/worldbank/stata-visual-library/tree/master/Library/do

[6]

stata-visual-library/Library/reg-het: https://github.com/worldbank/stata-visual-library/tree/master/Library/do

[7]

stata-visual-library/Library/reg-predicted: https://github.com/worldbank/stata-visual-library/tree/master/Library/do

[8]

回归分析: https://www.lianxh.cn/blogs/32.html

[9]

交乘项-调节: https://www.lianxh.cn/blogs/21.html 

本文转载自公众号：连享会（ID：lianxh_cn），作者：连享会。

星标⭐我们不迷路！想要文章及时到，文末“在看”少不了！

点击搜索你感兴趣的内容吧

往期推荐

数据Seminar

这里是大数据、分析技术与学术研究的三叉路口

推荐 | 彭绮荣

    欢迎扫描👇二维码添加关注    

点击下方“阅读全文”了解更多