[SciPy-dev] nbinom.ppf

[Pierre GM]

On Jul 27, 2009, at 4:58 PM, Robert Kern wrote:

> On Mon, Jul 27, 2009 at 15:29, <josef.pktd at gmail.com> wrote:
>> On Mon, Jul 27, 2009 at 3:55 PM, Robert Kern<robert.kern at gmail.com>  
>> wrote:
>>>
>>> ceil() is correct; round() is not. round() would be okay if the only
>>> inputs are expected to be outputs of the CDF, but one frequently  
>>> needs
>>> the PPF to take all values in [0,1], like for example doing random
>>> number generation via inversion.

Fair enough, but a bit frustrating nevertheless in this particular case.

>> Do you think it would be useful to round if we are epsilon (?) close
>> to the next integer? It is more likely that users have a case like
>> Pierre's where the answer might be the closest integer, instead of an
>> epsilon below. It would move the floating point error, but to an, in
>> actual usage, less likely location.



>> I think, it is in scipy.special where I saw some code that treats
>> anything that is within 1e-8 (?) of an integer as the integer.
>
> I think the code after the first line is an attempt to remedy this  
> situation:
>
>    def _ppf(self, q, n, pr):
>        vals = ceil(special.nbdtrik(q,n,pr))
>        vals1 = vals-1
>        temp = special.nbdtr(vals1,n,pr)
>        return where(temp >= q, vals1, vals)
>
> It is possible that "temp >= q" should be replaced with "temp >= q- 
> eps".

Doesn't work in the case I was presenting, as temp is here an array of  
NaNs. Using
 >>> vals = ceil(special.nbdtrik(q,n,pr)-eps)
seems to work well enough, provided that eps is around 1e-8 as Josef  
suggested.