【原创】用stata软件实现剂量反应meta分析之二(含word文档与视频)

自从上次发布用stata软件实现剂量反应关系meta分析之一，也就是单项研究的meta分析后，好久没有更新了，从今天开始，为大家陆续更新后面的系列内容，并且为了大家更好的学习，发布贴子的同时，也为大家录制了配套视频，欢迎学习！
文献来源：
Orsini N, Ruifeng L, Wolk A, Khudyakov P, Spiegelman D.Meta-analysis for linear and non-linear dose-response relationships: examples, an evaluation of approximations, and software.American Journal of Epidemiology. 2012; 175(1):66-73.
Summarized data about the relation between total alcohol intake and colorectal cancer risk in 8 prospective cohort studies participating in the Pooling Project of Prospective Studies of Diet and Cancer.
第一步 安装 GLST/ 读取数据：
ssc install glst
use http://nicolaorsini.altervista.org/data/ex_alcohol_crc,

这个就是数据， 8 个研究（计算的是发病率，也就是队列研究）， 7 个变量。

第二步，计算固定效应模型线性剂量反应关系 // Two-stage fixed-effect dose-response model assuming linearity

glst logrr dose , se(se) cov(peryears cases) pfirst(study type) ts(f)

命令ts(f)应该就是固定效应模型，从结果我们看出，模型有统计学意义，P=0.000,异质性较小，无统计学意义，P=0.6868应该说模型拟合是挺好的。

/* Relative risk for 12 grams/day incremental unit */
lincom dose*12 , eform

增加12个单位的情况是这样的。

第三步随机效应线性模型的结果
/ Two-stage random-effect dose-response model assuming linearity

glst logrr dose , se(se) cov(peryears cases) pfirst(study type) ts(r)

和固定效应模型是差不多的，原因是因为这个研究的异质性很小。后来我用R软件计算，才1%。

/* Relative risk for 12 grams/day incremental unit */
lincom dose*12 , eform

这是相对于参照（数据中最小的参照DOSE）增加12个单位的情况

第四步 非线性固定效应模型的结果
// Non-linearity using fixed-effect  

capture drop doses*
_pctile dose , percentile(5 35 65 95)
ret list
mkspline doses  = dose  , knots(`=r(r1)' `=r(r2)' `=r(r3)' `=r(r4)') cubic displayknots
glst logrr doses* , se(se)  cov(peryears cases)    pfirst(study type)  
testparm doses2 doses3

这里的命令很多，但大致可以这样分一下， mkspline 及以前是产生 4 个结点 5% ， 35% ， 65% ， 95% 的 DOSE, mkspline 命令就是进行样条回归。我们观察 mkpline 命令后数据里面就出现 dose1~dose3 。

GLST命令后的结果是这样的，也就是上面非线性mkspline 后产生的三个固定系数。模型是有意义的。

再看检验的结果：

第五步：接下来是画图了，先画一下之前的固定效应线性模型图。
// Figure 1 A of the paper  
glst  logrr dose , se(se)   cov(peryears cases)   pfirst(study type)  ts(f)
predictnl lrr_lin = _b[dose]*dose
gen rr_lin = exp(lrr_lin)

]
在数据框里出现了这两个数据LRR_LIN和RR_LIN,可以看出这是拟合的结果
下面是非线性的模型的结果图，使用 0 作为参考剂量
// using 0 as referent

glst logrr doses*, se(se) cov(peryears cases) pfirst(study type)
predictnl logrrwithref = _b[doses1]*doses1 + _b[doses2]*doses2 + _b[doses3]*doses3, ci(lo hi)
gen rrwithref  = exp(logrrwithref)
gen lbwithref  = exp(lo)
gen ubwithref  = exp(hi)

接下来是把画图的结果用表格的形式显示出来。在这之前，可能有的同学，没有安装 XBLC 命令，没有关系，在 STATA 中输入 help xblc 然后按照步骤一步一步在线安装就可以了。

// Tabulate result using xblc (findit xblc)

levelsof dose, local(level)
xblc doses*, c(dose) at(`r(levels)') ref(0) eform  

这就是最后的结果，实际上把这个结果放在EXCEL中了就可以绘图了，当然也可以用下面的STATA绘图的代码完成。

twoway  (line lbwithref ubwithref rrwithref dose, sort lp(longdash longdash l ) lc(black black black) ) (line rr_lin dose, sort lp(shortdash) lc(black) ) ,  scheme(s1mono) ylabel(.9 1 1.2 1.5 1.9, angle(horiz) format(%3.2fc)) xlabel(0(5)60) legend(off)  ytitle("Relative Risk", margin(right))  xtitle("Alcohol intake, grams/day" , margin(top_bottom) ) name(figure1B, replace) yscale(log)  plotregion(style(none))

这就是图表了，里面的控件是可以调的，大家自已更改。图中那条直线是线性固定效应模型的结果。
<br /