2

我有两个数组，看起来像这样：

X = np.array([ 157,  262,  368,  472,  577,  682,  786,  891,  996, 1100, 1204,
       1310, 1415, 1520, 1625, 1731, 1879])

Y = np.array([  30,  135,  240,  345,  450,  555,  660,  765,  870,  975, 1080,
       1185, 1290, 1395, 1500, 1605])

阵列将：

  • 使值从开始按升序排序。
  • 有时长度不等。

我想Z根据以下内容将这两个交错成一个新数组：

  • 每个元素只能使用一次
  • 不需要使用所有元素
  • 仅当存在元素in且其中没有其他元素的值差小于且不存在其值距离小于的元素时，Xi才可以包含元素。（同样的规则适用于 中的元素。）ZYjYYabs(Xi - Yj)XYjabs(Xi - Yj)Y

我看到我可以用一堆嵌套的 for 循环来做到这一点，但我想知道是否有一些更聪明、更整洁的方法来做到这一点？

（我意识到，我提出这个问题的方式，这听起来像是从教科书上剪下来的。它不是。但也许这是一个经典的排序功能，谁知道呢，但对于我作为一名生物学家来说......我只能说我不知道如何以有效，整洁的方式解决它。）

编辑：不是那么漂亮的例子

new_list = list()
for i in X:
    delta_i = np.abs(Y - i)
    delta_reciprocal = np.abs(X - Y[delta_i.argmin()])
    if delta_i.min() == delta_reciprocal.min():
        new_list += sorted([Y[delta_i.argmin()],
        X[delta_reciprocal.argmin()]])
Z = np.array(new_list)

我什至不完全确定它是否满足所有标准，但是在重写旧代码时，我只需要一个循环......仍然必须有一些更好的方法！

4

2 回答 2

5

让我们试着找出这个例子的解决方案：

In [1]: import numpy as np

In [5]: X = np.array([1879, 1731])

In [6]: Y = np.array([1481, 1691, 1586, 1796])

X我们可以像这样计算值 in和值 in之间的所有距离Y：

In [7]: dist = np.abs(np.subtract.outer(X,Y))

In [8]: dist
Out[8]: 
array([[398, 188, 293,  83],
       [250,  40, 145,  65]])

行对应于X值，列对应于Y值。

为了找到X最接近 中某个元素的值Y，我们正在寻找对应于矩阵列中的X最小值的值 。每列对应一个特定的，因此一列中的最小距离对应于一些和特定之间的最小值。distYXY

从视觉上讲，我们正在寻找的是它们所在行和dist所在列的最小值。我们称它们为“行列最小值”。

在上面的dist数组中，40 是行列最小值。65 是列最小值，但不是行列最小值。

对于每一列，我们可以通过这种方式找到最小化该列的 X-index：

In [6]: idx1 = np.argmin(dist, axis = 0)

In [7]: idx1
Out[7]: array([1, 1, 1, 1])

同样，对于每一行，我们可以通过以下方式找到 Y 索引：

In [8]: idx2 = np.argmin(dist, axis = 1)

In [9]: idx2
Out[9]: array([3, 1])

现在，让我们暂时忘记这个例子，假设idx1看起来像这样：

        0,1,2,3,4,5   # the index value 
idx1 = (_,_,_,_,_,2,...)

这就是说在第 5 列中，第 2 行具有最小值。

然后，如果第 2 行、第 5 列对应于行列最小值，则idx2 必须如下所示：

        0,1,2        # index value
idx2 = (_,_,5,...)

我们可以用 NumPy 来表达这种关系

idx1[idx2] == np.arange(len(X))
idx2[idx1] == np.arange(len(Y))

因此，对应于行列最小值的 X、Y 值是

X[idx1[idx2] == np.arange(len(X))]

和

Y[idx2[idx1] == np.arange(len(Y))]

import numpy as np
tests = [
    (np.array([1879, 1731]),
     np.array([1481, 1691, 1586, 1806])), 
    (np.array([1879, 1731]),
     np.array([1481, 1691, 1586, 1796])),
    (np.array([ 157,  262,  368,  472,  577,  682,  786,  891,  996, 1100, 1204]),
     np.array([  30,  135,  240,  345,  450,  555,  660,  765,  870,  975])),
    (np.array([ 157, 262, 368, 472, 577, 682, 786, 891, 996, 1100, 1204, 1310,
                1415, 1520, 1625, 1731, 1879]),
     np.array([ 221, 326, 431, 536, 641, 746, 851, 956, 1061, 1166, 1271, 1376,
                1481, 1586, 1691, 1796]))]

def find_close(X,Y):
    new_list = list()
    for i in X:
        delta_i = np.abs(Y - i)
        # print(delta_i)
        delta_reciprocal = np.abs(X - Y[delta_i.argmin()])
        if delta_i.min() == delta_reciprocal.min():
            new_list += sorted([Y[delta_i.argmin()],
                                X[delta_reciprocal.argmin()]])
    Z = np.array(new_list)
    return Z

def alt_find_close(X,Y):
    dist = np.abs(np.subtract.outer(X,Y))
    idx1 = np.argmin(dist, axis = 0)
    idx2 = np.argmin(dist, axis = 1)
    Z = np.r_[X[idx1[idx2] == np.arange(len(X))], Y[idx2[idx1] == np.arange(len(Y))]]
    return Z

for X, Y in tests:
    assert np.allclose(sorted(find_close(X,Y)), sorted(alt_find_close(X,Y)))

时间结果：

% python -mtimeit -s'import test' 'test.find_close(test.X, test.Y)'
1000 loops, best of 3: 454 usec per loop
% python -mtimeit -s'import test' 'test.alt_find_close(test.X, test.Y)'
10000 loops, best of 3: 40.6 usec per loop

所以alt_find_close明显快于find_close.

于 2012-12-11T15:57:59.690 回答
2

我想你可能想要使用scipy.spatial.cKDTree（你可以用 numpys 自己构建这样的东西，searchsorted但我看不出有什么意义，除非可能对等距离的问题有更多的控制）。

但是，通常您应该小心。您的示例是整数数组，根据发现的情况可能会出现问题（argmin 总是找到第一个，所以也许没关系，但是对于这个，如果距离相等，您可能会丢分）。

import numpy as np
from scipy.spatial import cKDTree

def find_close_fast(X, Y):
    kX = cKDTree(X[:,None]) # needs to be 2D
    kY = cKDTree(Y[:,None])

    nearest_X = kX.query(Y[:,None], p=1)[1] # might as well use 1-norm

    # Which Y corresponds the other way around?
    nearest_Y = kY.query(X[nearest_X][:,None], p=1)[1]

    w = nearest_Y == np.arange(len(Y))
    result = np.concatenate((X[nearest_X[w]], Y[w]))
    return result

如果您的阵列变大（每个阵列可能大约几百个），这将快得多。例如：

In [121]: X = np.random.random(5000)

In [122]: Y = np.random.random(5000)

In [123]: %timeit alt_find_close(X, Y)
1 loops, best of 3: 1.03 s per loop

In [124]: %timeit find_close_fast(X, Y)
10 loops, best of 3: 23.3 ms per loop

In [125]: np.all(np.sort(find_close_fast(X,Y)) == np.sort(alt_find_close(X, Y))) 
Out[125]: True

但是我在这里使用浮点数是有原因的，如果你可以有相等的距离，则不能保证结果。排序是不同的，并没有真正试图弄清楚为什么。

编辑：实际上，您也可以将两个数组归为一个（并记住哪个属于哪个类），然后从那里检查两个不同类的相邻位置。如果一个点有另一个类的两个邻居，手动选择更接近的一个。这可能更快，并且只使用 numpy。

于 2012-12-11T22:21:15.947 回答