R/rgarch-meanres.R

In rgarch: Flexible GARCH modelling in R

#################################################################################
##
##   R package rgarch by Alexios Ghalanos Copyright (C) 2008, 2009, 2010, 2011
##   This file is part of the R package rgarch.
##
##   The R package rgarch is free software: you can redistribute it and/or modify
##   it under the terms of the GNU General Public License as published by
##   the Free Software Foundation, either version 3 of the License, or
##   (at your option) any later version.
##
##   The R package rgarch is distributed in the hope that it will be useful,
##   but WITHOUT ANY WARRANTY; without even the implied warranty of
##   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
##   GNU General Public License for more details.
##
#################################################################################

# Gets the mean residuals for initializing the garch equation (initialize the starting
# variance) for every pass
.getresiduals = function(data, include.mean, mu, armaOrder, ar, ma, include.inmean, 
						inmean, darfima, mexdata, mxreg, h)
{
	# Case 1 garchInMean yes:
	armap = armaOrder[1]
	armaq = armaOrder[2]
	im = armaOrder[3]
	arfima = armaOrder[4]
	include.mex = armaOrder[5]
	zrf = 0
	if(include.mex!=0){
		mexdata=mexdata
		mxn = dim(as.matrix(mexdata))[2]
		mxreg=mxreg
		xm = as.numeric(mexdata%*%matrix(mxreg,nrow=mxn))
	} else{
		xm=0
	}
	if(include.inmean){
		ixm = inmean*(h^im)
	} else{
		ixm = 0
	}
	if(armap!=0) ar =ar else ar=0
	if(armaq!=0) ma = ma else ma=0
	if(include.mean) mu = mu else mu=0
	if(armap!=0 && armaq!=0) myfix=c(ar,ma)
	if(armap!=0 && armaq==0) myfix=ar
	if(armap==0 && armaq!=0) myfix=ma

	if((armap+armaq)>0){
		xm  = data - mu - xm - ixm
		res = arima(x = xm, order = c(armap, 0, armaq), include.mean = FALSE,
				fixed = myfix, transform.pars = FALSE, optim.control = list(maxit = 0))$residuals
		if(arfima)
		{
			zrf = .fracdiff(c(1,rep(0,length(data)-1)), darfima=darfima)
			#2 stage arfima (arma>arfima...maybe we can use the fracdiff package
			# but it does not yet return residuals/ fArma has an implementation...
			res = .fracdiff(res, darfima=darfima)
		}
		if(any(is.na(res))) res[which(is.na(res))]=0
	} else{
		res = data - mu - xm - ixm
		if(arfima)
		{
			zrf = .fracdiff(c(1,rep(0,length(data)-1)), darfima=darfima)
			
			#2 stage arfima (arma>arfima...maybe we can use the fracdiff package
			# but it does not yet return residuals/ fArma has an implementation...
			res = .fracdiff(res, darfima=darfima)
		}
	}
	ans = list(res = res, zrf = zrf)
	return(ans)
}

# we use this method (ignore the arma recursion for the residuals) instead as it 
# is much faster and has minimal impact on the likelihood (the impact on the mean
# for the res^2 is very very small)
.getresiduals2 = function(data, include.mean, mu, armaOrder, ar, ma, include.inmean, 
		inmean, darfima, mexdata, mxreg, h)
{
	armap = armaOrder[1]
	armaq = armaOrder[2]
	im = armaOrder[3]
	arfima = armaOrder[4]
	include.mex = armaOrder[5]
	zrf = 0
	if(include.mex!=0){
		mexdata=mexdata
		mxn = dim(as.matrix(mexdata))[2]
		mxreg=mxreg
		xm = as.numeric(mexdata%*%matrix(mxreg,nrow=mxn))
	} else{
		xm=0
	}
	if(include.inmean){
		ixm = inmean*(h^im)
	} else{
		ixm =0 
	}
	if(armap!=0) ar=ar else ar=0
	if(armaq!=0) ma=ma else ma=0
	if(include.mean) mu=mu else mu=0
	if(armap!=0 && armaq!=0) myfix=c(ar,ma)
	if(armap!=0 && armaq==0) myfix=ar
	if(armap==0 && armaq!=0) myfix=ma
	
	res  = data - mu - xm - ixm
	if(arfima)
	{
		zrf = .fracdiff(c(1,rep(0,length(data)-1)), darfima=darfima)
		res = .fracdiff(res, darfima=darfima)
	}

	ans = list(res = res, zrf = zrf)
	return(ans)
}