[R-SIG-Finance] efficient code for nonlinear garch model
Patrick Burns
patrick at burns-stat.com
Thu May 1 00:42:02 CEST 2014
It should be easy enough to use 'Rcpp' to
write the likelihood in C++ and gain a
bunch of speed.
Pat
On 30/04/2014 20:21, Milos Cipovic wrote:
> Hi,
> I'm trying to think of efficient code for NGARCH model,
> sigma[t]^2=omega + alpha*(R[t-1]-theta*sigma[t-1])^2 + beta*(sigma[t-1)]^2
> from book (Elements of Financial Risk Menagment-Peter F Christoffersen page
> 77).
>
> I already coded it the "usual way" but it's taking too long (around 30 s).
> Now I'm trying somehow to implement -filter- function like in paper -
> Parameter Estimation of ARMA Models with
> GARCH/APARCH Errors An R and SPlus Software Implementation - page 6.
> The problem I can't overcome is that in this model I can't figure out how
> to present the middle part ( 2 * Alpha * R[t-1]*theta*sigma[t-1]) with
> filter function or is it even possible to code it this way.
>
> Any help would appreciated.Here is my slower code:
>
> "Ngarch" <- function(rtn,value){
> # Estimation of a non-symmertic GARCH, NGARCH(1,1), model.
> # Assume normal innovations
> # rtn: return series
> #
> # The likelihood function "mxlk" can be modified to fit more general
> NGARCH
> # models.
> # obtain initial estimates
> rtn=as.matrix(rtn)
> mu=mean(rtn)
> par=c(0.0001,0.8,0.01,0.7)
> #
> mxlk <- function(par){
> mxlk=0
> ht=var(rtn)
> T=length(rtn)
> if(T > 40)ht=var(rtn[1:40])
> at=rtn-mu
> for (i in 2:T){
>
> sig2t=par[1]+par[2]*ht[i-1]+par[3]*ht[i-1]*(at[i-1]/sqrt(ht[i-1])-par[4])^2
> ht[i]=sig2t
> mxlk=mxlk+0.5*(log(sig2t) + (at[i])^2/ht[i])
> }
> mxlk
> }
>
> low=c(1e-6,1e-6,1e-6,1e-6)
> upp=c(100*var(rtn),1-1e-6,1-1e-6,5)
> mm=optim(par,mxlk,method="Nelder-Mead",hessian=T)
> #mm=nlminb(par,mxlk,hessian=T,scale=1,lower=low,upper=upp)
> #mm=optimx(par,mxlk,method="L-BFGS-B",hessian=T,lower=low,upper=upp)
> ## Print the results
> par=mm$par
> H=mm$hessian
> Hi = solve(H)
> cat(" ","\n")
> cat("Estimation results of NGARCH(1,1) model:","\n")
> cat("estimates: ",par,"\n")
> se=sqrt(diag(Hi))
> cat("std.errors: ",se,"\n")
> tra=par/se
> cat("t-ratio: ",tra,"\n")
> # compute the volatility series and residuals
>
> ht=var(rtn)
> T=length(rtn)
> if(T > 40)ht=var(rtn[1:40])
> at=rtn-mu
> for (i in 2:T){
> sig2t=par[1]+par[2]*ht[i-1]+par[3]*(at[i-1]-par[4]*sqrt(ht[i-1]))^2
> ht=c(ht,sig2t)
> }
> sigma.t=sqrt(ht)
> Ngarch <- as.matrix(cbind(as.numeric(at),as.numeric(sigma.t)))
> colnames(Ngarch)=c("residuals","volatility")
> Ngarch
> }
>
>
>
>
> I think that Alexios has this function in his package (called NAGARCH) but
> I need to code this myself for my master work.
> I already implemented (read ported) his ugarch function in mathematica
> wolfram and it's working great.
>
> [[alternative HTML version deleted]]
>
> _______________________________________________
> R-SIG-Finance at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-finance
> -- Subscriber-posting only. If you want to post, subscribe first.
> -- Also note that this is not the r-help list where general R questions should go.
>
--
Patrick Burns
patrick at burns-stat.com
http://www.burns-stat.com
http://www.portfolioprobe.com/blog
twitter: @burnsstat @portfolioprobe
More information about the R-SIG-Finance
mailing list