There is no error.  !x is not the gamma function but factorial,
and (gamma x) = factorial x-1.  The formula calls for

gamma 0.5*x+1
! _1+0.5*x+1
! 0.5*_2+x+1
! 0.5*x-1
! -: x-1

The last expression was what I had in my previous posts.

Results greater than 1 are due to accumulated numerical 
errors (I believe).  The requirement that (%:x)>|y is to
ensure that the argument to H. is within the theoretical
radius of convergence.  Large values of y, even when
bounded by %:x, likely cause the series to converge
slowly, thus requiring more terms, thus accumulates
more error.



----- Original Message -----
From: Robert Cyr <[EMAIL PROTECTED]>
Date: Tuesday, August 12, 2008 17:40
Subject: Re: [Jprogramming] general Gamma distribution
To: Programming forum <[email protected]>

> There appears to be a sign problem in the formula.
> 
> The formula either form Wolfram or wikipedia are:
> 
> tcdf=: 4 : 0
>  assert. (%:x)>|y
>  0.5 + y * (!-:x+1) * ((0.5,-:1+x) H. 1.5 x%~-*:y) % 
> (%:o.x) * !<:-:x
> )
> 
> As noted by others, the limit also needs to be adjusted.
> The limits probably are too high, as tcdf cannot exceed 1.  
> It seems to be
> in that region that problem results occur.
> 
> On Tue, Aug 12, 2008 at 5:38 PM, Roger Hui <[EMAIL PROTECTED]> wrote:
> 
> > For large x one may want to use the N(0,1) CDF
> > in a pinch.  From
> > http://www.jsoftware.com/jwiki/Essays/Normal_CDF
> >
> > erf   =: (1 H. 1.5)@*: * 2p_0.5&* % ^@:*:
> > n01cdf=: -: @ >: @ erf @ %&(%:2)
> >
> >   10 tcdf 3
> > 0.993328
> >    n01cdf 3
> > 0.99865
> >
> > The bad  100 tcdf 9  result is likely due to extreme
> > accumulation of numerical errors.  On my machine I get:
> >
> >   100 tcdf 9
> > 1.03437e17
> >
> >   n01cdf 9
> > 1
> >   1 - n01cdf 9
> > _2.22045e_16
> >
> >
> >
> > ----- Original Message -----
> > From: John Randall <[EMAIL PROTECTED]>
> > Date: Tuesday, August 12, 2008 14:07
> > Subject: Re: [Jprogramming] general Gamma distribution
> > To: Programming forum <[email protected]>
> >
> > > Roger Hui wrote:
> > > > It's the absolute value of y that should be bounded,
> > > > not just y itself.  Thus:
> > > >
> > > > tcdf=: 4 : 0
> > > >  assert. (%:x)>|y
> > > >  0.5 + y * (!-:x-1) * ((0.5,-:1+x) H. 1.5 x%~-*:y) %
> > > (%:o.x) * !<:-:x
> > > > )
> > > >
> > > >
> > > Mathworld <http://mathworld.wolfram.com/Studentst-
> > > Distribution.html> has
> > > different formulations (equations 7 and  8), which rely on
> > > the incomplete
> > > beta function and are not subject to the hypergeometric bounds.
> > >
> > > I don't think the restriction is bad when x (degrees of freedom)
> > > is above
> > > about 10.  For example
> > >
> > >    10 tcdf 3
> > > 0.993328
> > >
> > > but it is not good if x is smaller.   High values 
> of y
> > > give strange results:
> > >
> > >    100 tcdf 9
> > > _7.15204e15
----------------------------------------------------------------------
For information about J forums see http://www.jsoftware.com/forums.htm

Reply via email to