我认为指出我的问题最简单,一般情况很难解释.

说我有一个矩阵

a with dimensions NxMxT,

人们可以将T视为时间维度(使问题更容易).设(n,m)为通过NxM的指数.我可以调用(n,m)状态空间标识符.然后我需要找到python / scipy等价的

for each (n,m):

find a*(n,m) = min(a(n,m,:) s.t. a*(n,m) > a(n,m,T)

也就是说,找到仍然高于最后一个(在时间维度中)观察的最小状态空间值 – 对于整个状态空间.

我的第一个尝试是首先解决内部问题(找到一个高于[…, – 1]的)：

aHigherThanLast = a[ a > a[...,-1][...,newaxis] ]

然后我想在每个(n,m)中找到所有这些中最小的.不幸的是,aHigherThanLast现在包含所有这些值的一维数组,因此我不再具有(n,m)对应关系.什么是更好的方法呢？

作为一个额外的问题：状态空间是可变的,它也可能是3维或更多维度(NxMxKx …),我不能硬编码.任何一种

for (n,m,t) in nditer(a):

不可行.

非常感谢！

/编辑：

a = array([[[[[[[[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.]]]],

[[[[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.]]]]],

[[[[[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.]]]],

[[[[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.],

[ 0., 2., 1.]]]]]]]])

# a.shape = (1L, 1L, 2L, 2L, 1L, 1L, 10L, 3L). so in this case, T = 3.

# expected output would be the sort of

# b.shape = (1L, 1L, 2L, 2L, 1L, 1L, 10L), which solves

> b [a,b,c,d,e,f,g]> a [a,b,c,d,e,f,g,-1](b高于最新观察值)

> i中没有任何元素可以满足两者

– a [a,b,c,d,e,f,g,t]>一个[A,B,C,d,E,F,G,-1]

– a [a,b,c,d,e,f,g,t]< b [a,b,c,d,e,f,g](b是高于最新观测值的最小元素) 所以,假设前一个数组是一个简单的堆栈,如果[0,2,1]沿着最后一个观察,我会期待

b = ones((1,1,2,2,1,1,10))*2

然而,

– 如果在某些(a,b,c,d,e,f,g)中,不仅有{0,1,2}的值,还有{3}的值,那么我仍然希望2 (因为i = {2,3}中的较小者满足i> 1.

– 如果在某些(a,b,c,d,e,f,g)中只有值{0,1,3},我会想要3,因为i = 3将是满足i的最小数字> 1.

希望能把它清理一下吗？

/ EDIT2：

非常感谢答案,它有效.如果我想要相反的话,我会如何调整它,即较小的那些中最大的？我没有试图通过那个复杂的索引逻辑,所以我(仅)更改前三行的尝试没有成功：

b = sort(a[...,:-1], axis=-1)

b = b[...,::-1]

mask = b < a[..., -1:]

index = argmax(mask, axis=-1)

indices = tuple([arange(j) for j in a.shape[:-1]])

indices = meshgrid(*indices, indexing='ij', sparse=True)

indices.append(index)

indices = tuple(indices)

a[indices]

另外,我的第二次尝试[…,:: – 1] [指数]也没有成果.

最佳答案 我认为E先生走在正确的轨道上.你绝对可以在没有上一次值的情况下对数组进行排序：

b = np.sort(a[..., :-1], axis=-1)

你现在理想地使用`np.searchsorted来查找第一个项大于结束值的位置,但不幸的是np.searchsorted仅适用于平顶数组,因此我们必须做更多的工作,比如创建一个布尔掩码,然后找到使用np.argmax的第一个True：

mask = b > a[..., -1:]

index = np.argmax(mask, axis=-1)

你现在有了索引,要提取实际值,你需要做一些索引魔术：

indices = tuple([np.arange(j) for j in b.shape[:-1]])

indices = np.meshgrid(*indices, indexing='ij', sparse=True)

indices.append(index)

indices = tuple(indices)

你现在终于可以做到了：

>>> b[indices]

array([[[[[[[ 2., 2., 2., 2., 2., 2., 2., 2., 2., 2.]]],

[[[ 2., 2., 2., 2., 2., 2., 2., 2., 2., 2.]]]],

[[[[ 2., 2., 2., 2., 2., 2., 2., 2., 2., 2.]]],

[[[ 2., 2., 2., 2., 2., 2., 2., 2., 2., 2.]]]]]]])

>>> b[indices].shape

(1L, 1L, 2L, 2L, 1L, 1L, 10L)

要在较小的那些中获得最大的,你可以做类似的事情：

mask = b >= a[..., -1:]

index = np.argmax(mask, axis=-1) - 1

即,在较小的那些中最大的项目是在相等或更大的项目中最小的项目.第二种情况更清楚地表明,如果没有满足条件的项目,这种方法会产生垃圾结果.在第二种情况下,当发生这种情况时,您将获得一个-1的索引,因此您可以检查结果是否有效执行np.any(index == -1).

如果通过执行第一种情况无法满足条件,则可以将索引设置为-1

mask = b > a[..., -1:]

wrong = np.all(~mask, axis=-1)

index = np.argmax(mask, axis=-1)

index[wrong] = -1