library(bmeta)
## Loading required package: R2jags
## Loading required package: rjags
## Loading required package: coda
## Linked to JAGS 4.2.0
## Loaded modules: basemod,bugs
##
## Attaching package: 'R2jags'
## The following object is masked from 'package:coda':
##
## traceplot
## Loading required package: forestplot
## Loading required package: grid
## Loading required package: magrittr
## Loading required package: checkmate
### Read and format the data (binary)
data = read.csv(url("http://www.statistica.it/gianluca/bmeta/Data-bin.csv"))
### List data for binary outcome
data.list <- list(y0=data$y0,y1=data$y1,n0=data$n0,n1=data$n1)
### generate output from bmeta
x <- bmeta(data=data.list,outcome="bin",model="std.dt",type="ran")
## module glm loaded
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 18
## Unobserved stochastic nodes: 21
## Total graph size: 123
##
## Initializing model
### generate autocorrelation function plot
acf.plot(x,"alpha[1]")
### generate autocorrelation function plot and specify the title
acf.plot(x,"alpha[1]",title="Autocorrelation plot")
### Read and format the data (binary)
data = read.csv(url("http://www.statistica.it/gianluca/bmeta/Data-bin.csv"))
### List data for binary outcome (for meta-analysis)
d1 <- data.list <- list(y0=data$y0,y1=data$y1,n0=data$n0,n1=data$n1)
### List data for binary outcome when there is a covariate (for meta-regression)
d1 <- data.list <- list(y0=data$y0,y1=data$y1,n0=data$n0,n1=data$n1,X=cbind(data$X0))
### Select fixed-effects meta-analysis with normal prior for binary data
m1 <- bmeta(d1, outcome="bin", model="std.norm", type="fix",n.iter=100)
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 18
## Unobserved stochastic nodes: 10
## Total graph size: 95
##
## Initializing model
m1
## Inference for Bugs model at "model.txt", fit using jags,
## 2 chains, each with 100 iterations (first 50 discarded)
## n.sims = 100 iterations saved
## mu.vect sd.vect 2.5% 97.5% Rhat n.eff
## alpha[1] -1.688 0.452 -2.199 -0.451 1.144 15
## alpha[2] -2.008 0.291 -2.407 -1.382 1.085 100
## alpha[3] -0.469 0.225 -0.837 0.035 1.119 29
## alpha[4] -0.820 0.270 -1.531 -0.493 1.023 75
## alpha[5] -2.296 0.786 -3.098 -0.158 0.994 100
## alpha[6] -1.702 0.288 -1.944 -0.824 0.992 100
## alpha[7] -2.250 0.677 -2.901 -0.365 0.997 100
## alpha[8] -2.270 0.514 -2.822 -1.043 1.001 100
## alpha[9] -0.164 0.164 -0.430 0.282 1.016 94
## delta -0.223 0.176 -0.401 0.253 1.004 100
## rho 0.814 0.165 0.670 1.289 1.004 100
## deviance 595.691 1054.209 133.362 4116.464 0.992 100
##
## For each parameter, n.eff is a crude measure of effective sample size,
## and Rhat is the potential scale reduction factor (at convergence, Rhat=1).
##
## DIC info (using the rule, pD = var(deviance)/2)
## pD = 561312.9 and DIC = 561908.6
## DIC is an estimate of expected predictive error (lower deviance is better).
### Select random-effects meta-regression with t-distribution prior for binary
### data
m2 <- bmeta(data.list, outcome="bin", model="reg.dt", type="ran",n.iter=100)
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 18
## Unobserved stochastic nodes: 22
## Total graph size: 147
##
## Initializing model
m2
## Inference for Bugs model at "model.txt", fit using jags,
## 2 chains, each with 100 iterations (first 50 discarded)
## n.sims = 100 iterations saved
## mu.vect sd.vect 2.5% 97.5% Rhat n.eff
## alpha[1] -0.832 26.376 -48.076 47.225 1.022 63
## alpha[2] -3.527 32.924 -63.443 55.218 1.024 58
## alpha[3] -1.603 32.911 -61.529 57.203 1.022 62
## alpha[4] 0.655 26.358 -46.588 48.591 1.026 54
## alpha[5] -4.392 32.916 -64.394 54.305 1.022 61
## alpha[6] -0.803 26.346 -47.910 47.167 1.022 61
## alpha[7] -2.103 26.311 -49.009 46.009 1.020 66
## alpha[8] -4.179 32.988 -64.257 55.105 1.025 55
## alpha[9] 0.797 26.343 -46.182 48.578 1.022 61
## delta[1] -0.348 0.184 -0.707 -0.032 1.023 57
## delta[2] -0.587 0.198 -0.979 -0.265 2.327 3
## delta[3] -0.620 0.289 -1.186 -0.192 2.823 3
## delta[4] -0.610 0.311 -1.319 -0.223 1.637 5
## delta[5] -0.711 0.378 -1.556 -0.249 2.606 3
## delta[6] -0.246 0.147 -0.445 0.055 1.853 4
## delta[7] -0.247 0.241 -0.633 0.289 3.106 3
## delta[8] -0.298 0.226 -0.684 0.323 1.175 100
## delta[9] -0.106 0.247 -0.514 0.441 2.508 3
## gamma[1] 0.718 0.132 0.493 0.968 1.023 57
## gamma[2] 0.567 0.111 0.376 0.768 2.327 3
## gamma[3] 0.560 0.154 0.306 0.826 2.823 3
## gamma[4] 0.567 0.152 0.268 0.800 1.637 5
## gamma[5] 0.523 0.163 0.211 0.780 2.606 3
## gamma[6] 0.791 0.122 0.641 1.056 1.853 4
## gamma[7] 0.804 0.205 0.531 1.335 3.106 3
## gamma[8] 0.762 0.191 0.505 1.381 1.175 100
## gamma[9] 0.928 0.247 0.598 1.555 2.508 3
## mu -0.413 0.117 -0.650 -0.249 1.169 20
## rho 0.666 0.079 0.522 0.780 1.169 20
## sigma 0.314 0.204 0.054 0.785 3.339 3
## tau 44.915 99.479 1.630 346.268 3.339 3
## deviance 114.949 8.035 101.393 130.700 3.137 3
##
## For each parameter, n.eff is a crude measure of effective sample size,
## and Rhat is the potential scale reduction factor (at convergence, Rhat=1).
##
## DIC info (using the rule, pD = var(deviance)/2)
## pD = 13.0 and DIC = 127.9
## DIC is an estimate of expected predictive error (lower deviance is better).
### Read and format the data (continuous)
data = read.csv(url("http://www.statistica.it/gianluca/bmeta/Data-ctns.csv"))
### List data for continuous outcome for studies reporting two arms separately
### (for meta-analysis)
d1 <- data.list <- list(y0=data$y0,y1=data$y1,se0=data$se0,se1=data$se1)
### List data for continuous outcome for studies reporting mean difference and
### variance with a covariate (for meta-regression)
d2 <- data.list2 <- list(y=data$y,prec=data$prec,X=cbind(data$X0))
### Select fixed-effects meta-analysis with studies reporting information of
### both arm for continuous data
m1 <- bmeta(data.list, outcome="ctns", model="std.ta", type="fix",n.iter=100)
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 12
## Unobserved stochastic nodes: 7
## Total graph size: 77
##
## Initializing model
m1
## Inference for Bugs model at "model.txt", fit using jags,
## 2 chains, each with 100 iterations (first 50 discarded)
## n.sims = 100 iterations saved
## mu.vect sd.vect 2.5% 97.5% Rhat n.eff
## alpha0[1] 17.620 0.979 15.585 19.318 1.021 100
## alpha0[2] 8.306 0.910 6.544 10.244 1.003 100
## alpha0[3] 12.689 0.866 11.091 14.524 1.023 100
## alpha0[4] 10.358 0.598 9.274 11.314 1.003 100
## alpha0[5] 15.153 0.389 14.328 15.819 1.035 41
## alpha0[6] 18.408 0.628 17.158 19.618 1.000 100
## alpha1[1] 16.185 0.959 14.295 17.718 1.036 100
## alpha1[2] 6.871 0.833 5.202 8.696 1.002 100
## alpha1[3] 11.255 0.901 9.429 12.802 0.998 100
## alpha1[4] 8.923 0.483 8.025 9.803 1.001 100
## alpha1[5] 13.718 0.424 12.843 14.464 0.994 100
## alpha1[6] 16.974 0.591 15.759 17.955 1.033 100
## delta -1.435 0.477 -2.322 -0.569 1.050 36
## deviance 36.837 3.290 32.038 43.780 1.041 38
##
## For each parameter, n.eff is a crude measure of effective sample size,
## and Rhat is the potential scale reduction factor (at convergence, Rhat=1).
##
## DIC info (using the rule, pD = var(deviance)/2)
## pD = 5.3 and DIC = 42.2
## DIC is an estimate of expected predictive error (lower deviance is better).
### Select random-effects meta-regression with studies reporting mean difference and
### variance only for continuous data
m2 <- bmeta(data.list2, outcome="ctns", model="reg.mv", type="ran",n.iter=100)
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 6
## Unobserved stochastic nodes: 9
## Total graph size: 46
##
## Initializing model
m2
## Inference for Bugs model at "model.txt", fit using jags,
## 2 chains, each with 100 iterations (first 50 discarded)
## n.sims = 100 iterations saved
## mu.vect sd.vect 2.5% 97.5% Rhat n.eff
## alpha[1] -0.929 1.471 -3.699 2.269 1.316 9
## alpha[2] -2.578 1.965 -6.640 0.785 1.211 18
## alpha[3] -1.255 1.935 -5.468 3.119 1.086 100
## alpha[4] -3.274 1.301 -6.365 -1.032 1.007 100
## alpha[5] -1.391 1.137 -3.218 1.447 1.155 21
## alpha[6] -1.976 1.356 -4.453 0.039 1.062 100
## mu -1.876 1.211 -4.314 0.706 1.020 100
## sigma 2.205 1.649 0.441 5.954 1.375 7
## tau 0.883 1.392 0.028 5.175 1.375 7
## deviance 21.145 3.724 16.609 29.087 1.016 100
##
## For each parameter, n.eff is a crude measure of effective sample size,
## and Rhat is the potential scale reduction factor (at convergence, Rhat=1).
##
## DIC info (using the rule, pD = var(deviance)/2)
## pD = 7.0 and DIC = 28.1
## DIC is an estimate of expected predictive error (lower deviance is better).
### Read and format the data (count)
data = read.csv(url("http://www.statistica.it/gianluca/bmeta/Data-count.csv"))
### List data for count outcome (for meta-analysis)
d1 <- data.list <- list(y0=data$y0,y1=data$y1,p0=data[,6],p1=data[,10])
### List data for count outcome when there is a covariate (for meta-regression)
d2 <- data.list <- list(y0=data$y0,y1=data$y1,p0=data[,6],p1=data[,10],X=cbind(data$X0))
### Select fixed-effects meta-analysis for count data
m1 <- bmeta(d1, outcome="count", model="std", type="fix",n.iter=100)
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 30
## Unobserved stochastic nodes: 16
## Total graph size: 215
##
## Initializing model
m1
## Inference for Bugs model at "model.txt", fit using jags,
## 2 chains, each with 100 iterations (first 50 discarded)
## n.sims = 100 iterations saved
## mu.vect sd.vect 2.5% 97.5% Rhat n.eff
## IRR 0.797 0.009 0.779 0.816 0.992 100
## delta -0.227 0.011 -0.250 -0.204 0.992 100
## lambda0[1] 24.689 3.640 18.345 32.786 1.050 100
## lambda0[2] 1029.544 25.989 986.234 1089.534 1.089 36
## lambda0[3] 1566.310 30.418 1508.889 1616.678 1.008 91
## lambda0[4] 284.141 10.588 263.081 306.689 1.027 51
## lambda0[5] 352.101 12.030 330.197 373.815 0.991 100
## lambda0[6] 290.572 10.046 273.553 308.090 1.045 73
## lambda0[7] 216.595 9.250 201.809 233.810 1.013 94
## lambda0[8] 90.492 6.623 79.264 103.044 0.997 100
## lambda0[9] 145.415 8.854 128.772 160.769 1.029 100
## lambda0[10] 792.967 21.016 750.139 828.522 1.017 100
## lambda0[11] 36.246 4.823 28.646 46.982 0.998 100
## lambda0[12] 483.162 14.542 450.556 506.476 1.014 100
## lambda0[13] 110.940 8.319 95.432 126.082 1.037 58
## lambda0[14] 10565.667 93.501 10401.174 10743.397 0.992 100
## lambda0[15] 347.512 19.743 307.360 390.587 1.086 100
## lambda1[1] 19.684 2.931 14.503 26.580 1.050 100
## lambda1[2] 842.106 21.859 802.772 890.130 1.085 33
## lambda1[3] 1239.691 25.225 1191.947 1285.130 1.005 100
## lambda1[4] 227.636 8.320 210.835 244.124 1.029 44
## lambda1[5] 416.077 13.229 393.846 441.450 0.991 100
## lambda1[6] 461.764 16.289 434.199 492.699 1.033 87
## lambda1[7] 183.400 7.664 169.810 198.516 1.015 80
## lambda1[8] 66.587 4.852 58.468 75.657 0.996 100
## lambda1[9] 115.102 7.206 101.724 129.847 1.029 100
## lambda1[10] 619.730 16.172 585.931 649.266 1.025 100
## lambda1[11] 32.450 4.222 25.635 42.078 0.998 100
## lambda1[12] 384.753 11.965 360.484 409.255 1.032 100
## lambda1[13] 83.599 6.372 71.786 95.604 1.031 63
## lambda1[14] 8366.730 86.099 8183.100 8526.024 0.997 100
## lambda1[15] 255.864 14.508 226.986 286.446 1.075 100
## deviance 441.316 5.818 432.505 451.875 1.023 60
##
## For each parameter, n.eff is a crude measure of effective sample size,
## and Rhat is the potential scale reduction factor (at convergence, Rhat=1).
##
## DIC info (using the rule, pD = var(deviance)/2)
## pD = 16.8 and DIC = 458.1
## DIC is an estimate of expected predictive error (lower deviance is better).
### Select random-effects meta-analysis with half-Cauchy prior for count data
m2 <- bmeta(d1, outcome="count", model="std.hc", type="ran",n.iter=100)
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 30
## Unobserved stochastic nodes: 34
## Total graph size: 262
##
## Initializing model
m2
## Inference for Bugs model at "model.txt", fit using jags,
## 2 chains, each with 100 iterations (first 50 discarded)
## n.sims = 100 iterations saved
## mu.vect sd.vect 2.5% 97.5% Rhat n.eff
## IRR 0.694 0.050 0.607 0.791 1.097 19
## delta[1] -0.363 0.220 -0.752 0.023 0.998 100
## delta[2] -0.712 0.047 -0.822 -0.628 1.072 52
## delta[3] -0.238 0.030 -0.303 -0.183 1.011 100
## delta[4] -0.317 0.081 -0.469 -0.157 1.001 100
## delta[5] -0.237 0.071 -0.391 -0.124 1.050 47
## delta[6] -0.076 0.074 -0.214 0.048 1.091 28
## delta[7] -0.027 0.090 -0.197 0.142 1.071 24
## delta[8] -0.613 0.149 -0.889 -0.262 1.154 15
## delta[9] -0.310 0.104 -0.483 -0.129 1.039 39
## delta[10] -0.430 0.055 -0.547 -0.326 1.120 16
## delta[11] -0.539 0.162 -0.910 -0.271 1.141 16
## delta[12] -0.227 0.056 -0.349 -0.117 1.114 19
## delta[13] -0.626 0.138 -0.920 -0.395 1.082 100
## delta[14] -0.152 0.011 -0.169 -0.132 0.994 100
## delta[15] -0.532 0.071 -0.671 -0.399 1.012 93
## gamma[1] 0.712 0.153 0.472 1.024 0.998 100
## gamma[2] 0.491 0.023 0.440 0.534 1.072 52
## gamma[3] 0.788 0.024 0.738 0.833 1.011 100
## gamma[4] 0.731 0.059 0.626 0.855 1.001 100
## gamma[5] 0.791 0.056 0.676 0.883 1.050 47
## gamma[6] 0.929 0.069 0.807 1.049 1.091 28
## gamma[7] 0.978 0.088 0.821 1.152 1.071 24
## gamma[8] 0.548 0.085 0.411 0.770 1.154 15
## gamma[9] 0.737 0.077 0.617 0.879 1.039 39
## gamma[10] 0.651 0.036 0.579 0.722 1.120 16
## gamma[11] 0.591 0.093 0.402 0.763 1.141 16
## gamma[12] 0.798 0.045 0.705 0.889 1.114 19
## gamma[13] 0.540 0.074 0.399 0.674 1.082 100
## gamma[14] 0.859 0.009 0.844 0.877 0.994 100
## gamma[15] 0.589 0.042 0.511 0.671 1.012 93
## lambda0[1] 26.580 4.460 18.048 34.718 1.012 100
## lambda0[2] 1246.092 34.389 1170.775 1316.559 0.995 100
## lambda0[3] 1576.053 33.117 1517.010 1631.247 1.034 60
## lambda0[4] 295.348 14.621 266.124 320.030 1.001 100
## lambda0[5] 354.903 18.318 316.802 388.817 1.053 100
## lambda0[6] 264.710 16.577 238.871 298.764 1.117 18
## lambda0[7] 196.299 12.043 176.078 223.059 1.152 14
## lambda0[8] 103.426 9.963 87.076 129.544 1.073 28
## lambda0[9] 149.575 11.821 128.663 172.866 1.056 76
## lambda0[10] 860.902 31.427 802.118 922.848 1.074 28
## lambda0[11] 40.428 4.845 32.077 50.364 1.233 11
## lambda0[12] 483.795 20.383 445.355 526.423 1.022 61
## lambda0[13] 128.248 10.938 109.710 148.218 1.022 100
## lambda0[14] 10211.011 78.806 10055.981 10367.531 1.078 100
## lambda0[15] 391.655 19.667 353.197 428.019 1.005 100
## lambda1[1] 18.560 3.527 12.652 24.702 1.018 100
## lambda1[2] 627.691 25.471 581.483 682.902 1.040 50
## lambda1[3] 1232.923 31.102 1171.399 1287.055 1.002 100
## lambda1[4] 216.418 14.270 193.805 244.040 1.000 100
## lambda1[5] 415.185 19.959 374.302 452.979 1.013 71
## lambda1[6] 488.576 21.769 448.097 530.675 1.013 100
## lambda1[7] 203.169 14.137 179.353 227.572 0.992 100
## lambda1[8] 51.819 6.271 43.640 65.672 1.094 30
## lambda1[9] 109.001 9.943 92.701 131.671 1.003 100
## lambda1[10] 549.066 21.525 511.410 592.957 1.066 32
## lambda1[11] 26.594 3.916 19.755 35.178 1.004 100
## lambda1[12] 385.284 17.918 353.216 414.739 1.037 38
## lambda1[13] 64.984 7.043 51.864 76.660 1.032 100
## lambda1[14] 8714.643 87.627 8554.025 8878.344 1.020 100
## lambda1[15] 212.567 13.788 188.698 239.193 1.011 76
## mu -0.368 0.072 -0.499 -0.235 1.097 19
## sigma 0.238 0.046 0.166 0.355 0.996 100
## tau 19.681 7.705 7.953 36.096 0.996 100
## deviance 252.280 6.896 239.140 265.915 1.257 10
##
## For each parameter, n.eff is a crude measure of effective sample size,
## and Rhat is the potential scale reduction factor (at convergence, Rhat=1).
##
## DIC info (using the rule, pD = var(deviance)/2)
## pD = 21.2 and DIC = 273.5
## DIC is an estimate of expected predictive error (lower deviance is better).
### Select random-effects meta-regression with uniform prior for count data
m3 <- bmeta(d2, outcome="count", model="reg.unif", type="ran",n.iter=100)
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 30
## Unobserved stochastic nodes: 33
## Total graph size: 287
##
## Initializing model
m3
## Inference for Bugs model at "model.txt", fit using jags,
## 2 chains, each with 100 iterations (first 50 discarded)
## n.sims = 100 iterations saved
## mu.vect sd.vect 2.5% 97.5% Rhat n.eff
## IRR 0.702 0.058 0.611 0.855 1.046 50
## delta[1] -0.353 0.229 -0.858 0.046 0.995 100
## delta[2] -0.704 0.051 -0.795 -0.613 1.014 71
## delta[3] -0.238 0.043 -0.315 -0.162 1.019 100
## delta[4] -0.345 0.095 -0.514 -0.159 0.998 100
## delta[5] -0.229 0.073 -0.367 -0.102 0.992 100
## delta[6] -0.083 0.062 -0.215 0.019 1.060 100
## delta[7] -0.037 0.094 -0.227 0.131 1.056 29
## delta[8] -0.657 0.128 -0.889 -0.352 0.998 100
## delta[9] -0.312 0.109 -0.509 -0.089 0.994 100
## delta[10] -0.433 0.044 -0.518 -0.359 1.009 100
## delta[11] -0.548 0.211 -0.994 -0.193 1.129 16
## delta[12] -0.221 0.071 -0.341 -0.074 1.003 100
## delta[13] -0.680 0.128 -0.914 -0.446 1.035 39
## delta[14] -0.154 0.014 -0.180 -0.128 0.992 100
## delta[15] -0.529 0.069 -0.658 -0.397 0.992 100
## gamma[1] 0.720 0.160 0.424 1.047 0.995 100
## gamma[2] 0.495 0.025 0.452 0.542 1.014 71
## gamma[3] 0.789 0.034 0.730 0.851 1.019 100
## gamma[4] 0.711 0.068 0.598 0.853 0.998 100
## gamma[5] 0.798 0.058 0.693 0.903 0.992 100
## gamma[6] 0.922 0.056 0.806 1.019 1.060 100
## gamma[7] 0.968 0.090 0.797 1.140 1.056 29
## gamma[8] 0.523 0.070 0.411 0.704 0.998 100
## gamma[9] 0.736 0.081 0.601 0.915 0.994 100
## gamma[10] 0.649 0.029 0.596 0.698 1.009 100
## gamma[11] 0.591 0.121 0.370 0.824 1.129 16
## gamma[12] 0.804 0.058 0.711 0.928 1.003 100
## gamma[13] 0.511 0.066 0.401 0.641 1.035 39
## gamma[14] 0.857 0.012 0.836 0.880 0.992 100
## gamma[15] 0.591 0.041 0.518 0.673 0.992 100
## lambda0[1] 26.694 5.001 19.521 37.125 1.003 100
## lambda0[2] 1240.718 33.865 1177.320 1298.207 1.162 13
## lambda0[3] 1574.808 41.183 1495.213 1650.101 1.034 55
## lambda0[4] 299.153 19.255 260.959 331.961 0.995 100
## lambda0[5] 354.619 18.015 323.913 387.707 0.993 100
## lambda0[6] 267.419 12.302 243.965 289.862 1.025 100
## lambda0[7] 198.518 12.875 172.751 223.064 1.000 100
## lambda0[8] 106.573 9.800 85.932 121.577 0.994 100
## lambda0[9] 150.250 11.026 127.294 170.806 1.022 61
## lambda0[10] 861.271 22.949 816.955 912.516 1.033 41
## lambda0[11] 40.768 5.592 31.760 51.295 1.024 56
## lambda0[12] 483.759 23.575 441.504 521.464 0.997 100
## lambda0[13] 132.469 11.072 113.116 152.983 1.009 87
## lambda0[14] 10218.774 102.444 10004.462 10403.188 0.991 100
## lambda0[15] 388.709 15.708 363.336 417.894 1.001 100
## lambda1[1] 18.720 3.298 13.508 25.711 1.028 45
## lambda1[2] 629.706 24.742 584.233 677.168 0.994 100
## lambda1[3] 1232.682 41.005 1158.316 1307.613 1.006 100
## lambda1[4] 212.972 13.995 187.100 239.683 1.008 100
## lambda1[5] 418.168 19.636 382.863 455.211 0.992 100
## lambda1[6] 490.480 22.700 441.113 526.949 1.041 100
## lambda1[7] 203.426 14.763 174.729 229.883 1.051 32
## lambda1[8] 51.144 6.058 40.982 63.221 0.997 100
## lambda1[9] 109.279 8.954 92.866 124.030 1.000 100
## lambda1[10] 547.814 20.671 507.267 585.322 0.992 100
## lambda1[11] 26.601 4.460 18.109 36.750 1.082 22
## lambda1[12] 387.657 22.039 346.869 432.659 0.996 100
## lambda1[13] 63.556 6.405 52.470 76.279 1.012 75
## lambda1[14] 8700.385 91.335 8524.475 8872.235 0.994 100
## lambda1[15] 211.786 11.792 188.213 233.996 1.003 100
## mu -0.358 0.081 -0.492 -0.158 1.046 50
## sigma 0.290 0.086 0.182 0.556 0.997 100
## tau 14.630 7.359 3.233 30.065 0.997 100
## deviance 252.614 7.190 240.456 268.760 1.160 18
##
## For each parameter, n.eff is a crude measure of effective sample size,
## and Rhat is the potential scale reduction factor (at convergence, Rhat=1).
##
## DIC info (using the rule, pD = var(deviance)/2)
## pD = 24.6 and DIC = 277.2
## DIC is an estimate of expected predictive error (lower deviance is better).
### Read and format the data (binary)
data = read.csv(url("http://www.statistica.it/gianluca/bmeta/Data-bin.csv"))
### List data for binary outcome
data.list <- list(y0=data$y0,y1=data$y1,n0=data$n0,n1=data$n1)
### generate output using bmeta
x <- bmeta(data=data.list,outcome="bin",model="std.norm",type="fix")
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 18
## Unobserved stochastic nodes: 10
## Total graph size: 95
##
## Initializing model
### run the diagnostic plot to examine the Gelman-Rubin statistic
diag.plot(x)
### run the diagnostic plot to examine the effective sample size
diag.plot(x,diag="n.eff")
### Read and format the data (binary)
data = read.csv(url("http://www.statistica.it/gianluca/bmeta/Data-bin.csv"))
### List data for binary outcome
data.list <- list(y0=data$y0,y1=data$y1,n0=data$n0,n1=data$n1)
### Select fixed-effects meta-analysis with normal prior for binary data
x <- bmeta(data.list, outcome="bin", model="std.norm", type="fix")
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 18
## Unobserved stochastic nodes: 10
## Total graph size: 95
##
## Initializing model
### Plot forest plot
forest.plot(x)
### Plot forest plot on log scale
forest.plot(x,log=TRUE)
### Select random-effects meta-analysis with t-distribution prior for binary
### data
x <- bmeta(data.list, outcome="bin", model="std.dt", type="ran")
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 18
## Unobserved stochastic nodes: 21
## Total graph size: 123
##
## Initializing model
### Plot 'two-line' forest plot showing estimates from both randome-effects
### model and no-pooling effects model for comparison
forest.plot(x,add.null=TRUE,title="Two-line forestplot for comparison")
### Read and format the data (continuous)
data = read.csv(url("http://www.statistica.it/gianluca/bmeta/Data-ctns.csv"))
### List data for continuous outcome
data.list <- list(y0=data$y0,y1=data$y1,se0=data$se0,se1=data$se1)
### Select fix-effects meta-analysis for studies reporting two arms separately
x <- bmeta(data=data.list,outcome="ctns",model="std.ta",type="fix")
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 12
## Unobserved stochastic nodes: 7
## Total graph size: 77
##
## Initializing model
### Define for individual studies
study.label <- c(paste0(data$study,", ",data$year),"Summary estimate")
### Produce forest plot with label for each study and control the lower and upper
### limits for clipping credible intervals to arrows
forest.plot(x,study.label=study.label,clip=c(-7,4))
### Read and format the data (binary)
data = read.csv(url("http://www.statistica.it/gianluca/bmeta/Data-bin.csv"))
### List data for binary outcome
data.list <- list(y0=data$y0,y1=data$y1,n0=data$n0,n1=data$n1)
### Select random-effects meta-analysis with t-distribution prior for binary
### data
x <- bmeta(data.list, outcome="bin", model="std.dt", type="ran")
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 18
## Unobserved stochastic nodes: 21
## Total graph size: 123
##
## Initializing model
### using output from bmeta to produce funnel plot
funnel.plot(x)
### using output from bmeta and specify title of the plot
funnel.plot(x,title="funnel plot")
### using output from bmeta and specify the limit of x-axis and title
funnel.plot(x,title="funnel plot",xlim=c(-2,1))
### Read and format the data (binary)
data = read.csv(url("http://www.statistica.it/gianluca/bmeta/Data-bin.csv"))
### List data for binary outcome
data.list <- list(y0=data$y0,y1=data$y1,n0=data$n0,n1=data$n1)
### Select random-effects meta-analysis with t-distribution prior for binary
### data
x <- bmeta(data.list, outcome="bin", model="std.dt", type="ran")
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 18
## Unobserved stochastic nodes: 21
## Total graph size: 123
##
## Initializing model
### using output from bmeta to produce posterior plot
posterior.plot(x)
### using output from bmeta and specify the horizontal limits
posterior.plot(x,xlim=c(-2,1))
### using output from bmeta on natural scale and specify more options
posterior.plot(x,xlim=c(-0.5,2.5),xlab="odds ratio",main="Posterior distribution
of pooled odds ratio", scale="exp")
### examine heterogeneity by producing posterior plot for between-study standard
### deviation
posterior.plot(x,heterogeneity=TRUE,xlim=c(0,3),xlab="between-study standard
deviation")
### Read and format the data (binary)
data = read.csv(url("http://www.statistica.it/gianluca/bmeta/Data-bin.csv"))
### List data for binary outcome
data.list <- list(y0=data$y0,y1=data$y1,n0=data$n0,n1=data$n1)
### Select random-effects meta-analysis with t-distribution prior for binary
### data
x <- bmeta(data.list, outcome="bin", model="std.dt", type="ran")
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 18
## Unobserved stochastic nodes: 21
## Total graph size: 123
##
## Initializing model
### using output from bmeta to produce traceplot for a specific node
traceplot.bmeta(x,"mu")
### using output from bmeta to produce traceplot and specify the node used
traceplot.bmeta(x,"mu",lab="mu")
