[SciPy-Dev] scipy.stats: algorithm to for ticket 1493

[josef.pktd at gmail.com]

On Mon, May 14, 2012 at 2:45 PM, nicky van foreest <vanforeest at gmail.com> wrote:
>>> Nice example. The answer is negative, while it should be positive, but
>>> the answer is within numerical accuracy I would say.
>>
>> oops, didn't we have a case with negative sign already ?
>> maybe a check self.a <= p <= self.b  ?
>
> I included this. I also think that a check on whether left and right
> stay within  self.a and self.b should be included, perhaps just for
> safety reasons.
>
>>
>>>
>>>> I don't see anything yet to criticize in your latest version :(
>>>
>>> Ok. I just checked the tests in scipy/stats/tests.
>>
>> If you are curious, you could temporarily go closer to q=0 and q=1 in
>> the tests for ppf, and see whether it breaks for any distribution.
>
> Good idea. Just to see what would happen I changed the following code
> in test_continuous_basic.py:
>
> @_silence_fp_errors
> def check_cdf_ppf(distfn,arg,msg):
>    values = [-1.e-5, 0.,0.001,0.5,0.999,1.]
>    npt.assert_almost_equal(distfn.cdf(distfn.ppf(values, *arg), *arg),
>                            values, decimal=DECIMAL, err_msg= msg + \
>                            ' - cdf-ppf roundtrip')

roundtrip: looks like ppf should be ok, but cdf is not

>>> stats.norm.ppf(-1e-5)
nan
>>> stats.norm.cdf(np.nan)
0.0
>>> stats.norm.cdf(stats.norm.ppf(-1e-5))
0.0

I'm using scipy 0.9. but I don't think this has changed, not that I know of

I'm trying to track down when this got changed.
(github doesn't show changes in a file that has too many changes, need
to dig out git)

>
>
> Thus, I changed the values into an array. It should fail on the first
> value, as it is negative, but I get a pass. Specifically, I ran:
>
> nicky at chuck:~/prog/scipy/scipy/stats/tests$ python test_continuous_basic.py
> ..............................................................................................................................
> ----------------------------------------------------------------------
> Ran 126 tests in 93.990s
>
> OK
>
>>
>
> Weird result. If I add a q  = 1.0000001 I get a fail on the fourth
> test, as expected.
>
>>> - repair for the cases q =  0 and q = 1 by means of an explicit test.
>>
>> isn't ppf (generic part) taking care of this, if not then it should, I think
>
> Actually, from the code in lines:
>
> https://github.com/scipy/scipy/blob/master/scipy/stats/distributions.py#L1529
>
> I am inclined to believe you. However, in view of the above test ...
> Might it be that the conditions on L1529 have been added quite
> recently, and did not yet make it to my machine? I'll check this right
> now....As a matter of fact, my distributions.py contains the same
> check, i.e.,         cond1 = (q > 0) & (q < 1) . Hmmm.
>
> Now I admit that I do not understand in all nitty-gritty detail the
> entire implementation of ppf(), but I suspect that this is a bug.
>
>>
>> ppf(0) = self.a
>> ppf(1) = self.b
>
> Good idea.

this already looks correct in the generic ppf code

>>> stats.beta.ppf(0, 0.5)
0.0
>>> stats.beta.a
0.0

Josef
>
> I'll implement the code in my branch, and do a pull request.
> _______________________________________________
> SciPy-Dev mailing list
> SciPy-Dev at scipy.org
> http://mail.scipy.org/mailman/listinfo/scipy-dev