Conventional economic theory says that long term bonds are more sensitive to interest rate changes than short term bonds. Intuitively, we should see higher realized volatility in funds that track these types of bonds. Here, I want to test that if as you move from short duration bonds to higher and higher duration bonds, then the difference of mean volatility grows larger and larger. To do this, I selected Blackrock ETF’s of various lengths listed below:
SHY, 1-3 years. IEI, 3-7 years. IEF, 7-10 years. TLH, 10-20 years. TLT, 20+ years.
Comparing difference of means is fairly simple to do from a frequentist standpoint; the t-test is a straightforward procedure that does this. However, to better estimate difference of means I want to use a Bayesian approach and follow the BEST procedure, which was developed by Dr. John Krushke at Indiana University. A quick overview is that rather testing and rejecting null hypotheses, the Bayesian approach estimates distributions for group means and their differences. This is more robust and provides a richer picture of the underlying data.To learn more, I suggest reading the PyMC article as well as the original Krushke paper. In this article, I will use the Krushke package that he wrote himself.
To estimate realized volatility, there are a few ways to do it but I followed the “Realized Volatility Puzzle” method from Harel Jacobsen (article linked below) and decided on the Garman-Klass Volatility estimator over a 10 day window.
Here are the results:
SHY v IEI
library(tidyr)
library(TTR) ## Warning: package 'TTR' was built under R version 4.1.3library(BEST)## Warning: package 'BEST' was built under R version 4.1.3## Loading required package: HDInterval## Warning: package 'HDInterval' was built under R version 4.1.3shy1 = riingo_prices("shy", start_date = "01-01-2010", end_date = "07-04-2022")
iei = riingo_prices("iei", start_date = "01-01-2010", end_date = "07-04-2022")
# measuring 10 day realized volatility
# using garman klaus volatility
# https://medium.com/swlh/the-realized-volatility-puzzle-588a74ab3896
n = 10
shy1_vol = volatility(data.frame(shy1[6],shy1[,3:5]), n = n, N = 252, mean0 = FALSE)
shy1_vol = drop_na(as.data.frame(shy1_vol))
iei_vol = volatility(data.frame(iei[6],iei[,3:5]), n = n, N = 252, mean0 = FALSE)
iei_vol = drop_na(as.data.frame(iei_vol))
shy_v_iei = BESTmcmc(shy1_vol$shy1_vol*100, iei_vol$iei_vol*100,
parallel = TRUE, numSavedSteps = 1000)## Waiting for parallel processing to complete...## done.summary(shy_v_iei)## mean median mode HDI% HDIlo HDIup compVal %>compVal
## mu1 0.717 0.717 0.718 95 0.700 0.732
## mu2 2.627 2.627 2.628 95 2.589 2.669
## muDiff -1.910 -1.910 -1.917 95 -1.947 -1.866 0 0
## sigma1 0.347 0.347 0.346 95 0.334 0.361
## sigma2 0.942 0.943 0.944 95 0.909 0.977
## sigmaDiff -0.595 -0.596 -0.598 95 -0.629 -0.563 0 0
## nu 3.410 3.408 3.423 95 3.140 3.723
## log10nu 0.532 0.533 0.535 95 0.497 0.571
## effSz -2.690 -2.691 -2.699 95 -2.787 -2.596 0 0SHY v IEF
library(BEST)
library(TTR)
library(tidyr)
shy2 = riingo_prices("shy", start_date = "01-01-2010", end_date = "07-04-2022")
ief = riingo_prices("ief", start_date = "01-01-2010", end_date = "07-04-2022")
# measuring 10 day realized volatility
# using garman klaus volatility
# https://medium.com/swlh/the-realized-volatility-puzzle-588a74ab3896
n = 10
shy2_vol = volatility(data.frame(shy2[6],shy2[,3:5]), n = n, N = 252, mean0 = FALSE)
shy2_vol = drop_na(as.data.frame(shy2_vol))
ief_vol = volatility(data.frame(ief[6],ief[,3:5]), n = n, N = 252, mean0 = FALSE)
ief_vol = drop_na(as.data.frame(ief_vol))
shy_v_ief = BESTmcmc(shy2_vol$shy2_vol*100, ief_vol$ief_vol*100,
parallel = TRUE, numSavedSteps = 1000)## Waiting for parallel processing to complete...done.summary(shy_v_ief)## mean median mode HDI% HDIlo HDIup compVal %>compVal
## mu1 0.717 0.717 0.718 95 0.701 0.731
## mu2 4.786 4.786 4.785 95 4.720 4.849
## muDiff -4.069 -4.068 -4.066 95 -4.131 -4.000 0 0
## sigma1 0.348 0.348 0.348 95 0.335 0.361
## sigma2 1.518 1.518 1.519 95 1.459 1.570
## sigmaDiff -1.169 -1.170 -1.167 95 -1.221 -1.110 0 0
## nu 3.434 3.434 3.447 95 3.155 3.742
## log10nu 0.535 0.536 0.538 95 0.501 0.575
## effSz -3.697 -3.695 -3.689 95 -3.846 -3.565 0 0SHY v TLH
library(BEST)
library(TTR)
library(tidyr)
shy3 = riingo_prices("shy", start_date = "01-01-2010", end_date = "07-04-2022")
tlh = riingo_prices("tlh", start_date = "01-01-2010", end_date = "07-04-2022")
# measuring 10 day realized volatility
# using garman klaus volatility
# https://medium.com/swlh/the-realized-volatility-puzzle-588a74ab3896
n = 10
shy3_vol = volatility(data.frame(shy3[6],shy3[,3:5]), n = n, N = 252, mean0 = FALSE)
shy3_vol = drop_na(as.data.frame(shy3_vol))
tlh_vol = volatility(data.frame(tlh[6],tlh[,3:5]), n = n, N = 252, mean0 = FALSE)
tlh_vol = drop_na(as.data.frame(tlh_vol))
shy_v_tlh = BESTmcmc(shy3_vol$shy3_vol*100, tlh_vol$tlh_vol*100,
parallel = TRUE, numSavedSteps = 1000)## Waiting for parallel processing to complete...done.summary(shy_v_tlh)## mean median mode HDI% HDIlo HDIup compVal %>compVal
## mu1 0.713 0.713 0.712 95 0.698 0.730
## mu2 7.127 7.128 7.130 95 7.031 7.237
## muDiff -6.414 -6.414 -6.413 95 -6.521 -6.316 0 0
## sigma1 0.341 0.341 0.341 95 0.328 0.354
## sigma2 2.315 2.316 2.312 95 2.229 2.396
## sigmaDiff -1.974 -1.976 -1.986 95 -2.050 -1.886 0 0
## nu 3.109 3.100 3.077 95 2.889 3.376
## log10nu 0.492 0.491 0.488 95 0.461 0.528
## effSz -3.878 -3.878 -3.887 95 -4.022 -3.728 0 0SHY v TLT
library(BEST)
library(TTR)
library(tidyr)
shy4 = riingo_prices("shy", start_date = "01-01-2010", end_date = "07-04-2022")
tlt = riingo_prices("tlt", start_date = "01-01-2010", end_date = "07-04-2022")
# measuring 10 day realized volatility
# using garman klaus volatility
# https://medium.com/swlh/the-realized-volatility-puzzle-588a74ab3896
n = 10
shy4_vol = volatility(data.frame(shy4[6],shy4[,3:5]), n = n, N = 252, mean0 = FALSE)
shy4_vol = drop_na(as.data.frame(shy4_vol))
tlt_vol = volatility(data.frame(tlt[6],tlt[,3:5]), n = n, N = 252, mean0 = FALSE)
tlt_vol = drop_na(as.data.frame(tlt_vol))
shy_v_tlt = BESTmcmc(shy4_vol$shy4_vol*100, tlt_vol$tlt_vol*100,
parallel = TRUE, numSavedSteps = 1000)## Waiting for parallel processing to complete...done.summary(shy_v_tlt)## mean median mode HDI% HDIlo HDIup compVal %>compVal
## mu1 0.716 0.716 0.715 95 0.701 0.731
## mu2 11.120 11.119 11.109 95 10.983 11.274
## muDiff -10.405 -10.404 -10.406 95 -10.550 -10.260 0 0
## sigma1 0.345 0.345 0.342 95 0.332 0.357
## sigma2 3.394 3.393 3.397 95 3.276 3.519
## sigmaDiff -3.049 -3.048 -3.046 95 -3.171 -2.931 0 0
## nu 3.319 3.308 3.296 95 3.090 3.605
## log10nu 0.521 0.520 0.519 95 0.490 0.557
## effSz -4.315 -4.312 -4.296 95 -4.489 -4.168 0 0As you can see, it is clear that each ETF has different amounts of volatility, and that the longer duration ETF has the higher mean volatility.This confirms the economic tuition that long term bonds are more sensitive to interest rates, and thus leads to higher volatility.
Hopefully this article gave insight as to what the BEST process is and why it is more informative than a standard t-test. While the t-test would have worked for a simple example like the one done here, for more nuanced questions I would default to using the BEST process, with more informative priors.
Sources:
https://www.pymc.io/projects/examples/en/latest/case_studies/BEST.html https://jkkweb.sitehost.iu.edu/articles/Kruschke2013JEPG.pdf https://medium.com/swlh/the-realized-volatility-puzzle-588a74ab3896
Data:
Tiingo