Моделирование GARCH в R - PullRequest
Новая
       15

Моделирование GARCH в R

5 голосов
/ 02 апреля 2012

Я делаю имитацию модели GARCH.Сама модель не слишком актуальна, я хотел бы спросить вас об оптимизации симуляции в R. Больше всего, если вы видите место для векторизации, я думал об этом, но не вижу его.Итак, что у меня есть, это:

Позвольте: 

#    ht=cond.variance in t
#    zt= random number 
#    et = error term
#    ret= return
#    Horizon= n periods ahead

Итак, это код: 

randhelp= function(horizon=horizon){
    ret <- zt <- et <- rep(NA,horizon)#initialize ret and zt et
    for( j in 1:horizon){
      zt[j]= rnorm(1,0,1)
      et[j] = zt[j]*sqrt(ht[j])
      ret[j]=mu + et[j]

      ht[j+1]= omega+ alpha1*et[j]^2 + beta1*ht[j]
    }
    return(sum(ret))
  }

Я хочу сделать симуляцию возвращений 5периоды с этого момента, так что я буду запускать это, скажем, 10000. 

#initial values of the simulation
ndraws=10000
horizon=5 #5 periods ahead
ht=rep(NA,horizon) #initialize ht
ht[1] = 0.0002
alpha1=0.027
beta1 =0.963
mu=0.001
omega=0


sumret=sapply(1:ndraws,function(x) randhelp(horizon))

Я думаю, что это работает достаточно быстро, но я хотел бы спросить вас, есть ли какой-нибудь способ подойти к этой проблеме лучше.

Спасибо!

Ответы [ 2 ]

5 голосов
/ 02 апреля 2012

основываясь на ответе Винсента, все, что я изменил, - это определение zt сразу и переключение apply(ret, 1, sum) на rowSums(ret), и оно немного ускорилось.Я попробовал оба скомпилированных, но без больших различий: 

randhelp2 <- function(horizon = 5, N = 1e4, h0 = 2e-4, 
                       mu = 0, omega = 0, alpha1 = 0.027, 
                       beta1  = 0.963 ){
    ret <- et <- ht <- matrix(NA, nc = horizon, nr = N)
    zt  <- matrix(rnorm(N * horizon, 0, 1), nc = horizon)
    ht[, 1] <- h0
    for(j in 1:horizon){
        et[, j]  <- zt[, j] * sqrt(ht[, j])
        ret[,j]  <- mu + et[, j]
        if( j < horizon )
            ht[, j + 1] <- omega + alpha1 * et[, j] ^ 2 + beta1 * ht[, j]
    }
    rowSums(ret)
}

system.time(replicate(10,randhelp(N=1e5)))
   user  system elapsed 
  7.413   0.044   7.468 

system.time(replicate(10,randhelp2(N=1e5)))   
   user  system elapsed 
  2.096   0.012   2.112

, вероятно, все еще есть возможности для улучшения: -)

4 голосов
/ 02 апреля 2012

Вместо использования чисел в цикле, вы можете использовать векторы размера N: удаляет петлю, скрытую в sapply. 

randhelp <- function(
  horizon=5, N=1e4, 
  h0 = 2e-4, 
  mu = 0, omega=0,
  alpha1 = 0.027,
  beta1  = 0.963
){
  ret <- zt <- et <- ht <- matrix(NA, nc=horizon, nr=N)
  ht[,1] <- h0
  for(j in 1:horizon){
    zt[,j]   <- rnorm(N,0,1)
    et[,j]   <- zt[,j]*sqrt(ht[,j])
    ret[,j]  <- mu + et[,j]
    if( j < horizon )
      ht[,j+1] <- omega+ alpha1*et[,j]^2 + beta1*ht[,j]
  }
  apply(ret, 1, sum)
}
x <- randhelp(N=1e5)
...