Hi R users

I'm trying to replicate what the fGarch model does when any of the its
submodels incorporated are specified. I tried to look for help on the web
but not succeeded.

Up to now, i am able to get almost all but i got problems to replicate the
result of the aparch submodel.

I'm following the specification of the model set out in
"Introduction_to_the_rugarch_package.pdf"  (equation 26, to be more
precise) . It seems like i am not correctly setting the constraints of the
parameters needed. Based on the document, to get the aparch submodel from
the general equation of the fGarch model, the following constraints have to
be set:

- delta = lambda
- eta2 = 0
- abs(eta1) <= 1

It seems like i'm missing someone of them.

This is what i did.

 library(rugarch)

 data(sp500ret)

 seed=1

 ##################################################
 # SIMULATION WITH THE RUGARCH BUILT-IN FUNCTIONS

 spec=ugarchspec(variance.model = list(model = "fGARCH", submodel =
"APARCH",
                 garchOrder = c(1,1)),
                 mean.model = list(armaOrder = c(0,0), include.mean = TRUE),
                 distribution.model = "norm")

 model=ugarchfit(spec,as.xts(sp500ret),solver = "hybrid")

 sim.ruGarch = ugarchsim(model, n.sim = 100,n.start = 0, m.sim = 1,
                         startMethod = 'sample',rseed = seed)


##################################################
 ###  MY ATTEMPT TO MANUALLY REPLICATE IT.


 # STORE PARAMETERS TO USE THEM LATER ON

 ParIndx=model@model$pidx
 Param=model@model$pars

 mu=Param[ParIndx["mu",1]:ParIndx["mu",2],1]
 ar=Param[ParIndx["ar",1]:ParIndx["ar",2],1]
 orden.ar=length(ar)
 ma=Param[ParIndx["ma",1]:ParIndx["ma",2],1]
 orden.ma=length(ma)
 arfima=Param[ParIndx["arfima",1]:ParIndx["arfima",2],1]
 archm=Param[ParIndx["archm",1]:ParIndx["archm",2],1]
 mxreg=Param[ParIndx["mxreg",1]:ParIndx["mxreg",2],1]
 omega=Param[ParIndx["omega",1]:ParIndx["omega",2],1]
 alpha=Param[ParIndx["alpha",1]:ParIndx["alpha",2],1]
 orden.alpha=length(alpha)
 beta=Param[ParIndx["beta",1]:ParIndx["beta",2],1]
 orden.beta=length(beta)
 gamma=Param[ParIndx["gamma",1]:ParIndx["gamma",2],1]
 eta1=Param[ParIndx["eta1",1]:ParIndx["eta1",2],1]
 eta2=Param[ParIndx["eta2",1]:ParIndx["eta2",2],1]
 delta=Param[ParIndx["delta",1]:ParIndx["delta",2],1]
 lambda=Param[ParIndx["lambda",1]:ParIndx["lambda",2],1]
 vxreg=Param[ParIndx["vxreg",1]:ParIndx["vxreg",2],1]
 skew=Param[ParIndx["skew",1]:ParIndx["skew",2],1]
 shape=Param[ParIndx["shape",1]:ParIndx["shape",2],1]
 ghlambda=Param[ParIndx["ghlambda",1]:ParIndx["ghlambda",2],1]

 # SET INITIAL VALUES.

 h=array(0,c(101))
 epsilon=array(0,c(101))
 z=array(0,c(101))
 y=array(0,c(101))

 y[1]=tail(model@model$modeldata$data,1)
 h[1]=tail(model@fit$sigma,1)
 epsilon[1]=tail(model@fit$residuals,1)
 z[1]=epsilon[1]/h[1]


 set.seed(seed)

 z[2:101]=rdist(distribution = "norm", n=100, mu = 0, sigma = 1,
                        lambda = ghlambda, skew = skew,shape = shape)

 for(i in 2:101)
 {

    h[i] = omega+

 sum(alpha*h[i-(1:orden.alpha)]^lambda*((abs(z[i-(1:orden.alpha)]-eta2)-eta1*(z[i-(1:orden.alpha)]-eta2))^delta))+
           sum(beta*h[i-(1:orden.beta)]^lambda)

    h[i] = h[i]^(1/lambda)

    epsilon[i] = h[i]*z[i]

    y[i] = mu+sum(ar*(y[i-(1:orden.ar)]-mu))+sum(ma*epsilon[i-(1:orden.ma
)])+epsilon[i]

 }

 h=h[-1]
 epsilon=epsilon[-1]
 y=y[-1]


After running this, this is what i get:

> head(data.frame(round(h,7),round(sigma(sim.ruGarch),7)))
    round.h..7. round.sigma.sim.ruGarch...7.
T+1   0.0247601                    0.0260287
T+2   0.0247232                    0.0262280
T+3   0.0246867                    0.0246775
T+4   0.0246505                    0.0255587
T+5   0.0246145                    0.0243801
T+6   0.0245789                    0.0229765


Could anyone help me to figure out what i am doing wrong?

Any help would be strongly appreciated.

Cheers!!
-- 
*Jaimie.*

	[[alternative HTML version deleted]]

