[R-meta] Antw: RE: Network Meta-analysis metafor

[Caroline Woehl]

Dear Wolfgang,

it works great, thanks a lot for that!

I've already created a forest plot and calculated some subgroups, for
example based on intervention duration or cognitive function at
baseline. Could you give me some advice on how to add two subgroups in
one forest plot? Is there any possibility?

Thanks again for your help.

Best regards,
Caroline

>>> "Viechtbauer Wolfgang (SP)"
<wolfgang.viechtbauer at maastrichtuniversity.nl> 14.02.18 13.16 Uhr >>>
Dear Caroline,

The example in help(dat.senn2013) uses mean differences to illustrate a
network MA with metafor but is easy to adjust to standardized mean
differences. For SMDs, you have to use a contrast-based model. The
dat.senn2013 dataset is given in an arm-based format, so it first needs
to be restructured. Note that Willms (1999) is the only multiarm study
in this dataset, but the code should work regardless of how many
multiarm studies there are. Just make sure that the 'reference'
treatment for multiarm studies always goes under 'trt1'.

library(metafor)

############################################################################

### load data
dat <- get(data(dat.senn2013, package="metafor"))

### restructure dataset to a contrast-based format
dat <- dat.senn2013[,c(1,4,3,2,5,6)]
dat.c <- lapply(split(dat, dat$study),
                function(x) cbind(x[rep(1,nrow(x)-1),], x[-1,c(3:6)]))
dat.c <- do.call(rbind, dat.c)
names(dat.c)[3:10] <- c("trt1", "n1i", "m1i", "sd1i", "trt2", "n2i",
"m2i", "sd2i")
rownames(dat.c) <- 1:nrow(dat.c)
dat.c$id <- 1:nrow(dat.c)
dat <- dat.c
rm(dat.c)

############################################################################

### make sure dataset is sorted by 'study'
dat <- dat[order(dat$study),]

### compute the total sample sizes of the studies
dat$Ni <- unlist(lapply(split(dat, dat$study), function(x) rep(x$n1i[1]
+ sum(x$n2i), each=nrow(x))))

### compute the pooled SDs of the studies
dat$sdpi <- unlist(lapply(split(dat, dat$study), function(x)
sqrt(((x$n1i[1]-1)*x$sd1i^2 + sum((x$n2i-1)*x$sd2i^2)) / ((x$n1i[1] - 1)
+ sum(x$n2i - 1)))))

### calculate SMDs and corresponding sampling variances
dat$yi <- with(dat, (m1i-m2i)/sdpi)
dat$vi <- with(dat, 1/n1i + 1/n2i + yi^2/(2*Ni))

### calculate the variance-covariance matrix of the SMDs for
multitreatment studies
calc.v <- function(x) {
   v <- matrix(1/x$n1i[1] + outer(x$yi, x$yi, "*")/(2*x$Ni[1]),
nrow=nrow(x), ncol=nrow(x))
   diag(v) <- x$vi
   v
}
V <- bldiag(lapply(split(dat, dat$study), calc.v))

### convert trt1 and trt2 variables to factors and set levels
lvls <- levels(factor(c(dat$trt1, dat$trt2)))
dat$trt1 <- factor(dat$trt1, levels=lvls)
dat$trt2 <- factor(dat$trt2, levels=lvls)

### create variables to indicate the contrasts examined
dat <- cbind(dat, model.matrix(~ dat$trt1 - 1) - model.matrix(~ dat$trt2
- 1))
names(dat)[(ncol(dat)-9):ncol(dat)] <- lvls

############################################################################

### network meta-analysis using a contrast-based random-effects model
### by setting rho=1/2, tau^2 reflects the amount of heterogeneity for
all treatment comparisons
### the treatment left out (placebo) becomes the reference level for the
treatment comparisons
res <- rma.mv(yi, V, mods = ~ acarbose + benfluorex + metformin +
miglitol + pioglitazone +
                              rosiglitazone + sitagliptin + sulfonylurea
+ vildagliptin - 1,
              random = ~ factor(id) | study, rho=1/2, data=dat)
res

### estimate/test all pairwise differences between treatments
### need 'multcomp' package for contrMat() function
library(multcomp)
contr <- contrMat(setNames(rep(1,res$p), colnames(res$X)), type="Tukey")
sav <- predict(res, newmods=contr)
sav[["slab"]] <- rownames(contr)
sav$pval <- anova(res, L=contr)$pval
sav

############################################################################

Two points:

1) Note that I am using the pooled SD of all arms within each trial to
compute the SMDs. That way, results within triIf you just use the pooled SD of each pairwise comparison, then this
won't be the case.

2) I see that you have the reference treatment listed under treat2. If
you want to use the code above, just rename your variables to: study
trt2 trt1 n2i m2i sd2i n1i m1i sd1i.

Best,
Wolfgang

-- 
Wolfgang Viechtbauer, Ph.D., Statistician | Department of Psychiatry and

Neuropsychology | Maastricht University | P.O. Box 616 (VIJV1) | 6200 MD

Maastricht, The Netherlands | +31 (43) 388-4170 | http://www.wvbauer.com


>-----Original Message-----
>From: R-sig-meta-analysis [mailto:r-sig-meta-analysis-bounces at r-
>project.org] On Behalf Of Caroline Woehl
>Sent: Wednesday, 14 February, 2018 8:45
>To: r-sig-meta-analysis at r-project.org
>Subject: [R-meta] Network Meta-analysis metafor
>
>Dear all,
>
>I've only calucated network meta-analyses with the package netmeta yet.
>Due to the possibility to calculate subgroup analyses, I'd like to have
>a try at a network meta-analysis with metafor.
>
>The analysis is based on 20 randomized controlled trials. They compare
>physical activity or cognitive activity with usual care. One of them is
>a multi-arm trial comparing physical and cognitive activity with usual
>care. I could extract sample size, mean und standard deviation
>(posttreatment) from each study arm. Thus, the file contains the study,
>treat1, treat2, n1i, m1i, sd1i, n2i, m2i and sd2i.
>
>Study                                                treat1      
treat2
>        n1i          m1i      sd1i        n2i          m2i       sd2i
>       1      physical activity    usual care     37   17,16    4,33
>36   15,17    4,48
>       2                    physical activity    usual care     39
>19,40    8,68    35   18,48    8,40
>       2                    cognitive activity   usual care    36
>19,98    8,16    35   18,48    8,40
>       3      cognitive activity usual care    35    42,34   9,30    33
>  44,49   9,53
>
>So far, I’ve calculated the standardized mean differences (yi) and
>variances (vi) for each pairwise comparison.
>nma <- escalc(measure="SMD", m1i=m1i, sd1i=sd1i, n1i=n1i, m2i=m2i,
>sd2i=sd2i, n2i=n2i, data=nma)
>
>How can I define usual care as the reference group in the next step?
>
>Many thanks for your help.
>
>Best regards,
>Caroline