Analysis of Financial Time Series in Stan - Chapter 3

Chapter_3_1.utf8

In this series we are trying to reproduce the models and examples listed in the book “Analysis of Financial Time Series, 3rd Edition”, by Ruey S. Tsay, using Stan https://mc-stan.org/ and the package RStan https://cran.r-project.org/web/packages/rstan/index.html. The main repository for the presented models and data can be found at https://github.com/marcomarconi/AFTS_with_Stan.

Most of the models listed in chapter 3 were already translated into Stan as part of the official documentation and can be found at https://mc-stan.org/docs/2_26/stan-users-guide/time-series-chapter.html. Still, we will go through them in the current document.

{
  require(tidyverse)
  require(rstan)
  require(bayesplot)
  require(loo)
  require(LaplacesDemon)
  library(xts)
  library(quantmod)
  library(FKF)
  library(MARSS)
  library(forecast)
  rstan_options(auto_write = TRUE)
  rstan_options(javascript = FALSE)
  theme_set(theme_classic(base_size = 24))
}

The ARCH Model

The typical ARCH(K) model can be coded in Stan as follow:

data {
  int<lower=0> N;
  vector[N] y;
  int Kar;
  int family;
}

parameters {
  real mu;                       
  real<lower=0> nu[family==1];
  real<lower=0> alpha0;          
  real<lower=0,upper=1> alpha[Kar];  
}
transformed parameters {
  vector[N] sigma;
  for(i in 1:Kar)
    sigma[i] = alpha0;
  for (i in (Kar+1):N) {
    sigma[i] = alpha0;
    for(k in 1:Kar)
      sigma[i] += alpha[k] * pow(y[i-k] - mu, 2);
  }
  sigma = sqrt(sigma);
}
model {
  mu ~ normal(0, 1);
  nu ~ normal(0, 1);
  alpha0 ~ normal(1, 1);
  alpha ~ normal(0.5, 0.5);
  if(family == 0)
    y ~ normal(mu, sigma);
  else if(family == 1)  
    y ~ student_t(nu[1], mu, sigma);
}

generated quantities {
   vector[N] log_lik;
   {
     for(i in 1:N) 
      if(family == 0)
        log_lik[i] = normal_lpdf(y[i] | mu, sigma[i] );  
      else if(family == 1)    
        log_lik[i] = student_t_lpdf(y[i] | nu[1], mu, sigma[i] );  
   }
}

We will show the ARCH model on the monthly log returns of Intel stock:

intel <- read.table("data/m-intc7308.txt", header=T)
intel <- xts(intel$rtn, as.Date(as.character(intel$date), format = "%Y%m%d")) 
plot(intel)

We will fit both an ARCH(3) and an ARCH(1) model to the data:

ARCH <- stan_model("models/Chapter_3/ARCH.stan")
fit_intel_3 <- sampling(ARCH, data = list(N = length(intel), y = intel %>% as.vector(), Kar = 3, family = 0), chains = 4, refresh=0)
fit_intel_1 <- sampling(ARCH, data = list(N = length(intel), y = intel %>% as.vector(), Kar = 1, family = 0), chains = 4, refresh=0)
print(fit_intel_3, pars=names(fit_intel_3)[1:5])

## Inference for Stan model: ARCH.
## 4 chains, each with iter=2000; warmup=1000; thin=1; 
## post-warmup draws per chain=1000, total post-warmup draws=4000.
## 
##          mean se_mean   sd 2.5%  25%  50%  75% 97.5% n_eff Rhat
## mu       0.02       0 0.01 0.01 0.01 0.02 0.02  0.03  4724    1
## alpha0   0.01       0 0.00 0.01 0.01 0.01 0.01  0.01  2538    1
## alpha[1] 0.30       0 0.12 0.11 0.22 0.29 0.37  0.56  2407    1
## alpha[2] 0.10       0 0.06 0.01 0.05 0.09 0.13  0.25  3222    1
## alpha[3] 0.08       0 0.05 0.01 0.04 0.07 0.10  0.18  3870    1
## 
## Samples were drawn using NUTS(diag_e) at Sat May 15 20:53:27 2021.
## For each parameter, n_eff is a crude measure of effective sample size,
## and Rhat is the potential scale reduction factor on split chains (at 
## convergence, Rhat=1).

print(fit_intel_1, pars=names(fit_intel_1)[1:3])

## Inference for Stan model: ARCH.
## 4 chains, each with iter=2000; warmup=1000; thin=1; 
## post-warmup draws per chain=1000, total post-warmup draws=4000.
## 
##          mean se_mean   sd 2.5%  25%  50%  75% 97.5% n_eff Rhat
## mu       0.02       0 0.01 0.01 0.02 0.02 0.02  0.03  2865    1
## alpha0   0.01       0 0.00 0.01 0.01 0.01 0.01  0.01  2340    1
## alpha[1] 0.41       0 0.11 0.22 0.33 0.40 0.48  0.64  2113    1
## 
## Samples were drawn using NUTS(diag_e) at Sat May 15 20:53:34 2021.
## For each parameter, n_eff is a crude measure of effective sample size,
## and Rhat is the potential scale reduction factor on split chains (at 
## convergence, Rhat=1).

It can be seen that the parameter estimation are in reasonable accordance with the results obtained in book with the maximum likelihhod approach.

We further perform some model checking by testing the independence of the standardized residuals:

sq <- extract(fit_intel_1, pars="sigma")[[1]] %>% colMeans()
mu <- extract(fit_intel_1, pars="mu")[[1]] %>% mean()
residuals <- ((intel - mu) / sq)
plot(residuals)

Box.test(residuals, lag=10, type="Ljung-Box")

## 
##  Box-Ljung test
## 
## data:  residuals
## X-squared = 11.822, df = 10, p-value = 0.2971

A quick comparison with the loo package show little difference between the ARCH(3) and ARCH(1), as already discussed in the book:

loo_compare(loo(fit_intel_3), loo(fit_intel_1))

##        elpd_diff se_diff
## model1  0.0       0.0   
## model2 -1.1       2.7

a <- extract_log_lik(fit_intel_3, merge_chains = FALSE)[,,4:(length(intel))]
b <- extract_log_lik(fit_intel_1, merge_chains = FALSE)[,,4:(length(intel))]
r_eff_a <- relative_eff(exp(a))
r_eff_b <- relative_eff(exp(b))
loo_compare(loo(a, r_eff = r_eff_a, moment_match = TRUE), loo(b, r_eff = r_eff_b, moment_match = TRUE))

##        elpd_diff se_diff
## model1  0.0       0.0   
## model2 -0.8       2.6

For comparison, we also fit an ARCH(1) model with Student-t innovations to the series:

fit_intel_3_t <- sampling(ARCH, data = list(N = length(intel), y = intel %>% as.vector(), Kar = 3, family = 1), chains = 4, refresh=0)
print(fit_intel_3_t, pars=names(fit_intel_3_t)[1:5])

## Inference for Stan model: ARCH.
## 4 chains, each with iter=2000; warmup=1000; thin=1; 
## post-warmup draws per chain=1000, total post-warmup draws=4000.
## 
##          mean se_mean   sd 2.5%  25%  50%  75% 97.5% n_eff Rhat
## mu       0.02    0.00 0.01 0.01 0.02 0.02 0.02  0.03  6502    1
## nu[1]    3.78    0.01 0.49 2.92 3.42 3.75 4.10  4.81  3953    1
## alpha0   0.01    0.00 0.00 0.00 0.01 0.01 0.01  0.01  3910    1
## alpha[1] 0.13    0.00 0.07 0.02 0.08 0.12 0.17  0.28  3890    1
## alpha[2] 0.08    0.00 0.05 0.01 0.04 0.07 0.11  0.18  3938    1
## 
## Samples were drawn using NUTS(diag_e) at Sat May 15 20:53:58 2021.
## For each parameter, n_eff is a crude measure of effective sample size,
## and Rhat is the potential scale reduction factor on split chains (at 
## convergence, Rhat=1).

Let’s also check the independence assumption for standardized residuals:

sq <- extract(fit_intel_3_t, pars="sigma")[[1]] %>% colMeans()
mu <- extract(fit_intel_3_t, pars="mu")[[1]] %>% mean()
residuals <- ((intel - mu) / sq)
plot(residuals)

Box.test(residuals, lag=12, type="Ljung-Box")

## 
##  Box-Ljung test
## 
## data:  residuals
## X-squared = 10.648, df = 12, p-value = 0.5593

The Student-t ARCH(3) seems to perform slightly between in comparison with the normal ARCH(3) model, however, it must be kept in mind that with this kind of analysis we are allowing data from the future to influence the current estimations, and therefore the result from loo will likely over-estimate the difference between models:

loo_compare(loo(fit_intel_3), loo(fit_intel_3_t))

##        elpd_diff se_diff
## model2  0.0       0.0   
## model1 -8.3       9.1

a <- extract_log_lik(fit_intel_3, merge_chains = FALSE)[,,4:(length(intel))]
b <- extract_log_lik(fit_intel_3_t, merge_chains = FALSE)[,,4:(length(intel))]
r_eff_a <- relative_eff(exp(a))
r_eff_b <- relative_eff(exp(b))
loo_compare(loo(a, r_eff = r_eff_a, moment_match = TRUE), loo(b, r_eff = r_eff_b, moment_match = TRUE))

##        elpd_diff se_diff
## model2  0.0       0.0   
## model1 -8.4       9.1

Search This Blog

No Time for Bayes

Analysis of Financial Time Series in Stan - Chapter 3 - ARCH

The ARCH Model

Comments

Post a Comment

Popular posts from this blog

Dynamic Linear Models with Switching in Stan

Analysis of Financial Time Series in Stan - Chapter 10 - Multivariate Volatility Models and Their Applications

Analysis of Financial Time Series in Stan - Chapter 4 - Nonlinear Models