有没有一种快速的方法可以将numpy数组中的一个元素与该数组中其余元素进行比较？

Question

我有一个数组，我想查看该数组中的任何元素是否大于或等于该数组中的其他任何元素。我可以使用两个for循环来实现，但是我的数组长度为10,000或更长，因此这会创建一个非常缓慢的程序。有没有什么方法可以更快地完成这个任务？

[编辑] 我只需要查看它是否大于或等于我正在查看的元素后面的元素，如果是，则需要知道它的索引。

[编辑] 我将重新解释我的问题，因为当前的解决方案不适用于我所需要的内容。首先，这是一些代码。

x=linspace(-10, 10, 10000)
t=linspace(0,5,10000)

u=np.exp(-x**2)

k=u*t+x

所以，我有一个x数组，通过将其放入高斯函数来获取其高度，然后基于该高度，这是x值通过空间传播的速度，在此我使用k来找到。我的问题是，我需要找到高斯函数变为双值函数的时候（或换句话说，发生激波时）。如果我使用argmax解决方案，我总是会得到k中的最后一个值，因为它非常接近零，我需要在元素之后的第一个值，这将为我提供双值函数。

[编辑]小例子

x=[0,1,2,3,4,5,6,7,8,9,10] #Input 
k=[0,1,2,3,4,5,6,5,4,10] #adjusted for speed

output I want
in this case, 5 is the first number that goes above a number that comes after it.
So I need to know the index of where 5 is located and possibly the index 
of the number that it is greater than

- Wesley Bowman

---

你能更精确地描述一下你的问题吗？按照目前的情况，你可以简单地打印“是的”，因为数组中的某个元素实际上大于或等于所有其他元素。你是在试图找到最大的元素吗？ - Robᵩ

我想最重要的是一旦知道某个元素大于前面的元素，就要找到它的索引。 - Wesley Bowman

让我尝试重新陈述您的问题：给定列表 A 和索引 I，确定 A[I] 是否大于 A 中所有后续值。如果不是，则确定 A 中最大后续值的索引。这样对吗？ - Robᵩ

是的，那就是它。 - Wesley Bowman

Answer 1

第一个大于后面值的值必然对应于局部最小值中的最小值：

k = np.array([0,1,2,3,4,5,6,5,4,10])
lm_i = np.where(np.diff(np.sign(np.diff(k))) > 0)[0] + 1
mlm = np.min(k[lm_i])
mlm_i = lm_i[np.argmin(k[lm_i])]

首先，找到第一个比后面数值大的数的索引，然后该最小局部最小值之后第一个索引就是所要寻找的索引：

i = np.where(k > mlm)[0][0]

请忽略图像似乎在切线处穿过水平线的显示问题，那只是一个显示问题。

一句话总结：

np.where(k > np.min(k[np.where(np.diff(np.sign(np.diff(k))) > 0)[0] + 1]))[0][0]

请注意，这种方法比root的解决方案快大约1000倍，因为它完全向量化。

%timeit np.where(k > np.min(k[np.where(np.diff(np.sign(np.diff(k))) > 0)[0] + 1]))[0][0]
1000 loops, best of 3: 228 us per loop

- ecatmur

---

太棒了，这绝对是最快的解决方案。由于我一直在遇到值错误，所以我不得不添加一个try except代码块，但是之后它表现得非常好。谢谢！ - Wesley Bowman

你可以使用 np.where(k > np.min(k[np.where(np.diff(k) < 0)[0][0]:]))[0][0] 再次减少25%的时间，而且你实际上不需要 np.min() 的调用。 - root

@root 如果可能存在多个局部最小值，则np.min是必要的。不过，我喜欢你的解决方案；在第一个局部最大值之后找到全局最小值比考虑多个局部最小值要简单得多。 - ecatmur

Answer 2

矢量化解决方案，比ecatmur的解决方案快约25％：

np.where(k > np.min(k[np.where(np.diff(k) < 0)[0][0]:]))[0][0]

一个幼稚的方法：

next(i for i in np.arange(len(arr)) if arr[i:].argmin() != 0)

- root

---

仍然不确定 OP 在问什么，但是 arr[i:].argmax() 显然是解决方案。 - ecatmur

请看我编辑后的问题，就会明白为什么这对我行不通。不过我确实学到了一些东西，所以感谢你。 - Wesley Bowman

@NightHallow -- 即使您进行了当前的编辑，理解您想要实现的内容仍然代价高昂。请提供一个小的输入示例（长度为10），并附上您希望得到的输出和解释... - root

@NightHallow - 这是您想要的输出吗？还是您只想要大于其直接后继的索引 - 如果是这样，那么它甚至更简单... - root

我已经有了直接后继的代码，所以这就是我要找的。 - Wesley Bowman

Answer 3

编辑 实际上，使用一个包含10,000个元素的Python循环比操作一个包含100,000,000个元素的数组更便宜：

In [14]: np.where(np.array([True if np.all(k[:j] <= k[j]) else
                            False for j in xrange(len(k))]) == 0)
Out[14]: (array([5129, 5130, 5131, ..., 6324, 6325, 6326]),)

In [15]: %timeit np.where(np.array([True if np.all(k[:j] <= k[j]) else
                                    False for j in xrange(len(k))]) == 0)
1 loops, best of 3: 201 ms per loop

就内存而言，这会很昂贵，但您可以使用广播来矢量化搜索。 如果您这样做：

>>> k <= k[:, None]
array([[ True, False, False, ..., False, False, False],
       [ True,  True, False, ..., False, False, False],
       [ True,  True,  True, ..., False, False, False],
       ..., 
       [ True,  True,  True, ...,  True, False, False],
       [ True,  True,  True, ...,  True,  True, False],
       [ True,  True,  True, ...,  True,  True,  True]], dtype=bool)

返回的是布尔值数组，其中位置上的元素告诉你k[j]是否小于或等于k[i]。我们可以使用以下方法：np.cumprod。

>>> np.cumprod(k <= k[:, None], axis=1)
array([[1, 0, 0, ..., 0, 0, 0],
       [1, 1, 0, ..., 0, 0, 0],
       [1, 1, 1, ..., 0, 0, 0],
       ..., 
       [1, 1, 1, ..., 1, 0, 0],
       [1, 1, 1, ..., 1, 1, 0],
       [1, 1, 1, ..., 1, 1, 1]])

在位置 [i, j] 的项目告诉您是否 k[j] 小于或等于 k[:i] 中的所有项目。 如果您取该矩阵的对角线：

>>> np.cumprod(k <= k[:, None], axis=1)[np.diag_indices(k.shape[0])]
array([1, 1, 1, ..., 1, 1, 1])

在位置[i]的项目告诉您是否k[i]小于或等于其前面所有项目。找到该数组为零的位置：

>>> np.where(np.cumprod(k <= k[:, None],
...                     axis=1)[np.diag_indices(k.shape[0])] == 0)
(array([5129, 5130, 5131, ..., 6324, 6325, 6326]),)

你将得到所有满足所需条件的值的索引。

如果您只对第一个感兴趣：

>>> np.argmax(np.cumprod(k <= k[:, None],
...                      axis=1)[np.diag_indices(k.shape[0])] == 0)
5129

这不是一个轻松的操作，但如果您有足够的内存来容纳所有的布尔数组，它不会让您等待太长时间：

In [3]: %timeit np.argmax(np.cumprod(k <= k[:, None],
                                     axis=1)[np.diag_indices(k.shape[0])] == 0)
1 loops, best of 3: 948 ms per loop

- Jaime

---

你的修改导致出现错误，似乎是由于一个]放错了位置？我无法让它工作。 - Wesley Bowman

现在可以用了！但是我接受的解决方案对我需要的东西更好一些，但我很可能也会将此解决方案用于我计划中的其他事情。谢谢！ - Wesley Bowman

@NightHallow ecatmur的回答显然是正确的方法：好的数学总是胜过编程。 - Jaime