用r绘制不同尺度上的脉冲响应

  • 本文关键字:脉冲 响应 绘制
  • 更新时间 :
  • 英文 :


目前我正在使用r中的vars包。

library(vars)
data(Canada)
var.2c <- VAR(Canada, p = 2, type = "const")
plot(irf(var.2c, impulse = "e", response = c("prod", "rw", "U"), boot = 
  T))

可以看到,每个图都有相同的y尺度。我如何用单独的缩放来代替这种统一的缩放?在我自己的数据集中,我遇到的问题是,一些脉冲响应的范围是[0.80:-0.80],而另一些脉冲响应的范围是从0.001到-0.001。

提前感谢!

取决于你想显示什么。注意在图表上混合缩放或缩放数据,因为这很容易误导读者。那么,您自担风险,:-),这里有两种可能性。

首先,手动重新缩放数据,即规范化y-data的每个向量。

其次,如果结果属于两个基本的缩放类别,则将一组绘制在左y轴上,另一组绘制在次(右)y轴上。查看?par?axis

你可以试试我的代码:plot。varirfGM (x,阴谋。类型="多个",数控= 2,3月。Multi = c(0,4,2,1), oma。

= c(1, 2, 4, 2)

结果示例

    plot.varirfGM=
function (x, plot.type = c("multiple", "single"), names = NULL, 
          main = NULL, sub = NULL, lty = NULL, lwd = NULL, col = NULL, 
          ylim = NULL, ylab = NULL, xlab = NULL, nc, 
          mar.multi = c(0, 4, 0, 4), oma.multi = c(6, 4, 6, 4), adj.mtext = NA, 
          padj.mtext = NA, col.mtext = NA, ...) 
{
  op <- par(no.readonly = TRUE)
  on.exit(par(op))
  plot.type <- match.arg(plot.type)
  inames <- x$impulse
  rnames <- x$response
  if (is.null(names)) {
    names <- inames
  }
  else {
    names <- as.character(names)
    if (!(all(names %in% inames))) {
      warning("nInvalid variable name(s) supplied, using first variable.n")
      inames <- inames[1]
    }
    else {
      inames <- names
    }
  }
  nvi <- length(inames)
  nvr <- length(rnames)
  ifelse(is.null(lty), lty <- c(1, 1, 2, 2), lty <- rep(lty, 
                                                        4)[1:4])
  ifelse(is.null(lwd), lwd <- c(1, 1, 1, 1), lwd <- rep(lwd, 
                                                        4)[1:4])
  ifelse(is.null(col), col <- c("black", "gray", "red", "red"), 
         col <- rep(col, 4)[1:4])
###1. dataplot function ########
   dataplot <- function(x, iname) {
    impulses <- x$irf[[iname]]
    range <- range(impulses)
    upper <- NULL
    lower <- NULL
    if (x$boot) {
      upper <- x$Upper[[iname]]
      lower <- x$Lower[[iname]]
      range <- range(cbind(impulses, upper, lower))
      rangeList=getRespRangeList(x,iname) ############## rangelist
      #print(rangeList)
    }
    if ((x$model == "varest") || (x$model == "vec2var")) {
      if (x$ortho) {
        text1 <- paste("Orthogonal Impulse Response from", 
                       iname, sep = " ")
      }
      else {
        text1 <- paste("Impulse Response from", iname, 
                       sep = " ")
      }
    }
    else if (x$model == "svarest") {
      text1 <- paste("SVAR Impulse Response from", iname, 
                     sep = " ")
    }
    else if (x$model == "svecest") {
      text1 <- paste("SVECM Impulse Response from", iname, 
                     sep = " ")
    }
    if (x$cumulative) 
      text1 <- paste(text1, "(cumulative)", sep = " ")
    text2 <- ""
    if (x$boot) 
      text2 <- paste((1 - x$ci) * 100, "% Bootstrap CI, ", 
                     x$runs, "runs")
    result <- list(impulses = impulses, upper = upper, lower = lower, 
                   range = range, text1 = text1, text2 = text2, rangeList=rangeList)
    return(result)
   }
###2. plot.single function ########
  plot.single <- function(x, iname, rname, ...) {
    ifelse(is.null(main), main <- x$text1, main <- main)
    ifelse(is.null(sub), sub <- x$text2, sub <- sub)
    xy <- xy.coords(x$impulse[, rname])
    ifelse(is.null(ylab), ylabel <- rname, ylabel <- ylab)
    ifelse(is.null(xlab), xlabel <- "", xlabel <- xlab)
    ifelse(is.null(ylim), ylim <- x$range, ylim <- ylim) #ylim assignment
    plot(xy, type = "l", ylim = ylim, col = col[1], lty = lty[1], 
         lwd = lwd[1], axes = FALSE, ylab = paste(ylabel), 
         xlab = paste(xlab), ...)
    title(main = main, sub = sub, ...)
    axis(1, at = xy$x, labels = c(0:(length(xy$x) - 1)))
    axis(2, ...)
    box()
    if (!is.null(x$upper)) 
      lines(x$upper[, rname], col = col[3], lty = lty[3], 
            lwd = lwd[3])
    if (!is.null(x$lower)) 
      lines(x$lower[, rname], col = col[3], lty = lty[3], 
            lwd = lwd[3])
    abline(h = 0, col = col[2], lty = lty[2], lwd = lwd[2])
  }
### 3. plot.multiple function ## 
  plot.multiple <- function(dp, nc = nc, ...) {
    #print(dp$rangeList)
    x <- dp$impulses
    y <- dp$upper
    z <- dp$lower
    ifelse(is.null(main), main <- dp$text1, main <- main)
    ifelse(is.null(sub), sub <- dp$text2, sub <- sub)
    ifelse(is.null(ylim), ylim <- dp$range, ylim <- ylim) #ylim assignment
    range <- range(c(x, y, z)) #unique range for each response in a dp$impulse
    rangeList=dp$rangeList
    nvr <- ncol(x)
    if (missing(nc)) {
      nc <- ifelse(nvr > 4, 2, 1)
    }
    nr <- ceiling(nvr/nc)
    par(mfrow = c(nr, nc), mar = mar.multi, oma = oma.multi)
    if (nr > 1) {
      for (i in 1:(nvr - nc)) {
        ylimL=dp$rangeList[i,]################## ylim rangeList
        #print(ylimL) #print("ylimL i")
        ifelse(is.null(ylab), ylabel <- colnames(x)[i], 
               ylabel <- ylab)
        xy <- xy.coords(x[, i])
        plot(xy, axes = FALSE, type = "l", ylab = ylabel, 
             ylim = ylimL, col = col[1], lty = lty[1], lwd = lwd[1], 
             ...)
        axis(2, at = pretty(rangeList[i,]))
        abline(h = 0, col = "red")
        if (!is.null(y)) 
          lines(y[, i], col = col[3], lty = lty[3], lwd = lwd[3])
        if (!is.null(z)) 
          lines(z[, i], col = col[3], lty = lty[3], lwd = lwd[3])
        box()
      }
      for (j in (nvr - nc + 1):nvr) {
        ylimL=dp$rangeList[j,]
        #print(ylimL) #;print("ylimL j")
        ifelse(is.null(ylab), ylabel <- colnames(x)[j], 
               ylabel <- ylab)
        xy <- xy.coords(x[, j])
        plot(xy, axes = FALSE, type = "l", ylab = ylabel, 
             ylim = ylimL, col = col[1], lty = lty[1], lwd = lwd[1], 
             ...)
        axis(2, at = pretty(rangeList[j,]))
        axis(1, at = pretty(1:(nrow(x))) ) #labels = c(0:(nrow(x) - 1))
        box()
        abline(h = 0, col = "red")
        if (!is.null(y)) 
          lines(y[, j], col = col[3], lty = lty[3], lwd = lwd[3])
        if (!is.null(z)) 
          lines(z[, j], col = col[3], lty = lty[3], lwd = lwd[3])
      }
      mtext(main, 3, line = 2, outer = TRUE, adj = adj.mtext, 
            padj = padj.mtext, col = col.mtext, ...)
      mtext(sub, 1, line = 4, outer = TRUE, adj = adj.mtext, 
            padj = padj.mtext, col = col.mtext, ...)
    }
    else {
      for (j in 1:nvr) {
        ylimL=dp$rangeList[j,] #ylimL
        ifelse(is.null(ylab), ylabel <- colnames(x)[j], 
               ylabel <- ylab)
        xy <- xy.coords(x[, j])
        plot(xy, type = "l", ylab = ylabel, ylim = ylimL, ##ylimL
             col = col[1], lty = lty[1], lwd = lwd[1], ...)
        if (!is.null(y)) 
          lines(y[, j], col = col[3], lty = lty[3], lwd = lwd[3])
        if (!is.null(z)) 
          lines(z[, j], col = col[3], lty = lty[3], lwd = lwd[3])
        abline(h = 0, col = "red")
      }
      mtext(main, 3, line = 2, outer = TRUE, adj = adj.mtext, 
            padj = padj.mtext, col = col.mtext, ...)
      mtext(sub, 1, line = 4, outer = TRUE, adj = adj.mtext, 
            padj = padj.mtext, col = col.mtext, ...)
    }
  }
  if (plot.type == "single") {
    for (i in 1:nvi) {
      dp <- dataplot(x, iname = inames[i])
      for (j in 1:nvr) {
        plot.single(dp, iname = inames[i], rname = rnames[j], 
                    ...)
        if (nvr > 1) 
          par(ask = FALSE) #edit
      }
    }
  }
  if (plot.type == "multiple") {
    for (i in 1:nvi) {
      dp <- dataplot(x, iname = inames[i])
      plot.multiple(dp, nc = nc, ...)
      if (nvi > 1) 
        par(ask = FALSE) #edit
    }
  }
}

getRespRangeList=function(irfResults,iname){
  impulse=irfResults$irf[[iname]] #employment as impulses
  upper <- irfResults$Upper[[iname]]
  lower=irfResults$Lower[[iname]]
  rangeList=matrix(ncol = 2)
  for (r in colnames(impulse)){
    row=range(c(lower[,c(r)], impulse[,c(r)], upper[,c(r)]))
    rangeList=rbind(rangeList, row)
  }
  rangeList=rangeList[-1,]
  return(rangeList)
}

最新更新