广播 lstsq（最小二乘）

Question

我有一堆 3x2 矩阵，比如说 777 个，以及大小为 3 的许多右侧。对于每个矩阵，我想知道最小二乘解决方案，所以我在做

import numpy

A = numpy.random.rand(3, 2, 777)
b = numpy.random.rand(3, 777)

for k in range(777):
    numpy.linalg.lstsq(A[..., k], b[..., k])

这可行，但速度很慢。我宁愿一次性计算出所有的解决方案，但是在

numpy.linalg.lstsq(A, b)

我来了

numpy.linalg.linalg.LinAlgError: 3-dimensional array given. Array must be two-dimensional

关于如何广播numpy.linalg.lstsq的任何提示？

Answer 1

如果A = U \Sigma V^T是A的奇异值分解，可以利用这一事实，

x = V \Sigma^+ U^T b

是Ax = b 的最小二乘解。 SVD is broadcasted in numpy。现在只需要稍微摆弄einsums 就可以了：

A = numpy.random.rand(7, 3, 2)
b = numpy.random.rand(7, 3)
for k in range(7):
    x, res, rank, sigma = numpy.linalg.lstsq(A[k], b[k])
    print(x)

print

u, s, v = numpy.linalg.svd(A, full_matrices=False)
uTb = numpy.einsum('ijk,ij->ik', u, b)
xx = numpy.einsum('ijk, ij->ik', v, uTb / s)
print(xx)