计算矩阵乘法的子集

当我有两个非稀疏矩阵A和B时，有什么方法可以有效地计算C=A.T.dot(B)当我只想要C的元素的子集？我有指定here CSC格式存储所需的指数C。计算矩阵乘法的子集

2012-12-05 bluecat

尽管这是编程社区，但您可能需要在Math SO上检查这一点。 – jdotjdot

@jdotjdot：对我来说，这纯粹是一个编程问题，几乎没有数学内容。 – NPE

我检查mathexchange，同时为帖子制作标签，但它没有像scipy和numpy或甚至稀疏似乎相关的。 – bluecat

忽略CSC业务，也许回答一个比你问的更简单的问题。给出C索引值的元组列表，我将如何计算C的元素子集。

既然你正在评估C = ATdot（B）你是B.所以列乘以列，

for i, j in indexList: 
    C[i, j] = np.dot(A[:,i], B[:,j])

我猜，是不是你在找什么，但我发现简单的答案有时有助于澄清问题。

来源

2012-12-07 16:50:32 GaryBishop

是的，这或多或少是我现在使用的方法，但不幸的是对np.dot（）的循环调用非常缓慢。 – bluecat

如果您事先知道您需要哪些C部分以及这些部分中的某些部分是连续的和矩形区域*，那么您可以使用与分区矩阵乘法（1）或块矩阵乘法相关的矩阵代数规则2）来加速其中一些计算。因此，例如，您可以使用@GaryBishop的相同基本逻辑，但不是具有“i”和“j”元素的列表，而是具有i_start，i_end和j_start的四元组的列表（或数组），j_end定义C的子矩阵，那么你可以使用这些指数（尽管在这些链接中建立的规则）来找出需要为所需的C块求解的A和B的子矩阵。

对于简单的例子，假设你只想要C的中间块，所以我们将C分割成C1，C2和C3，我们关心的只有C2。如果A^{T}同样被划分成三组A1，A2，A3，那么C2 = A2 * B。这个想法推广到任何形状的矩形，它只需要计算A和B的不同分区。这个想法是一样的。

- 这很平凡，但如果区域大于单个元素，则只会节省时间。

来源

2012-12-21 12:54:36 BKay

这和我想的一样。我怀疑构成C2和A2的行甚至可能是非连续的（不连续的切片） – heltonbiker

我认为你是正确的，它们不必是连续的，尽管代数可能有点棘手。非连续的A2可能被认为是更精细分区的元素，例如， A1（输出），A2a（输入），A2b（输出），A2C（输入），A3（输出），但仔细索引并很好地理解您关心它的哪些元素可能会更快地使用非连续子集直接使用一块连续的分区进行处理。 – BKay

相反在使用Python（GaryBishop的答案）坐标迭代的，你可以有numpy的做循环，从而构成实质性加速（以下计时）：

def sparse_mult(a, b, coords) : 
    rows, cols = zip(*coords) 
    rows, r_idx = np.unique(rows, return_inverse=True) 
    cols, c_idx = np.unique(cols, return_inverse=True) 
    C = np.dot(a[rows, :], b[:, cols]) 
    return C[r_idx, c_idx] 

>>> A = np.arange(12).reshape(3, 4) 
>>> B = np.arange(15).reshape(3, 5) 
>>> np.dot(A.T, B) 
array([[100, 112, 124, 136, 148], 
     [115, 130, 145, 160, 175], 
     [130, 148, 166, 184, 202], 
     [145, 166, 187, 208, 229]]) 
>>> sparse_mult(A.T, B, [(0, 0), (1, 2), (2, 4), (3, 3)]) 
array([100, 145, 202, 208])

sparse_mult返回扁平您提供的坐标值的数组作为元组列表。我不是很熟悉，稀疏矩阵格式，所以我不知道如何从上面的数据定义CSC，但以下工作：

>>> coords = [(0, 0), (1, 2), (2, 4), (3, 3)] 
>>> sparse.coo_matrix((sparse_mult(A.T, B, coords), zip(*coords))).tocsc() 
<4x5 sparse matrix of type '<type 'numpy.int32'>' 
    with 4 stored elements in Compressed Sparse Column format>

这是各种替代方案的时机：

>>> import timeit 
>>> a = np.random.rand(2000, 3000) 
>>> b = np.random.rand(3000, 5000) 
>>> timeit.timeit('np.dot(a,b)[[0, 0, 1999, 1999], [0, 4999, 0, 4999]]', 'from __main__ import np, a, b', number=1) 
5.848562187263569 
>>> timeit.timeit('sparse_mult(a, b, [(0, 0), (0, 4999), (1999, 0), (1999, 4999)])', 'from __main__ import np, a, b, sparse_mult', number=1) 
0.0018596387374678613 
>>> np.dot(a,b)[[0, 0, 1999, 1999], [0, 4999, 0, 4999]] 
array([ 758.76351111, 750.32613815, 751.4614542 , 758.8989648 ]) 
>>> sparse_mult(a, b, [(0, 0), (0, 4999), (1999, 0), (1999, 4999)]) 
array([ 758.76351111, 750.32613815, 751.4614542 , 758.8989648 ])

来源

2013-01-07 23:23:15 Jaime

我想弄明白这是如何工作，但我无法得到它的工作。为什么np.dot（A.T，B）是5乘4矩阵？它不应该是4x2吗？使用你的例子，我在函数的行“C = np.dot（a [rows，：]，b [：，cols]）”，“index error：index（2）out of range（0 <索引<1）“ – BKay

@BKay B在示例中的大小不正确，将其替换为'B = np.arange（15）.reshape（3，5）'（编辑提交） – Maxim

计算矩阵乘法的子集

回答

相关问题