0

Well hello there, could anyone please tell me if there is a random number generator for Poisson distributed random variables implemented in QuantLib?If yes, where do I find the code for this?I'm trying to simulate a Jump-Diffusion process and need the number of jumps between steps in time ( i.e. for each time interval [t_(i-1);t_i[.Is there a way to do this in QuantLib directly or do I need to use the boost library? Thanks in advance!

p.s. Or would you recommend using the actual jump arrival times by generating exponentially distributed numbers instead?

sunshine
  • 1

2 Answers2

0

Whether or not you simulate the jump times or the jump density depends on how you're writing your diffusion loop. IMHO, simulating density is cleaner, because it requires less state to be carried forward through the simulation.

I don't know whether you'll find something already written in Boost or QuantLib. But sampling the Poisson distribution is actually very simple, if you already have a uniform RNG. For instance (pseudo-code):

p = exp(-lambda);
F = p;   % cumulative distribution function
N = 0;
U = rand();
while (U > F)
    N = N + 1;
    p = p*lambda/N;
    F = F + p;
end
return N;

This is based on inverse transform sampling. There are several other techniques out there.

Oliver Charlesworth
  • 267,707
  • 33
  • 569
  • 680
  • Thanks, I might need to use this after all although that was I was trying to avoid. I stumbled onto some previously asked questions here on Stack Overflow and read that there might be numerical issues with this approach. I am looking for sth like BigInteger seed =12324; MersenneTwisterUniformRng unifMt(seed); BoxMullerGaussianRng bmGauss(unifMt); but for a poisson distributed r.v. – sunshine Mar 16 '12 at 10:37
  • @sunshine: If there are numerical issues, I'm not sure what they are. If you find out, I'd be grateful to hear about it. – Oliver Charlesworth Mar 16 '12 at 10:43
  • funnily enough I can't seem to find the thread! If I come across it again I'll post it here. – sunshine Mar 16 '12 at 15:32
0

The closest you can currently get in QuantLib is the InverseCumulativeRng class template together with the InverseCumulativePoisson class; something like

MersenneTwisterUniformRng unifMt(seed);
InverseCumulativePoisson f(lambda);
InverseCumulativeRng<MersenneTwisterUniformRng, InverseCumulativePoisson> rng(unifMt, f);

will give you a Poisson generator. Be aware that it will return the samples as doubles, not ints: they will be whole numbers, but expressed in the wrong type.

Also, it looks like for some reason InverseCumulativeRng doesn't provide a constructor taking the function. Strange that we've overlooked that... Anyway, you'll have to edit <ql/math/randomnumbers/inversecumulativerng.hpp> and add it; when you're done, please send the patch to the QuantLib mailing list and I'll add it to the repository.

Luigi Ballabio
  • 4,128
  • 21
  • 29
  • hehe I'd love to, unfortunately I have absolutely no idea how to do that..I think I'll stick to the boost rng then.thanks! – sunshine Mar 16 '12 at 15:31
  • Just add another constructor that takes a second arguments of type IC and copies it into the ICND_ data member. – Luigi Ballabio Mar 16 '12 at 15:54