您可以使用在numpy.lib.stride_tricks 中找到的as_strided() 来操纵数组的步幅以更快地迭代它。
您的计算可以被视为在阵列的 (1,1,3) 窗口上运行。我喜欢使用一个通用函数(sliding_window() 使用as_strided() 在 n 个窗口中创建 n。我在这里找到了它 - Efficient Overlapping Windows with Numpy;该函数的功劳显然归功于 johnvinyard。那个博客页面很好地描述了正在发生的事情。
制作一些 1x1x3 的窗口
import numpy as np
x = np.array([[[1,2,3],[6,np.nan,4]],[[0.5,2,1],[9,3,np.nan]]])
for thing in sliding_window(x, (1,1,3)):
print thing
# [ 1. 2. 3.]
# [ 6. nan 4.]
# [ 0.5 2. 1. ]
# [ 9. 3. nan]
应用 ```np.percentile()'' - 忽略 NaN
for thing in sliding_window(x, (1,1,3)):
print np.percentile(thing[np.isfinite(thing)], 40)
# 1.8
# 4.8
# 0.9
# 5.4
制作一个结果数组:
per_s = [np.percentile(thing[np.isfinite(thing)], 40)
for thing in sliding_window(x, (1,1,3))]
print per_s
# [1.8, 4.8000000000000007, 0.90000000000000002, 5.4000000000000004]
per_s = np.array(per_s)
print per_s
# array([ 1.8, 4.8, 0.9, 5.4])
让它恢复到您期望的形状
print per_s.reshape((2,2))
# array([[ 1.8, 4.8],
# [ 0.9, 5.4]])
print per_s.reshape(x.shape[:-1])
# array([[ 1.8, 4.8],
# [ 0.9, 5.4]])
这应该更快。我很好奇它是否会 - 我没有任何现实世界问题来测试它。
对 numpy as_strided 的谷歌搜索出现了一些不错的结果:我有这个书签,http://scipy-lectures.github.io/advanced/advanced_numpy/
sliding_window() 来自Efficient Overlapping Windows with Numpy
from numpy.lib.stride_tricks import as_strided as ast
from itertools import product
def norm_shape(shape):
'''
Normalize numpy array shapes so they're always expressed as a tuple,
even for one-dimensional shapes.
Parameters
shape - an int, or a tuple of ints
Returns
a shape tuple
'''
try:
i = int(shape)
return (i,)
except TypeError:
# shape was not a number
pass
try:
t = tuple(shape)
return t
except TypeError:
# shape was not iterable
pass
raise TypeError('shape must be an int, or a tuple of ints')
def sliding_window(a,ws,ss = None,flatten = True):
'''
Return a sliding window over a in any number of dimensions
Parameters:
a - an n-dimensional numpy array
ws - an int (a is 1D) or tuple (a is 2D or greater) representing the size
of each dimension of the window
ss - an int (a is 1D) or tuple (a is 2D or greater) representing the
amount to slide the window in each dimension. If not specified, it
defaults to ws.
flatten - if True, all slices are flattened, otherwise, there is an
extra dimension for each dimension of the input.
Returns
an array containing each n-dimensional window from a
'''
if None is ss:
# ss was not provided. the windows will not overlap in any direction.
ss = ws
ws = norm_shape(ws)
ss = norm_shape(ss)
# convert ws, ss, and a.shape to numpy arrays so that we can do math in every
# dimension at once.
ws = np.array(ws)
ss = np.array(ss)
shape = np.array(a.shape)
# ensure that ws, ss, and a.shape all have the same number of dimensions
ls = [len(shape),len(ws),len(ss)]
if 1 != len(set(ls)):
raise ValueError(\
'a.shape, ws and ss must all have the same length. They were %s' % str(ls))
# ensure that ws is smaller than a in every dimension
if np.any(ws > shape):
raise ValueError('ws cannot be larger than a in any dimension. a.shape was %s and ws was %s' % (str(a.shape),str(ws)))
# how many slices will there be in each dimension?
newshape = norm_shape(((shape - ws) // ss) + 1)
# the shape of the strided array will be the number of slices in each dimension
# plus the shape of the window (tuple addition)
newshape += norm_shape(ws)
# the strides tuple will be the array's strides multiplied by step size, plus
# the array's strides (tuple addition)
newstrides = norm_shape(np.array(a.strides) * ss) + a.strides
strided = ast(a,shape = newshape,strides = newstrides)
if not flatten:
return strided
# Collapse strided so that it has one more dimension than the window. I.e.,
# the new array is a flat list of slices.
meat = len(ws) if ws.shape else 0
firstdim = (np.product(newshape[:-meat]),) if ws.shape else ()
dim = firstdim + (newshape[-meat:])
# remove any dimensions with size 1
#dim = filter(lambda i : i != 1,dim)
dim = tuple(thing for thing in dim if thing != 1)
return strided.reshape(dim)