Speedup a code using apply_along_axis

Hi, I'm sure I reinventing the wheel with the following code: from numpy import * from scipy import polyfit,stats def f(x,y,z): return x+y+z M=fromfunction(f,(2000,2000,10)) def foo(M): ramp=where(M<1000)[0] l=len(ramp) t=arange(l) if(l>1): return polyfit(t,ramp,1)[0] else: return 0 print apply_along_axis(foo,2,M) In real life M is not the result of one fromfunction call but it does not matter. The basic idea is to compute the slope (and only the slope) along one axis of 3D array. Only the values below a given threshold should be taken into account. The current code is ugly and slow. How to remove the len and the if statement? How to rewrite the code in a numpy oriented way? Xavier

On Sun, Feb 28, 2010 at 1:51 PM, Xavier Gnata <xavier.gnata@gmail.com> wrote:
Hi,
I'm sure I reinventing the wheel with the following code: from numpy import * from scipy import polyfit,stats
def f(x,y,z): return x+y+z M=fromfunction(f,(2000,2000,10))
def foo(M): ramp=where(M<1000)[0]
is this really what you want? I think this returns the indices not the values
l=len(ramp) t=arange(l) if(l>1): return polyfit(t,ramp,1)[0] else: return 0
print apply_along_axis(foo,2,M)
In real life M is not the result of one fromfunction call but it does not matter. The basic idea is to compute the slope (and only the slope) along one axis of 3D array. Only the values below a given threshold should be taken into account.
The current code is ugly and slow. How to remove the len and the if statement? How to rewrite the code in a numpy oriented way?
Getting the slope or the linear fit can be done completely vectorized see numpy-discussion threads last April with titles "polyfit on multiple data points" "polyfit performance" Josef
Xavier
_______________________________________________ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion

On 02/28/2010 08:17 PM, josef.pktd@gmail.com wrote:
On Sun, Feb 28, 2010 at 1:51 PM, Xavier Gnata <xavier.gnata@gmail.com> wrote:
Hi,
I'm sure I reinventing the wheel with the following code: from numpy import * from scipy import polyfit,stats
def f(x,y,z): return x+y+z M=fromfunction(f,(2000,2000,10))
def foo(M): ramp=where(M<1000)[0]
is this really what you want? I think this returns the indices not the values
Correct! It should be M[where(M<1000)]
l=len(ramp) t=arange(l) if(l>1): return polyfit(t,ramp,1)[0] else: return 0
print apply_along_axis(foo,2,M)
In real life M is not the result of one fromfunction call but it does not matter. The basic idea is to compute the slope (and only the slope) along one axis of 3D array. Only the values below a given threshold should be taken into account.
The current code is ugly and slow. How to remove the len and the if statement? How to rewrite the code in a numpy oriented way?
Getting the slope or the linear fit can be done completely vectorized see numpy-discussion threads last April with titles "polyfit on multiple data points" "polyfit performance"
Josef
Ok but the problem is that I also want to apply a threshold. In some cases, I end up less than 2 values below the threshold: There is nothing to fit and it should return 0. Hum....sounds like masked arrays could help...but I'm not familiar with masked arrays... Xavier
participants (2)
-
josef.pktd@gmail.com
-
Xavier Gnata