Below we address implementation of volatility models as discussed in chapter three of the book chapter3

GARCH models have to be solved using maximum likelihood methods, how the GARCH packages, like `rmgarch`, fit their models, discussed in Chapter 10. For the mathematical details, see the book or slides, see package documentation for how to use `rmgarch`.

We use maximum likelihood methods for estimating volatility models, and use the R package `rmgarch` for the actual implementation. See the package documentation for more details, manual and a more detailed vignette. It is developed by Alexios Galanos and is constantly maintained and updated. The development code can be found on github.

## 12.1 Libraries

``suppressPackageStartupMessages(library(rmgarch))``

## 12.2 Data

The data we use is `stocks.csv`.

``source('common.r')``
``````load('data/data.RData')
Return=data\$Return
Ticker=data\$Ticker``````

## 12.3 EWMA

``````y=as.matrix(Return[,c("JPM","INTC")])
EWMA = matrix(nrow=dim(y)[1],ncol=3)
lambda = 0.94
S = cov(y)
EWMA[1,] = c(S)[c(1,4,2)]
for (i in 2:dim(y)[1]){
S = lambda*S+(1-lambda)*
y[i-1,] %*% t(y[i-1,])
EWMA[i,] = c(S)[c(1,4,2)]
}
rhoEWMA = EWMA[,3]/sqrt(EWMA[,1]*EWMA[,2])
plot(rhoEWMA,type='l')``````

## 12.4 DCC

### 12.4.1`rmgarch`

The package `rmgarch` allows us easily estimate models with multiple assets in the same fashion as `rugarch`. The procedure is analogous, we first need to specify the model we want to use, and then fit it to the data. However, specifying the model has extra steps. To explain this, we will focus on the estimation of a DCC model.

In a DCC model, we assume that each individual asset follows some type of univariate model, usually a GARCH. And then we model the correlation between the assets using an ARMA-like process.

The process for fitting a DCC model using `rmgarch` is then:

1. Specify the univariate volatility model each asset follows using `ugarchspec()`;
2. Use the `multispec()` function to create a multivariate specification. This is essentially a list of univariate specifications. If we are going to use the same for every asset, we can use `replicate()`;
3. Then we need to create a `dccspec()` object, which takes in the list of univariate specifications for every asset, and the additional DCC joint specifications, like `dccOrder` and `distribution`;
4. Fit the specification to the data.

We will use the returns for JPM and Intel:

``y=as.matrix(Return[,c("JPM","INTC")])``

We will assume a simple GARCH(1,1) with mean zero for each stock. We will create a single univariate specification, and then replicate it into `multispec()`:

``````# Create the univariate specification
uni_spec = ugarchspec(
variance.model = list(
garchOrder = c(1,1)),
mean.model = list(
armaOrder = c(0,0),
include.mean = FALSE)
)
# Replicate it into a multispec element
mspec = multispec(replicate(2, uni_spec))``````

Let’s take a look inside the `mspec` object. We will see there are two univariate specifications, one per asset, and they are equal:

``mspec``
``````
*-----------------------------*
*     GARCH Multi-Spec        *
*-----------------------------*
Multiple Specifications : 2
Multi-Spec Type         : equal``````

You can check the specifications with `mspec@spec`

Now we proceed to create the specification for the DCC model using `dccspec()`:

``````spec = dccspec(
# Univariate specifications - Needs to be multispec
uspec = mspec,

# DCC specification. We will assume an ARMA(1,1)-like process
dccOrder = c(1,1),

# Distribution, here multivariate normal
distribution = "mvnorm"
)``````

We can call `spec` to see what is inside:

``spec``
``````
*------------------------------*
*       DCC GARCH Spec         *
*------------------------------*
Model          :  DCC(1,1)
Estimation     :  2-step
Distribution   :  mvnorm
No. Parameters :  9
No. Series     :  2``````

We can again see more details in `spec@model` and `spec@umodel`.

Now we can proceed to fit the specification to the data:

``````res = dccfit(spec, data = y)
res``````
``````
*---------------------------------*
*          DCC GARCH Fit          *
*---------------------------------*

Distribution         :  mvnorm
Model                :  DCC(1,1)
No. Parameters       :  9
[VAR GARCH DCC UncQ] : [0+6+2+1]
No. Series           :  2
No. Obs.             :  5034
Log-Likelihood       :  27227.76
Av.Log-Likelihood    :  5.41

Optimal Parameters
-----------------------------------
Estimate  Std. Error  t value Pr(>|t|)
[JPM].omega    0.000004    0.000003   1.5613  0.11846
[JPM].alpha1   0.092326    0.017139   5.3868  0.00000
[JPM].beta1    0.896338    0.020491  43.7438  0.00000
[INTC].omega   0.000004    0.000003   1.3599  0.17385
[INTC].alpha1  0.030623    0.004388   6.9782  0.00000
[INTC].beta1   0.959936    0.001737 552.5567  0.00000
[Joint]dcca1   0.009528    0.005876   1.6215  0.10491
[Joint]dccb1   0.973320    0.020506  47.4651  0.00000

Information Criteria
---------------------

Akaike       -10.814
Bayes        -10.802
Shibata      -10.814
Hannan-Quinn -10.810

Elapsed time : 1.237428 ``````

We can check what slots are inside:

``````names(res@model)
names(res@mfit)``````
1. 'modelinc'
2. 'modeldesc'
3. 'modeldata'
4. 'varmodel'
5. 'pars'
6. 'start.pars'
7. 'fixed.pars'
8. 'maxgarchOrder'
9. 'maxdccOrder'
10. 'pos.matrix'
11. 'pidx'
12. 'DCC'
13. 'mu'
14. 'residuals'
15. 'sigma'
16. 'mpars'
17. 'ipars'
18. 'midx'
19. 'eidx'
20. 'umodel'
1. 'coef'
2. 'matcoef'
3. 'garchnames'
4. 'dccnames'
5. 'cvar'
6. 'scores'
7. 'R'
8. 'H'
9. 'Q'
10. 'stdresid'
11. 'llh'
12. 'log.likelihoods'
13. 'timer'
14. 'convergence'
15. 'Nbar'
16. 'Qbar'
17. 'plik'
``````# Coefficient matrix
res@mfit\$matcoef``````
A matrix: 8 × 4 of type dbl
Estimate Std. Error t valuePr(>|t|)
[JPM].omega4.320478e-062.767254e-06 1.5612871.184561e-01
[JPM].alpha19.232608e-021.713926e-02 5.3868197.171575e-08
[JPM].beta18.963376e-012.049062e-02 43.7438010.000000e+00
[INTC].omega3.501196e-062.574546e-06 1.3599271.738529e-01
[INTC].alpha13.062290e-024.388355e-03 6.9782192.989387e-12
[INTC].beta19.599362e-011.737263e-03552.5567270.000000e+00
[Joint]dcca19.528089e-035.876088e-03 1.6215021.049100e-01
[Joint]dccb19.733201e-012.050600e-02 47.4651360.000000e+00
``````# Log likelihood
res@mfit\$llh``````
27227.7569568891

The matrix H inside `res@mfit` includes the covariances. It is 3-dimensional, since it includes the 2x2 covariance matrix for each of the T time periods:

``````H = res@mfit\$H
dim(H)``````
1. 2
2. 2
3. 5034
``````# First period's covariances
H[,,1]``````
 0.000522062 0.000194688 0.000194688 0.00039073

We can extract the conditional correlation in two ways. One is computing it from H:

``````# Initializing the vector
rhoDCC = vector(length = dim(y)[1])

# Populate with the correlations
rhoDCC = H[1,2,] / sqrt(H[1,1,]*H[2,2,])``````

## 12.5 Compare EWMA and DCC

``````par(mar=c(2,4,2,0))
matplot(cbind(rhoEWMA,rhoDCC),
type='l',
bty='l',
lty=1,
col=c("green","blue"),
main="EWMA and DCC correlations for JPM and C",
ylab="Correlations",
las=1
)
legend("bottomright",
legend=c("EWMA","DCC"),
lty=1,
col=c("green","blue"),
bty='n'
)``````