Octave mrdivide代码转Numpy结果不一致问题求助

Octave的mrdivide操作（/）在执行x = b'/A时，本质是求解线性方程组 x * A = b'，而不是简单地计算b'乘以A的逆。直接使用np.dot(b.T, np.linalg.inv(A))会因为矩阵求逆的数值不稳定性（尤其是当A条件数较大时）导致结果与Octave不一致。

下面是两种更准确的Numpy替代方案：

1. 当A是方阵且非奇异时：使用np.linalg.solve #

x * A = b'可以转置为 A.T @ x.T = b，我们可以直接用solve求解这个线性方程组，这和Octave内部的求解逻辑更一致，数值稳定性更好：

import numpy as np

# 假设A是方阵，b是列向量
A = np.array([[1, 2], [3, 4]])
b = np.array([[5], [6]])

# 求解 A.T @ x_t = b，其中x_t = x.T
x_t = np.linalg.solve(A.T, b)
# 转置得到最终的x（行向量）
x = x_t.T

2. 当A不是方阵或接近奇异时：使用np.linalg.lstsq #

如果A是长方形矩阵，或者存在数值上的奇异性，Octave的mrdivide会返回最小二乘解。此时用lstsq（最小二乘求解）能得到一致的结果：

import numpy as np

A = np.array([[1, 2, 3], [4, 5, 6]])  # 非方阵示例
b = np.array([[7], [8]])

# 求解 A.T @ x_t = b 的最小二乘解
x_t, residuals, rank, s = np.linalg.lstsq(A.T, b, rcond=None)
x = x_t.T

为什么直接用inv会有差异？ #

矩阵求逆（np.linalg.inv）的数值稳定性远不如直接求解线性方程组。当A的条件数较大（即矩阵接近奇异）时，求逆操作会放大计算中的浮点误差，导致结果偏离Octave的计算结果。而Octave的mrdivide和Numpy的solve/lstsq都会使用更稳定的数值分解算法（如LU分解、QR分解）来求解，避免了直接求逆的误差问题。

你可以用自己的实际数据测试上述两种方法，结果应该会和Octave的b'/A高度一致。

内容的提问来源于stack exchange，提问作者Merakli