Beta

[Martin Frankel]

Hi.

Does anyone know of a good, fast algorithm for sampling from a beta distribution? Something that can generate, in VB, several hundred thousand samples per minute?

Thanks.

[Actuarybert]

If you're using Excel Visual Basic, you could probably use something like this:

For i=1 to 10,000
Randomize
p = Rnd(1)
Observation=Application.WorksheetFunction.BetaInv( p,alpha,beta)
Next i

I'm not sure how fast this would be, but I think it would get the job done. Hope this helps.

[Martin Frankel]

Thanks bert.

Unfortunately, I'm working in VB, not VBA. I'd prefer to leave Excel out of the picture. You're right, that algorithm will work (although the BetaInv function crashes about once per five hundred attempts, so additional error-handling code is required) but it's very slow.

I appreciate the reply.

[Mopus]

I've never used VB so I don't know whether you can call external functions written in C, but if you can, you might want to try getting your hands on the GNU Scientific Library. This is basically a toolbox of functions which include a whole whack of stuff useful for simulation. I don't know whether it'll give you the performance that you want but it might be worth a shot.

Added - GSL uses Knuth's method, which relies on generating 2 gamma RV's. Other packages like NAG use algorithms that depend on the parameters of the Beta. The so called Chen BB and BC algorithms appear common. You could find some really ugly source code (lots of goto's) in the package randlib 1.3


<font size=-1>[ This Message was edited by: samiam on 2002-03-11 17:25 ]</font>

[General Kenobi (ret.)]

How accurate do you have to be? Especially if alpha and beta are fixed, you can fill a 100 (or 1000, or whatever) element array with the points of the inverse distribution, and then do a binary search on it. Not real accurate, but if you code your binary search well, it oughta go pretty quickly.

[Cho Da]

Numerical Recipes (http://www.nr.com) has an algorithm for the incomplete beta function, which should be the CDF&hellip;
It relies on a continued fraction approximation which converges fairly rapidly for <blockquote><large>x &lt; (&alpha; + 1) / (&alpha; + &beta; + 2)</large></blockquote>

[Martin Frankel]

Thanks for the replies.

sam, what is randlib 1.3 and how can I check it out?

[Mopus]

Randlib is a library of random number generating routines. You can find a description at
http://odin.mdacc.tmc.edu/anonftp/page_2.html

I've never used it so can't vouch for it's correctness. Again, the source is very ugly. Nag's selection of algorithms can be found in
http://www.nag.com/numeric/cl/manual/pdf/G05/g05fec_cl05.pdf

[Packet_Storm]

another good thing for random fuction is to see with the clock time if you don't you run the chance of repeating your random number every time you call the function, which is not good.

Also, for good random generators, look under encryption. Encryption programs normally have the best random generators, also check out SourceCode.com or the GNU lib mentioned before. Check out this page for more info on the random number function, because I think you will need to reseed the function with that amount of number needed, for you not to repeat the pattern.
http://www.vbexplorer.com/random/random_numbers_1.asp

If you actualy need random number, use this function to call certain cells out of a large set of data, and being an actuary that should not be hard to come by. Example if you get 51 and 12 then get the first digit of the number in row 12 col 51; however, but this will ensure more "random" numbers.

If you need any more, just ask. I could send you one in C++ that I used to test encryption and network preformance.

[Martin Frankel]

samiam,

Knuth's double-gamma method generates, on my PC, over 2.3 million samples per minute. So it's roughly 6 times as fast as the Excel BetaInv function. And it converges every time!

Thanks for your assistance.