使用多元回归最小二乘法拟合cos(x)（二次多项式）时损失偏高不下降，求排查

嘿，我来帮你排查下问题所在！你的代码里有几个关键问题导致损失没下降甚至居高不下，咱们一个个拆解：

首先最核心的问题是：你的循环里完全没有更新模型参数theta的步骤——你只计算了预测值和损失，但没有根据损失的梯度去调整a、b、c对应的theta[0]、theta[1]、theta[2]，所以theta一直是初始的随机值，损失当然不会有任何变化，甚至因为初始值太差而损失极高。

除此之外，还有几个细节问题会影响拟合效果，我整理成了以下几点：

1. 缺失梯度更新逻辑 #

梯度下降的核心是通过计算损失对每个参数的偏导数（梯度），然后沿着梯度反方向更新参数。对于你的二次模型y_hat = a*x² + b*x + c，损失函数是均方误差cost = (1/(2m))*sum((y_hat - y)^2)（乘1/2是为了求导后抵消系数，不影响更新方向），各参数的梯度公式如下：

• 对a（theta[0]）的梯度：(1/m) * sum( (y_hat - y) * x² )
• 对b（theta[1]）的梯度：(1/m) * sum( (y_hat - y) * x )
• 对c（theta[2]）的梯度：(1/m) * sum( y_hat - y )

你需要在循环中计算这些梯度，然后用theta = theta - learning_rate * gradients来更新参数（注意选择合适的学习率，比如0.1或0.01）。

2. 特征尺度不统一导致收敛困难 #

你的特征x²的范围是[0,4]，x的范围是[-2,2]，常数项是1，不同特征的尺度差异很大，这会让梯度下降的步长对不同参数的影响不均衡，导致收敛极慢甚至无法收敛。建议对特征做标准化处理，比如将每个特征缩放到均值为0、方差为1的范围：

x_scaled = (x - x.mean()) / x.std()
x_squared_scaled = (x**2 - (x**2).mean()) / (x**2).std()

3. 低效的预测计算（可选优化） #

你用列表推导式计算y_hat的方式既低效又容易出错，建议把特征整理成矩阵形式，用矩阵乘法来计算预测值，代码更简洁且运行更快：

# 构造特征矩阵：每一行是 [x², x, 1]
X = np.column_stack([x**2, x, np.ones_like(x)])
y_hat = X @ theta  # 矩阵乘法，等价于逐个计算每个x的预测值

修正后的完整代码示例 #

import numpy as np
import matplotlib.pyplot as plt

# 生成数据
x = np.linspace(-2, 2, 100)
y = np.cos(x)

# 特征标准化（关键步骤）
x_scaled = (x - x.mean()) / x.std()
x_squared_scaled = (x**2 - (x**2).mean()) / (x**2).std()
# 构造标准化后的特征矩阵
X = np.column_stack([x_squared_scaled, x_scaled, np.ones_like(x)])
y = y.reshape(-1, 1)  # 转为列向量，和theta形状匹配

# 初始化参数
theta = np.random.random((3, 1)) * 0.1  # 小范围初始化，避免初始损失过大
m = len(y)
learning_rate = 0.1
epochs = 10000

# 梯度下降循环
for i in range(epochs):
    # 计算预测值
    y_hat = X @ theta
    # 计算损失
    cost = np.sum((y_hat - y) ** 2) / (2 * m)
    # 每1000轮打印一次损失，观察变化
    if i % 1000 == 0:
        print(f"Epoch {i}, Cost: {cost:.4f}")
    # 计算梯度
    gradients = (1/m) * X.T @ (y_hat - y)
    # 更新参数
    theta = theta - learning_rate * gradients

# 还原参数到原始特征尺度（因为我们标准化了特征）
a = theta[0] / (x**2).std()
b = theta[1] / x.std()
c = theta[2] - theta[0]*(x**2).mean()/(x**2).std() - theta[1]*x.mean()/x.std()

# 绘制拟合结果
y_pred = a * x**2 + b * x + c
plt.scatter(x, y, label='cos(x)')
plt.plot(x, y_pred, color='red', label=f'Fit: {a:.4f}x² + {b:.4f}x + {c:.4f}')
plt.legend()
plt.show()

运行这段代码你会看到损失随着迭代逐步下降，最终得到一个不错的二次拟合曲线。

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