fev() returns tidy df

R/fev.R

@@ -1,19 +1,25 @@
 #' Calculate the forecast error variance for a VAR out to a specified horizon
 #'
 #' @param var vars::VAR or svars object
-#' @param h number of steps out to calculate fev
+#' @param n_ahead number of steps out to calculate fev
 #'
 #' @return forecast error variance
 #'
-fev <- function(var, h) {
+fev <- function(var, n_ahead) {
   n <- colnames(var$y)
   k <- length(n)
 
+  impulse_names <- n
+  if (inherits(var, "fevdvar")) {
+    impulse_names <- c("Main", paste0("Orth_", 2:k))
+  }
+  response_names <- n
+
   ## Calculate IRFs out to horizon (then adj to 3-dim matrix from DF)
   if (inherits(var, "fevdvar")) {
-    irf <- vars::irf(var, n.ahead = h, as_matrix = TRUE)
+    irf <- vars::irf(var, n.ahead = n_ahead, as_matrix = TRUE)
   } else {
-    irf <- vars::irf(var, n.ahead = h)[[1]][, -1] |>
+    irf <- vars::irf(var, n.ahead = n_ahead)[[1]][, -1] |>
       apply(1, matrix, simplify = FALSE, nrow = k, ncol = k, byrow = TRUE) |>
       simplify2array()
   }
@@ -29,10 +35,19 @@ fev <- function(var, h) {
   fe2 <- fe
 
   for (i in seq_along(fe)) {
-    for (j in 1:h) {
+    for (j in 1:n_ahead) {
       fe2[[i]][j, ] <- (colSums(fe[[i]][j:1, ]^2))
     }
   }
 
-  return(fe2)
+  ## Tidy to a DF
+  df <- data.frame(
+    h = rep(1:n_ahead, times = k * k),
+    impulse = rep(impulse_names, each = n_ahead, times = k),
+    response = rep(response_names, each = k * n_ahead),
+    fev = unlist(lapply(fe2, c)),
+    row.names = NULL
+  )
+
+  return(df)
 }

tests/testthat/test-test_id_fevdtd.R

@@ -16,14 +16,22 @@ test_that("Check that fevdtd > id.chol shocks", {
 
     mvar <- id_fevdtd(v, ti, hi)
 
-    fevdm <- fev(mvar, h = max(hi))
-    fevds <- fev(svar, h = max(hi))
-
-    m <- sum(fevdm[[ti]][, 1][hi])
-    m2 <- colSums(fevds[[ti]][hi, ])
+    fevm <- fev(mvar, max(hi))
+    fevs <- fev(svar, max(hi))
+
+    m <-
+      fevm[
+        fevm$impulse == "Main" &
+        fevm$response == ti &
+        fevm$h %in% hi,
+        "fev"] |>
+      sum()
+    m2 <-
+      fevm[fevm$impulse == ti & fevm$response == ti & fevm$h %in% hi, "fev"] |>
+      sum()
 
     expect_true(
-      all(m >= m2),
+      m >= m2,
       label = paste0("Iter: ", i, "\n", m, "\n", paste(m2, collapse = ", "))
       )
   }