遍历1D Numpy阵列进行群集

Question

我有一个1维的numpy数组,其中每个元素的值i指向数组中的另一个索引.每个集群中心都有一个唯一的负(整数)值.目标是将每个元素分配给一个集群.

例

# Generated/pre-computed elsewhere
a = np.array([6, 8, 1, -1, 0, 3, -2, 4, -3, 10, 5])

因此,聚类中心是元素3,6和8(因为它们具有负值)并且它们分别被标记为聚类-1,-2和-3.所以,从那以后

a[0] = 6 --> a[6] = -2, 

然后a [0]可以指定为-2.同样,因为

a[5] = 3 --> a[3] = -1,

然后a [5]可以指定为-1.遵循此逻辑,然后可以将所有元素分配给集群中心.结果数组将是:

[-2, -3, -3, -1, -2, -1, -2, -2, -3, -1, -1]

我知道如何在纸上实现这一点,但我不知道如何在代码中实现这一点或在numpy中使用矢量化代码.

更新:根据下面unutbu的答案,我用for循环替换了while循环,以避免无限循环:

a = np.array([6, 8, 1, -1, 0, 3, -2, 4, -3, 10, 5])
for i in range(len(a)):
    mask = a >= 0
    if not mask.any(): break
    a[mask] = a[a[mask]]

Answer 1

我们不知道先验需要多少"跳"找到一个集群中心.因此,我们必须进行一些迭代并检查我们是否已登陆集群中心:

import numpy as np
a = np.array([6, 8, 1, -1, 0, 3, -2, 4, -3, 10, 5])

for i in a:
    mask = a>=0
    # We can stop when all the values in `a` are negative
    if not mask.any(): break
    # perform a hop
    a[mask] = a[a[mask]]

print(a)

产量

[-2 -3 -3 -1 -2 -1 -2 -2 -3 -1 -1]

要了解正在发生的事情,可以更清楚地看一下更简单的值a:

a = np.arange(-1, 10)
print(a)
for i in a:
    mask = a>=0
    # We can stop when all the values in `a` are negative
    if not mask.any(): break
    # perform a hop
    a[mask] = a[a[mask]]
    print(a)

版画

[-1  0  1  2  3  4  5  6  7  8  9]
[-1 -1  0  1  2  3  4  5  6  7  8]
[-1 -1 -1 -1  0  1  2  3  4  5  6]
[-1 -1 -1 -1 -1 -1 -1 -1  0  1  2]
[-1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1]

在这里我们可以清楚地看到每个索引(即输出中的一列)的值正在减少.您可能希望值​​减少1,但实际上由于先前跃点的影响,它们在后续迭代中会减少多个值.

更一般地,让G一系列索引跳转到相同的簇值. 循环的不变量是G中的值映射到G中的值.此外,在每个索引处,对于每次迭代,到簇值的距离(即,跳数)减小(或者为零).这两个事实共同意味着Banach不动点定理,该算法收敛于一个固定点.

此外,每次迭代时,填充簇值的位置数量呈指数增长.如果在迭代之前存在n填充了簇值的位置,则在跳跃之后将存在2n填充有簇值的位置.所以收敛是指数级的.

在OP领先后,我将while True循环更改为a for-loop.固定点将在少于len(a)步骤中找到,但仍然for i in a足够.

请注意,上面的代码假定每个索引都收敛到一个簇值.如果不是这样,那么原始while True循环可能永远循环.该for i in a会结束,但也许不是在所有群集值.一个说明问题的简单例子是a = np.array([1,0]).

有趣的是,就简单np.arange(-1, 10)长度为11的数组所需的最大迭代而言,"简单示例" 也是最坏情况.由于簇值的数量增长为2**n,因此循环最多需要log2(len(a))迭代次序.因此,我们可以编写一个函数,它将返回簇值或在a包含无限循环时引发ValueError :

def converge(a):
    for i in range(int(np.log2(len(a)))+1):
        mask = a>=0
        # We can stop when all the values in `a` are negative
        if not mask.any(): break
        # perform a hop
        a[mask] = a[a[mask]]
    else:
        # for-loop completed without reaching break
        mask = a>=0
        if mask.any(): 
            raise ValueError('{!r} contains infinite cycles'.format(a))
    return a