In the script 'statdist.ijs' the verb normalrand is
provided for generating standard normal variates and afaik
is a widely used method based on Box and Muller [1] . Their
original method uses n independent uniform and n independent
exponential variates to generate **2n** (NB. just **n** in
the next paragraph) standard normal variates. Half of the
uniform variates are used as the argument for sine and half
are used as the argument for a cosine, and standard normal
variates are produced in such pairs.
In contrast, the verb normalrand uses n independent
uniform and n independent exponential random variates to
generate **n** standard normal variates; all uniform
variates are used as the argument for cosine. The normalrand
revision should work correctly.
In the interest of finding a more efficient
generator of standard normal variates, I have tried the
following two which follow more closely the original Box
and Muller approach.
normalrand2=: 3 : 0
n =. >. -: y
r1 =. (2 o. +: o. b =.rand01 n) * %: - +: ^. a =. rand01 n
r2 =. (1 o. +: o. b ) * %: - +: ^. a
NB. y{.,r1,.r2 NB. unnecessary imo
y{.r1,r2
)
normalrand3=: 3 : 0
n =. >. -: y
y{.,((2&o.,.1&o.) +: o. rand01 n) *"1 0 %: - +: ^. rand01 n
)
Using ts on my computer I have found only slight
improvements with normalrand2 and normalrand3.
normalrand2 takes about the same time as normalrand
and uses about half the space.
normalrand3 takes about 80 per cent of the time and
space of normalrand.
In all of these verbs the exponential variate is
generated by the following phrase.
- ^. rand01 n (or y in normalrand)
[1] Box, GEP and GM Jenkins, "A Note on the Generation of
Random Normal Deviates," Annals of Mathematical Statistics,
Vol 29, 1958, pp 610-611.
----------------------------------------------------------------------
For information about J forums see http://www.jsoftware.com/forums.htm
