[SciPy-user] Strange discontinuity in noncentral chisquare
Robert Kern
robert.kern@gmail....
Thu May 28 14:03:25 CDT 2009
On Thu, May 28, 2009 at 13:53, <josef.pktd@gmail.com> wrote:
> On Thu, May 28, 2009 at 1:56 PM, Neal Becker <ndbecker2@gmail.com> wrote:
>> def pmiss2 (x, esnodB, N):
>> esno = 10**(0.1 * esnodB) * N
>> var = 1/esno
>> _lambda = 1/(0.5*var)
>>
>> return ncx2.cdf (x, 2, _lambda)
>>
>> x = np.arange (0, 50, 0.1)
>> p1 = [pmiss2 (e, 3.5, 24) for e in x]
>>
>> What's with this strange discontinuity?
>> print p1:
>> ...
>> 3.475382846574262e-21,
>> 4.2226227741362447e-21,
>> 5.1248671653198949e-21,
>> 6.2130949241675783e-21,
>> 7.5242411687161146e-21,
>> 9.1022970542215721e-21,
>> 5.8787514615651664e-09,
>> 6.2565279721932619e-09,
>> 6.656924144742753e-09,
>> 7.0811937641411923e-09,
>> 7.5306544171300622e-09,
>> 8.0066904400027596e-09,
>> 8.5107559880675483e-09,
>> 9.044378231221098e-09,
>> 9.6091606801490766e-09,
>> ...
>>
>
> This must be a bug in scipy.special.chndtr
I notice the following snippets of code, which appear guilty:
C .. Statement Functions ..
LOGICAL qsmall
C ..
C .. Statement Function definitions ..
qsmall(xx) = sum .LT. 1.0D-20 .OR. xx .LT. eps*sum
# That is a feature of Fortran I knew nothing about.
60 IF (qsmall(term)) GO TO 80
80 cum = sum
RETURN
--
Robert Kern
"I have come to believe that the whole world is an enigma, a harmless
enigma that is made terrible by our own mad attempt to interpret it as
though it had an underlying truth."
-- Umberto Eco
More information about the SciPy-user
mailing list