1

我有一个非常大的 scipy 稀疏 csr 矩阵。它是一个 100,000x2,000,000 维矩阵。让我们称之为X。每行是 2,000,000 维空间中的样本向量。

我需要非常有效地计算每对样本之间的余弦距离。我一直在使用sklearn pairwise_distances带有向量子集的函数，X其中给了我一个密集矩阵 D：包含冗余​​条目的成对距离的平方形式。如何使用sklearn pairwise_distances直接获取压缩形式？请参阅http://docs.scipy.org/doc/scipy/reference/generated/scipy.spatial.distance.pdist.html以了解压缩形式是什么。它是scipy pdist函数的输出。

我有内存限制，我无法计算平方形式，然后得到压缩形式。由于内存限制，我也不能使用scipy pdist它，因为它需要一个密集的矩阵X，它不再适合内存。我想过循环遍历不同的块X并计算每个块的压缩形式并将它们连接在一起以获得完整的压缩形式，但这相对繁琐。有更好的想法吗？

非常感谢任何帮助。提前致谢。

下面是一个可重现的例子（当然为了演示目的X要小得多）：

from scipy.sparse import rand
from scipy.spatial.distance import pdist
from sklearn.metrics.pairwise import pairwise_distances
X = rand(1000, 10000, density=0.01, format='csr')
dist1 = pairwise_distances(X, metric='cosine')
dist2 = pdist(X.A, 'cosine')

如您所见dist2，它是压缩形式，是一个 499500 维向量。但是dist1是对称的正方形，是一个 1000x1000 的矩阵。

4

1 回答 1

4

我深入研究了这两个版本的代码，并认为我了解两者都在做什么。

从一个小的简单X（密集）开始：

X = np.arange(9.).reshape(3,3)

pdist余弦会：

norms = _row_norms(X)
_distance_wrap.pdist_cosine_wrap(_convert_to_double(X), dm, norms)

行点在哪里_row_norms- 使用einsum：

norms = np.sqrt(np.einsum('ij,ij->i', X,X)

所以这是第一个X必须是数组的地方。

我还没有深入研究 cosine_wrap，但它似乎可以（可能在 cython 中）

xy = np.dot(X, X.T)
# or xy = np.einsum('ij,kj',X,X)

d = np.zeros((3,3),float)   # square receiver
d2 = []                     # condensed receiver
for i in range(3):
    for j in range(i+1,3):
         val=1-xy[i,j]/(norms[i]*norms[j])
         d2.append(val)
         d[j,i]=d[i,j]=val

print('array')
print(d)
print('condensed',np.array(d2))

from scipy.spatial import distance
d1=distance.pdist(X,'cosine')
print('    pdist',d1)

生产：

array
[[ 0.          0.11456226  0.1573452 ]
 [ 0.11456226  0.          0.00363075]
 [ 0.1573452   0.00363075  0.        ]]

condensed [ 0.11456226  0.1573452   0.00363075]
    pdist [ 0.11456226  0.1573452   0.00363075]

distance.squareform(d1)d产生与我的数组相同的东西。

xy我可以通过将点积与适当的norm外积相除来生成相同的方阵：

dd=1-xy/(norms[:,None]*norms)
dd[range(dd.shape[0]),range(dd.shape[1])]=0 # clean up 0s

X或者在采用点积之前进行标准化。这似乎是scikit版本所做的。

Xnorm = X/norms[:,None]
1-np.einsum('ij,kj',Xnorm,Xnorm)

scikit添加了一些 cython 代码来进行更快的稀疏计算（超出 提供的那些sparse.sparse，但使用相同的csr格式）：

from scipy import sparse
Xc=sparse.csr_matrix(X)

# csr_row_norm - pyx of following
cnorm = Xc.multiply(Xc).sum(axis=1)
cnorm = np.sqrt(cnorm)
X1 = Xc.multiply(1/cnorm)  # dense matrix
dd = 1-X1*X1.T

要获得具有稀疏矩阵的快速压缩形式，我认为您需要实现X1*X1.T. 这意味着您需要了解如何在c代码中实现稀疏矩阵乘法。cython scikit'fast sparse' 代码也可能提供一些想法。

numpy有一些tri...函数是直接的 Python 代码。它不会试图通过直接实施三计算来节省时间或空间。迭代 nd 数组的矩形布局（具有形状和步幅）比执行更复杂的三角形数组的可变长度步骤更容易。精简的形式只是将空间和计算步骤减少了一半。

=============

c这是函数 的主要部分pdist_cosine，它迭代i和上层j，计算dot(x[i],y[j])/(norm[i]*norm[j])。

for (i = 0; i < m; i++) {
    for (j = i + 1; j < m; j++, dm++) {
        u = X + (n * i);
        v = X + (n * j);
        cosine = dot_product(u, v, n) / (norms[i] * norms[j]);
        if (fabs(cosine) > 1.) {
            /* Clip to correct rounding error. */
            cosine = npy_copysign(1, cosine);
        }
        *dm = 1. - cosine;
    }
}

https://github.com/scipy/scipy/blob/master/scipy/spatial/src/distance_impl.h

于 2016-07-14T01:41:55.777 回答