拿到一份带噪声的模拟数据点在平面图上呈现出明显的曲线趋势用普通线性回归直接拟合效果很差。我当时的处理方式是转成多项式回归建模并同时用最小二乘法闭式解和梯度下降法迭代解分别求解最后对比两者的系数、收敛速度和评估指标。这篇内容就把整个思路和代码完整拆开覆盖从建模、特征构造、两种求解方式到过拟合控制和常见坑的完整过程。适合正在学习回归算法、想把数学推导和Python代码对上的朋友也适合做特征工程时想弄明白“什么时候用闭式解、什么时候用迭代解”的开发者。全程只用numpy从零实现不依赖机器学习框架理解起来会更通透。1. 多项式回归的建模思路为什么两种方法能解同一个模型1.1 从“直线不够用”说起普通线性回归假设的是 y 与 x 之间呈直线关系也就是 y b0 b1·x。但现实里很多关系并不是线性的比如销量随广告投入增长会出现边际递减体温随时间变化会有波峰和波谷。这时候如果硬套一条直线就会出现系统性偏差中间预测偏高、两端预测偏低之类的。多项式回归的思路非常简单把原始特征 x 扩展成 x²、x³ 甚至更高次幂然后继续用线性回归的方式拟合。模型变成y b0 b1·x b2·x² b3·x³ ...这里的关键点是虽然对 x 是非线性的但对参数 b0、b1、b2 来说依然是线性的。这个性质非常重要它意味着我们不需要发明一套全新的求解算法而是可以直接复用线性回归里现成的数学工具——最小二乘法和梯度下降法都能用来求这组最优参数。1.2 同一个损失函数两条求解路径两种方法的起点完全一样都是最小化残差平方和L Σ(yi - f(xi))²写成矩阵形式更简洁L ||y - X·w||²其中 X 是设计矩阵每一行对应一个样本列是 [1, x, x², x³, ...]w 是待求参数向量。所谓“求解多项式回归模型”本质上就是找一组 w让这个 L 最小。最小二乘法走的是解析路线直接对 L 求导并令导数等于零一步算出最优解。梯度下降法走的是迭代路线先随机给一组 w然后反复沿着梯度方向更新参数慢慢逼近最优解。两条路殊途同归但是各自的适用场景、计算开销和调参要求差异非常大。1.3 求解路线选型的基本原则我自己的经验法则很简单数据量小比如几千条且特征维度不高时闭式解一步到位省心省力数据量大几十万条以上或特征维度很高时闭式解需要计算 (XᵀX)⁻¹矩阵求逆的复杂度高达 O(d³)而且内存里存不下整个矩阵这时候梯度下降的优势就出来了——它每轮只需要算 Xᵀ 乘以残差还可以分批次处理。另外还有一个容易被忽略的点如果用了正则化闭式解公式会稍微变一下但梯度下降只需要在梯度里加一个惩罚项就行改动成本很低。所以很多实际项目里即便数据量不大团队也会默认用梯度下降因为后续扩展更灵活。2. 最小二乘法从正规方程到Python实现2.1 正规方程的推导与矩阵视角最小二乘法的核心是正规方程。我们从矩阵形式出发L ||y - X·w||² yᵀy - 2wᵀXᵀy wᵀXᵀXw对 w 求导∂L/∂w -2Xᵀy 2XᵀXw令导数等于零得到XᵀXw Xᵀy所以w (XᵀX)⁻¹Xᵀy这个公式就是很多教材里写的“闭式解”。它的几何意义也很直观预测值 ŷ Xw 是 y 在 X 列空间上的投影残差向量 y - ŷ 与 X 的所有列都正交。正规方程本质上就是在找这个正交投影。聊一点实际经验代码里我一般不会真的用 np.linalg.inv(X.T X) 去算逆矩阵而是用 np.linalg.solve(X.T X, X.T y) 或者 np.linalg.lstsq(X, y, rcondNone)。前者是解线性方程组数值稳定性比先求逆再乘好很多后者直接基于奇异值分解能自动处理 XᵀX 接近奇异的情况。原理上都是在求同一个解但工程上的稳定性差别真的挺大。2.2 最小二乘法的Python实现步骤我用一组带噪声的模拟数据来演示。数据由真实规律 y 2 0.5x - 0.3x² 加上高斯噪声生成这样我们心里有底知道求出来的系数应该接近什么。import numpy as np import matplotlib.pyplot as plt rng np.random.default_rng(42) x np.linspace(0, 10, 300) y 2.0 0.5 * x - 0.3 * x ** 2 rng.normal(0, 1.5, sizex.size) def design_matrix(x, degree): return np.vander(x, degree 1, increasingTrue) degree 2 X design_matrix(x, degree) # 方式一正规方程 solve A X.T X b X.T y coef np.linalg.solve(A, b) # 方式二直接调用 lstsq更稳定 coef2, residuals, rank, singular_values np.linalg.lstsq(X, y, rcondNone) y_pred X coef ss_res np.sum((y - y_pred) ** 2) ss_tot np.sum((y - np.mean(y)) ** 2) r2 1 - ss_res / ss_tot print(solve得到的系数:, coef) print(lstsq得到的系数:, coef2) print(R²:, r2)两个方式求出来的系数基本一致degree2 时 R² 一般能在 0.9 以上。这里有个小细节值得注意np.vander 是构造多项式设计矩阵的一个快捷函数increasingTrue 表示按 x 的幂次从小到大排列也就是 [1, x, x², ...]跟模型表达式对应得上。预测新样本时必须用同样的特征变换。比如新样本 x_new 5.5就要构造 [1, 5.5, 5.5²]再与系数做点积。这个“训练和预测保持一致”的原则很多人会漏一旦漏了结果会错得毫无意义。2.3 数值稳定性特征缩放与条件数最小二乘看起来简单真正跑起来坑也不少。最典型的问题出在 XᵀX 的条件数上。当 degree 逐渐增大比如 d 10X 里的列包括 x 的 10 次方这一列和其他列的数值尺度差异会大到离谱。x 本来在 0~10 之间x^10 直接到 10 的 10 次方量级另一列却是常数 1。矩阵 XᵀX 会变得高度病态条件数巨大求得的 w 对数据的微小扰动极其敏感稍微换一组噪声结果就飘了。解决办法是在构造多项式特征之前先做标准化z (x - mean(x)) / std(x)标准化之后 x 的幂次依然存在但每一列的量级被压缩到可接受的范围内条件数会降好几个数量级。我自己测过degree15 时不做标准化np.linalg.solve 跑完的结果基本是废的做完标准化闭式解依然能求出一个结构合理的系数向量。这个习惯建议从一开始就养成后面讲到梯度下降时它更是决定成败的关键。3. 梯度下降法从梯度推导到从零手写3.1 损失函数对参数的梯度推导梯度下降的本质是反复执行“算出梯度方向沿着它迈一小步”。所以第一步就是把损失函数对每个参数的导数搞清楚。还是用均方误差版本的损失函数L (1/n)||y - X·w||²这里除以 n 是为了让损失和梯度不随样本量线性增长这样学习率可以设得更稳定。展开L (1/n)(yᵀy - 2wᵀXᵀy wᵀXᵀXw)对 w 求导注意 wᵀXᵀXw 的导数是 2XᵀXw所以∂L/∂w (-2Xᵀy 2XᵀXw) / n (2/n)·Xᵀ(Xw - y)这个式子非常干净残差向量 (Xw - y) 先算出来然后左乘 Xᵀ再乘以 2/n。代码实现里这就是两行矩阵运算的事。我建议你手推一次这个梯度花不了几分钟但对后续调参和诊断模型会非常有帮助。3.2 三种更新策略批量、随机、小批量梯度下降不是一个算法而是一族算法。差别在于每一步用多少样本去计算梯度。批量梯度下降Batch GD每轮用全部 n 个样本计算梯度方向稳定损失曲线平滑下降但每轮的计算成本高内存压力大。随机梯度下降SGD每步只随机抽一个样本算梯度方向噪声非常大损失曲线会震荡但偶尔能跳出局部极小点而且计算量小适合非常大的数据集。小批量梯度下降Mini-batch GD每次用一小堆样本计算梯度通常 32、64、128 条兼顾稳定性和计算效率是实际项目中最常用的。这里补一个容易误解的点我做拟合实验时经常有人问“为什么梯度下降求出的系数和最小二乘不一样”。其实如果损失函数是凸函数理论上两者最终会收敛到同一个全局最优解。不一样的原因要么是还没收敛完要么是学习率不合适导致在最优解附近反复震荡要么是特征没缩放导致收敛路径曲折、早早被停止条件截断。3.3 学习率与停止准则的选法学习率是梯度下降最重要的超参数没有之一。设得太大参数会越过最优点损失可能震荡甚至发散设得太小收敛慢到让人怀疑人生。我的做法是先跑 100 轮打印每轮损失观察曲线形状损失持续下降学习率太小可以适当加大。损失先下降后反弹学习率偏大需要调小。损失直接变成 NaN学习率过大或者特征尺度差异过大。除了固定学习率还可以用衰减策略每迭代一定轮次就让学习率乘以一个小于 1 的系数比如 0.9。这样前期步子大、快速接近后期步子小、精细收敛。停止准则一般用三个条件组合达到最大迭代次数、损失变化量小于阈值、梯度范数小于阈值。注意阈值不要设太小否则会在最优点附近空转很多轮。我自己常用 1e-6 作为损失变化的容忍度再配合几千到几万轮的上限。3.4 Python实现从零手写批量梯度下降先对 x 做标准化再构造多项式特征这是让梯度下降正常工作的大前提。标准化后的公式预测也很容易还原预测样本先做同样的变换。x_mean x.mean() x_std x.std() xs (x - x_mean) / x_std X_poly design_matrix(xs, degree) n len(y) theta np.zeros(degree 1) alpha 0.1 max_iter 5000 tol 1e-8 loss_history [] for i in range(max_iter): grad (2.0 / n) * X_poly.T (X_poly theta - y) theta - alpha * grad loss np.mean((X_poly theta - y) ** 2) loss_history.append(loss) if len(loss_history) 1 and abs(loss_history[-2] - loss) tol: break print(迭代轮数:, len(loss_history)) print(theta:, theta)在这组模拟数据上degree2、alpha0.1 通常 100 轮左右就能收敛到和最小二乘几乎一样的系数。如果你画损失历史曲线会看到一条快速下降后趋于平稳的典型收敛曲线。有一点必须提醒这里 theta 对应的是标准化后的特征空间。给新样本预测时也要先按 x_mean、x_std 做变换再构造多项式特征。千万不要拿原始 x 直接构造特征和 theta 做点积那样结果一定不对。3.5 特征缩放为什么是梯度下降的关键前提没有特征缩放时x 和 x² 的取值范围相差巨大损失函数的等高线会拉成非常扁的椭圆。梯度下降在这样的地形上会走“之”字形每一步在某一方向上跨得很大在另一个方向上却几乎不动收敛速度慢到难以接受。标准化之后特征尺度接近等高线趋向圆形梯度方向直接指向中心收敛速度会快很多。我在这类实验里见过太多类似现象不标准化的时候跑 5000 轮还不收敛标准化之后 50 轮就已经到达同样的损失水平。所以只要见到梯度下降收敛异常优先检查特征缩放而不是先去调学习率。4. 两种方法对比与模型效果评估4.1 最小二乘与梯度下降的结果对比把两种方法放在同一份数据上对比是最直观的学习方式。我用一张表整理两者的核心差异对比项最小二乘闭式解小批量梯度下降核心计算解正规方程 XᵀXw Xᵀy迭代更新 w ← w - α·∇L时间复杂度约 O(nd² d³)每轮 O(nd)收敛前重复多轮内存占用需要构造并存储 XᵀX每次仅需一个批次调参需求几乎无学习率、批量大小、停止阈值特征缩放影响影响矩阵条件数与数值稳定性影响收敛速度强烈依赖大数据集表现差矩阵求逆是瓶颈好支持流式数据扩展正则化公式略有改动梯度增加一个惩罚项即可实际跑出来的结果很有意思数据量只有几百时最小二乘几乎瞬间完成梯度下降反而要先调学习率还要迭代几百轮——所以小数据集上用闭式解是明智的。但一旦数据量涨到百万级闭式解光是构建设计矩阵就可能把内存吃满而小批量梯度下降可以一直稳定运行。4.2 多项式阶数选择与过拟合控制多项式回归最让人头疼的选择是 degree 到底取多少。degree1 就是直线通常欠拟合degree 太高比如 15模型会把噪声也拟合进去曲线在数据点之间剧烈扭曲训练集 R² 能到 0.99但一到测试集就露馅。一个可靠的评估流程是先把数据集划分成训练集和测试集在训练集上拟合在测试集上计算 RMSE 或 MAE然后尝试不同的 degree画出“degree-RMSE”曲线。曲线最低点对应的是测试误差最小的阶数。比单纯看训练 R² 可靠得多。还有一种做法是交叉验证把数据切成 K 份轮流拿其中一份做验证其余训练最后平均验证误差。虽然计算量大一些但在样本量不大时能得出更稳定的结论。我的习惯是 degree 从 1 试到 8 左右找到测试 RMSE 最低的位置再小幅上下试探。4.3 从多项式回归到正则化岭回归如果业务上必须用较高的阶数或者特征之间存在明显共线性直接再用普通最小二乘会很不稳。这时可以在损失函数后面加一个正则项L (1/n)||y - Xw||² λ·||w||²这就是岭回归。闭式解会变成w (XᵀX λI)⁻¹Xᵀy注意这里 λ 是惩罚系数I 是单位矩阵。这个加项的效果相当于让参数向量尽量小避免模型过分依赖某个特征的极大系数。用梯度下降来做也只需要在梯度后面追加一个 2λw一般不对截距项做惩罚实际实现时可以把截距单独处理。所以当你想控制过拟合时加正则化比生硬地降低多项式阶数更平滑尤其适合特征间关联较强的场景。5. 常见问题与避坑实录5.1 梯度下降不收敛先检查学习率还是特征直接给结论先检查特征是否标准化再检查学习率。特征尺度差异过大会带来两个恶果第一是梯度方向失真参数更新会被某个量级巨大的特征主导第二是学习率不好设小了收敛极慢大了直接发散。标准化之后学习率 0.1 或 0.01 通常就能跑得很好。症状可能原因排查方法损失一直不降学习率太小 / 停止条件过早打印损失历史曲线调大学习率损失出现 NaN学习率过大 / 梯度爆炸调小学习率检查特征最大值损失反复震荡学习率偏大 / 数据没打乱调小学习率或用小批量替换全量收敛很慢特征尺度差异大先做标准化 / z-score结果与最小二乘差距大迭代不足 / 停止阈值太宽松加大 max_iter收紧 tol5.2 矩阵求逆失败或警告XᵀX 奇异怎么处理当 degree 不低且没有做特征缩放时XᵀX 的条件数会非常高。表面上 np.linalg.solve 不会报错但结果可能已经严重失真。真正报错的情况通常有两种一是特征之间存在严格线性关系比如不小心把 x 的某次幂重复放进去了二是样本量少于特征数XᵀX 必然奇异。遇到这类问题我的处理顺序是先检查设计矩阵的构造逻辑排除重复列再用 np.linalg.matrix_rank 看一眼 XᵀX 的秩最后换成 np.linalg.lstsq它对秩亏矩阵依然能给出一个最小二乘范数意义下的解。更根本的方案是采样增加数据量或者降低 degree或者在 XᵀX 上做一点正则化。5.3 训练 R² 很高但测试效果差过拟合的典型症状我先说一个区分方法论训练集 R² 高不代表模型好。如果训练集误差小、测试集误差明显变大模型就是在记忆训练数据而不是学习背后的规律。高次多项式特别容易出现这种状况因为多项式的表达能力强到可以穿过每一个噪声点。遇到这种情况不要急着加数据或加特征。先画一下拟合曲线如果曲线形如锯齿多半是过拟合了。处理手段有三个方向降低多项式阶数增加正则化项增加训练样本量或在数据层面做平滑。最稳的做法是交叉验证选阶数而不是拍脑袋定 degree。5.4 数据处理顺序泄露先划分数据集再缩放一个特别容易被忽视的问题在构造多项式特征或计算标准化参数时如果把全部数据包括测试集都用来 fit 标准化函数测试集的信息就“泄露”进训练过程了。这样的评估结果会偏乐观但放到真实环境里立刻露馅。正确的流程是先划分训练集和测试集然后在训练集上计算 x_mean、x_std再对测试集应用同一组参数做变换。多项式特征也是一样degree 的选择应该基于训练集上的交叉验证结果而不是基于测试集反复试探。5.5 一个常被忽略的小细节标准化还原预测用标准化后的特征训练模型时得到的参数是相对于 z-score 空间的。预测时如果直接输入新样本的原始 x忘了先减去均值、除以标准差结果会错得非常离谱而且这种错不太容易看出来——预测出的 y 量级可能看起来正常只是完全不在合理范围内。所以我在实操中习惯写一个简单的 predict 函数把“缩放—构造特征—点积”三个步骤封装在一起避免每次手工计算时出错。代码不长但能节省大量排查时间。def predict(x_new, theta, x_mean, x_std, degree): xs (x_new - x_mean) / x_std features np.vander(np.array([xs]), degree 1, increasingTrue) return features theta这个函数只处理单个样本批量预测时改成对整个数组操作即可。最后再分享一点个人体会做这类实验时最值得养成的习惯是“先把原理公式看懂再写代码”。我不止一次看到有人拿着 sklearn 直接一把梭跳过数据标准化、跳过训练测试划分最后只能靠直觉调参出了问题也无从下手。相反用这段从零实现的代码跑一遍你会对那些抽象概念产生很直观的感受——梯度下降到底是怎么沿着等高线走的、标准化对收敛有多大影响、闭式解为什么在小规模数据上如此优雅。这些经验在之后看其他模型时同样用得上。
