R/gfevd.r

In joergrieger/bayesianConnectedness: Calculates spillover measures as in Diebold-Yilmaz (2009)

# Function to calculate generalized forecast error variance decomposition
# Inputs:
# bvObj - An object of the class bvar
# p     - order of the moving average
#

gfevd2 <- function(bvObj,h=10){
  # test if bvObj is of class bvar

  #
  constant <- 0
  if(bvObj$intercept){constant <- 1}
  testdraw <- bvObj$betadraws[,,1]
  if(constant){testdraw <- testdraw[-c(1),]}
  nVar <- dim(testdraw)[2]


  nLength <- dim(bvObj$betadraws)[2]
  delta <- array(0,dim=c(nVar,nVar,nLength))

  for(ii in 1:nLength){ # loop over all sample draws

    # transform coefficients into a threedimensional array

    bdraw <- bvObj$betadraws[,,ii]
    if(constant){
      bdraw <- bdraw[-c(1),]
    }
    varcoef <- array(0,dim=c(nVar,nVar,bvObj$NoLags))
    for( jj in 1:(bvObj$NoLags) ){

      varcoef[,,jj] <- bdraw[((jj-1)*nVar+1):(jj*nVar),]

    }

    # calculate moving averages
    vmacoef <- phi(varcoef,h)

    # Calculate the generalized fevd
    deltatemp <- array(0,dim=c(nVar,nVar))
    for(i in 1:nVar){
      ei <- array(0,dim=c(nVar))
      ei[i] <- 1
      tsum2 <- 0
      for(hh in 1:h){
        tmp <- t(ei) %*% vmacoef[,,hh] %*% bvObj$sigmadraws[,,ii] %*% t(vmacoef[,,hh]) %*% ei
        tsum2 <- tsum2 + tmp
      }
      for(j in 1:nVar){
        ej <- array(0, dim=c(nVar))
        ej[j] <- 1
        tsum <- 0
        for(hh in 1:h){
          tmp <- (ei %*% vmacoef[,,hh] %*% ej)^2
          tsum <- tsum + tmp
        }
        deltatemp[i,j] <- bvObj$sigmadraws[i,i,ii]^(-1)*tsum/tsum2
      }
    }

    delta[,,ii] <- deltatemp

  }
  deltafinal <- array(0,dim=c(nVar,nVar))
  for(ii in 1:nVar){
    for(jj in 1:nVar){
      deltafinal[ii,jj] <- mean(delta[ii,jj,])
    }
  }
  return(deltafinal)

}