R使用BUGSnet程序包实现二分类数据的贝叶斯网状meta分析

本文主要是使用BUGSnet程序包对2007年发表在Lancet上的《Incident diabetes in clinical trials of antihypertensive drugs: a network meta-analysis》一文结果进行复现。

公众号回:diabetes 获取效应数据。

1. 程序包加载及数据导入

install.packages(c("remotes", "knitr"))
remotes::install_github("audrey-b/BUGSnet@v1.1.0", upgrade = TRUE,
                        build_vignettes = TRUE, dependencies = TRUE)
library(BUGSnet)
library(readr)
data <- read_csv("diabetes.csv")

2. 干预措施信息比较

在进行分析时首先使用data.prep函数对数据进行规整,主要用于干预措施(treatment)和研究(study)。

data <- data.prep(arm.data = data,
                 varname.t = "treatment", 
                 varname.s = "study")

BUGSnet程序包的 net.tab()函数可输出网状、干预和比较特征的统计数据,其中:

  • network 输出的结果包含了网状的连通性、网状中的处理数量及事件数量等一些网状关系的基本特征

  • intervention 输出的结果包含了不同干预措施的总发生事件数,总样本量等信息

  • comparison 输出的结果包含了按不同治疗比较的事件数和总样本量等信息

network.char <- net.tab(data = data,
                        outcome = "responders",
                        N = "sampleSize", 
                        type.outcome = "binomial",
                        time = NULL)
> network.char$network
# A tibble: 12 × 2
   Characteristic                                        Value 
   <chr>                                                 <chr> 
 1 Number of Interventions                               6     
 2 Number of Studies                                     22    
 3 Total Number of Patients in Network                   154176
 4 Total Possible Pairwise Comparisons                   15    
 5 Total Number of Pairwise Comparisons With Direct Data 14    
 6 Is the network connected?                             TRUE  
 7 Number of Two-arm Studies                             18    
 8 Number of Multi-Arms Studies                          4     
 9 Total Number of Events in Network                     10962 
10 Number of Studies With  No Zero Events                22    
11 Number of Studies With At Least One Zero Event        0     
12 Number of Studies with All Zero Events                0  

> network.char$intervention 
# A tibble: 6 × 7
  treatment n.studies n.events n.patients min.outcome max.outcome av.outcome
  <chr>         <int>    <int>      <int>       <dbl>       <dbl>      <dbl>
1 ACE               8     1618      23351     0.0291       0.171      0.0693
2 ARB               5     1189      14185     0.00510      0.136      0.0838
3 BBlocker          9     2705      36150     0.0261       0.173      0.0748
4 CCB               9     2791      38809     0.0366       0.167      0.0719
5 Diuretic          8      973      18699     0.0229       0.0858     0.0520
6 Placebo           9     1686      22982     0.0154       0.185      0.0734

> network.char$comparison
# A tibble: 14 × 5
   comparison            n.studies n.patients n.outcomes proportion
   <chr>                     <int>      <int>      <int>      <dbl>
 1 ACE vs. BBlocker              3      15158       1022     0.0674
 2 ACE vs. CCB                   3      12597        538     0.0427
 3 ACE vs. Diuretic              2      16488        759     0.0460
 4 ACE vs. Placebo               3      17893       1929     0.108 
 5 ARB vs. BBlocker              1       7999        562     0.0703
 6 ARB vs. CCB                   1      10161       1535     0.151 
 7 ARB vs. Diuretic              1        392          9     0.0230
 8 ARB vs. Placebo               2       9778        573     0.0586
 9 BBlocker vs. CCB              5      44974       3361     0.0747
10 BBlocker vs. Diuretic         2       8752        241     0.0275
11 BBlocker vs. Placebo          1       3315         71     0.0214
12 CCB vs. Diuretic              2      15739        768     0.0488
13 CCB vs. Placebo               1       9711        331     0.0341
14 Diuretic vs. Placebo          3       7343        384     0.0523

3. 网络证据图( Network plot)

网络证据图使用net.plot()函数进行绘制,不需进行过多设置,便可以绘制比较好的图片。

net.plot(data)

图片

也可以进行复杂的图形绘制:

net.plot(data,
         node.scale=4, 
         edge.scale=1.5, 
         flag="Placebo", 
         label.offset1=0.5, 
         label.offset2=0.5,
         study.counts=T,
         graph.scale = TRUE, 
         node.lab.cex = 1, 
         edge.lab.cex = 1,
         node.colour = "steelblue", 
         edge.colour = "black",
         edge.lab.colour = "blue", 
         flag.edge.colour = "black",
         layout = "layout_in_circle", 
         layout.params = NULL)

图片

4. 模型选择(固定vs随机,一致vs不一致)

nma.model()函数是该程序包的核心函数,用来构建NMA模型可以运行固定效应模型、随机效应模型及一致性模型、不一致性模型;并提供了模型比较图,来进行模型选择。

4.1 固定效应模型和随机效应模型选择

fixmodel <- nma.model(data=data,#固定效应模型
                      outcome = "responders",
                      N = "sampleSize", 
                      reference="Placebo",
                      family="binomial",
                      link="log",
                      effects="fixed")

ranmodel <- nma.model(data = data,#随机效应模型
                      outcome = "responders",
                      N = "sampleSize", 
                      reference="Placebo",
                      family="binomial",
                      link="log",
                      effects="random")
fixresults <- nma.run(fixmodel,
                      n.adapt=1000,
                      n.burnin=10000, 
                      n.iter=50000)
ranresults <- nma.run(ranmodel,
                      n.adapt=1000,
                      n.burnin=10000,
                      n.iter=50000)

上述代码中,

  • 参数“outcome”用于指定结局变量;

  • 参数reference用于设置对照组,通常是某种安慰剂或对照药物等;

  • 参数family用于指定结果的分布类型,通常有“binomial”(二分类数据),“normal”(连续性数据),“poisson”(计数数据)可供选择;

  • 参数 link用于指 定NMA 模型使用的函数,其中“logit”用于二分类数据的比值比(Odds Ratio),“log”用于二分类数据的风险比(Risk Ratio)或计数数据的比率(RateRatio);“cloglog”用于二分类数据的危险比(Hazard Ratio);identity用于连续性数据;

  • effects用于设置效应模型的类型。

使用图形展示模型结果:

par(mfrow=c(1,2))
nma.fit(fixresults, main ="fixed model" )
nma.fit(ranresults, main= "random model")

图片

可使用nma.diag()函数绘制模型诊断图:

nma.diag(ranresults)

图片

4.2 一致性模型和不一性模型选择

根据DIC最小原则,上述模型结果显示随机效应模型拟合更佳,可进一步选择一致和不一致模型:

inconmodel <- nma.model(data = data,
                        outcome = "responders",
                        N = "sampleSize", 
                        reference="Placebo",
                        family="binomial",
                        link="log",
                        type = "inconsistency",
                        effects="random")

inconresults<- nma.run(inconmodel,
                        n.adapt=5000,
                        n.burnin=3,
                        n.iter=20000)
conresults<- nma.model(data = data,
                       outcome = "responders",
                       N = "sampleSize", 
                       reference="Placebo",
                           family = "binomial",
                           link = "log",
                           effects = "random",
                           type="consistency",)

conresults <- nma.run(conresults,
                      n.adapt=5000,
                      n.burnin=3,
                      n.iter=20000)

对结果进行可视化:

par(mfrow=c(1,2))
nma.fit(conresults,main = "consistency" )
nma.fit(inconresults,main = "inconsistency")

图片

根据DIC最小原则,最终选择随机效应模型和一致性模型。

5. 森林图 (Forest plot)

5.1 网状meta分析森林图

使用nma.forest() 函数输出森林图,对药物在不同研究和不同效应量水平的合并结果进行比较,并显示了有效的干预措施。

nma.forest(nma = ranresults,
           log.scale=FALSE,
           central.tdcy="mean",
           comparator = "Placebo")

图片

5.2 成对比较森林图

除了绘制网状meta结果的森林图外,该程序包还可以调用metafor程序包使用pma对两两干预措施进行直接比较,并绘制森林图,具体代码如下:

pma(data = data,
    type.outcome="binomial",
    sm="OR",
    name.trt1 = "CCB", 
    name.trt2 = "BBlocker", 
    outcome = "responders",
    N = "sampleSize") 
$summary
        comparison n.studies i.squared re.estimate   re.lci   re.uci fe.estimate   fe.lci   fe.uci
1 CCB vs. BBlocker         5 0.6249776    1.230353 1.075073 1.408062    1.277449 1.189803 1.371551

$forest
$forest$xlim
[1] -0.5615787  0.5615787

$forest$addrows.below.overall
[1] 0

$forest$colgap
[1] 2mm

$forest$colgap.left
[1] 2mm

$forest$colgap.right
[1] 2mm

$forest$colgap.studlab
[1] 2mm

$forest$colgap.forest
[1] 2mm

$forest$colgap.forest.left
[1] 2mm

$forest$colgap.forest.right
[1] 2mm

$forest$studlab
[1] "1"  "5"  "14" "17" "21"

$forest$TE.format
[1] ""     ""     ""     "0.10" "0.38" "0.17" "0.16" "0.02"

$forest$seTE.format
[1] ""       ""       ""       "0.2332" "0.0577" "0.0594" "0.0950" "0.1480"

$forest$cluster.format
[1] "" "" ""

$forest$effect.format
[1] "1.28" "1.23" ""     "1.11" "1.47" "1.19" "1.18" "1.02"

$forest$ci.format
[1] "[1.19; 1.37]" "[1.08; 1.41]" ""             "[0.70; 1.75]" "[1.31; 1.64]" "[1.06; 1.33]"
[7] "[0.98; 1.42]" "[0.77; 1.37]"


$raw
Number of studies combined: k = 5
Number of observations: o = 44974
Number of events: e = 3361

                         OR           95%-CI    z  p-value
Common effect model  1.2774 [1.1898; 1.3716] 6.75 < 0.0001
Random effects model 1.2304 [1.0751; 1.4081] 3.01   0.0026

Quantifying heterogeneity:
 tau^2 = 0.0130; tau = 0.1138; I^2 = 62.5% [0.6%; 85.8%]; H = 1.63 [1.00; 2.66]

Test of heterogeneity:
     Q d.f. p-value
 10.67    4  0.0306

Details on meta-analytical method:
- Mantel-Haenszel method
- DerSimonian-Laird estimator for tau^2
- Mantel-Haenszel estimator used in calculation of Q and tau^2 (like RevMan 5)

图片

6. 累计疗效排序图(SUCRA plot)

使用nma.rank()函数进行干预措施的排序,结果可产生rankogramlongtablesucraplot三种形式。

sucra.out <- nma.rank(ranresults,
                      largerbetter=T, 
                      sucra.palette= "Set1")
sucra.out$rankogram

图片

sucra.out$sucraplot

图片

league.out$table
write.csv(league.out$table,"effectR.csv"

7. 赛联表及赛联热图( league heat plot)

两两比较使用nma.league() 函数可实现,结果可产出赛联表和赛联热图。

league.out <- nma.league(ranresults,  
                         central.tdcy="mean",
                         high.colour = "cornflowerblue")      
league.out$heatplot

图片

也可以产出赛联表,并导出EXCEL。

league.out$table
write.csv(league.out$table,"effectR.csv")

参考文献

1. Béliveau A, Boyne DJ, Slater J, Brenner D, Arora P. BUGSnet: an R package to facilitate the conduct and reporting of Bayesian network Meta-analyses. BMC Med Res Methodol. 2019;19(1):196. Published 2019 Oct 22. doi:10.1186/s12874-019-0829-2.

2. .https://bugsnetsoftware.github.io/instructions.html

评论 1
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

皮肤小白生

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值