\documentclass{chapman} %%% copy Sweave.sty definitions %%% keeps `sweave' from adding `\usepackage{Sweave}': DO NOT REMOVE %\usepackage{Sweave} \RequirePackage[T1]{fontenc} \RequirePackage{graphicx,ae,fancyvrb} \IfFileExists{upquote.sty}{\RequirePackage{upquote}}{} \usepackage{relsize} \DefineVerbatimEnvironment{Sinput}{Verbatim}{} \DefineVerbatimEnvironment{Soutput}{Verbatim}{fontfamily=courier, fontshape=it, fontsize=\relsize{-1}} \DefineVerbatimEnvironment{Scode}{Verbatim}{} \newenvironment{Schunk}{}{} %%% environment for raw output \newcommand{\SchunkRaw}{\renewenvironment{Schunk}{}{} \DefineVerbatimEnvironment{Soutput}{Verbatim}{fontfamily=courier, fontshape=it, fontsize=\small} \rawSinput } %%% environment for labeled output \newcommand{\nextcaption}{} \newcommand{\SchunkLabel}{ \renewenvironment{Schunk}{\begin{figure}[ht] }{\caption{\nextcaption} \end{figure} } \DefineVerbatimEnvironment{Sinput}{Verbatim}{frame = topline} \DefineVerbatimEnvironment{Soutput}{Verbatim}{frame = bottomline, samepage = true, fontfamily=courier, fontshape=it, fontsize=\relsize{-1}} } %%% S code with line numbers \DefineVerbatimEnvironment{Sinput} {Verbatim} { %% numbers=left } \newcommand{\numberSinput}{ \DefineVerbatimEnvironment{Sinput}{Verbatim}{numbers=left} } \newcommand{\rawSinput}{ \DefineVerbatimEnvironment{Sinput}{Verbatim}{} } %%% R / System symbols \newcommand{\R}{\textsf{R}} \newcommand{\rR}{{R}} \renewcommand{\S}{\textsf{S}} \newcommand{\SPLUS}{\textsf{S-PLUS}} \newcommand{\rSPLUS}{{S-PLUS}} \newcommand{\SPSS}{\textsf{SPSS}} \newcommand{\EXCEL}{\textsf{Excel}} \newcommand{\ACCESS}{\textsf{Access}} \newcommand{\SQL}{\textsf{SQL}} %%\newcommand{\Rpackage}[1]{\hbox{\rm\textit{#1}}} %%\newcommand{\Robject}[1]{\hbox{\rm\texttt{#1}}} %%\newcommand{\Rclass}[1]{\hbox{\rm\textit{#1}}} %%\newcommand{\Rcmd}[1]{\hbox{\rm\texttt{#1}}} \newcommand{\Rpackage}[1]{\index{#1 package@{\fontseries{b}\selectfont #1} package} {\fontseries{b}\selectfont #1}} \newcommand{\rpackage}[1]{{\fontseries{b}\selectfont #1}} \newcommand{\Robject}[1]{\texttt{#1}} \newcommand{\Rclass}[1]{\index{#1 class@\textit{#1} class}\textit{#1}} \newcommand{\Rcmd}[1]{\index{#1 function@\texttt{#1} function}\texttt{#1}} \newcommand{\Roperator}[1]{\texttt{#1}} \newcommand{\Rarg}[1]{\texttt{#1}} \newcommand{\Rlevel}[1]{\texttt{#1}} %%% other symbols \newcommand{\file}[1]{\hbox{\rm\texttt{#1}}} %%\newcommand{\stress}[1]{\index{#1}\textit{#1}} \newcommand{\stress}[1]{\textit{#1}} \newcommand{\booktitle}[1]{\textit{#1}} %%' %%% Math symbols \usepackage{amstext} \usepackage{amsmath} \newcommand{\E}{\mathsf{E}} \newcommand{\Var}{\mathsf{Var}} \newcommand{\Cov}{\mathsf{Cov}} \newcommand{\Cor}{\mathsf{Cor}} \newcommand{\x}{\mathbf{x}} \newcommand{\y}{\mathbf{y}} \renewcommand{\a}{\mathbf{a}} \newcommand{\W}{\mathbf{W}} \newcommand{\C}{\mathbf{C}} \renewcommand{\H}{\mathbf{H}} \newcommand{\X}{\mathbf{X}} \newcommand{\B}{\mathbf{B}} \newcommand{\V}{\mathbf{V}} \newcommand{\I}{\mathbf{I}} \newcommand{\D}{\mathbf{D}} \newcommand{\bS}{\mathbf{S}} \newcommand{\N}{\mathcal{N}} \renewcommand{\L}{L} \renewcommand{\P}{\mathsf{P}} \newcommand{\K}{\mathbf{K}} \newcommand{\m}{\mathbf{m}} \newcommand{\argmin}{\operatorname{argmin}\displaylimits} \newcommand{\argmax}{\operatorname{argmax}\displaylimits} \newcommand{\bx}{\mathbf{x}} \newcommand{\bbeta}{\mathbf{\beta}} %%% links \usepackage{hyperref} \hypersetup{% pdftitle = {A Handbook of Statistical Analyses Using R (3rd Edition)}, pdfsubject = {Book}, pdfauthor = {Torsten Hothorn and Brian S. Everitt}, colorlinks = {black}, linkcolor = {black}, citecolor = {black}, urlcolor = {black}, hyperindex = {true}, linktocpage = {true}, } %%% captions & tables %% : conflics with figure definition in chapman.cls %%\usepackage[format=hang,margin=10pt,labelfont=bf]{caption} %% \usepackage{longtable} \usepackage[figuresright]{rotating} %%% R symbol in chapter 1 \usepackage{wrapfig} %%% Bibliography \usepackage[round,comma]{natbib} \renewcommand{\refname}{References \addcontentsline{toc}{chapter}{References}} \citeindexfalse %%% texi2dvi complains that \newblock is undefined, hm... \def\newblock{\hskip .11em plus .33em minus .07em} %%% Example sections \newcounter{exercise}[chapter] \setcounter{exercise}{0} \newcommand{\exercise}{\stepcounter{exercise} \item{Ex.~\arabic{chapter}.\arabic{exercise} }} %% URLs \newcommand{\curl}[1]{\begin{center} \url{#1} \end{center}} %%% for manual corrections %\renewcommand{\baselinestretch}{2} %%% plot sizes \setkeys{Gin}{width=0.95\textwidth} %%% color \usepackage{color} %%% hyphenations \hyphenation{drop-out} \hyphenation{mar-gi-nal} %%% new bidirectional quotes need \usepackage[utf8]{inputenc} %\usepackage{setspace} \definecolor{sidebox_todo}{rgb}{1,1,0.2} \newcommand{\todo}[1]{ \hspace{0pt}% \marginpar{% \fcolorbox{black}{sidebox_todo}{% \parbox{\marginparwidth} { \raggedright\sffamily\footnotesize{TODO: #1}% } }% } } \begin{document} %% Title page \title{A Handbook of Statistical Analyses Using \R{} --- 3rd Edition} \author{Torsten Hothorn and Brian S. Everitt} \maketitle %%\VignetteIndexEntry{Chapter Meta-Analysis} %%\VignetteDepends{rmeta} \setcounter{chapter}{16} \SweaveOpts{prefix.string=figures/HSAUR,eps=FALSE,keep.source=TRUE} <>= rm(list = ls()) s <- search()[-1] s <- s[-match(c("package:base", "package:stats", "package:graphics", "package:grDevices", "package:utils", "package:datasets", "package:methods", "Autoloads"), s)] if (length(s) > 0) sapply(s, detach, character.only = TRUE) if (!file.exists("tables")) dir.create("tables") if (!file.exists("figures")) dir.create("figures") set.seed(290875) options(prompt = "R> ", continue = "+ ", width = 63, # digits = 4, show.signif.stars = FALSE, SweaveHooks = list(leftpar = function() par(mai = par("mai") * c(1, 1.05, 1, 1)), bigleftpar = function() par(mai = par("mai") * c(1, 1.7, 1, 1)))) HSAURpkg <- require("HSAUR3") if (!HSAURpkg) stop("cannot load package ", sQuote("HSAUR3")) rm(HSAURpkg) ### hm, R-2.4.0 --vanilla seems to need this a <- Sys.setlocale("LC_ALL", "C") ### book <- TRUE refs <- cbind(c("AItR", "DAGD", "SI", "CI", "ANOVA", "MLR", "GLM", "DE", "RP", "GAM", "SA", "ALDI", "ALDII", "SIMC", "MA", "PCA", "MDS", "CA"), 1:18) ch <- function(x) { ch <- refs[which(refs[,1] == x),] if (book) { return(paste("Chapter~\\\\ref{", ch[1], "}", sep = "")) } else { return(paste("Chapter~", ch[2], sep = "")) } } if (file.exists("deparse.R")) source("deparse.R") setHook(packageEvent("lattice", "attach"), function(...) { lattice.options(default.theme = function() standard.theme("pdf", color = FALSE)) }) @ \pagestyle{headings} <>= book <- FALSE @ \chapter[Meta-Analysis]{Meta-Analysis: Nicotine Gum and Smoking Cessation and the Efficacy of BCG Vaccine in the Treatment of Tuberculosis \label{MA}} \section{Introduction} \section{Systematic Reviews and Meta-Analysis} \section{Analysis Using \R{}} The aim in collecting the results from the randomized trials of using nicotine gum to help smokers quit was to estimate the overall \stress{odds ratio}, the odds of quitting smoking for those given the gum, divided by the odds of quitting for those not receiving the gum. Following formula (\ref{MA:barY}), we can compute the pooled odds ratio as follows: <>= data("smoking", package = "HSAUR3") odds <- function(x) (x[1] * (x[4] - x[3])) / ((x[2] - x[1]) * x[3]) weight <- function(x) ((x[2] - x[1]) * x[3]) / sum(x) W <- apply(smoking, 1, weight) Y <- apply(smoking, 1, odds) sum(W * Y) / sum(W) @ Of course, the computations are more conveniently done using the functionality provided in package \Rpackage{rmeta}. The odds ratios and corresponding confidence intervals are computed by means of the \Rcmd{meta.MH} function for fixed effects meta-analysis as shown here <>= library("rmeta") smokingOR <- meta.MH(smoking[["tt"]], smoking[["tc"]], smoking[["qt"]], smoking[["qc"]], names = rownames(smoking)) @ and the results can be inspected via a \Rcmd{summary} method -- see Figure~\ref{MA-smoking-summary}. \renewcommand{\nextcaption}{\R{} output of the \Rcmd{summary} method for \Robject{smokingOR}. \label{MA-smoking-summary}} \SchunkLabel <>= summary(smokingOR) @ \SchunkRaw \begin{figure} \begin{center} <>= plot(smokingOR, ylab = "") @ \caption{Forest plot of observed effect sizes and $95\%$ confidence intervals for the nicotine gum studies. \label{MA:smokingplot}} \end{center} \end{figure} We shall use both the fixed effects and random effects approaches here so that we can compare results. For the fixed effects model (see Figure~\ref{MA-smoking-summary}) the estimated overall log-odds ratio is \Sexpr{round(smokingOR$logMH, 3)} with a standard error of \Sexpr{round(smokingOR$selogMH, 3)}. This leads to an estimate of the overall odds ratio of \Sexpr{round(exp(smokingOR$logMH), 3)}, with a 95\% confidence interval as given above. For the random effects model, which is fitted by applying function \Rcmd{meta.DSL} to the \Robject{smoking} data as follows \vspace{1cm} <>= (smokingDSL <- meta.DSL(smoking[["tt"]], smoking[["tc"]], smoking[["qt"]], smoking[["qc"]], names = rownames(smoking))) @ the corresponding estimate is \Sexpr{round(exp(smokingDSL$logDSL), 3)}. Both models suggest that there is clear evidence that nicotine gum increases the odds of quitting. The random effects confidence interval is considerably wider than that from the fixed effects model; here the test of homogeneity of the studies is not significant implying that we might use the fixed effects results. But the test is not particularly powerful and it is more sensible to assume a priori that heterogeneity is present and so we use the results from the random effects model. \section{Meta-Regression} The examination of heterogeneity of the effect sizes from the studies in a meta-analysis begins with the formal test for its presence, although in most meta-analyses such heterogeneity can almost be assumed to be present. There will be many possible sources of such heterogeneity and estimating how these various factors affect the observed effect sizes in the studies chosen is often of considerable interest and importance, indeed usually more important than the relatively simplistic use of meta-analysis to determine a single summary estimate of overall effect size. We can illustrate the process using the BCG vaccine data. We first find the estimate of the overall effect size from applying the fixed effects and the random effects models described previously: <>= data("BCG", package = "HSAUR3") BCG_OR <- meta.MH(BCG[["BCGVacc"]], BCG[["NoVacc"]], BCG[["BCGTB"]], BCG[["NoVaccTB"]], names = BCG$Study) BCG_DSL <- meta.DSL(BCG[["BCGVacc"]], BCG[["NoVacc"]], BCG[["BCGTB"]], BCG[["NoVaccTB"]], names = BCG$Study) @ The results are inspected using the \Rcmd{summary} method as shown in Figures~\ref{MA-BCGOR-summary} and \ref{MA-BCGDSL-summary}. \renewcommand{\nextcaption}{\R{} output of the \Rcmd{summary} method for \Robject{BCG\_OR}. \label{MA-BCGOR-summary}} \SchunkLabel <>= summary(BCG_OR) @ \SchunkRaw \renewcommand{\nextcaption}{\R{} output of the \Rcmd{summary} method for \Robject{BCG\_DSL}. \label{MA-BCGDSL-summary}} \SchunkLabel <>= summary(BCG_DSL) @ \SchunkRaw To assess how the two covariates, latitude and year, relate to the observed effect sizes we shall use multiple linear regression but will weight each observation by $W_i = (\hat{\sigma}^2 + V_i^2)^{-1}, i = 1, \dots, 13$, where $\hat{\sigma}^2$ is the estimated between-study variance and $V_i^2$ is the estimated variance from the $i$th study. The required \R{} code to fit the linear model via weighted least squares is: \index{Meta-Analysis!weighted least squares regression} <>= studyweights <- 1 / (BCG_DSL$tau2 + BCG_DSL$selogs^2) y <- BCG_DSL$logs BCG_mod <- lm(y ~ Latitude + Year, data = BCG, weights = studyweights) @ and the results of the \Rcmd{summary} method are shown in Figure~\ref{MA-mod-summary}. There is some evidence that latitude is associated with observed effect size, the log-odds ratio becoming increasingly negative as latitude increases, as we can see from a scatterplot of the two variables with the added weighted regression fit seen in Figure~\ref{MA-BCG}. \renewcommand{\nextcaption}{\R{} output of the \Rcmd{summary} method for \Robject{BCG\_mod}. \label{MA-mod-summary}} \SchunkLabel <>= summary(BCG_mod) @ \SchunkRaw \begin{figure} \begin{center} <>= plot(y ~ Latitude, data = BCG, ylab = "Estimated log-OR") abline(lm(y ~ Latitude, data = BCG, weights = studyweights)) @ \caption{Plot of observed effect size for the \Robject{BCG} vaccine data against latitude, with a weighted least squares regression fit shown in addition. \label{MA-BCG}} \end{center} \end{figure} \section{Publication Bias} \begin{figure} \begin{center} <>= set.seed(290875) sigma <- seq(from = 1/10, to = 1, length.out = 35) y <- rnorm(35) * sigma gr <- (y > -0.5) layout(matrix(1:2, ncol = 1)) plot(y, 1/sigma, xlab = "Effect size", ylab = "1 / standard error") plot(y[gr], 1/(sigma[gr]), xlim = range(y), xlab = "Effect size", ylab = "1 / standard error") @ \caption{Example funnel plots from simulated data. The asymmetry in the lower plot is a hint that a publication bias might be a problem. \label{MA-funnel}} \end{center} \end{figure} We can construct a funnel plot for the nicotine gum data using the \R{} code depicted with Figure~\ref{MA:funnel}. There does not appear to be any strong evidence of publication bias here. \begin{figure} \begin{center} <>= funnelplot(smokingDSL$logs, smokingDSL$selogs, summ = smokingDSL$logDSL, xlim = c(-1.7, 1.7)) abline(v = 0, lty = 2) @ \caption{Funnel plot for nicotine gum data. \label{MA:funnel}} \end{center} \end{figure} \index{Meta-analysis!funnel plots|)} %%\bibliographystyle{LaTeXBibTeX/refstyle} %%\bibliography{LaTeXBibTeX/HSAUR} \end{document}