我有三个numpy数组:
X: a 3073 x 49000矩阵W:一个10 × 3073矩阵y:一个49000 × 1向量
y包含0到9之间的值，每个值代表W中的一行。
我想将X的第一列添加到W中由y中的第一个元素给定的行。例如，如果y中的第一个元素是3，则将X的第一列添加到W的第四行。然后将X的第二列添加到y中第二个元素给定的W中的行，以此类推，直到X的所有列都被添加到y指定的W中的行，即总共添加了49000行。
W[y] += X.T不适合我，因为这不会将多个向量添加到W中的一行。
请注意:我只寻找矢量化解决方案。即没有for循环。
EDIT:为了澄清，我将添加一个小矩阵大小的例子，改编自Salvador Dali下面的例子。
In [1]: import numpy as np
In [2]: a, b, c = 3, 4, 5
In [3]: np.random.seed(0)
In [4]: X = np.random.randint(10, size=(b,c))
In [5]: W = np.random.randint(10, size=(a,b))
In [6]: y = np.random.randint(a, size=(c,1))
In [7]: X
Out[7]:
array([[5, 0, 3, 3, 7],
[9, 3, 5, 2, 4],
[7, 6, 8, 8, 1],
[6, 7, 7, 8, 1]])
In [8]: W
Out[8]:
array([[5, 9, 8, 9],
[4, 3, 0, 3],
[5, 0, 2, 3]])
In [9]: y
Out[9]:
array([[0],
[1],
[1],
[2],
[0]])
In [10]: W[y.ravel()] + X.T
Out[10]:
array([[10, 18, 15, 15],
[ 4, 6, 6, 10],
[ 7, 8, 8, 10],
[ 8, 2, 10, 11],
[12, 13, 9, 10]])
In [11]: W[y.ravel()] = W[y.ravel()] + X.T
In [12]: W
Out[12]:
array([[12, 13, 9, 10],
[ 7, 8, 8, 10],
[ 8, 2, 10, 11]])
问题是将 X中的第0列和第4列都添加到W中的第0行，以及将X中的第1列和第2列都添加到W中的第1行。
期望的结果是:
W = [[17, 22, 16, 16],
[ 7, 11, 14, 17],
[ 8, 2, 10, 11]]
首先是直接循环解决方案作为参考:
In [65]: for i,j in enumerate(y):
W[j]+=X[:,i]
....:
In [66]: W
Out[66]:
array([[17, 22, 16, 16],
[ 7, 11, 14, 17],
[ 8, 2, 10, 11]])
add.at解决方案:
In [67]: W=W1.copy()
In [68]: np.add.at(W,(y.ravel()),X.T)
In [69]: W
Out[69]:
array([[17, 22, 16, 16],
[ 7, 11, 14, 17],
[ 8, 2, 10, 11]])
add.at执行非缓冲计算，绕过阻止W[y.ravel()] += X.T工作的缓冲。它仍然是迭代的，但是循环已经移动到编译的代码中。这不是真正的矢量化，因为应用程序的顺序很重要。X.T的一行的加法取决于前一行的结果。
https://stackoverflow.com/a/20811014/901925是几年前我对一个类似问题的答案(对于1d数组)。
但是当处理大型数组时:
X: a 3073 x 49000 matrix
W: a 10 x 3073 matrix
y: a 49000 x 1 vector
可能会遇到速度问题。请注意，W[y.ravel()]与X.T的大小相同(为什么选择这些需要转置的大小?)这是拷贝，不是视图。所以已经有时间惩罚了。
bincount在之前的问题中已经建议过了，我认为它更快。使用索引数组使for循环更快(包括bincount和add.at解决方案)
迭代较小的3073维度也可能具有速度优势。或者像Divakar演示的那样，更好的是在尺寸为10的维度上。
对于小测试用例a,b,c=3,4,5, add.at解决方案最快，其次是Divakar's bincount和einseum。对于较大的a,b,c=10,1000,20000, add.at变得非常慢，而bincount是最快的。
相关SO答案
https://stackoverflow.com/a/28205888/901925(注意bincount要求完全覆盖y)。
https://stackoverflow.com/a/30041823/901925(这里Divakar再次显示了bincount的规则!)
矢量化方法
方法# 1
基于this answer，这里是一个使用np.bincount -
N = y.max()+1
id = y.ravel() + np.arange(X.shape[0])[:,None]*N
W[:N] += np.bincount(id.ravel(), weights=X.ravel()).reshape(-1,N).T
方法# 2
您可以很好地利用boolean indexing和np.einsum以简洁的矢量化方式完成工作-
N = y.max()+1
W[:N] += np.einsum('ijk,lk->il',(np.arange(N)[:,None,None] == y.ravel()),X)
呆头呆脑的方法
方法# 3
由于您正在选择和添加每个唯一y的X中的大量列，因此在性能方面运行一个循环，complexity等于这种唯一y's的数量，这似乎在max等于W中的行数，而在您的情况下只是10。因此，循环只有10次迭代，还不错!以下是实现这些愿望的实现-
for k in range(W.shape[0]):
W[k] += X[:,(y==k).ravel()].sum(1)
方法# 4
您可以引入np.einsum来执行按列求和，最终输出如下-
for k in range(W.shape[0]):
W[k] += np.einsum('ij->i',X[:,(y==k).ravel()])
这会达到你想要的效果:X + W[y.ravel()].T
要看看这是否真的起作用，这里有一个可重复的例子:
import numpy as np
np.random.seed(0)
a, b, c = 3, 5, 4 # you can use your 3073, 49000, 10 later
X = np.random.rand(a, b)
W = np.random.rand(c, a)
y = np.random.randint(c, size=(b, 1))
现在你的矩阵是:
[[ 0.0871293 0.0202184 0.83261985]
[ 0.77815675 0.87001215 0.97861834]
[ 0.79915856 0.46147936 0.78052918]
[ 0.11827443 0.63992102 0.14335329]]
[[3]
[0]
[3]
[2]
[0]]
[[ 0.5488135 0.71518937 0.60276338 0.54488318 0.4236548 ]
[ 0.64589411 0.43758721 0.891773 0.96366276 0.38344152]
[ 0.79172504 0.52889492 0.56804456 0.92559664 0.07103606]]
和W[y.ravel()]给你"W由y中的第一个元素给出"。通过对它的转置，你将得到一个准备添加到X的矩阵:
[[ 0.11827443 0.0871293 0.11827443 0.79915856 0.0871293 ]
[ 0.63992102 0.0202184 0.63992102 0.46147936 0.0202184 ]
[ 0.14335329 0.83261985 0.14335329 0.78052918 0.83261985]]
虽然我不能说这是非常python的，但这是一个解决方案(我认为):
for column in range(x.shape[1]):
w[y[column]] = x[:,column].T