右 |在内部合并使用 par(mfrow = ...) 的图

Question

以plot.acf 为例。 acf 和 pacf 都在内部调用此函数。如何将它们并排绘制？

示例：

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

我尝试使用par(mfrow = c(2,4)) 和layout 来组合它们，但stats:::plot.acf 覆盖了它。预期的输出是：

Answer 1

与我的其他答案不同的方法：使用 ggplot2 绘制 ACF。

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
}

这旨在创建与plot.acf() 类似的界面。然后，您可以使用 gridExtra 包中的 ggplot2 绘图可用的所有强大功能。

library(ggplot2)
library(gridExtra)

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

然后你得到这个：

很遗憾，grid.arrange() 不适用于基本图形，因此建议使用ggplot2。

Answer 2

这不是一个理想的解决方案，但您可以通过定义 plot.acf() 来重新定义绘制 ACF/PACF 的含义。

先存储现有版本。

old.plot.acf <- plot.acf

现在您可以使用stats:::plot.acf 获取源代码并将其复制/粘贴到编辑器中。移除重置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()
}

现在这是在本地定义的，您可以根据需要设置mfrow，进行绘图，然后重置函数或从命名空间中清除它。

plot.acf <- old.plot.acf

为避免同时更改plot.pacf()，您可以只使用acf(..., type="partial")，它会获取PACF。

Answer 3

您可以使用PerformanceAnalytics 包：

library(PerformanceAnalytics)
chart.ACFplus(TS)