[Numpy-discussion] custom accumlators
Matt Knox
mattknox_ca at hotmail.com
Fri Jan 5 22:01:44 EST 2007
> >>
> > You might want to look at frompyfunc:
> >
> > def expmave2(x, k):
> > def expmave_sub(a, b):
> > return a + k * (b - a)
> > return np.frompyfunc(expmave_sub, 2, 1).accumulate(x)
> >
> >
> > It's amazing what you find when you dig around.
Thanks a lot everyone. This has been informative. For what it's worth, I did
some performance comparisons...
import numpy as np
import profile
def expmave1(x, k):
def expmave_sub(a, b):
return b + k * (a - b)
return np.frompyfunc(expmave_sub, 2, 1).accumulate(x).astype(x.dtype)
def expmave2(x, k):
result = np.array(x, copy=True)
for i in range(1, result.size):
result[i] = result[i-1] + k * (result[i] - result[i-1])
return result
testArray = np.cumprod(1 + np.random.normal(size=10000)/100)
profile.run('expmave1(testArray, 0.2)')
profile.run('expmave2(testArray, 0.2)')
and the second function is faster, which I guess makes sense if frompyfunc is
pure python, although the first one does have a nice elegance to it I think.
And actually, the speed of expmave2 is quite sufficient for my use cases when
working with standard arrays... it's faster than I thought because initially I
was doing this with masked arrays directly and the __getitem__ and __setitem__
calls are very expensive for masked arrays (expensive to the point where they
can't really be used in looping functions like this), but I can get around that
by just filling the array first and then computing the appropriate resulting
mask separately. Obviously it would be faster in pure C, but the effort isn't
worth it for my needs.
> That said, it might also be worth looking at whether numexpr can do
> this sort of thing - it's supposed to take a numpy expression and
> evaluate it efficiently in C (avoiding intermediate arrays and so on).
numexpr does look interesting. If numexpr was adapted to be able to do
accumulations with expressions passed to it, that would be super slick... but
that is probably a few leaps beyond my current knowledge level to code that.
Thanks again everybody,
- Matt
More information about the NumPy-Discussion
mailing list