我们仅仅是代码的编辑者、整合者、搬运工,仅免费传授方法,下文数据和代码取自于网络和免费软件“R语言说明书”,如果您觉得我们侵犯了您的版权,请通知我们撤稿。请大家谅解,谢谢!
# stata软件 R软件请查看本系列之前推文,点击链接查看#
#第一部分 感谢“赵坤同学”的代码和说明 #
#1 代码来源####
#第二部分 meta分析的实现#
1.1 软件安装及程序包加载
R 软件当前最新版本为R-3.1.2,获取网址为http://www.r-project.org/;GeMTC软件最新版本为0.14.3,获取网址为http://drugis.org/software/addis1/gemtc。GeMTC软件在运行时需要Java环境,Java可从 http://www.java.com/zh_CN/下载。本文采用OpenBUGS软件执行数据分析,因此还需从http://www.openbugs.net/w/Downloads下载该软件。所有软件下载完成后,双击安装包进行安装。
软件安装完毕后,双击桌面的R软件图标启动R 的交互式窗口(R—GUI);接着在菜单栏处选择程序包下的安装程序包→选择CRAN镜像→gemtc,安装完成后可由library(gemtc)命令完成gemtc程序包加载。此外还需加载coda程序包、lattice程序包,具体方法如上。
1.2. 数据处理及录入
本文以Multiple-treatments Meta-analysis of a network of interventions[12]的数据为例进行演示,处理后的数据如表1所示。
表1 本文使用数据
Study ID | treatment | mean.change | sd | n |
1 | 1 | -10.26 | 13.43 | 70 |
1 | 2 | -4.88 | 11.64 | 66 |
2 | 1 | -14.78 | 12.49 | 54 |
2 | 2 | -8.13 | 12.72 | 56 |
3 | 1 | -13.36 | 8.78 | 125 |
3 | 3 | -10.39 | 10.45 | 123 |
4 | 1 | -13.11 | 8.53 | 220 |
4 | 2 | -9.1 | 9.36 | 114 |
5 | 1 | -16.6 | 12.49 | 55 |
5 | 3 | -14.9 | 11.39 | 60 |
6 | 1 | -15.49 | 8.148877223 | 58 |
6 | 2 | -15.25 | 8.372448865 | 59 |
7 | 3 | -11.5 | 10.9 | 187 |
7 | 2 | -9 | 10.9 | 177 |
8 | 3 | -10.1 | 10.9 | 144 |
8 | 2 | -8.7 | 10.9 | 78 |
9 | 1 | -11.71 | 7 | 100 |
9 | 2 | -8.97 | 6.934414179 | 101 |
10 | 1 | -9.4 | 8.506468127 | 201 |
10 | 3 | -8.2 | 8.455672652 | 186 |
10 | 2 | -7.4 | 7.959899497 | 99 |
注:设Treatment中“1”代表文拉法辛,“2”代表西酞普兰,“3”代表帕罗西汀
数据处理完成后,采用GeMTC软件录入(操作简便且不易出错)。
1.2.1录入干预措施组
选择Data下的Treatment进行干预措施的数据输入,可点击红色框区的Add teartment进行添加、编辑与删除,ID编号与表格一致,输入1,具体描述为文法拉辛(图1)。西酞普兰与帕罗西汀的添加方法与此相同。
图1 干预措施组数据录入
图2 研究组数据录入
1.2.2录入研究组
选择Data下的Study进行研究数据输入,可点击Adds tudy进行添加,ID编号与表格一致,可知1号研究对应的干预措施为1和2,因此ID 为1,选择干预措施1 2并打勾。在红色框区选择数据类型并对该数据集命名,本例选择CONTINOUS(连续性变量),如图2所示。其余研究同理进行添加。
1.2.3录入结局指标
连续型变量录入的是mean、sd和n。图3为数据录入完成后的界面,需注意的是ID为自动排序因此编号为10的排列第二,且数据列表中干预措施的排列顺序与在研究组录入时的点击顺序有关。录入数据完成后,点击Save进行保存及文件命名,本例命名为depression。
图3 数据完成界面
1.3.数据分析
采用R软件gemtc程序包调用OpenBUGS软件进行NMA以及相关图形的绘制。
1.3.1 读取录入数据
读取数据的命令为:
network<-read.mtc.network("C:/Users/Administrator/Desktop/depression.gemtc")
注意括号内容应于1.2.3的文件命名及保存位置一致,命令前的network<-相当于赋值,不能省略,命令符号应用英文输入状态且不应有空格,以下同理。如需查看数据可输入“print(network)”命令。
1.3.2 设置model
设置model的命令为:
model <- mtc.model(network,type = "consistency",n.chain = 3)
其中network 为network 数据,type 为是否选取一致性模型(注:在模型的选择上,通常认为由一致性模型得出的结果更加接近真实情况),n.chain为迭代运算中链的条数(一般来说,保证n≥3的链条数目使结果更加可信)。
1.3.3 分析数据
在加载程序包之前,应确保OpenBUGS安装在C盘,否则不能调用。执行数据分析的命令为:
results <-mtc.run(model, sampler ="BRugs",n.adapt = 5000, n.iter = 20000, thin = 1)
命令中,“n.adapt”为预迭代次数,“n.iter”为迭代运算次数,“thin”为步长。
1.4 汇总结果
输入命令summary(results),查看均数差与百分数,结果如图4所示
图4 Meta分析图
图5收敛诊断图
1.5 图形绘制
1.5.1收敛诊断图
收敛诊断图(convergence diagnostics)的绘制方法为通过输入命令“gelman.plot(results)”实现,图形见图5所示。保存格式及位置可以自行选择。
1.5.2森林图
森林图(forest)的绘制命令为forest(results),命令执行后结果见图6。
图6 森林图
1.5.3轨迹密度图
轨迹密度图(trace and density)的实现命令为“plot(results)”,结果如图7。
图7 轨迹密度图
1.5.4排序概率图
排序概率图(rank probability)分别有一下三种表示方式(图8,图9),分别对应的命令为1.ranks<- rank.probability(results),输入命令print(ranks)查看。
2.plot(ranks)
3.plot(ranks, beside=TRUE)
图8 排序概率表
图9排序概率图
1.5.5内部关系总结图
输入命令summary(network),可查看每种干预措施的研究数量以及试验臂数对应的研究数量。如图10所示,本研究中研究文法拉辛的研究有8篇,2臂试验9篇,3臂实验1篇。
图10.内部关系总结图
1.5.6网状关系图
因为利用plot(network)制作出的网状关系图的命令仅反映该研究之间的直接比较关系,无法通过顶点反映研究数量的变化,因此可选择R软件或STATA绘制网状关系图(network plot)[9],因STATA还可制作贡献图(contribution plot)及发表偏移检测图(publication bias test),因此本文使用STATA进行连续性变量的NMA的详细介绍。
2.1 STATA下载与操作
用户可从http://www.stata.com/stata13/购买最新版stata13.0,本文使用的为该版本[10]。
2.2模块加载
加载metan模块的命令为:ssc install metan。
加载图形分析network graphs模块的命令为:net from http://www.mtm.uoi.gr,点击 network_graphs后按照提示进行安装。
加载metareg
2.3数据录入
首先需要将表1中的数据整理成2臂实验的格式[7],整理后的数据如表2所示。整理完成后,输入到STATA中。
表2 STATA应用数据格式
STUDY | t1 | t2 | x1 | s1 | n1 | x2 | s2 | n2 |
1 | 1 | 2 | -10.26 | 13.43 | 70 | -4.88 | 11.64 | 66 |
2 | 1 | 2 | -14.78 | 12.49 | 54 | -8.13 | 12.72 | 56 |
3 | 1 | 3 | -13.36 | 8.78 | 125 | -10.39 | 10.45 | 123 |
4 | 1 | 2 | -13.11 | 8.53 | 220 | -9.1 | 9.36 | 114 |
5 | 1 | 3 | -16.6 | 12.49 | 55 | -14.9 | 11.39 | 60 |
6 | 1 | 2 | -15.49 | 8.148877 | 58 | -15.25 | 8.372449 | 59 |
7 | 3 | 2 | -11.5 | 10.9 | 187 | -9 | 10.9 | 177 |
8 | 3 | 2 | -10.1 | 10.9 | 144 | -8.7 | 10.9 | 78 |
9 | 1 | 2 | -11.71 | 7 | 100 | -8.97 | 6.934414 | 101 |
10 | 1 | 2 | -9.4 | 8.506468 | 201 | -7.4 | 7.959899 | 99 |
10 | 1 | 3 | -9.4 | 8.506468 | 201 | -8.2 | 8.455673 | 186 |
10 | 2 | 3 | -7.4 | 7.959899 | 99 | -8.2 | 8.455673 | 186 |
2.4 图形绘制
2.4.1 网状关系图
在命令输入框输入networkplot t1 t2[13],如图11,顶点与线条的粗细均表示研究数量的多少,与样本量无关,这应与R软件相区别[10]。用户可以使用Graph Editor进行网状控制图的编辑,调整字体的大小、线条的颜色等等。
图11 网状关系图
图12 贡献图
2.4.2 贡献图
输入命令metan n1 x1 s1 n2 x2 s2,可得出各个研究的SMD以及seSMD,在STATA中表示为_ES, _seES。
输入命令netweight _ES _seES t1 t2[13],结果如图12。
2.4.3 发表偏倚检测图
输入命令netfunnel _ES _seES t1 t2, bycomparison[13]。结果如图13,
图13 发表偏移检测图
拓展知识:
SNPmeta
网状emta分析必备技能
网状meta分析必备技能6_R+Rstudio运用meta包做简单meta分析
网状meta分析必备技能7_使用R、GeMTC和STATA软件实现连续变量的网状Meta分析
网状meta分析stata简易教程
NMA(网状meta分析)stata简易教程(4)漏斗图的制作
NMA(网状meta分析)stata简易教程(5)排序图的制作
诊断性meta分析简单实现
其他
本文由 GCBI学院 作者:其明技术专家 发表,转载请注明来源!