Here is an alternative approach: numerically integrating the pdf. It is not subject to the bound in the beta function, and gives sensible values for large x and y. The downside is figuring how to get the required accuracy. For the normal distribution, you can get 4-digit accuracy by using h=0.05: I have used h=0.01.
gamma=:!@:<: beta =:[EMAIL PROTECTED] * [EMAIL PROTECTED] % [EMAIL PROTECTED] NB. Probability density function tpdf=:4 : '((1+ x %~ *:y)^- -: >: x) % (0.5 beta -:x)*%: x' ppr=: +//.@(*/) NB. polynomial multiplication dp=:+/ .* NB. dot product NB. Simson's rule integrator int=:1 : 0 n=.2>.+:>.y%0.01 h=.y%n c=.1 4 1 ppr 0=2 | i. <: n (h%3)*c mp u h*i.>:n ) NB. Find cdf by numerically integrating pdf tcdf2=:4 : 0 if. y<0 do. 1-x tcdf -y return. end. if. y=0 do. 0.5 return. end. 0.5+x&tpdf int y ) 5 tcdf2 1.47588 0.899999 5 tcdf2 2.01505 0.95 100 tcdf2 9 1 Best wishes, John ---------------------------------------------------------------------- For information about J forums see http://www.jsoftware.com/forums.htm
