定义一个类似于numpy.isclose
的isclose
函数,但速度更快(主要是因为不检查任何输入并且不支持相对和绝对容差):
import numpy as np
array1 = np.array([[(1.22, 5.64)],
[(2.31, 7.63)],
[(4.94, 4.15)]])
array2 = np.array([[(1.23, 5.63)],
[(6.31, 10.63)],
[(2.32, 7.65)]])
def isclose(x, y, atol):
return np.abs(x - y) < atol
现在来到了难点。我们需要计算内部维度中是否有任何两个值是接近的。为此,我会重新塑造数组,使第一个数组的值沿着第二个维度复制到第一维和第二维上,而第二个数组的值沿着第一维复制到第二维上(注意1, 3
和3, 1
):
In [92]: isclose(array1.reshape(1,3,2), array2.reshape(3,1,2), 0.03)
Out[92]:
array([[[ True, True],
[False, False],
[False, False]],
[[False, False],
[False, False],
[False, False]],
[[False, False],
[ True, True],
[False, False]]], dtype=bool)
现在我们想要所有数值接近其他数值的条目(在同一维度上):
In [93]: isclose(array1.reshape(1,3,2), array2.reshape(3,1,2), 0.03).any(axis=0)
Out[93]:
array([[ True, True],
[ True, True],
[False, False]], dtype=bool)
然后我们只想要那些元组中两个值都接近的:
In [111]: isclose(array1.reshape(1,3,2), array2.reshape(3,1,2), 0.03).any(axis=0).all(axis=-1)
Out[111]: array([ True, True, False], dtype=bool)
最后,我们可以使用这个来索引
array1
:
In [112]: array1[isclose(array1.reshape(1,3,2), array2.reshape(3,1,2), 0.03).any(axis=0).all(axis=-1)]
Out[112]:
array([[[ 1.22, 5.64]],
[[ 2.31, 7.63]]])
如果您想的话,可以交换
any
和
all
的调用。在您的情况下,其中一个可能比另一个更快。
reshape
调用中的
3
需要替换为您数据的实际长度。
该算法将具有与使用
itertools.product
的其他答案相同的糟糕运行时,但至少实际循环是由
numpy
隐式执行的并且是用C实现的。这在计时中可见。
In [122]:
11.6 µs ± 493 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)
In [126]:
267 µs ± 8.72 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)
函数pares
是由@Ferran Parés在另一个答案中定义的代码,数组已经被重塑。
对于更大的数组,这变得更加明显:
array1 = np.random.normal(0, 0.1, size=(1000, 1, 2))
array2 = np.random.normal(0, 0.1, size=(1000, 1, 2))
array1_pares = array1.reshape(1000, 2)
array2_pares = arra2.reshape(1000, 2)
In [149]: %timeit array1[isclose(array1.reshape(1,len(array1),2), array2.reshape(len(array2),1,2), 0.03).any(axis=0).all(axis=-1)]
135 µs ± 5.34 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)
In [157]: %timeit pares(array1_pares, array2_pares)
1min 36s ± 6.85 s per loop (mean ± std. dev. of 7 runs, 1 loop each)
最终,这受可用系统内存的限制。我的机器(16GB RAM)仍然可以处理长度为20000的数组,但这几乎将其推到了100%。这也需要大约12秒:
In [14]: array1 = np.random.normal(0, 0.1, size=(20000, 1, 2))
In [15]: array2 = np.random.normal(0, 0.1, size=(20000, 1, 2))
In [16]:
12.3 s ± 514 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)
np.around(array1, 1)
或向上取整值np.ceil(array1)
来近似你的结果。 - J. Doe