R | combine plots that use par(mfrow = ...) internally

Question

Take plot.acf as an example. Both acf and pacf call this function internally. How can i plot them side by side?

Example:

TS <- ts.union(mdeaths, fdeaths)
acf(TS)
pacf(TS)

I tried to use par(mfrow = c(2,4)) and layout to combine them, but stats:::plot.acf overwrites this. The expected output would be:

Answer 1

You can use the PerformanceAnalytics package:

library(PerformanceAnalytics)
chart.ACFplus(TS)

Answer 2

A different approach than my other answer: Plot the ACF using ggplot2.

ggacf <- function(x, ci=0.95, type="correlation", xlab="Lag", ylab=NULL,
                  ylim=NULL, main=NULL, ci.col="blue", lag.max=NULL) {

    x <- as.data.frame(x)

    x.acf <- acf(x, plot=F, lag.max=lag.max, type=type)

    ci.line <- qnorm((1 - ci) / 2) / sqrt(x.acf$n.used)

    d.acf <- data.frame(lag=x.acf$lag, acf=x.acf$acf)

    g <- ggplot(d.acf, aes(x=lag, y=acf)) +
        geom_hline(yintercept=0) +
        geom_segment(aes(xend=lag, yend=0)) +
        geom_hline(yintercept=ci.line, color=ci.col, linetype="dashed") +
        geom_hline(yintercept=-ci.line, color=ci.col, linetype="dashed") +
        theme_bw() +
        xlab("Lag") +
        ggtitle(ifelse(is.null(main), "", main)) +
        if (is.null(ylab))
            ylab(ifelse(type=="partial", "PACF", "ACF"))
        else
            ylab(ylab)

    g
}

This seeks to create a similar interface to plot.acf(). Then you can use all of the great features available to ggplot2 plots from the gridExtra package.

library(ggplot2)
library(gridExtra)

grid.arrange(ggacf(lh), ggacf(lh, type="partial"), ncol=2)

Then you get this:

Unfortunately grid.arrange() doesn't work with base graphics, hence the ggplot2 suggestion.

Answer 3

This isn't an ideal solution, but you can redefine what it means to plot an ACF/PACF by defining plot.acf().

First store the existing version.

old.plot.acf <- plot.acf

Now you can use stats:::plot.acf to get the source and copy/paste into the editor. Remove the part that resets mfrow.

plot.acf <- function(x, ci = 0.95, type = "h", xlab = "Lag", ylab = NULL,
                     ylim = NULL, main = NULL, ci.col = "blue",
                     ci.type = c("white", "ma"), max.mfrow = 6,
                     ask = Npgs > 1 && dev.interactive(), 
                     mar = if (nser > 2) c(3, 2, 2, 0.8) else par("mar"),
                     oma = if (nser > 2) c(1, 1.2, 1, 1) else par("oma"),
                     mgp = if (nser > 2) c(1.5, 0.6, 0) else par("mgp"),
                     xpd = par("xpd"), cex.main = if (nser > 2) 1 else
                     par("cex.main"), verbose = getOption("verbose"), ...) 
{
    ci.type <- match.arg(ci.type)
    if ((nser <- ncol(x$lag)) < 1L) 
        stop("x$lag must have at least 1 column")
    if (is.null(ylab)) 
        ylab <- switch(x$type, correlation = "ACF", covariance = "ACF (cov)", 
                       partial = "Partial ACF")
    if (is.null(snames <- x$snames)) 
        snames <- paste("Series ", if (nser == 1L) 
            x$series
            else 1L:nser)
    with.ci <- ci > 0 && x$type != "covariance"
    with.ci.ma <- with.ci && ci.type == "ma" && x$type == "correlation"
    if (with.ci.ma && x$lag[1L, 1L, 1L] != 0L) {
        warning("can use ci.type=\"ma\" only if first lag is 0")
        with.ci.ma <- FALSE
    }
    clim0 <- if (with.ci) 
        qnorm((1 + ci)/2)/sqrt(x$n.used)
    else c(0, 0)
    Npgs <- 1L
    nr <- nser
    if (nser > 1L) {
        sn.abbr <- if (nser > 2L) 
            abbreviate(snames)
        else snames
        if (nser > max.mfrow) {
            Npgs <- ceiling(nser/max.mfrow)
            nr <- ceiling(nser/Npgs)
        }

        ### NOT INCLUDED: mfrow = rep(nr, 2L)

        opar <- par(mar = mar, oma = oma, 
                    mgp = mgp, ask = ask, xpd = xpd, cex.main = cex.main)
        on.exit(par(opar))
        if (verbose) {
            message("par(*) : ", appendLF = FALSE, domain = NA)
            str(par("mfrow", "cex", "cex.main", "cex.axis", "cex.lab", 
                    "cex.sub"))
        }
    }
    if (is.null(ylim)) {
        ylim <- range(x$acf[, 1L:nser, 1L:nser], na.rm = TRUE)
        if (with.ci) 
            ylim <- range(c(-clim0, clim0, ylim))
        if (with.ci.ma) {
            for (i in 1L:nser) {
                clim <- clim0 * sqrt(cumsum(c(1, 2 * x$acf[-1, 
                                                           i, i]^2)))
                ylim <- range(c(-clim, clim, ylim))
            }
        }
    }
    for (I in 1L:Npgs) for (J in 1L:Npgs) {
        dev.hold()
        iind <- (I - 1) * nr + 1L:nr
        jind <- (J - 1) * nr + 1L:nr
        if (verbose) 
            message("Page [", I, ",", J, "]: i =", paste(iind, 
                                                         collapse = ","), "; j =", paste(jind, collapse = ","), 
                    domain = NA)
        for (i in iind) for (j in jind) if (max(i, j) > nser) {
            frame()
            box(col = "light gray")
        }
        else {
            clim <- if (with.ci.ma && i == j) 
                clim0 * sqrt(cumsum(c(1, 2 * x$acf[-1, i, j]^2)))
            else clim0
            plot(x$lag[, i, j], x$acf[, i, j], type = type, xlab = xlab, 
                 ylab = if (j == 1) 
                     ylab
                 else "", ylim = ylim, ...)
            abline(h = 0)
            if (with.ci && ci.type == "white") 
                abline(h = c(clim, -clim), col = ci.col, lty = 2)
            else if (with.ci.ma && i == j) {
                clim <- clim[-length(clim)]
                lines(x$lag[-1, i, j], clim, col = ci.col, lty = 2)
                lines(x$lag[-1, i, j], -clim, col = ci.col, lty = 2)
            }
            title(if (!is.null(main)) 
                main
                else if (i == j) 
                    snames[i]
                else paste(sn.abbr[i], "&", sn.abbr[j]), line = if (nser > 
                                                                        2) 
                    1
                else 2)
        }
        if (Npgs > 1) {
            mtext(paste("[", I, ",", J, "]"), side = 1, line = -0.2, 
                  adj = 1, col = "dark gray", cex = 1, outer = TRUE)
        }
        dev.flush()
    }
    invisible()
}

Now that this is defined locally, you can set mfrow as needed, do your plotting, then reset the function or clear it from the namespace.

plot.acf <- old.plot.acf

To avoid also having to change plot.pacf() as well, you can just use acf(..., type="partial"), which gets the PACF.