并行构造距离矩阵

Question

我在大量多维向量上进行分层凝聚聚类，我注意到最大的瓶颈是构造距离矩阵。这项任务的天真实现如下（在Python中）：

''' v = an array (N,d), where rows are the observations
and columns the dimensions'''
def create_dist_matrix(v):
   N = v.shape[0]
   D = np.zeros((N,N))
   for i in range(N):
      for j in range(i+1):
          D[i,j] = cosine(v[i,:],v[j,:]) # scipy.spatial.distance.cosine()
   return D

我想知道哪个是为这个例程添加一些并行性的最佳方法。一种简单的方法是打破并将外循环分配给许多作业，例如如果你有10个处理器，为不同范围的i创建10个不同的作业，然后连接结果。然而，这种“横向”解决方案似乎并不合适。是否有任何其他并行算法（或现有库）用于此任务？任何帮助将受到高度赞赏。

Answer 1

看起来像scikit-learn有一个名为pairwise_distances的并行版本的pdist

from sklearn.metrics.pairwise import pairwise_distances

D = pairwise_distances(X = v, metric = 'cosine', n_jobs = -1)

其中n_jobs = -1指定将使用所有CPU。

Answer 2

请参阅@agartland答案 - 您可以在n_jobs中指定sklearn.metrics.pairwise.pairwise_distances，或者使用sklearn.cluster参数在n_jobs查找聚类算法。 E. g。 sklearn.cluster.KMeans。

不过，如果你喜欢冒险，你可以实现自己的计算。例如，如果您需要scipy.cluster.hierarchy.linkage的1D距离矩阵，您可以使用：

#!/usr/bin/env python3
from multiprocessing import Pool
import numpy as np
from time import time as ts


data = np.zeros((100,10)) # YOUR data: np.array[n_samples x m_features]
n_processes = 4           # YOUR number of processors
def metric(a, b):         # YOUR dist function
    return np.sum(np.abs(a-b)) 


n = data.shape[0]
k_max = n * (n - 1) // 2  # maximum elements in 1D dist array
k_step = n ** 2 // 500    # ~500 bulks
dist = np.zeros(k_max)    # resulting 1D dist array


def proc(start):
    dist = []
    k1 = start
    k2 = min(start + k_step, k_max)
    for k in range(k1, k2):
        # get (i, j) for 2D distance matrix knowing (k) for 1D distance matrix
        i = int(n - 2 - int(np.sqrt(-8 * k + 4 * n * (n - 1) - 7) / 2.0 - 0.5))
        j = int(k + i + 1 - n * (n - 1) / 2 + (n - i) * ((n - i) - 1) / 2)
        # store distance
        a = data[i, :]
        b = data[j, :]
        d = metric(a, b)
        dist.append(d)
    return k1, k2, dist


ts_start = ts()
with Pool(n_processes) as pool:
    for k1, k2, res in pool.imap_unordered(proc, range(0, k_max, k_step)):
        dist[k1:k2] = res
        print("{:.0f} minutes, {:,}..{:,} out of {:,}".format(
            (ts() - ts_start)/60, k1, k2, k_max))


print("Elapsed %.0f minutes" % ((ts() - ts_start) / 60))
print("Saving...")
np.savez("dist.npz", dist=dist)
print("DONE")

你知道，scipy.cluster.hierarchy.linkage的实现并不是平行的，它的复杂性至少是O（N * N）。我不确定scipy是否具有此功能的并行实现。

Answer 3

我怀疑你会在pdist模块中比scipy更快地得到它。可能这就是它说的原因

请注意，应避免将引用传递给此库中定义的距离函数之一。例如，：

dm = pdist(X, sokalsneath)

将使用Python函数sokalsneath计算X中向量之间的成对距离。这将导致sokalsneath被称为n选择2次，这是低效的。相反，优化的C版本更有效，我们使用以下语法调用它：

dm = pdist(X, 'sokalsneath')

So no Python function is used, if you use pdist(X, 'cosine'). When I run it, to me it seems, that it does only use one core, so if you have a lot of cores, you might get it faster. But bear in mind, that to achieve this, your native implementation has to be as fast as SciPy's. That won't be trivial. You'd rather be patient or go for a different clustering method, e. g. an algorithm which supports a spatial index.

Answer 4

除了@agartland提出的我喜欢用pairwise_distances或pairwise_disances_chunked和numpy.triu_indices来获得浓缩的距离向量。这是scipy.spatial.distance.pdist提供的确切输出

重要的是要注意k kwarg为triu_indices控制对角线的偏移量。默认值k=0将返回零的对角线以及实际距离值，并应设置为k=1以避免这种情况。

对于大型数据集，我遇到了一个问题，pairwise_distances在从工作线程返回值时从ValueError引发struct.unpack。因此我在下面使用pairwise_distances_chunked。

gen = pairwise_distances_chunked(X, method='cosine', n_jobs=-1)
Z = np.concatenate(list(gen), axis=0)
Z_cond = Z[np.triu_indices(Z.shape[0], k=1)

对我来说，这比使用pdist要快得多，并且可以很好地扩展可用内核的数量。

注：我认为值得指出的是，过去对scipy.cluster.hierarchy.linkage的论点存在一些混淆，因为文档在某一点上表明用户可以传递浓缩或方形距离矢量/矩阵（linkage() function mistakes distance matrix as observation vectors #2614）。实际情况并非如此，传递给连接的值应该是压缩距离矢量或原始观测值的m×n数组。

Answer 5

如果您决定自己编排多处理，则可能需要在CPU之间均匀地分配计算次数，以便最大限度地缩短计算时间。然后回复this question on equally splitting the diagonal matrix可能会派上用场。