我试图弄清楚如何最好地执行稀疏矩阵和稀疏向量的元素加法(和减法)。我在 SO 上找到了 this trick:
mat = sp.csc_matrix([[1,0,0],[0,1,0],[0,0,1]])
vec = sp.csr_matrix([[1,2,1]])
mat.data += np.repeat(vec.toarray()[0], np.diff(mat.indptr))
但不幸的是,它只更新非零值:
print(mat.todense())
[[2 0 0]
[0 3 0]
[0 0 2]]
SO线程上实际接受的答案:
def sum(X,v):
rows, cols = X.shape
row_start_stop = as_strided(X.indptr, shape=(rows, 2),
strides=2*X.indptr.strides)
for row, (start, stop) in enumerate(row_start_stop):
data = X.data[start:stop]
data -= v[row]
sum(mat,vec.A[0])
做同样的事情。不幸的是,我现在已经没有想法了,所以我希望你能帮助我找出解决这个问题的最佳方法。
编辑:
我期待它的功能与密集版本相同:
np.eye(3) + np.asarray([[1,2,1]])
array([[ 2., 2., 1.],
[ 1., 3., 1.],
[ 1., 2., 2.]])
谢谢
最佳答案
使用 10x10 稀疏垫和 vec 进行的一些测试:
In [375]: mat=sparse.rand(10,10,.1)
In [376]: mat
Out[376]:
<10x10 sparse matrix of type '<class 'numpy.float64'>'
with 10 stored elements in COOrdinate format>
In [377]: mat.A
Out[377]:
array([[ 0. , 0. , 0. , 0. , 0. ,
0. , 0. , 0. , 0. , 0. ],
[ 0. , 0. , 0. , 0. , 0. ,
0.15568621, 0.59916335, 0. , 0. , 0. ],
...
[ 0. , 0. , 0.15552687, 0. , 0. ,
0.47483064, 0. , 0. , 0. , 0. ]])
In [378]: vec=sparse.coo_matrix([0,1,0,2,0,0,0,3,0,0]).tocsr()
<1x10 sparse matrix of type '<class 'numpy.int32'>'
with 3 stored elements in Compressed Sparse Row format>
maxymoo 的解决方案:
def addvec(mat,vec):
Mc = mat.tocsc()
for i in vec.nonzero()[1]:
Mc[:,i]=sparse.csc_matrix(Mc[:,i].todense()+vec[0,i])
return Mc
以及使用
lil
格式的变体,在更改稀疏结构时应该更有效:def addvec2(mat,vec):
Ml=mat.tolil()
vec=vec.tocoo()
for i,v in zip(vec.col, vec.data):
Ml[:,i]=sparse.coo_matrix(Ml[:,i].A+v)
return Ml
总和有 38 个非零项,高于原始
mat
中的 10 个。它添加了 vec
中的 3 列。这是稀疏性的一个巨大变化。In [382]: addvec(mat,vec)
Out[382]:
<10x10 sparse matrix of type '<class 'numpy.float64'>'
with 38 stored elements in Compressed Sparse Column format>
In [383]: _.A
Out[383]:
array([[ 0. , 1. , 0. , 2. , 0. ,
0. , 0. , 3. , 0. , 0. ],
[ 0. , 1. , 0. , 2. , 0. ,
0.15568621, 0.59916335, 3. , 0. , 0. ],
...
[ 0. , 1. , 0.15552687, 2. , 0. ,
0.47483064, 0. , 3. , 0. , 0. ]])
与 addvec2 相同的输出:
In [384]: addvec2(mat,vec)
Out[384]:
<10x10 sparse matrix of type '<class 'numpy.float64'>'
with 38 stored elements in LInked List format>
在时间上,
addvec2
比 2x 好In [385]: timeit addvec(mat,vec)
100 loops, best of 3: 6.51 ms per loop
In [386]: timeit addvec2(mat,vec)
100 loops, best of 3: 2.54 ms per loop
和密集的等价物:
In [388]: sparse.coo_matrix(mat+vec.A)
Out[388]:
<10x10 sparse matrix of type '<class 'numpy.float64'>'
with 38 stored elements in COOrdinate format>
In [389]: timeit sparse.coo_matrix(mat+vec.A)
1000 loops, best of 3: 716 µs per loop
In [390]: timeit sparse.coo_matrix(mat.A+vec.A)
1000 loops, best of 3: 338 µs per loop
一个可以节省临时密集矩阵空间的版本,同时运行:
In [393]: timeit temp=mat.A; temp+=vec.A; sparse.coo_matrix(temp)
1000 loops, best of 3: 334 µs per loop
所以密集版本比我的稀疏版本好 5-7 倍。
对于非常大的
mat
,内存问题可能会影响密集性能,但迭代稀疏解决方案也不会发光。我可以通过更有效地索引
addvec2
来从 Ml
中榨取更多性能。 Ml.data[3],Ml.rows[3]
比 Ml[3,:]
或 Ml[:,3]
快得多。def addvec3(mat,vec):
Mtl=mat.T.tolil()
vec=vec.tocoo()
n = mat.shape[0]
for i,v in zip(vec.col, vec.data):
t = np.zeros((n,))+v
t[Mtl.rows[i]] += Mtl.data[i]
t = sparse.coo_matrix(t)
Mtl.rows[i] = t.col
Mtl.data[i] = t.data
return Mtl.T
In [468]: timeit addvec3(mat,vec)
1000 loops, best of 3: 1.8 ms per loop
适度的改进,但没有我希望的那么大。再挤一点:
def addvec3(mat,vec):
Mtl = mat.T.tolil()
vec = vec.tocoo();
t0 = np.zeros((mat.shape[0],))
r0 = np.arange(mat.shape[0])
for i,v in zip(vec.col, vec.data):
t = t0+v
t[Mtl.rows[i]] += Mtl.data[i]
Mtl.rows[i] = r0
Mtl.data[i] = t
return Mtl.T
In [531]: timeit mm=addvec3(mat,vec)
1000 loops, best of 3: 1.37 ms per loop
关于python - 具有广播的稀疏 scipy 矩阵向量的元素相加,我们在Stack Overflow上找到一个类似的问题:https://stackoverflow.com/questions/30741461/