R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

Time:2022-8-5

Original link:http://tecdat.cn/?p=20015 

This article will illustrate different models for univariate and multivariate financial time series, in particular models for conditional means and conditional covariance matrices, volatility.

mean model

This section explores conditional mean models.

iid model

We start with a simple iid model. The iid model assumes that the log-return xt is an N-dimensional Gaussian time series:

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

The sample estimators for the mean and covariance matrix are the sample mean, respectively

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

and the sample covariance matrix

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

We start by generating the data, become familiar with the process and ensure that the estimation process gives correct results (i.e. sanity checks). Then use real market data and fit different models.

Let's generate synthetic iid data and estimate mean and covariance matrices:

# Generate comprehensive income data X <- rmvnorm(n = T, mean = mu, sigma = Sigma)# sample estimates (sample mean and sample covariance matrix) mu\_sm <- colMeans(X)Sigma\_scm <- cov( X)# Error norm(mu\_sm - mu, "2")#> \[1\] 2.44norm(Sigma\_scm - Sigma, "F")#> \[1\] 70.79

Now, let's do it again for a different number of observations T:

# First generate all data X <- rmvnorm(n = T\_max, mean = mu, sigma = Sigma) # Now iterate over a subset of samples for (T\_ in T\_sweep) { # Sample estimation mu\_sm <- colMeans(X_) Sigma\_scm <- cov(X\_) # Calculate error error\_mu\_vs\_T <- c(error\_mu\_vs\_T, norm(mu\_sm - mu, "2")) error\_Sigma\_vs\_T <- c(error\_Sigma\_vs\_T, norm(Sigma\_scm - Sigma, "F"))# plot(T\_sweep, error\_mu\_vs\_T, main = "mu estimation error",

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

plot(T\_sweep, error\_Sigma\_vs\_T main = "Error in Sigma estimation", ylab = "Error"

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

Univariate ARMA Model

The ARMA(p,q) model on the logarithmic rate of return xt is

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

where wt is a white noise sequence with zero mean and variance σ2. The parameters of the model are the coefficients ϕi, θi and the noise variance σ2.

Note that the ARIMA(p,d,q) model is an ARMA(p,q) model with a time difference of order d. So if we replace the log price with xt, the previous log return model is actually an ARIMA(p, 1, q) model because once the log price is differentiated, we get the log return.

rugarchGenerate data

we will userugarchThe package generates univariate ARMA data, estimates parameters and makes predictions.

First, we need to define the model:

# Specify AR(1) model with given coefficients and parameters #> #> *------------------------------- ---*#> * ARFIMA Model Spec *#> *------------------------------------------------*# > Conditional Mean Dynamics#> ------------------------------------#> Mean Model : ARFIMA(1 ,0,0)#> Include Mean : TRUE #> #> Conditional Distribution#> ----------------------------------------- -----#> Distribution : norm #> Includes Skew : FALSE #> Includes Shape : FALSE #> Includes Lambda : FALSE#> Level Fixed Include Estimate LB UB#> mu 0.01 1 1 0 NA NA#> ar1 -0.90 1 1 0 NA NA#> ma 0.00 0 0 0 NA NA#> arfima 0.00 0 0 0 NA NA#> archm 0.00 0 0 0 NA NA#> mxreg 0.00 0 0 0 NA NA#> sigma 0.20 1 1 0 NA NA #> alpha 0.00 0 0 0 NA NA#> beta 0.00 0 0 0 NA NA#> gamma 0.00 0 0 0 NA NA#> eta1 0.00 0 0 0 NA NA#> eta2 0.00 0 0 0 NA NA#> delta 0.00 0 0 0 NA NA#> lambda 0.00 0 0 0 NA NA#> vxreg 0.00 0 0 0 NA NA#> skew 0.00 0 0 0 NA NA#> shape 0.00 0 0 0 NA NA#> ghlambda 0.00 0 0 0 NA NA#> xi 0.00 0 0 0 NA NA NAfixed.pars#> $mu#> \[1\] 0.01#> #> $ar1#> \[1\ ] -0.9#> #> $sigma#> \[1\] 0.2true_params#> mu ar1 sigma #> 0.01 -0.90 0.20

Then, we can generate the time series:

# Simulate a path apath(spec, n.sim = T)# Convert to xts and plot plot(synth\_log\_returns, main = "Log returns for ARMA model" plot(synth\_log\_prices, main = " Logarithmic price of ARMA model"

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

ARMA model

Now, we can estimate the parameters (which we already know):

# Specify AR(1) model arfimaspec(mean.model = list(armaOrder = c(1,0), include.mean = TRUE))# Estimation model #> mu ar1 sigma #> 0.0083 -0.8887 0.1987#> mu ar1 sigma #> 0.01 -0.90 0.20

We can also study the effect of the sample size T on the parameter estimation error:

# loop for (T_ in T\_sweep) { estim\_coeffs\_vs\_T <- rbind(estim\_coeffs\_vs\_T, coef(arma\_fit)) error\_coeffs\_vs\_T <- rbind(error\ _coeffs\_vs\_T, abs(coef(arma\_fit) - true\_params)/true\_params)# plot matplot(T\_sweep, estim\_coeffs\_vs_T, main = "Estimated ARMA coefficients", xlab = " T", ylab = "value",

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

matplot(T\_sweep, 100*error\_coeffs\_vs\_T, main = "Relative error in estimating ARMA coefficients", xlab = "T", ylab = "Error (%)",

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

First, the true μ is almost zero, so the relative error may appear unstable. After T = 800 samples, the other coefficients are well estimated.

ARMApredict

For a sanity check, we will now compare the two packagesForecastandrugarch results

# specify AR(1) model with given coefficients and parameters spec(mean.model = list(armaOrder = c(1,0), include.mean = TRUE), fixed.pars = list(mu = 0.005, ar1 = -0.9, sigma = 0.1))# Generate a sequence of length 1000 arfima(arma\_fixed\_spec, n.sim = 1000)@path$seriesSim # Specify and fit a model using the rugarch package spec(mean.model = list( armaOrder = c(1,0), include.mean = TRUE))# Fit a model using package "forecast" #> ARIMA(1,0,0) with non-zero mean #> #> Coefficients:#> ar1 mean #> -0.8982 0.0036#> se 0.0139 0.0017#> #> sigma^2 estimated as 0.01004: log likelihood=881.6#> AIC=-1757.2 AICc=-1757.17 BIC=-1742.47# Compare model coefficients#> ar1 intercept sigma #> -0.898181148 0.003574781 0.100222964#> mu ar1 sigma #> 0.003605805 -0.898750138 0.100199956

Indeed, both packages gave the same results.

ARMA model selection

In previous experiments, we assumed that we knew the order of the ARMA model, i.e. p = 1 and q = 0. In reality, the order is unknown, so different order combinations have to be tried. The higher the order, the better the fit, but this will inevitably lead to overfitting. Many methods have been developed to penalize the increase in complexity to avoid overfitting, such as AIC, BIC, SIC, HQIC, etc.

# Try different combinations# View ranking#> AR MA Mean ARFIMA BIC converged#> 1 1 0 1 0 -0.38249098 1#> 2 1 1 1 0 -0.37883157 1#> 3 2 0 1 0 -0.37736340 1#> 4 1 2 1 0 -0.37503980 1#> 5 2 1 1 0 -0.37459177 1#> 6 3 0 1 0 -0.37164609 1#> 7 1 3 1 0 -0.37143480 1#> 8 2 2 1 0 -0.37107841 1#> 9 3 1 1 0 -0.36795491 1#> 10 2 3 1 0 -0.36732669 1#> 11 3 2 1 0 -0.36379209 1#> 12 3 3 1 0 -0.36058264 1#> 13 0 3 1 0 -0.11875575 1#> 2 1 0 0.02957266 1#> 15 0 1 1 0 0.39326050 1#> 16 0 0 1 0 1.17294875 1# choose the best armaOrder#> AR MA #> 1 0

In this case, since the number of observations T = 1000 is large enough, the order is correctly detected. Conversely, if you try to use T=200, the detected order is p=1, q=3.

ARMA forecast

Once the ARMA model parameters ϕi^i and θ^j have been estimated, the model can be used to predict future values. For example, the prediction of xt based on past information is

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

And the prediction error will be xt-x^t=wt (assuming the parameters have been estimated) with a variance of σ2. packagerugarchMakes predictions on out-of-sample data simple:

# Estimated model (excluding out-of-sample) coef(arma\_fit)#> mu ar1 sigma #> 0.007212069 -0.898745183 0.200400119# Predicted log-returns out-of-sample forecast\_log\_returns <- xts(arma\[email protected]$ seriesFor\[1, \], dates\_out\_of\_sample)# restore log price prev\_log\_price <- head(tail(synth\_log\_prices, out\_of\_sample+1), out\_of \_sample)# log-return graph plot(cbind("fitted" = fitted(arma\_fit),# log-price graph plot(cbind("forecast" = forecast\_log\_prices, main = "log price forecast" , legend.loc = "topleft")

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

Multivariate VARMA Model

The VARMA(p,q) model on the logarithmic rate of return xt is

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

where wt is a white noise sequence with zero mean and covariance matrix Σw. The parameters of this model are the vector/matrix coefficients ϕ0, Φi, Θj and the noise covariance matrix Σw.

Compare

Let's load the S&P500 first:

# Load S&P 500 data head(SP500\_index\_prices)#> SP500#> 2012-01-03 1277.06#> 2012-01-04 1277.30#> 2012-01-05 1281.06#> 2012-01-06 1277.81# > 2012-01-09 1280.70#> 2012-01-10 1292.08# Prepare training and test data logreturns\_trn <- logreturns\[1:T\_trn\]logreturns\_tst <- logreturns\[-c(1:T \_trn)\]# Plot { plot(logreturns, addEventLines(xts("training"

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

Now, we use the training data (i.e., for t=1,…,Ttrnt=1,…,Ttrn) to fit different models (note that out-of-sample data is excluded by indicatingout.sample = T_tst). In particular, we will consider iid models, AR models, ARMA models, and some ARCH and GARCH models (variance modeling will be studied in more detail later).

# Fit the iid model coef(iid\_fit)#> mu sigma #> 0.0005712982 0.0073516993mean(logreturns\_trn)#> \[1\] 0.0005681388sd(logreturns\_trn)#> \[1\] 0.007360208# Fitting AR(1) model coef(ar\_fit)#> mu ar1 sigma #> 0.0005678014 -0.0220185181 0.0073532716# Fit ARMA(2,2) model coef(arma\_fit)#> mu ar1 ar2 ma1 ma2 sigma #> 0.0007223304 0.0268612 0.9095552008 -0.0832923604 -0.9328475211 0.0072573570# Fit ARMA(1,1) + ARCH(1) model coef(arch\_fit)#> mu ar1 ma1 omega alpha1 #> 6.321441e-04 8.720929e-02 -9.398810 e-05 9.986975e-02#Fit ARMA(0,0)+ARCH(10) model coef(long\_arch\_fit)#> mu omega alpha1 alpha2 alpha3 alpha4 alpha5 #> 7.490786e-04 2.452099e-05 6.888561 e-02 7.207551e-02 1.419938e-01 1.909541e-02 3.082806e-02 #> alpha6 alpha7 alpha8 alpha9 alpha10 #> 4.026539e-02 3.050040e-07 9.260183e -02 1.150128e-01 1.068426e-06# Fit ARMA(1,1)+GARCH(1,1) model coef(garch_fit)#> mu ar1 ma1 omega alpha1 beta1 #> 6.660346e-04 9.664597e-01 - 1.000000e+00 7.066506e-06 1.257786e-01 7.470725e-01

We use different models to predict log returns:

# Prepare to forecast logarithmic returns for out-of-sample periods # iid model forecast forecast(iid\_fit, n.ahead = 1, n.roll = T\_tst - 1) dates\_out\_of\_sample)# AR(1) model Make prediction forecast(ar\_fit, n.ahead = 1, n.roll = T\_tst - 1) dates\_out\_of\_sample)# ARMA(2,2) model make prediction forecast(arma\_fit, n. ahead = 1, n.roll = T\_tst - 1) dates\_out\_of\_sample)# Use ARMA(1,1) + ARCH(1) model for prediction forecast(arch\_fit, n.ahead = 1, n.roll = T\_tst - 1) dates\_out\_of\_sample)# ARMA(0,0)+ARCH(10) model forecast forecast(long\_arch\_fit, n.ahead = 1, n.roll = T\_tst - 1) dates\_out\_of\_sample)# ARMA(1,1)+GARCH(1,1) model forecast forecast(garch\_fit, n.ahead = 1, n.roll = T\_tst - 1) dates\_out\_of_sample)

We can calculate the prediction errors (in-sample and out-of-sample) for different models:

print(error_var)#>                           in-sample out-of-sample#> iid                    5.417266e-05  8.975710e-05#> AR(1)                  5.414645e-05  9.006139e-05#> ARMA(2,2)              5.265204e-05  1.353213e-04#> ARMA(1,1) + ARCH(1)    5.415836e-05  8.983266e-05#> ARCH(10)               5.417266e-05  8.975710e-05#> ARMA(1,1) + GARCH(1,1) 5.339071e-05  9.244012e-05

We can observe that as the model complexity increases, the in-sample error tends to be smaller (due to higher degrees of freedom in fitting the data), although the difference is negligible. What matters is actually the out-of-sample error: we can see that increasing the model complexity may give poorer results. In terms of error in predicting returns, it seems that the simplest iid model is sufficient.

Finally, let's show some graphs of out-of-sample error:

plot(error, main = "Out-of-Sample Error in Return Forecasts for Different Models",

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

Note that since we did not refit the model, the error increases over time (especially for ARCH modeling).

Scrolling window comparison

Let's first compare the concepts of static forecasting versus rolling forecasting with a simple example:

#ARMA(2,2) model spec <- spec(mean.model = list(armaOrder = c(2,2), include.mean = TRUE))# static fit and prediction ar\_static\_fit <- fit( spec = spec, data = logreturns, out.sample = T\_tst)#Roll fit and prediction modelroll <- aroll(spec = spec, data = logreturns, n.ahead = 1, #Prediction plot plot(cbind("static forecast" = ar\_static\_fore\_logreturns, main = "Prediction using ARMA(2,2) model", legend.loc = "topleft")# forecast error plot plot(error_logreturns, col = c("black", "red"), lwd = 2, main = "Prediction error of ARMA(2,2) model", legend.loc = "topleft")

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

We can clearly observe the effect of the rolling window process on the time series.

Now, we can redo all predictions for all models on a rolling window basis:

# Rolling forecast based on iid model roll(iid\_spec, data = logreturns, n.ahead = 1, forecast.length = T\_t# Rolling forecast of AR(1) model roll(ar\_spec, data = logreturns, n .ahead = 1, forecast.length = T\_tst, # Rolling forecast for ARMA(2,2) model roll(arma\_spec, data = logreturns, n.ahead = 1, forecast.length = T\_tst, # ARMA (1,1) + ARCH(1) model rolling forecast roll(arch\_spec, data = logreturns, n.ahead = 1, forecast.length = T\_tst, refit.every = 50, refit.win# ARMA( 0,0) + ARCH(10) model rolling forecast roll(long\_arch\_spec, data = logreturns, n.ahead = 1, forecast.length = T\_tst, refit.every = 50, # ARMA(1, 1) + GARCH(1,1) model rolling forecast roll(garch\_spec, data = logreturns, n.ahead = 1, forecast.length = T_tst, refit.every = 50, refit.window

Let's look at the forecast error in the rolling benchmark case:

print(rolling\_error\_var)#>                           in-sample out-of-sample#> iid                    5.417266e-05  8.974166e-05#> AR(1)                  5.414645e-05  9.038057e-05#> ARMA(2,2)              5.265204e-05  8.924223e-05#> ARMA(1,1) + ARCH(1)    5.415836e-05  8.991902e-05#> ARCH(10)               5.417266e-05  8.976736e-05#> ARMA(1,1) + GARCH(1,1) 5.339071e-05  8.895682e-05

and some graphs:

plot(error_logreturns, main = "Rolling prediction errors for different models", legend.loc = "topleft"

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

We see that now all models fit the time series. Furthermore, we did not find any significant differences between the models.

We can finally compare static and rolling errors:

barplot(rbind(error\_var\[, "out-of-sample"\], rolling\_error_var\[, "out-of-sample"\]) col = c("darkblue", "darkgoldenrod"), legend = c("Static Forecast", "Rolling Forecast"),

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

We can see that rolling forecasts are necessary in some cases. So, in practice, we need rolling forecast improvements on a regular basis.

variance model

ARCH and GARCH models

The ARCH(m) model for the log-return residual wt is

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

where zt is a white noise sequence with zero mean and constant variance, and the conditional variance σ2t is modeled as

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

Among them, m is the model order, ω>0, αi≥0 are parameters.

The GARCH(m,s) model extends the ARCH model with a recursive term on σ2t:

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

The parameters ω\> 0, αi≥0, βj≥0 need to satisfy the stability of ∑mi = 1αi+ ∑sj = 1βj≤1.

rugarchGenerate data

First, we need to define the model:

# Specify a GARCH model with given coefficients and parameters #> #> *---------------------------------* #> * GARCH Model Spec *#> *---------------------------------*#> #> Conditional Variance Dynamics #> ------------------------------------#> GARCH Model : sGARCH(1,1) #> Variance Targeting : FALSE #> #> Conditional Mean Dynamics#> ----------------------------------- -#> Mean Model : ARFIMA(1,0,0)#> Include Mean : TRUE #> GARCH-in-Mean : FALSE #> #> Conditional Distribution#> ------------- -----------------------#> Distribution : norm #> Includes Skew : FALSE #> Includes Shape : FALSE #> Includes Lambda : FALSE#> Level Fixed Include Estimate LB UB#> mu 0.005 1 1 0 NA NA#> ar1 -0.900 1 1 0 NA NA#> ma 0.000 0 0 0 NA NA#> arfima 0.000 0 0 0 NA NA#> archm 0.000 0 0 0 NA NA #> mxreg 0.000 0 0 0 NA NA#> omega 0.001 1 1 0 NA NA #> alpha1 0.300 1 1 0 NA NA#> beta1 0.650 1 1 0 NA NA#> gamma 0.000 0 0 0 NA NA#> eta1 0.000 0 0 0 NA NA#> eta2 0.000 0 0 0 NA NA#> delta 0.000 0 0 0 NA NA#> lambda 0.000 0 0 0 NA NA#> vxreg 0.000 0 0 0 NA NA#> skew 0.000 0 0 0 NA NA#> shape 0.000 0 0 0 NA NA#> ghlambda 0.000 0 0 0 NA NA# > xi 0.000 0 0 0 NA NA#> $mu#> \[1\] 0.005#> #> $ar1#> \[1\] -0.9#> #> $omega#> \[1\] 0.001# > #> $alpha1#> \[1\] 0.3#> #> $beta1#> \[1\] 0.65true_params#> mu ar1 omega alpha1 beta1 #> 0.005 -0.900 0.001 0.300 0.650

Then, we can generate the returns time series:

# Simulate a path hpath(garch\_spec, n.sim = T)#> num \[1:2000, 1\] 0.167 -0.217 # Plot log returns { plot(synth\_log\_returns, main = "GARCH model ", lwd = 1.5) lines(synth\_volatility

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

GARCH

Now, we can estimate the parameters:

# Specify a GARCH model ugarchspec(mean.model = list(armaOrder = c(1,0)# Estimate the model coef(garch_fit)#> mu ar1 omega alpha1 beta1 #> 0.0036510100 -0.8902333595 0.00081811434 0.2810460728 0.6717486402 omega ar beta1 o #> 0.005 -0.900 0.001 0.300 0.650# Coefficient error#> mu ar1 omega alpha1 beta1 #> 0.0013489900 0.0097666405 0.0001188566 0.0189539272 0.0217486402

We can also study the effect of the sample size T on the parameter estimation error:

# 循环for (T_ in T\_sweep) {  garch\_fit   error\_coeffs\_vs\_T <- rbind(error\_coeffs\_vs\_T, abs((coef(garch\_fit) - true\_params)/true\_params))  estim\_coeffs\_vs\_T <- rbind(estim\_coeffs\_vs\_T, coef(garch\_fit))# 绘图matplot(T\_sweep, 100*error\_coeffs\_vs\_T,         main = "估计GARCH系数的相对误差", xlab = "T", ylab = "误差 (%)",

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

The true ω is almost zero, so the error is very unstable. As for the other coefficients, as in the ARMA case, the estimates of μ are really poor (more than 50% relative error), while the other coefficients seem to be well estimated after T=800 samples.

Comparison of GARCH results

As a sanity check, we will now compare the two packagesfGarchandrugarch results

# Specify ARMA(0,0)-GARCH(1,1) with specific parameter values ​​as data generation process garch\_spec #Generate data with length 1000 path(garch\_fixed\_spec, n.sim = 1000)@path $# Specify and fit model using 'rugarch' package rugarch\_fit &lt;- ugarchfit(spec = garch\_spec, data = x) # Fit model using 'fGarch' package garchFit(formula = ~ garch(1, 1), data = x, trace = FALSE)# compare model coefficients #&gt; mu omega alpha1 beta1 #&gt; 0.09749904 0.01395109 0.13510445 0.73938595#&gt; mu omega alpha1 beta1 #&gt; 0.09750394 0.0139G2648 0.13527024 compare standard deviation of fit f5_839716 #&gt; \[1\] 0.3513549 0.3254788 0.3037747 0.2869034 0.2735266 0.2708994

Indeed, both packages gave the same results.

userugarchpackage for GARCH prediction

Once the parameters of a GARCH model have been estimated, the model can be used to predict future values. For example, a one-step prediction of conditional variance based on past information is

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

Given ω^/(1-∑mi=1α^i-∑sj=1β^j). packagerugarchMakes predictions on out-of-sample data simple:

# Estimate the model, excluding out-of-sample garch\_fit coef(garch\_fit)#&gt; mu ar1 omega alpha1 beta1 #&gt; 0.0034964331 -0.8996287630 0.0006531088 0.3058756796 0.6815452241# Predict the log return for the entire sample garch\[email protected]$sigmaFor\[forecast$sigmaFor\ , \]# log-return plot plot(cbind(&quot;fitted&quot; = fitted(garch\_fit), main = &quot;synthetic log-return forecast&quot;, legend.loc = &quot;topleft&quot;)

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

#Volatility log-return plot plot(cbind(&quot;fitted volatility&quot; = sigma(garch_fit), main = &quot;Predicted volatility of synthetic log-return&quot;, legend.loc = &quot;topleft&quot;)

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

different methods

Let's load the S&amp;P500 first:

# Load S&amp;P 500 index data head(SP500\_index\_prices)#&gt; SP500#&gt; 2008-01-02 1447.16#&gt; 2008-01-03 1447.16#&gt; 2008-01-04 1411.63#&gt; 2008-01-07 1416.18 #&gt; 2008-01-08 1390.19#&gt; 2008-01-09 1409.13# Prepare training and test data x\_trn &lt;- x\[1:T\_trn\]x\_tst &lt;- x\[-c(1: T\_trn)\]# plot { plot(x, main = &quot;revenue&quot; addEventLines(xts(&quot;training&quot;, in

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

constant

Let's start with constants:

plot(cbind(sqrt(var\_constant), x\_trn) main = &quot;constant&quot;)

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

moving average

Now, let's use a moving average of squared returns:

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

plot(cbind(sqrt(var\_t), x\_trn), main = &quot;Envelope based on simple rolling square mean (time period=20)

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

EWMA

Exponentially Weighted Moving Average (EWMA):

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

Note that this can also be modeled as an ETS(A,N,N) state space model:

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

plot(cbind(std\_t, x\_trn), main = &quot;Envelope based on squared EWMA&quot;)

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

Multiplication ETS

We can also try different variants of the ETS model. For example, a multiplicative noise version ETS(M,N,N) with a state-space model:

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

plot(cbind(std\_t, x\_trn), col = c(&quot;red&quot;, &quot;black&quot;) main = &quot;Envelope of square-based ETS(M,N,N)&quot;

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

ARCH

Now, we can use more complex ARCH modeling:

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

plot(cbind(std\_t, x\_trn), col = c(&quot;red&quot;, &quot;black&quot;) main = &quot;ARCH(5) based envelope&quot;)

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

GARCH

We can lift the model to GARCH:

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

plot(cbind(std\_t, x\_trn), col = c(&quot;red&quot;, &quot;black&quot;) main = &quot;GARCH(1,1) based envelope&quot;)

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

SV Stochastic Volatility

Finally, we can use stochastic volatility modeling:

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

or, equivalently,

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

plot(cbind(std\_t, x\_trn), col = c(&quot;red&quot;, &quot;black&quot;), main = &quot;Envelope Analysis Based on Stochastic Volatility&quot;)

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

Compare

We can now compare the error in the variance estimates for each method over the out-of-sample period:

#&gt; MA EWMA ETS(M,N,N) ARCH(5) GARCH(1,1) SV #&gt; 2.204965e-05 7.226188e-06 3.284057e-06 7.879039e-05 6.496545e-06 6.705059e-06barplot( error_all, main = &quot;Error in out-of-sample variance estimation&quot;

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

Scrolling window comparison

Rolling window comparison of six methods: MA, EWMA, ETS(MNN), ARCH(5), GARCH(1,1) and SV.

#rolling window lookback &lt;- 200len\_tst &lt;- 40for (i in seq(lookback, T-len\_tst, by = len\_tst)) { # MA var\_t &lt;- roll\_meanr(x\_trn^2, n = 20, fill = NA) var\_fore &lt;- var(x\_trn/sqrt(var\_t), na.rm = TRUE) * tail(var\_t, 1) error\_ma &lt;- c(error\ _ma, abs(var\_fore - var\_tst)) # EWMA error\_ewma &lt;- c(error\_ewma, abs(var\_fore - var\_tst)) # ETS(M,N,N) error\_ets\ _mnn &lt;- c(error\_ets\_mnn, abs(var\_fore - var\_tst)) # ARCH error\_arch &lt;- c(error\_arch, abs(var\_fore - var\_tst)) # GARCH error\ _garch &lt;- c(error\_garch, abs(var\_fore - var\_tst)) # SV error\_sv &lt;- c(error\_sv, abs(var\_fore - var\_tst))}barplot(error_all, main = &quot;Variance Estimation Error&quot;,

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

Multivariate GARCH Model

For illustration purposes, we will only consider the Constant Conditioned Correlation (CCC) and Dynamic Conditioned Correlation (DCC) models as they are the most popular. The logarithmic rate of return residual wt is modeled as

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

where zt is an iid white noise sequence with zero mean and constant covariance matrix II. The conditional covariance matrix Σt is modeled as

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

where Dt = Diag(σ1,t,…,σN,t) is the normalized noise vector C, and the covariance matrix ηt=C-1wt (ie, it contains diagonal elements equal to 1).

Basically, with this model, the diagonal matrix Dt contains a set of univariate GARCH models, and then the matrix C contains some correlations between the sequences. The main disadvantage of this model is that the matrix C is constant. To overcome this problem, DCC is proposed as

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

where Ct contains diagonal elements equal to 1. To enforce a diagonal element equal to 1, Engle models it as

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

Qt has arbitrary diagonal elements and follows the model

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

We will generate data, estimate parameters and forecasts.

Start by loading the multivariate ETF data:

  • SPDR S&P 500 ETF
  • 20+ Year Treasury Bond ETFs
  • IEF: 7-10 Year Treasury Bond ETF
# Download data prices &lt;- xts()head(prices)#&gt; SPY TLT IEF#&gt; 2013-01-02 127.8779 99.85183 93.65224#&gt; 2013-01-03 127.5890 98.49886 93.17085#&gt; 2013-01-04 126# 9&gt;3.849 2013-01-07 127.7991 98.92480 93.26714#&gt; 2013-01-08 127.4314 99.57622 93.49468#&gt; 2013-01-09 127.7553 99.48438 93.54719# Plot three &quot;logarithmic price series of three ETFs&quot; Number price&quot;, legend.loc = &quot;topleft&quot;)

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

First, we define the model:

# Specify iid univariate time series model ugarch_spec # Specify DCC model spec( multispec(replicate(spec, n = 3))

Next, we fit the model:

# 估计模型#> #> *---------------------------------*#> *          DCC GARCH Fit          *#> *---------------------------------*#> #> Distribution         :  mvnorm#> Model                :  DCC(1,1)#> No. Parameters       :  44#> \[VAR GARCH DCC UncQ\] : \[30+9+2+3\]#> No. Series           :  3#> No. Obs.             :  1007#> Log-Likelihood       :  12198.4#> Av.Log-Likelihood    :  12.11 #> #> Optimal Parameters#> -----------------------------------#>               Estimate  Std. Error   t value Pr(>|t|)#> \[SPY\].omega   0.000004    0.000000  11.71585 0.000000#> \[SPY\].alpha1  0.050124    0.005307   9.44472 0.000000#> \[SPY\].beta1   0.870051    0.011160  77.96041 0.000000#> \[TLT\].omega   0.000001    0.000001   0.93156 0.351563#> \[TLT\].alpha1  0.019716    0.010126   1.94707 0.051527#> \[TLT\].beta1   0.963760    0.006434 149.79210 0.000000#> \[IEF\].omega   0.000000    0.000001   0.46913 0.638979#> \[IEF\].alpha1  0.031741    0.023152   1.37097 0.170385#> \[IEF\].beta1   0.937777    0.016498  56.84336 0.000000#> \[Joint\]dcca1  0.033573    0.014918   2.25044 0.024421#> \[Joint\]dccb1  0.859787    0.079589  10.80278 0.000000#> #> Information Criteria#> ---------------------#>                     #> Akaike       -24.140#> Bayes        -23.925#> Shibata      -24.143#> Hannan-Quinn -24.058#> #> #> Elapsed time : 0.8804049

We can plot the time-varying correlation:

# Extract time-varying covariance and correlation matrix dim(dcc\_cor)#&gt; \[1\] 3 3 1007#plot(corr\_t main = &quot;time-varying correlation&quot;, legend.loc = &quot;left&quot;)

R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

We see a very high and fairly stable correlation between the two income ETFs. The correlation with SPY is small and fluctuates in an interval less than 0.


R language uses multivariate ARMA, GARCH, EWMA, ETS, stochastic volatility SV model to model financial time series data

Most Popular Insights

1.HAR-RV-J and Recurrent Neural Network (RNN) Hybrid Model Predicting and Trading High-Frequency Volatility in Large Equity Indices

2.HAR-RV Model Based on Mixed Data Sampling (MIDAS) Regression in R Language to Predict GDP GrowthThe regression HAR-RV model predicts GDP growth&quot;)

3.Realization of Volatility: ARCH Model and HAR-RV Model

4.R language ARMA-EGARCH model and integrated forecasting algorithm to forecast the actual volatility of SPX

5.VaR comparison of GARCH(1,1), MA and historical simulation method

6.R language multivariate COPULA GARCH model time series forecast

7.VAR fitting and prediction based on ARMA-GARCH process in R language

8.Matlab predicts ARMA-GARCH conditional mean and variance model

9.ARIMA+GARCH trading strategy for S&amp;P500 stock index in R language