On Fri, Aug 8, 2008 at 9:29 AM, Robert Cyr <[EMAIL PROTECTED]> wrote:
> I am afraid I do not understand the basis of
> your proposed variation, sorry.
poissonrand generates a sequence of numbers such as would
be generated from [EMAIL PROTECTED]@? bind 0, and counts how many such
numbers are needed to exceed the desired total (100, in your case).
This approach:
ps=. >.l % ([EMAIL PROTECTED]@[EMAIL PROTECTED] 0:) m
generates just one number and then divides the desired total by
that number to find how many times it would need to be added
to reach that total.
This does not give a poisson distribution but in the context of your
program (which sums many such numbers) the distribution of
the sums is very like the distribution of the sums of numbers
obtained from poisson distributions.
This approach:
prd=: 4 :0
n=. >:&.(p:inv) 1+y+2*x
r=. n [EMAIL PROTECTED]@[EMAIL PROTECTED] 0
>. _0.5+y{.x}. +//. (x,_1+n)$r
)
does seem to give a poisson distribution, though I would want
someone with a better analysis background than mine to confirm
that that is actually the case. It's an approach similar to my first
simplification but instead of taking just one random number I
use a sequence of random numbers that's based on the desired
total (and I add them rather than multiplying just the one number).
You would use it in your program like this:
ps=. l prd m
In both cases the time savings come from recycling the rather expensive
calculations required to generate the logarithms of random numbers.
FYI,
--
Raul
----------------------------------------------------------------------
For information about J forums see http://www.jsoftware.com/forums.htm
