How to compute confidence intervals for this impulse response?

Question

I am trying to create bootstrapped confidence intervals for my impulse responses. Basically, I need a confidence intervall for each impulse response observation of each variable. It's easier for me to show you with an example where I am stuck.

library(vars)
library(varexternalinstrument)
data(GKdata)

gkvar <- VAR(GKdata[, c("logip", "logcpi", "gs1", "ebp")], p = 12, type = "const")
shockcol <- externalinstrument(gkvar, GKdata$ff4_tc, "gs1")
shockcol

ma_representation <- Phi(gkvar, 50)
irfs <- apply(ma_representation, 3, function(x) x %*% shockcol) #unclear
irfs <- as.data.frame(t(irfs))
colnames(irfs) <- names(shockcol)
irfs

So far so good, we get the impulse responses for the 4 variables in the dataset. Now, for each observation in each variable I want to compute bootstrapped confidence intervals. Yet, I got stuck probably on the statistics to insert in the function. I indicated with a question mark below the missing code that I can't write

? <- function() # I need to create a function that says bootstrap for each observation of each variable, I don't know how to do it

boot <- boot(irfs, ?, R = 1000, stype = "w", sim = "ordinary")

boot_ci <- boot.ci(boot, conf = c(0.90, 0.95), type = c("norm", "basic", "perc", "bca"))

Can anyone help me with that?

Thanks a lot!

Answer 1

So the data is collected in a period from 1979:7 to 2012:6, hence every 12 rows is a year and the bootstrap should sample the year in a block.

We first add this information to the data:

Data = GKdata
Data$year = (1:nrow(Data) - 1) %/% 12

And we create a function, that takes in a list, with sampled indices, recombines them into a data.frame and performs the same fitting function. For boot, you cannot return a data.frame so we return a matrix:

func = function(mylist,ind){

  mydf = do.call(rbind,mylist[ind])
  gkvar <- VAR(mydf[, c("logip", "logcpi", "gs1", "ebp")], p = 12, type = "const")
  shockcol <- externalinstrument(gkvar, mydf$ff4_tc, "gs1")

  ma_representation <- Phi(gkvar, 50)
  irfs <- t(apply(ma_representation, 3, function(x) x %*% shockcol)) #unclear
  colnames(irfs) <- names(shockcol)
  irfs
}

So we do the bootstrap by feeding it a list, that is split by the Year. The bootstrapping will sample the indices of this list, combine the sampled years into a new data frame and perform the fitting:

res = boot(split(Data,Data$year),func,100)

We can check whether the function is doing the right thing on the original observation:

table(res$t0 == irfs)

Your irfs has 51 x 4 values. If we look at the estimates inside boot results, this is flattened out:

dim(res$t)
[1] 100 204

To get back the 95% confidence interval for each value or variable, we do:

boot_est= sapply(1:ncol(res$t),function(i){
  boot.ci(res,index=i,type = "perc")$percent[4:5]
})

And then create the upper and lower bound by putting them back into a matrix:

lb = matrix(boot_est[1,],ncol=ncol(res$t0))
ub = matrix(boot_est[2,],ncol=ncol(res$t0))

You can plot lets say the first column of irfs, somehow the ci is really huge for the first few entries:

plot(irfs[,1],ylim=c(-8,8))
lines(ub[,1])
lines(lb[,1])