使用 numpy 线性代数计算二次响应

Question

我在 numpy 中操作矩阵时遇到问题。我有一个包含 N 个输入的向量 theta，我想使用二次模型 theta^T A theta + theta^T b + c 来计算标量响应，其中 ^T 表示转置，A 是 NxN 方阵，b 是N维向量，c是标量。当 theta 是一个 (NxM) 矩阵时，这意味着我有 M 个 theta 值要传播，我必须计算 theta^T A theta 以产生一个 M 维矩阵。在索引符号中，计算是 theta_{mj} A_{ji} theta_{im}，其中 m 不求和。

如果我只有一组 theta 值，numpy 的线性代数将按预期工作（此处，N=10 和 M=1）：

In [1]: import numpy as np                                            

In [2]: theta = np.ones(10)                                           

In [3]: A = np.ones((10, 10))                                         

In [4]: b = np.ones(10)                                               

In [5]: theta.T.dot(A).dot(theta) + theta.T.dot(b) + 1                
Out[5]: 111.0

我认为 theta^T 操作 theta 可能与此处的点积不同。当我将 theta 设为 NxM 矩阵时，我不明白有什么根本不同。我认为额外的维度自然会贯穿这段代码，就像它对 b 和 c 项所做的那样。

如何通过 theta^T A theta 运算使 numpy return 成为 M 维数组？

我只能return方阵。不幸的是，点积函数将此操作视为矩阵乘法（此处，N=10 和 M=5）：

In [6]: theta = np.ones((10, 5))
#       theta.T.dot(A).dot(theta) is equivalient to:
#             
#                 (M x N)           (N x N)  (N x M)
In [7]: np.matmul(theta.T, np.matmul(A,     theta))                      
Out[7]: 
array([[100., 100., 100., 100., 100.],
       [100., 100., 100., 100., 100.],
       [100., 100., 100., 100., 100.],
       [100., 100., 100., 100., 100.],
       [100., 100., 100., 100., 100.]])

相比之下，b和c项自然带有额外的theta项，提供我想要的M维输出：

In [8]: theta.T.dot(b) + 1                                           
Out[8]: array([11., 11., 11., 11., 11.])

Answer 1

两种可能性：

N,M = 10,5
A = np.random.randint(0,10,(N,N))
theta = np.random.randint(0,10,(N,M))
b = np.random.randint(0,10,N)

1) 使用 >2D 操作数被 matmul 视为堆栈的事实：

(theta.T[:,None]@A@theta.T[...,None])[...,0,0] + b@theta + 1
# array([ 8188, 14837,  7697,  9719,  7262])

或

2) 使用einsum

np.einsum("ik,ij,jk->k",theta,A,theta) + b@theta + 1
# array([ 8188, 14837,  7697,  9719,  7262])

与一对一评估进行比较以进行验证：

[t@A@t + b@t + 1 for t in theta.T]
# [8188, 14837, 7697, 9719, 7262]

使用 numpy 线性代数计算二次响应

Computing a quadratic response using numpy linear algebra

python

arrays

numpy

matrix

linear-algebra