[SciPy-User] scipy.stats.nanmedian

[josef.pktd at gmail.com]

On Thu, Jan 21, 2010 at 10:01 PM, Keith Goodman <kwgoodman at gmail.com> wrote:
> On Thu, Jan 21, 2010 at 6:41 PM, Pierre GM <pgmdevlist at gmail.com> wrote:
>> On Jan 21, 2010, at 9:28 PM, Keith Goodman wrote:
>>> That's the only was I was able to figure out how to pull 1.0 out of
>>> np.array(1.0). Is there a better way?
>>
>>
>> .item()
>
> Thanks. item() looks better than tolist().
>
> I simplified the function:
>
> def nanmedian(x, axis=0):
>    x, axis = _chk_asarray(x,axis)
>    if x.ndim == 0:
>        return float(x.item())
>    x = x.copy()
>    x = np.apply_along_axis(_nanmedian,axis,x)
>    if x.ndim == 0:
>        x = float(x.item())
>    return x
>
> and opened a ticket:
>
> http://projects.scipy.org/scipy/ticket/1098


How about getting rid of apply_along_axis?    see attachment

I don't know whether or how much faster it is, but there is a ticket
that the current version is slow.
No hidden bug or corner case guarantee yet.


Josef
-------------- next part --------------
# -*- coding: utf-8 -*-
"""
Created on Wed Jan 20 10:18:32 2010
Author: josef-pktd
"""

import numpy as np
from scipy import stats


def nanmedian(x, axis = 0):
    x, axis = stats.stats._chk_asarray(x, axis)
    if x.ndim == 0:
       return 1.0*x[()]
    x = np.sort(x, axis=axis)
    nall = x.shape[axis]
    notnancount = nall - np.isnan(x).sum(axis=axis)
    
    (idx, rmd) = divmod(notnancount, 2)
    indx = [np.arange(x.shape[ax]) for ax in range(x.ndim)]
    indxlo = indx[:]
    indxlo[axis] = idx
    indxhi = indx[:]
    indxhi[axis] = idx - (1-rmd)
    nanmed = (x[indxlo] + x[indxhi])/2.
    if nanmed.ndim == 0:
       return nanmed[()]
    return nanmed
    



for axis in [0,1]:
    for i in range(5):
        # for complex
        #x = 1j+np.arange(20).reshape(4,5)
        x = np.arange(20).reshape(4,5).astype(float)
        x[zip(np.random.randint(4, size=(2,5)))] = np.nan
        print nanmedian(x, axis=0)
        print stats.nanmedian(x, axis=0)
    
print nanmedian(1)
print nanmedian(np.array(1))
print nanmedian(np.array([1]))