Numpy:将矩阵与向量数组相乘

6

我在学习numpy方面遇到了一些困难。我的最终目标是通过一个矩阵变换得到简单的矢量箭头图。我已经多次听说要使用数组来表示矩阵,这很合理。同时,我已经创建了网格以表示x和y坐标。

X,Y = np.meshgrid( np.arange(0,10,2),np.arange(0,10,1) )
a = np.array([[1,0],[0,1.1]])

但是,即使经过2个多小时的谷歌搜索和尝试,我仍然无法从矩阵乘法和每个向量中获取结果向量。我知道quiver需要输入分量长度,因此进入quiver函数的结果向量应该类似于x分量的np.dot(a, [X[i,j], Y[i,j]]) - X[i,j],其中i和j在范围内迭代。

当然我可以用循环编程,但是numpy有太多内置工具来方便这些向量化的事情,以至于我相信一定有更好的方法。

编辑:好吧,这里是循环版本。

import numpy as np
import matplotlib.pyplot as plt

plt.figure(figsize=(10,10))

n=10
X,Y = np.meshgrid( np.arange(-5,5),np.arange(-5,5) )
print("val test", X[5,3])
a = np.array([[0.5,0],[0,1.3]])
U = np.zeros((n,n))
V = np.zeros((n,n))
for i in range(10):
    for j in range(10):
        product = np.dot(a, [X[i,j], Y[i,j]]) #matrix with vector
        U[i,j] = product[0]-X[i,j]  # have to substract the position since quiver accepts magnitudes
        V[i,j] = product[1]-Y[i,j]

Q = plt.quiver( X,Y, U, V)

你期望的结果是什么样子?a.shape = (2,2), X.shape = (10,5), Y.shape = (10,5)。我不明白你的意思... - jkalden
你能否展示一个简单的循环版本(以便更容易理解你想要做什么)? - Lee
完成。实际上这只是一个矩阵在向量场中的表现,结果通过箭头可视化。 - Basti
我在这里发布了一个关于如何进行点积的答案(http://stackoverflow.com/a/22081723/553404),可能对你有帮助。 - YXD
3个回答

6
您可以使用NumPy广播手动执行矩阵乘法,如下所示:
import numpy as np
import matplotlib.pyplot as plt

X,Y = np.meshgrid(np.arange(-5,5), np.arange(-5,5))
a = np.array([[0.5, 0], [0, 1.3]])

U = (a[0,0] - 1)*X + a[0,1]*Y
V = a[1,0]*X + (a[1,1] - 1)*Y

Q = plt.quiver(X, Y, U, V)

如果你想使用 np.dot,你需要将XY 数组扁平化,并将它们合并到适当的形状中,如下所示:

import numpy as np
import matplotlib.pyplot as plt

X,Y = np.meshgrid(np.arange(-5,5), np.arange(-5,5))
a = np.array([[0.5, 0], [0, 1.3]])

U,V = np.dot(a-np.eye(2), [X.ravel(), Y.ravel()])

Q = plt.quiver(X, Y, U, V)

嗯,我很感激。但这两种情况都非常有限,对吧?第一种是手写矩阵乘法。那相当的不专业,而且对于任何更高维度的情况都行不通。第二种是依赖于零。 - Basti
两种情况都是针对2D的。对于3D,应该是这样的:U,V,W = np.dot(a-np.eye(3), [X.ravel(), Y.ravel(), Z.ravel()]) - Ondro

5

根据文档所述,对于多维数据,np.mul(或@)的操作方式如下:

For N dimensions it is a sum product over the last axis of a and the second-to-last of b:

dot(a, b)[i,j,k,m] = sum(a[i,j,:] * b[k,:,m])
这不是我们想要的。然而,有一些简单的替代方案没有涉及到展平-还原或手动矩阵乘法:np.tensordotnp.einsum
第一个例子直接来自文档
要使用最左边的指数而不是右边来进行矩阵乘积,您可以执行np.einsum('ij...,jk...->ik...', a, b)
import numpy as np
import matplotlib.pyplot as plt

X,Y = np.meshgrid(np.arange(-5,5), np.arange(-5,5))
a = np.array([[0.5, 0], [0, 1.3]])

U, V = np.einsum('ij...,jk...->ik...', a - np.eye(2), np.array([X, Y]))
Q = plt.quiver(X, Y, U, V)

第二个方法是使用简单的np.tensordot函数。我们只需指定第一个参数对第二个轴(列)求和,对第二个参数对第一个轴(行)求和即可。
import numpy as np
import matplotlib.pyplot as plt

X,Y = np.meshgrid(np.arange(-5,5), np.arange(-5,5))
a = np.array([[0.5, 0], [0, 1.3]])

U, V = np.tensordot(a - np.eye(2), np.array([X, Y]), axes=(1, 0))
plt.quiver(X, Y, U, V)

0

我曾经为同样的问题苦苦挣扎,最终使用了numpy.matrix类。请考虑下面的例子。

import numpy as np

>>> transformation_matrix = np.array([(1, 0, 0, 1),
...                                   (0, 1, 0, 0),
...                                   (0, 0, 1, 0),
...                                   (0, 0, 0, 1)])
>>> coordinates = np.array([(0,0,0),
...                         (1,0,0)])
>>> coordinates = np.hstack((coordinates, np.ones((len(coordinates), 1))))
>>> coordinates
array([[ 0.,  0.,  0.,  0.],
       [ 1.,  0.,  0.,  0.]])

在这种情况下,numpy.matrix类会有所帮助。通过将坐标转置为列向量并使用numpy.matrix类的指定矩阵乘法重载,以下代码可得到预期结果。
>>> (np.asmatrix(transformation_matrix) * np.asmatrix(coordinates).T).T
matrix([[ 1.,  0.,  0.,  1.],
        [ 2.,  0.,  0.,  1.]])

1
你为什么先将它们定义为数组,然后再转换为矩阵呢?你本可以这样做的:transformation_matrix = np.matrix([(1,0,0,1), (0,1,0,0), (0,0,1,0), (0,0,0,1)]) - R. Navega
2
@R.Navega 这是真的,只是因为 OP 在 ndarray 上操作。请注意,asmatrix() 不会复制数据:“与矩阵不同,如果输入已经是矩阵或 ndarray,则 asmatrix 不会进行复制。等效于 matrix(data, copy=False)。”https://docs.scipy.org/doc/numpy/reference/generated/numpy.asmatrix.html - user1556435

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接