【问题标题】:Scipy sparse... arrays?Scipy稀疏...数组?
【发布时间】:2011-02-02 03:51:49
【问题描述】:

所以,我正在使用非常稀疏的 numpy 数组进行一些 Kmeans 分类 - 很多很多零。我想我会使用 scipy 的“稀疏”包来减少存储开销,但我对如何创建数组而不是矩阵有点困惑。

我已经阅读了有关如何创建稀疏矩阵的教程: http://www.scipy.org/SciPy_Tutorial#head-c60163f2fd2bab79edd94be43682414f18b90df7

为了模拟一个数组,我只创建了一个 1xN 矩阵,但正如您可能猜到的那样,Asp.dot(Bsp) 并不能很好地工作,因为您不能将两个 1xN 矩阵相乘。我必须将每个数组转置为 Nx1,这很糟糕,因为我会为每个点积计算都这样做。

接下来,我尝试创建一个 NxN 矩阵,其中第 1 列 == 第 1 行(这样您可以将两个矩阵相乘并将左上角作为点积),但结果证明效率非常低.

我很想使用 scipy 的 sparse 包作为 numpy 的 array() 的神奇替代品,但到目前为止,我还不确定该怎么做。

有什么建议吗?

【问题讨论】:

  • 见下面的 cmets,但我最终只是滚动了我自己的稀疏向量实现,使用类似于“dok”矩阵的东西。
  • 原来的问题链接好像失效了。 @spitzanator。

标签: python matrix numpy scipy sparse-matrix


【解决方案1】:

使用基于行或列的scipy.sparse 格式:csc_matrixcsr_matrix

这些在底层使用高效的 C 实现(包括乘法),并且转置是无操作的(尤其是如果您调用 transpose(copy=False)),就像使用 numpy 数组一样。

编辑:通过ipython进行一些时间安排:

import numpy, scipy.sparse
n = 100000
x = (numpy.random.rand(n) * 2).astype(int).astype(float) # 50% sparse vector
x_csr = scipy.sparse.csr_matrix(x)
x_dok = scipy.sparse.dok_matrix(x.reshape(x_csr.shape))

现在 x_csrx_dok 是 50% 稀疏的:

print repr(x_csr)
<1x100000 sparse matrix of type '<type 'numpy.float64'>'
        with 49757 stored elements in Compressed Sparse Row format>

还有时间:

timeit numpy.dot(x, x)
10000 loops, best of 3: 123 us per loop

timeit x_dok * x_dok.T
1 loops, best of 3: 1.73 s per loop

timeit x_csr.multiply(x_csr).sum()
1000 loops, best of 3: 1.64 ms per loop

timeit x_csr * x_csr.T
100 loops, best of 3: 3.62 ms per loop

看来我撒了谎。转置 非常便宜,但没有 csr * csc 的有效 C 实现(在最新的 scipy 0.9.0 中)。每次调用都会构造一个新的 csr 对象 :-(

作为一个hack(虽然scipy现在比较稳定),你可以直接在稀疏数据上做点积:

timeit numpy.dot(x_csr.data, x_csr.data)
10000 loops, best of 3: 62.9 us per loop

请注意,最后一种方法再次执行 numpy 密集乘法。稀疏度为 50%,因此它实际上比 dot(x, x) 快 2 倍。

【讨论】:

  • +1 表示普通的 numpy.dot。对于 kmeans,您需要 argmax( dot( k x N 个中心,每个 Nvec x ));无论如何,中心都会变得密集,所以不妨保持它们密集。 (不过,为新中心平均许多稀疏 x 非常慢。)
  • 好吧,如果我们把乘法速度放在一边,OP还不如使用scipy.cluster.kmeans...
  • 似是而非。我更喜欢(advt)this code,它可以使用 scipy.spatial.distance 中的 20 多个指标中的任何一个;对于高维度 kmeans,度量比算法更重要。
【解决方案2】:

您可以创建现有二维稀疏数组之一的子类

from scipy.sparse import dok_matrix

class sparse1d(dok_matrix):
    def __init__(self, v):
        dok_matrix.__init__(self, (v,))
    def dot(self, other):
        return dok_matrix.dot(self, other.transpose())[0,0]

a=sparse1d((1,2,3))
b=sparse1d((4,5,6))
print a.dot(b)

【讨论】:

  • 不幸的是,这样做的问题是您必须即时转置这些问题,当您进行数百万次比较时,这没有多大意义。我尝试缓存点积,但不幸的是,我们不经常做相同的点积,所以没有太大帮助。
【解决方案3】:

我不确定它是否真的更好或更快,但您可以这样做以避免使用转置:

Asp.multiply(Bsp).sum()

这只是取两个矩阵的逐元素乘积并将乘积相加。您可以创建您使用的任何矩阵格式的子类,将上述语句作为点积。

但是,转置它们可能更容易:

Asp*Bsp.T

这似乎没什么大不了的,但您也可以创建一个子类并修改 mul() 方法。

【讨论】:

  • 我也试过,对于一个向量 [1, 2, 3],创建一个矩阵: [1, 2, 3] [2, 0, 0] [3, 0, 0] 取两个这些和乘法(以任何顺序)在结果矩阵的左上角给出所需的点积。不幸的是,这严重影响了速度。
猜你喜欢
  • 2018-01-11
  • 2021-05-15
  • 1970-01-01
  • 2017-03-26
  • 2017-03-31
  • 2016-05-26
  • 2021-03-01
  • 1970-01-01
  • 2011-09-18
相关资源
最近更新 更多