[R] Simulating inhomogeneous Poisson process without loop

Ben Bolker bbolker at gmail.com
Tue Jul 5 13:53:45 CEST 2011


Tristan Linke <tristan.linke <at> gmail.com> writes:

> 
> Dear all
> 
> I want to simulate a stochastic jump variance process where N is Bernoulli
> with intensity lambda0 + lambda1*Vt. lambda0 is constant and lambda1 can be
> interpreted as a regression coefficient on the current variance level Vt. J
> is a scaling factor
> 
> How can I rewrite this avoiding the loop structure which is very
> time-consuming for long simulations?
> 
> for (i in 1:N){
> ...
> N <- rbinom(n=1, size=1, prob=(lambda0+lambda1*Vt))
> Vt <- ... + J*N
> ..
> }

  Is it a typo that you are using N both as your loop variable
and as part of the state of your simulation?

> P.S. This is going towards the Duffie, Pan, Singleton 2000 Transform Pricing
> paper, here stochastic volatility with state-dependent correlated jumps
> (Eraker 2004).

  I don't think there's any way to rearrange this completely without loops,
because each call of rbinom() depends on the previous updating.

  One thing that might speed things up a lot would be to pick
N random uniform variates as a single vector in advance:

  rvec <- runif(N)
  for (i in 1:N) {
    Vt <- ... + J*(rvec[i]<lambda0+lambda1*Vt)
  }   

  Depending on how badly you need it to be fast/how much effort you're
willing to invest, you could try compiling/using inline or Rcpp etc.

  The other point is that if the probability stays constant over several
steps when the jump doesn't happen, and if the probability is low,
you could skip a few loop steps by drawing a geometric RV to determine
the time of the next jump, rather than testing one step at a time.



More information about the R-help mailing list