R/fanchart.R

In vars: VAR Modelling

"fanchart" <-
function(x, colors = NULL, cis = NULL, names = NULL, main = NULL, ylab = NULL, xlab = NULL, col.y = NULL, nc, plot.type = c("multiple", "single"), mar = par("mar"), oma = par("oma"), ...){
  if(!is(x, "varprd")){
      stop("\nPlease provide an object of class 'varprd',\ngenerated by predict-method for objects of class 'varest'.\n")
      }
  if(is.null(colors))
    colors <- gray(sqrt(seq(from = 0.05, to = 1.0, length = 9))) 
  if(is.null(cis)){
    cis <- seq(0.1, 0.9, by = 0.1)
  } else {
    if((min(cis) <= 0) || (max(cis) >= 1))
      stop("\nValues of confidence intervals must be in(0, 1).\n")
    if(length(cis) > length(colors))
      stop("\nSize of 'colors' vector must be at least as long as\nsize of 'cis' vector\n")
  }
  n.regions <- length(cis)
  n.ahead <- nrow(x$fcst[[1]])
  K <- ncol(x$endog)
  e.sample <- nrow(x$endog)
  endog <- x$endog
  fcst <- NULL 
  for(j in 1:n.regions){
    fcst[[j]] <- predict(x$model, n.ahead = n.ahead, ci = cis[j], dumvar = x$exo.fcst)$fcst
  }
  xx <- seq(e.sample, length.out = n.ahead + 1)
  xx <- c(xx, rev(xx))
  op <- par(no.readonly = TRUE)
  plot.type <-  match.arg(plot.type)
  ynames <- colnames(endog)
  if (is.null(names)) {
    names <- ynames
  } else {
    names <- as.character(names)
    if (!(all(names %in% ynames))) {
      warning("\nInvalid variable name(s) supplied, using first variable.\n")
      names <- ynames[1]
    }
  }

  nv <- length(names)

  ifelse(is.null(main), main <- paste("Fanchart for variable", names), main <- rep(main, nv)[1:nv])
  ifelse(is.null(ylab), ylab <- "", ylab <- ylab)
  ifelse(is.null(xlab), xlab <- "", xlab <- xlab)
  ifelse(is.null(col.y), col.y <- "black", col.y <- col.y)

    if(plot.type == "single") {
      if(nv > 1) par(ask = TRUE)
      par(mar = mar, oma = oma)
    } else if(plot.type == "multiple"){
      if (missing(nc)) {
        nc <- ifelse(nv > 4, 2, 1)
      }
      nr <- ceiling(nv/nc)
      par(mfcol = c(nr, nc), mar = mar, oma = oma)
    }

  
  for(i in 1 : nv){
    ymax <- max(c(fcst[[n.regions]][names[i]][[1]][, 3]), endog[, names[i]])
    ymin <- min(c(fcst[[n.regions]][names[i]][[1]][, 2]), endog[, names[i]])
    yy1 <- c(endog[e.sample, names[i]], fcst[[1]][names[i]][[1]][, 2], rev(c(endog[e.sample, names[i]], fcst[[1]][names[i]][[1]][, 3])))
    plot.ts(c(endog[, names[i]], rep(NA, n.ahead)), main = main[i], ylim = c(ymin, ymax), ylab = ylab, xlab = xlab, col = col.y, ...)
    polygon(xx, yy1, col = colors[1], border = colors[1])
    if(n.regions > 1){
      for(l in 2:n.regions){
        yyu <- c(endog[e.sample, names[i]], fcst[[l]][names[i]][[1]][, 3], rev(c(endog[e.sample, names[i]], fcst[[l-1]][names[i]][[1]][, 3])))
        yyl <- c(endog[e.sample, names[i]], fcst[[l-1]][names[i]][[1]][, 2], rev(c(endog[e.sample, names[i]], fcst[[l]][names[i]][[1]][, 2])))
        polygon(xx, yyu, col = colors[l], border = colors[l])
        polygon(xx, yyl, col = colors[l], border = colors[l])
      }
    }
  }
  on.exit(par(op))
}