[R-SIG-Finance] efficient code for nonlinear garch model

[Patrick Burns]

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