如何在Python中简便计算回归系数(X'X)^{-1}X'Y？

嘿，这个问题问得好！确实，拆分三次np.dot加求逆的写法显得有点啰嗦，Python里有两种更简洁的方式——既能贴近你熟悉的MATLAB风格，还能兼顾数值稳定性：

1. 完全贴近数学公式的链式写法 #

Numpy从1.10版本开始支持@运算符作为矩阵乘法的简写，配合X.T表示转置，你可以直接把原始数学公式翻译成一行代码，和MATLAB的写法几乎一模一样：

import numpy as np

# 假设X是(100,3)的数组，Y是(100,1)的数组
betahat = np.linalg.inv(X.T @ X) @ X.T @ Y

这里X.T对应MATLAB里的X'，@对应矩阵乘法*，完美还原了公式结构，看起来清爽多了。

2. 更推荐的数值稳定方案：最小二乘求解 #

不过要提一句：直接计算(X'X)^{-1}在数值上可能存在风险（比如当X的列存在多重共线性时，求逆结果会不稳定）。Numpy提供了专门的最小二乘求解函数np.linalg.lstsq，它会用更稳定的QR分解或SVD算法计算回归系数，写法同样简洁：

# 调用最小二乘，只取第一个返回值就是回归系数
betahat = np.linalg.lstsq(X, Y, rcond=None)[0]

这个方法还能额外返回残差、矩阵秩、奇异值等信息，实际工程场景中更常用。

验证两种方法的一致性 #

你可以用一组测试数据确认结果一致：

# 构造随机测试数据
X = np.random.rand(100, 3)
Y = np.random.rand(100, 1)

# 方法1的结果
beta1 = np.linalg.inv(X.T @ X) @ X.T @ Y
# 方法2的结果
beta2 = np.linalg.lstsq(X, Y, rcond=None)[0]

# 检查浮点误差范围内是否一致
print(np.allclose(beta1, beta2))  # 输出True

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