[R] Overlaying density plot on forest plot

Viechtbauer Wolfgang (STAT) wolfgang.viechtbauer at maastrichtuniversity.nl
Tue Dec 13 01:37:55 CET 2011


Here is an example. It requires a bit of extra work, but that's the beauty of R -- you can essentially customize a graph to your taste with low-level functions. I also added the prediction interval using the addcred argument (as suggested by Michael), since it gives the same information as that fancy normal distribution at the bottom. Essentially, that distribution is just window-dressing. I am not sure if I would go as far as calling that "chart-junk", but it's certainly somewhat gratuitous. Well, I hope the code below gives you enough of an idea so that you can adjust this to your specific case.


library(metafor)

### load BCG vaccine data
data(dat.bcg)

### meta-analysis of the log relative risks using a random-effects model
res <- rma(ai=tpos, bi=tneg, ci=cpos, di=cneg, data=dat.bcg, measure="RR", 
           slab=paste(author, year, sep=", "), method="REML")

### forest plot with extra annotations
forest(res, slab=paste(dat.bcg$author, dat.bcg$year, sep=", "), 
       xlim=c(-16, 6), at=log(c(.05, .25, 1, 4)), atransf=exp,
       ilab=cbind(dat.bcg$tpos, dat.bcg$tneg, dat.bcg$cpos, dat.bcg$cneg), 
       ilab.xpos=c(-9.5,-8,-6,-4.5), cex=.7, ylim=c(-3,16), addcred=TRUE)
text(c(-9.5,-8,-6,-4.5), 15, c("TB+", "TB-", "TB+", "TB-"), font=2, cex=.70)
text(c(-8.75,-5.25),     16, c("Vaccinated", "Control"), font=2, cex=.70)
text(-16,                15, "Author(s) and Year",     pos=4, font=2, cex=.70)
text(6,                  15, "Relative Risk [95% CI]", pos=2, font=2, cex=.70)
text(-16,              -2.5, "Prediction Distribution", pos=4, cex=.70)

### add prediction distribution
res <- predict(res)
pred.m <- res$pred
pred.s <- (res$cr.ub - res$cr.lb) / (2*1.96)
xs <- seq(-2.5,1,.01)
cords.xs <- c(-2,xs,-2)
cords.ys <- c(0, 1.5*dnorm(xs, pred.m, sd=pred.s), 0) - 3
polygon(cords.xs, cords.ys, col="darkgray", border=NA)
ptext <- paste(round(exp(res$pred),2), " [ ", round(exp(res$cr.lb),2), " , ", round(exp(res$cr.ub),2), "]", sep="")
text(6, -2.5, ptext, pos=2, cex=.70)

--   
Wolfgang Viechtbauer, Ph.D., Statistician   
Department of Psychiatry and Psychology   
School for Mental Health and Neuroscience   
Faculty of Health, Medicine, and Life Sciences   
Maastricht University, P.O. Box 616 (VIJV1)   
6200 MD Maastricht, The Netherlands   
+31 (43) 368-5248 | http://www.wvbauer.com   

> -----Original Message-----
> From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org]
> On Behalf Of Michael Dewey
> Sent: Monday, December 12, 2011 18:50
> To: Frank Peter; r-help at r-project.org
> Subject: Re: [R] Overlaying density plot on forest plot
> 
> At 15:09 11/12/2011, you wrote:
> >Dear Michael,
> >
> >Thanks for the email. This is the kind of forest plot, I want to
> >replicate
> >
> >http://www.biomedcentral.com/1471-2334/11/122/figure/F2
> 
> It would be helpful to cc to r-help in case someone else there knows
> better.
> 
> If you follow my suggestion you get a similar plot to sub-plot A in your
> example but the credible interval is shown as a dashed line superimposed
> on the summary diamond and not as a normal density as shown in sub-plot A.
> Whether that is what you want is up to you I think.
> 
> >Regards
> >Frank Peter
> >
> >-------- Original Message --------
> >From: Michael Dewey <info at aghmed.fsnet.co.uk>
> >To: "Frank Peter" <frankpeter at Safe-mail.net>, r-help at r-project.org
> >Subject: Re: [R] Overlaying density plot on forest plot
> >Date: Sun, 11 Dec 2011 14:20:13 +0000
> >
> > > At 07:16 10/12/2011, Frank Peter wrote:
> > > >Dear R User,
> > > >
> > > >Please, I am new to R. I want to overlay density plot for
> > > >predictive interval pooled result in meta-analysis.
> > > >http://addictedtor.free.fr/graphiques/graphcode.php?graph=114
> > >
> > > It is hard to be sure from your rather brief question but does the
> > > addcred parameter to forest.rma in package metafor do what you want?
> > >
> > > >
> > > >
> > > >Regards
> > > >Frank Peter
> > >
> > > Michael Dewey
> > > info at aghmed.fsnet.co.uk
> > > http://www.aghmed.fsnet.co.uk/home.html



More information about the R-help mailing list