tvarGIRF is an R package that calculates generalised impulse response functions to reduced form shocks for threshold vector autoregressions estimated using the tsDyn package.
Please note: this package is something I put together quickly for an uncompleted project because the only other implementation of GIRFs in R that existed at the time looks like it implements the bootstrap incorrectly to me. According to the tsDyn wiki, GIRFs are now included in the tsDyn package, so you might have better luck with those.
Installation
Install the lastest released version of the package using the R remotes package:
You can also install the latest development version:
Basic usage
The library extends the R tsDyn package. The following example illustrates how to create a simple GIRF for a threshold VAR using the zeroyld dataset provided with the tsDyn package.
GIRF is given a reduced form shock - in the example below a shock to only the second variable c(0,1). It is up to you to identify shocks; tvarGIRF won’t do it for you.
## Registered S3 method overwritten by 'xts':
## method from
## as.zoo.xts zoo## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo## Registered S3 methods overwritten by 'forecast':
## method from
## fitted.fracdiff fracdiff
## residuals.fracdiff fracdifflibrary(tvarGIRF)
# Estimate an example TVAR using the zeroyld dataset included in the tsDyn package
data(zeroyld)
exampleTVAR <- TVAR(zeroyld, lag=2, nthresh=1, thDelay=1, mTh=1, plot=FALSE)## Best unique threshold 10.653# Calculate GIRFs for a reduced form shock to the second variable (long.run)
GIRF(exampleTVAR, c(0,1), H = 10, R = 10)## # GIRF of tvar exampleTVAR
## # A tibble: 20 x 2
## short.run long.run
## <dbl> <dbl>
## 1 -0.00390 1.00
## 2 0.0476 1.04
## 3 0.0744 0.988
## 4 0.0644 0.941
## 5 0.0789 0.888
## 6 0.125 0.913
## 7 0.104 0.803
## 8 0.0954 0.755
## 9 0.188 0.827
## 10 0.231 0.839
## 11 0.265 0.888
## 12 0.275 0.862
## 13 0.326 0.840
## 14 0.288 0.694
## 15 0.325 0.728
## 16 0.300 0.595
## 17 0.274 0.449
## 18 0.198 0.353
## 19 0.117 0.254
## 20 0.174 0.373NB: These GIRFs are meaningless as I have used very few repetitions to calculate them. I have done this so that the examples don’t take too long to run. But you should use much larger values for your actual estimation.
Shocks to more than one variable
Note that you can shock more than one variable. If you had a three equation TVAR, passing shock = c(1,0,0) would give you the GIRF to a reduced-form shock to the first equation, and shock = c(1,1,0) would be a reduced form shock to both the first and second. The values need not be 1 either. You can estimate the response to a 0.9 unit shock to one variable and a 0.645 to the other (just as an example!):
## # GIRF of tvar exampleTVAR
## # A tibble: 20 x 2
## short.run long.run
## <dbl> <dbl>
## 1 0.874 0.623
## 2 0.969 1.00
## 3 0.901 0.855
## 4 0.923 0.820
## 5 0.907 0.807
## 6 0.873 0.650
## 7 0.853 0.629
## 8 0.841 0.619
## 9 0.701 0.475
## 10 0.552 0.202
## 11 0.522 0.163
## 12 0.610 0.292
## 13 0.679 0.518
## 14 0.623 0.450
## 15 0.632 0.385
## 16 0.641 0.464
## 17 0.667 0.493
## 18 0.736 0.546
## 19 0.706 0.610
## 20 0.702 0.615This is helpful, as you will typically want to identify structural shocks. The reduced form shock will be some linear combination of those structural shocks. See the section “A note on shocks” for a bit more discussion.
Viewing the results
tvarGIRF has a couple helpful functions to help you view the results from GIRF.
The first is to simply print the results of an estimation:
saved_girfs <- GIRF(exampleTVAR, c(0,1), H = 10, R = 10)
saved_girfs # or, equivalently print(saved_girfs)## # GIRF of tvar exampleTVAR
## # A tibble: 20 x 2
## short.run long.run
## <dbl> <dbl>
## 1 0.00591 0.992
## 2 0.0961 1.04
## 3 0.203 1.11
## 4 0.231 1.14
## 5 0.233 1.11
## 6 0.197 0.988
## 7 0.201 0.894
## 8 0.182 0.780
## 9 0.141 0.604
## 10 0.174 0.654
## 11 0.118 0.557
## 12 0.0777 0.512
## 13 0.0681 0.346
## 14 0.0558 0.347
## 15 0.0885 0.377
## 16 0.0747 0.316
## 17 0.0976 0.315
## 18 0.126 0.279
## 19 0.143 0.313
## 20 0.0692 0.237There is also a tidier for GIRFs (a la broom):
## # A tibble: 40 x 3
## horizon variable response
## <int> <chr> <dbl>
## 1 1 short.run 0.00591
## 2 2 short.run 0.0961
## 3 3 short.run 0.203
## 4 4 short.run 0.231
## 5 5 short.run 0.233
## 6 6 short.run 0.197
## 7 7 short.run 0.201
## 8 8 short.run 0.182
## 9 9 short.run 0.141
## 10 10 short.run 0.174
## # … with 30 more rowsYou can View GIRFs in the spreadsheet viewer. You can also plot your GIRFs, using a built in plot function that uses ggplot2:
(Again, these GIRFs are meaningless as I have used very few repetitions to calculate them, which is why they look so wonky. I have done this so that the examples don’t take too long to run. But you should use much larger values for your actual estimation.)
A note on shocks
tvarGIRF does not identify shocks for you. The argument you supply to tvarGIRF is for a reduced form shock. This is rarely what you want. The reduced form shocks are extremely unlikely to be uncorrelated with one another, and so it makes little sense to consider one such isolated shock.
Instead, you will want to identify structural shocks that are uncorrelated with one another.
If you want to use structural shocks, you will have to identify the shocks yourself and supply the reduced form shock that corresponds to your chosen structural shock. It is up to you to identify shocks; tvarGIRF won’t do it for you. There are many strategies for identifying structural shocks, tvarGIRf can’t help with this.
As an example, let’s identify a structural shock for the TVAR we estimated above using timing restrictions / Cholesky decomposition. Our timing assumption is that shocks to short.run affect both variables contemparenously, but shocks to long.run affect short.run with a one period lag.1
Using our estimated TVAR, we can get the residuals by:
We calculate the covariance matrix by:
## short.run long.run
## short.run 42.70877 54.85373
## long.run 54.85373 124.82032And take an upper triangular Cholesky decomposition:
## short.run long.run
## short.run 6.535195 8.393587
## long.run 0.000000 7.373466I could then calculate the GIRF to the identified structural shock to the variable short.run by running:
Which is what my identification scheme tells me the structural shock looks like as a reduced form shock (the first row of the cholesky-decomposed covariance matrix).
GIRFs conditional on starting in a particular regime
The restrict.to argument takes an integer corresponding to which regime you wish to restrict the GIRF to. For instance, using the toy example we started with before, if I estimate the following TVAR:
data(zeroyld)
exampleTVAR <- TVAR(zeroyld, lag = 2, nthresh = 1, thDelay = 1, mTh = 1, plot = FALSE)## Best unique threshold 10.653I have a TVAR with two regimes. These are regimes are denoted by 1 when then threshold variable (in this case short.run) is below my estimated threshold (10.653 in this example) and denoted by 2 when it is above the threshold value. This scales to more regimes: 1 corresponds to the lowest regime and nthresh+1 to the highest regime.
You can see the names of your regimes by printing your TVAR:
## Model TVAR with 1 thresholds
##
## $Bdown
## Intercept short.run -1 long.run -1 short.run -2
## Equation short.run 0.048644738 0.9211813 0.02872384 0.05347747
## Equation long.run 0.007964404 0.1835720 0.97111008 -0.14157883
## long.run -2
## Equation short.run -0.006114365
## Equation long.run -0.015307980
##
## $Bup
## Intercept short.run -1 long.run -1 short.run -2
## Equation short.run 1.561918 1.1780904 -0.03076567 -0.3468208
## Equation long.run 2.037192 0.7989221 0.81823826 -0.8562527
## long.run -2
## Equation short.run 0.06986963
## Equation long.run 0.05771228
##
##
## Threshold value[1] "10.653"My regimes are called Bdown (the regime denoted by 1, as it appears first) and Bup (the regime denoted by 2, as it appears second).
You can also see which regime your estimated model is in at each point in your estimation sample by looking at:
## [1] NA NA 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [24] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [47] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [70] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [93] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [116] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [139] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [162] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [185] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [208] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [231] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [254] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [277] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [300] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [323] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [346] 1 1 1 1 2 2 2 2 1 1 2 2 2 2 2 2 2 2 2 2 2 2 2
## [369] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 1 2 1 1 1 2 2
## [392] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 1
## [415] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [438] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [461] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1And you could confirm that these correspond to short.run being above or below the estimated threshold value.
The restrict.to argument in tvarGIRF uses those regime numbers. Suppose I want to examine GIRFs conditional on starting in regime Bdown, when short.run is below my threshold. This corresponds to regime 1. So I would use restrict.to = 1 in my call to GIRF:
## # GIRF of tvar exampleTVAR
## # A tibble: 20 x 2
## short.run long.run
## <dbl> <dbl>
## 1 -0.0273 0.959
## 2 0.0107 0.934
## 3 0.0397 0.900
## 4 0.106 0.933
## 5 0.171 0.988
## 6 0.155 0.901
## 7 0.161 0.879
## 8 0.193 0.862
## 9 0.272 0.912
## 10 0.300 0.936
## 11 0.291 0.892
## 12 0.351 0.915
## 13 0.356 0.836
## 14 0.312 0.746
## 15 0.286 0.728
## 16 0.361 0.791
## 17 0.377 0.796
## 18 0.424 0.841
## 19 0.415 0.728
## 20 0.367 0.650(NB: Again, I am using far too few repetitions to calculate meaningful GIRFs, but I don’t want the examples to take too long.)
This is not a very sensible assumption in this context and is meant purely to illustrate how you might go about this identification strategy using
tvarGIRF. There are many problems with the identification strategy I use here, including that I assume that the impact of the structural shocks is the same in both regimes.↩