On Fri, Aug 8, 2008 at 2:55 PM, Oleg Kobchenko <[EMAIL PROTECTED]> wrote:
> dstat poissonrand 100 10000
Thanks, I was needing something like that.
Really, really needing it -- my last post included so
many bugs that the routine I posted was utterly broken.
Here's a fixed version:
pssr0=:4 :0"0
m=. >.0>. 1.3*x
n=: [EMAIL PROTECTED]@[EMAIL PROTECTED]&0 >:&.(p:inv) m (%:@<. * >.) y
'r t'=. ((,: {:"[EMAIL PROTECTED])~ >:&x i."1 1:) +/\"1 (y,m)$n
i=. I. x > t
t=. i { t
while. # i do.
m=. >.0>. 1.7 * (+/%#) v=. x - t
n=. [EMAIL PROTECTED]@[EMAIL PROTECTED]&0 >:&.(p:inv) m (%:@<. * >.) y
'r2 t2'=. ((,: {:"1)~ >:&v i."1 1:) +/\"1 ((#i),m)$n
b=. x > t2+t
r=. (r2+i{r) i} r
i=. b#i
t=. b#t+t2
end.
r
)
pssr=:4 :0"0
blksiz=. >. y <. 1e5 % x
r=. $0
while. y > #r do.
r=. r, x pssrnd0 blksiz <. y-#r
end.
)
ts '100 pssr 1e4'
0.403073 2.89376e6
dstat 100 pssr 1e4
sample size: 10000
minimum: 66
maximum: 142
median: 99
mean: 99.7684
std devn: 9.86424
skewness: 0.102645
kurtosis: 2.93427
That's not as good as I was hoping, but rather than spend any more
time tuning this routine I thought I should post to keep people from having
to waste too much time on my previous errors.
Thanks,
--
Raul
----------------------------------------------------------------------
For information about J forums see http://www.jsoftware.com/forums.htm