二进制 numpy 数组之间的快速汉明距离计算

Question

我有两个包含二进制值的相同长度的 numpy 数组

import numpy as np
a=np.array([1, 1, 1, 1, 1, 1, 0, 1, 1, 0, 1, 1, 1, 0, 0, 0, 0, 1, 1, 1, 0])
b=np.array([1, 1, 1, 1, 0, 1, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 0, 1, 1, 0, 1])

我想尽快计算它们之间的汉明距离，因为我要进行数百万次这样的距离计算。

这是一个简单但缓慢的选项（取自维基百科）：

%timeit sum(ch1 != ch2 for ch1, ch2 in zip(a, b))
10000 loops, best of 3: 79 us per loop

我想出了更快的选项，灵感来自这里关于堆栈溢出的一些答案。

%timeit np.sum(np.bitwise_xor(a,b))
100000 loops, best of 3: 6.94 us per loop

%timeit len(np.bitwise_xor(a,b).nonzero()[0])
100000 loops, best of 3: 2.43 us per loop

我想知道是否有更快的方法来计算这个，可能使用 cython？

Answer 1

有一个现成的 numpy 函数胜过 len((a != b).nonzero()[0]) ;)

np.count_nonzero(a!=b)

Answer 2

与我平台上 np.count_nonzero(a!=b) 的 1.07µs 相比，gmpy2.hamdist 在将每个数组转换为 mpz（多精度整数）后降至约 143ns :

import numpy as np
from gmpy2 import mpz, hamdist, pack

a = np.array([1, 1, 1, 1, 1, 1, 0, 1, 1, 0, 1, 1, 1, 0, 0, 0, 0, 1, 1, 1, 0])
b = np.array([1, 1, 1, 1, 0, 1, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 0, 1, 1, 0, 1])

根据@casevh 的提示，使用 gmpy2.pack(list(reversed(list(array))) 可以相当高效地完成从 1 和 0 的一维数组到 gmpy2 mpz 对象的转换， 1).

# gmpy2.pack reverses bit order but that does not affect
# hamdist since both its arguments are reversed
ampz = pack(list(a),1) # takes about 4.29µs
bmpz = pack(list(b),1)

hamdist(ampz,bmpz)
Out[8]: 7

%timeit hamdist(ampz,bmpz)
10000000 loops, best of 3: 143 ns per loop

相对比较，在我的平台上：

%timeit np.count_nonzero(a!=b)
1000000 loops, best of 3: 1.07 µs per loop

%timeit len((a != b).nonzero()[0])
1000000 loops, best of 3: 1.55 µs per loop

%timeit len(np.bitwise_xor(a,b).nonzero()[0])
1000000 loops, best of 3: 1.7 µs per loop

%timeit np.sum(np.bitwise_xor(a,b))
100000 loops, best of 3: 5.8 µs per loop

Answer 3

使用pythran可以在这里带来额外的好处：

$ cat hamm.py
#pythran export hamm(int[], int[])
from numpy import nonzero
def hamm(a,b):
    return len(nonzero(a != b)[0])

作为参考（没有pythran）：

$ python -m timeit -s 'import numpy as np; a = np.random.randint(0,2, 100); b = np.random.randint(0,2, 100); from hamm import hamm' 'hamm(a,b)'
100000 loops, best of 3: 4.66 usec per loop

在 pythran 编译之后：

$ python -m pythran.run hamm.py
$ python -m timeit -s 'import numpy as np; a = np.random.randint(0,2, 100); b = np.random.randint(0,2, 100); from hamm import hamm' 'hamm(a,b)'
1000000 loops, best of 3: 0.745 usec per loop

这大约是 numpy 实现的 6x 加速，因为 pythran 在评估元素明智比较时跳过了中间数组的创建。

我也测量了：

def hamm(a,b):
    return count_nonzero(a != b)

我得到 Python 版本的 3.11 usec per loop 和 Pythran 版本的 0.427 usec per loop。

免责声明：我是 Pythran 开发人员之一。

Answer 4

对于字符串，它工作得更快

def Hamm(a, b):
    c = 0
    for i in range(a.shape[0]):
        if a[i] != b[i]:
            c += 1
    return c

Answer 5

我建议您使用 np.packbits

将 numpy 位数组转换为 numpy uint8 数组

看看 scipy 的 spatial.distance.hamming： https://docs.scipy.org/doc/scipy/reference/generated/scipy.spatial.distance.hamming.html

否则，这是一个小片段，它只需要受 Fast way of counting non-zero bits in positive integer 启发的 numpy :

bit_counts = np.array([int(bin(x).count("1")) for x in range(256)]).astype(np.uint8)
def hamming_dist(a,b,axis=None):
    return np.sum(bit_counts[np.bitwise_xor(a,b)],axis=axis)

轴=-1，这允许在一个条目和一个大数组之间取hammig距离；例如：

inp = np.uint8(np.random.random((512,8))*255) #512 entries of 8 byte
hd = hamming_dist(inp, inp[123], axis=-1) #results in 512 hamming distances to entry 123
idx_best = np.argmin(hd)    # should point to identity 123
hd[123] = 255 #mask out identity
idx_nearest= np.argmin(hd)    # should point entry in list with shortest distance to entry 123
dist_hist = np.bincount(np.uint8(hd)) # distribution of hamming distances; for me this started at 18bits to 44bits with a maximum at 31

二进制 numpy 数组之间的快速汉明距离计算

Fast hamming distance computation between binary numpy arrays

python

arrays

numpy

cython

hamming-distance