heavyModelC: HEAVY Model estimation using C code
In highfrequency: Tools for Highfrequency Data Analysis

Description Usage Arguments Details Value Author(s) References Examples

This function is the same as heavyModel function, except for using C code to speed up the calculation process.

This function calculatest the High frEquency bAsed VolatilitY (HEAVY) model proposed in Shephard and Sheppard (2010). This function is used as a predictive volatility model built to exploit highfrequency data.

1
2
3

heavyModelC(data, p=matrix( c(0,0,1,1),ncol=2 ), q=matrix( c(1,0,0,1),ncol=2 ), 
          startingvalues = NULL, LB = NULL, UB = NULL, 
            backcast = NULL, compconst = FALSE);

`data`	a (T x K) matrix containing the data, with T the number of days. For the traditional HEAVY model: K = 2, the first column contains the squared daily demeaned returns, the second column contains the realized measures.
`p`	a (K x K) matrix containing the lag length for the model innovations. Position (i, j) in the matrix indicates the number of lags in equation i of the model for the innovations in data column j. For the traditional heavy model p is given by matrix( c(0,0,1,1),ncol=2 ) (default).
`q`	a (K x K) matrix containing the lag length for the conditional variances. Position (i, j) in the matrix indicates the number of lags in equation i of the model for conditional variances corresponding to series j. For the traditionalheavy model introduced above q is given by matrix( c(1,0,0,1),ncol=2 ) (default).
`startingvalues`	a vector containing the starting values to be used in the optimization to find the optimal parameters estimates.
`LB`	a vector of length K indicating the lower bounds to be used in the estimation. If NULL it is set to a vector of zeros by default.
`UB`	a vector of length K indicating the upper bounds to be used in the estimation. If NULL it is set to a vector of Inf by default.
`backcast`	a vector of length K used to initialize the estimation. If NULL the unconditional estimates are taken.
`compconst`	a boolean variable. In case TRUE, the omega values are estimated in the optimization. In case FALSE, volatility targeting is done and omega is just 1 minus the sum of all relevant alpha's and beta's multiplied by the unconditional variance.

Assume there are T daily returns and realized measures in the period t. Let r_i and RM_i be the i^{th} daily return and daily realized measure respectively (with i=1, …,T).

The most basic heavy model is the one with lag matrices p of ≤ft( \begin{array}{ccc} 0 & 1 \\ 0 & 1 \end{array} \right) and q of ≤ft( \begin{array}{ccc} 1 & 0 \\ 0 & 1 \end{array} \right). This can be reprensented by the following equations:

\mbox{var}{≤ft(r_t \right)} = h_t = w + α RM_{t-1} + β h_{t-1}; w,α ≥q 0, β \in [0,1]

\mbox{E}{≤ft(RM_t \right)} = μ_t = w_R + α_R RM_{t-1} + β_R μ_{t-1}; w_R,α_R, β_R ≥q 0, α_R+β_R \in [0,1]

Equivalently, they can be presented in terms of matrix notation as below:

≤ft( \begin{array}{ccc} h_t \\ μ_t \end{array} \right) = ≤ft( \begin{array}{ccc} w \\ w_R \end{array} \right) + ≤ft( \begin{array}{ccc} 0 & α \\ 0 & α_R \end{array} \right) ≤ft( \begin{array}{ccc} r^2_{t-1} \\ RM_{t-1} \end{array} \right) + ≤ft( \begin{array}{ccc} β & 0 \\ 0 & β_R \end{array} \right) ≤ft( \begin{array}{ccc} h_{t-1} \\ μ_{t-1} \end{array} \right)

In this version, the parameters vector to be estimated is ≤ft( w, w_R,α, α_R, β, β_R \right) .

In terms of startingvalues, Shephard and Sheppard recommend for this version of the Heavy model to set β be around 0.6 and sum of α+β to be close to but slightly less than one.

In general, the lag length for the model innovation and the conditional covariance can be greater than 1. Consider, for example, matrix p is ≤ft( \begin{array}{ccc} 0 & 2 \\ 0 & 1 \end{array} \right) and matrix q is the same as above. Matrix notation will be as below:

≤ft( \begin{array}{ccc} h_t \\ μ_t \end{array} \right) = ≤ft( \begin{array}{ccc} w \\ w_R \end{array} \right) + ≤ft( \begin{array}{ccc} 0 & α_1 \\ 0 & α_R \end{array} \right) ≤ft( \begin{array}{ccc} r^2_{t-1} \\ RM_{t-1} \end{array} \right) +≤ft( \begin{array}{ccc} 0 & α_2 \\ 0 & 0 \end{array} \right) ≤ft( \begin{array}{ccc} r^2_{t-2} \\ RM_{t-2} \end{array} \right) + ≤ft( \begin{array}{ccc} β & 0 \\ 0 & β_R \end{array} \right) ≤ft( \begin{array}{ccc} h_{t-1} \\ μ_{t-1} \end{array} \right)

In this version, the parameters vector to be estimated is ≤ft( w, w_R,α_1, α_R, α_2, β, β_R \right) .

A list with the following values: (i) loglikelihood: The log likelihood evaluated at the parameter estimates. (ii) likelihoods: an xts object of length T containing the log likelihoods per day. (iii) condvar: a (T x K) xts object containing the conditional variances (iv) estparams: a vector with the parameter estimates. The order in which the parameters are reported is as follows: First the estimates for omega then the estimates for the non-zero alpha's with the most recent lags first in case max(p) > 1, then the estimates for the non-zero beta's with the most recent lag first in case max(q) > 1. (v) convergence: an integer code indicating the successfulness of the optimization. See optim for more information.

Giang Nguyen, Jonathan Cornelissen and Kris Boudt

Shephard, N. and K. Sheppard (2010). Realising the future: forecasting with high frequency based volatility (heavy) models. Journal of Applied Econometrics 25, 197-231.

 
 # Implementation of the heavy model on DJI:
 data("realized_library");
 returns =  realized_library$Dow.Jones.Industrials.Returns; 
 rk      =  realized_library$Dow.Jones.Industrials.Realized.Kernel; 
 returns = returns[!is.na(rk)];  rk = rk[!is.na(rk)]; # Remove NA's 
 data = cbind( returns^2, rk ); # Make data matrix with returns and realized measures
 backcast = matrix( c(var(returns),mean(rk)) ,ncol=1);
 
 #For traditional (default) version:
 startvalues = c(0.004,0.02,0.44,0.41,0.74,0.56); # Initial values;
 output = heavyModelC( data = as.matrix(data,ncol=2), compconst=FALSE, 
                    startingvalues = startvalues, backcast=backcast); 
                    
 #For general version: 
 startvalues = c(0.004,0.02,0.44,0.4,0.41,0.74,0.56); # Initial values;
 p = matrix(c(2, 0,0 , 1), ncol = 2);
 q = matrix(c(1,0, 0, 1), ncol = 2);
 
 output = heavyModelC( data = as.matrix(data,ncol=2), p=p, q=q, compconst=FALSE, 
                    startingvalues = startvalues, backcast=backcast);