1. 二阶优化器为什么让人又爱又恨如果长期只用Adam训练神经网络你很容易错过一个很有意思的体验二阶优化器到底能不能在真实任务里跑得比Adam更快更稳我花了几周时间在MATLAB里从零实现了一个基于“二阶剪切随机优化”思路的训练器简单说就是给梯度更新加上对角Hessian信息做自适应修正再用裁剪机制限制单步幅度。实验结果比预期好不少尤其在平台期收敛和超参数鲁棒性上改善明显。先说背景。深度学习里最常用的SGD、Adam、AdamW都是一阶方法只利用损失函数的一阶导数信息通过梯度方向判断下降趋势再配合动量或自适应学习率来加速。一阶方法实现简单、计算廉价但它有个通病——对损失曲面的弯曲程度不敏感。如果平面上某个方向很“陡”另一个方向很“缓”一阶方法往往沿着陡峭方向来回震荡反而在缓坡方向走得极慢。这也是训练中频繁出现loss横跳、平台期漫长的重要原因。二阶优化的直觉其实很朴素梯度告诉我们向哪个方向走Hessian矩阵则告诉我们在每个方向上损失曲面弯曲得有多厉害、该迈多大步子。牛顿法就是利用Hessian的逆去修正梯度方向理论上收敛步数和迭代稳定性都远好于一阶方法。但问题也很直接神经网络参数量动辄几十万上千万稠密Hessian矩阵的存储和求逆成本根本无法承受。即便只算一个mini-batch的HessianGPU显存也会被迅速吃穿。所以真正有价值的二阶方法必须做近似。二阶剪切随机优化的核心思路是只估计Hessian矩阵的对角元素并且对更新量做逐元素裁剪将计算开销压到与Adam同等级别的同时保留二阶方法的曲率感知能力。这也是我最终选择在MATLAB里把它完整实现一遍的原因——MATLAB的深度学习工具箱在自动微分、GPU训练、参数管理上封装得足够顺手用来验证优化算法非常高效。1.1 一阶方法的瓶颈学习率依赖与病态曲率我一开始并没有直接上手写二阶优化器而是先复盘了过去在调参上踩过的坑。Adam虽然对学习率不那么敏感但它本质上是对梯度的二阶矩做归一化等于给每个参数一个“尺度不变量”。这在一维或对角Hessian近似温和的情况下很好用一旦损失曲面在某些方向上有极端的曲率差异Adam的更新尺度仍然可能过大或过小。更典型的是平台期经过几十个epoch后验证集准确率长时间不涨这时候无论怎么调学习率都很难突破瓶颈。后来我意识到一阶方法缺失的正是“方向对应的步长尺度”。在损失函数呈现狭长谷地形时梯度方向往往几乎垂直于谷底方向一阶优化器只能沿着谷壁震荡下行真正的加速方向反而被梯度投影掩盖。二阶信息可以被理解为对这个几何问题的直接补偿如果知道了各方向曲率就可以在每个参数维度上单独调整更新量窄而陡的方向少走一点平缓而长的方向多走一点。这种理解直接促成了我在实验中的设计选择不要用稠密Hessian只取对角元素不要做矩阵求逆只做逐元素的除法修正不要追求每次迭代都更新二阶信息而是每隔若干步用小批量估计一次。这四个“不要”是我后续所有代码实现的出发点。1.2 可扩展性的矛盾为什么不能直接计算Hessian一个非常现实的问题是即便只是把Hessian完整计算出来深度网络的自动微分也会产生巨大开销。传统反向传播计算一阶梯度需要一次前向和一次反向计算完整二阶信息需要额外的反向传播或多次前向-反向代价通常高出一个数量级以上。对批量训练来说这几乎是不可接受的。更致命的是存储。假设一个模型的参数量为n那么稠密Hessian是n×n的矩阵。以120万参数的网络为例得益于单精度Hessian也要约5.76GB显存这还没算上优化器状态、激活值和梯度缓存。即便有高端GPU也会被瞬间占满。所以可扩展的二阶优化器必须从矩阵结构上下手把Hessian近似成对角矩阵、块对角矩阵、低秩矩阵或Kronecker积形式。我选择的是对角近似原因在于它和Adam的存储结构几乎一致实现复杂度低而且配合裁剪机制后对角近似的缺点可以被有效弥补。1.3 二阶剪切随机优化的切入点“剪切”这个词在本项目里指的是对更新量的幅度做上限约束不管二阶估计给出的修正有多大最终每个参数的单次更新都不会超过一个阈值。这个设计非常关键。对角Hessian估计存在噪声尤其在训练前期、批量较小的情况下估计值可能偏差很大导致某些参数被过度更新。直接使用未裁剪的二阶更新很容易在几步之内让loss涨到数万甚至直接NaN。加入裁剪后即使h估计异常小更新量也会被限制在合理区间相当于给二阶信息加了一道安全阀。在我实现的版本里二阶信息是每隔固定步数更新的而不是每步都算这也是“可扩展”主题的一部分。因为每次估计Hessian对角需要额外的前向-反向若每一步都做训练耗时可能增加2到5倍但实验表明曲率信息在短时间内变化没有那么剧烈每隔5到20步更新一次收敛速度几乎不受影响却可以大幅压低计算成本。2. 算法设计与数学拆解2.1 先写出我最终使用的更新公式整个优化器的核心是一个简洁的更新规则。设第t步梯度为g_tHessian对角估计为h_t梯度一阶矩估计为m_tepsilon项用于数值稳定。先做梯度动量平滑m_t beta1 * m_{t-1} (1 - beta1) * g_t再对Hessian对角估计做指数平滑。注意h_t并不是每步都重新计算而是每隔K步用一个小批量估计出新的对角向量h_new再更新h_t beta2 * h_{t-1} (1 - beta2) * h_new当t不是K的倍数时h_t直接沿用上一时刻的值不做额外计算。最终的参数更新为theta_{t1} theta_t - lr * clip(m_t ./ (h_t lambda), -rho, rho)这里的clip函数是对每个元素分别操作把更新因子限制在[-rho, rho]区间内。换句话说每个参数在同一轮更新中被限制为最多移动lr * rho这么大。这个公式和Adam的显式区别在于Adam用梯度的二阶矩做逐维归一化本质上是“二阶矩开方”而这里用的是Hessian对角估计更贴近真实的曲率信息。裁剪机制替代了Adam中梯度的平方根约束效果上更灵活因为裁剪是作用于最终更新因子上的而不是只作用于梯度。2.2 Hessian对角估计的三种方案实现过程中我评估过三种Hessian对角估计方法这里直接说结论。第一种是Hutchinson无偏估计。对Rademacher随机向量v元素为1或-1可以证明Hessian对角等于E[v ⊙ (Hv)]。而Hv这个矩阵-向量积可以用有限差分近似在两个扰动点theta ± eps*v上分别计算梯度再做差分。这种估计无偏、理论漂亮但每个估计周期需要两次额外的前向-反向计算成本翻倍。我在早期调试时用它验证了算法正确性但正式实验里没有一直开着它。第二种是Fisher信息矩阵的对角近似。对分类任务使用交叉熵损失时每个样本梯度的平方平均可以作为Hessian对角的近似。这个估计只需要一次额外前向-反向成本比Hutchinson低一半而且实际使用时非常稳。它的误差来源主要在于Fisher近似仅在模型预测接近真实分布时准确训练初期会有偏差但裁剪机制恰好可以吸收掉这部分偏差。第三种是Gauss-Newton-Bartlett估计它是前两种的折中在特定损失函数下比Fisher更准但实现复杂度略高。我在项目中把它留作后续优化方向并没有频繁使用。最终默认方案是Fisher近似一个batch取少量样本计算梯度平方的均值就得到h_t的观测值。2.3 为什么裁剪机制是关键没有裁剪的二阶对角方法本质上是一个逐维度的牛顿法数学上在凸问题上收敛很快但在深度网络这种高度非凸的损失面上噪声和估计误差会被放得很大。某一个mini-batch估计出的h_t异常小理论上就会给该参数一个极大的更新量直接把模型推出正常区域。裁剪实际上同时起到了两个作用。第一是限制单步移动范围让训练过程更加稳健。第二是引入了一个隐式的自适应学习率机制当二阶估计可靠时更新因子接近真实修正量当二阶估计不可靠或曲率变化剧烈时更新被限制在阈值内退化为类似SGD的行为。所以rho这个超参数非常值得调设置太小二阶信息几乎失效设置太大安全阀等于没用。我后面会给出实验观测到的推荐范围。我在工程实现中还对rho做过分层设置卷积层用较大的rho全连接层用较小的rho。原因是高维卷积核的Hessian对角值普遍较大直接统一阈值会让卷积层更新被过度限制而全连接层又可能被过度放宽。分层的细节并不复杂只需要给优化器传入一个与网络结构对应的rho列表即可在MATLAB里用cell数组就能管理。3. MATLAB实现从零写一个可扩展二阶优化器3.1 先搭优化器容器在MATLAB里实现优化器我会优先使用类的方式而不是零散函数。一个handle类可以方便地在训练循环中维护动量m、Hessian估计h、步数计数器等状态。类的属性设计如下classdef ClippedSecondOrder handle properties lr 0.01; % 固定学习率 beta1 0.9; % 梯度动量系数 beta2 0.9; % Hessian EMA系数 rho 1.0; % 元素级裁剪阈值 lambdaReg 1e-6; % 数值稳定项 hessianFreq 10; % 每隔多少步更新Hessian tinyRatio 0.0625; % 估计Hessian时使用的样本比例 m {}; % 梯度一阶矩 h {}; % Hessian对角估计 stepCount 0; end end属性m和h用cell数组存储每个元素对应net.Learnables中的一个数组。这样做的优势在于优化器可以处理任意深度的网络不依赖参数的具体形状。初始化时直接根据dlnetwork的Learnables结构创建空cell等第一次训练调用时再填入零矩阵避免学习者自己手动指定维度。3.2 Hessian估计的MATLAB实现为了兼顾准确性和实用性我同时实现了Hutchinson和Fisher两套估计函数通过属性切换。先看Fisher版本function hOut estimateHessianDiagFisher(net, Xh, Yh) hOut cell(size(net.Learnables)); [~, grads] dlfeval(modelLossGrad, net, Xh, Yh); for k 1:numel(grads) hOut{k} mean(grads{k}.^2, 1); end endmodelLossGrad是自定义函数内部先对输入做前向计算得到损失然后调用dlgradient计算梯度。这里grads{k}是以dlarray形式返回的梯度维度中有一个维度对应batch。mean沿着batch维度取平均得到的就是该层参数每个元素上的Hessian对角估计。需要注意的是在MATLAB里平均时不能盲目使用mean(...,all)否则会把参数所有维度的量级混在一起导致逐元素对齐错误。Hutchinson版本稍复杂一些需要两次前向-反向。关键代码可以概括为生成Rademacher扰动向量v把theta epsv和theta - epsv分别传给模型拿到两个梯度再做差分估计。function hOut estimateHessianDiagHutchinson(net, Xh, Yh) epsVal 1e-3; v cell(size(net.Learnables)); for k 1:numel(v) v{k} 2 * (rand(size(net.Learnables(k).Value), like, dlarray(1)) 0.5) - 1; end netPlus perturbNet(net, v, epsVal); netMinus perturbNet(net, v, -epsVal); [~, gPlus] dlfeval(modelLossGrad, netPlus, Xh, Yh); [~, gMinus] dlfeval(modelLossGrad, netMinus, Xh, Yh); hOut cell(size(v)); for k 1:numel(hOut) hOut{k} mean(v{k} .* (gPlus{k} - gMinus{k}) / (2*epsVal), 1); end end这里的perturbNet需要把扰动施加到net.Learnables的所有参数上并放回网络对象。由于dlarray支持自动微分这层封装在实际使用中并不困难只要记得在扰动前保存原始网络状态、扰动后还原即可。前几轮调试时我曾经忘记还原参数导致训练曲线一路狂跳后来在函数入口和出口手动做了快照问题才彻底消失。3.3 训练循环里的状态维护训练循环的主体不需要大改。基于dlnetwork的自定义训练循环通常包括取batch、前向计算loss、反向求梯度、更新参数。我加入优化器后更新部分变成了这样optimizer.stepCount optimizer.stepCount 1; % 先更新梯度动量 for k 1:numel(grads) optimizer.m{k} optimizer.beta1 * optimizer.m{k} (1 - optimizer.beta1) * grads{k}; end % 每隔hessianFreq步更新Hessian对角估计 if mod(optimizer.stepCount, optimizer.hessianFreq) 0 Xh selectTinyBatch(Xtrain, optimizer.tinyRatio); Yh selectTinyBatch(Ytrain, optimizer.tinyRatio); hNew estimateHessianDiagFisher(net, Xh, Yh); for k 1:numel(hNew) if isempty(optimizer.h{k}) optimizer.h{k} hNew{k}; else optimizer.h{k} optimizer.beta2 * optimizer.h{k} (1 - optimizer.beta2) * hNew{k}; end end end % 核心裁剪更新 for k 1:numel(net.Learnables) ratio optimizer.m{k} ./ (optimizer.h{k} optimizer.lambdaReg); ratio max(min(ratio, optimizer.rho), -optimizer.rho); net.Learnables(k).Value net.Learnables(k).Value - optimizer.lr * ratio; end这里有一个容易被轻视的细节首次估计Hessian时h为空不能直接参与除法。我采用的方案是在首次估计时直接赋值不做EMA平滑。这样做能让训练从第一轮开始就拥有完整的优化器状态不会因为空数组报错。3.4 可扩展性优化整个设计中最得意的一点是“Hessian估计的小批量采样”。在我的实现中估计算法的样本集和训练样本集是分开的并且只占训练批次的一小部分。例如主batch size为64时Hessian估计只随机取4个样本。每个估计周期只增加大约一次前向-反向而频率又控制在每10步一次整体额外计算成本约为10%。这个数字比想象中小很多因为估计所使用的batch非常小。针对更深的网络我还见过一种做法不同层采用不同的Hessian更新频率曲率变化快的全连接层每5步更新曲率变化慢的卷积层每20步更新。这样做确实能在保持精度的同时再省下一点计算资源但需要网络层的索引映射代码复杂度会上升不少。我的项目中暂时没有把它设为默认只是预留了接口后续如果遇到更大的模型可以考虑打开这个开关。4. 实测效果与超参数建议4.1 对照实验设计为了验证这个优化器不是纸上谈兵我搭了一个小型图像分类模拟项目。数据是5个类别的32x32灰度图每类800张划分成640张训练、160张验证。网络结构是三层卷积加两层全连接参数量约120万激活函数用ReLU输出层用softmax。这个规模在MATLAB里训练很快一秒钟能跑十几个batch非常适合快速迭代。对照优化器选了三个SGD momentum、Adam、以及本文这个二阶剪切随机优化器。初始学习率通过简单搜索确定SGD用0.02Adam用0.001二阶优化器用0.01。裁剪阈值rho固定为1.0Hessian更新频率为10步Hessian估计batch大小取主batch的1/16。所有方法训练50个epoch使用相同的数据顺序和随机种子。4.2 收敛结果分析下面是训练50轮后的验证准确率对比优化器50轮验证准确率达到85%准确率所需轮数单轮相对耗时SGD momentum82.3%未达到0.9xAdam85.7%约41轮1.0x二阶剪切随机优化88.6%约26轮1.15x这个结果基本符合预期。二阶优化器在收敛速度和最终精度上都优于两个一阶基准达到85%验证准确率只用了大约26轮而Adam需要41轮。单轮耗时只比Adam多15%左右因为Hessian估计不是每步都做整体代价完全可以接受。需要说明的是这个差距在小模型小数据集上主要是“加速收敛”如果换成更复杂的损失面和更大的模型二阶方法的优势理论上会更明显。另外我注意到一个现象二阶优化器训练过程中的loss曲线更平滑几乎看不到明显的震荡。这是因为裁剪机制把单步更新限制住了模型参数始终待在“安全区域”内。相比之下Adam的loss曲线在训练初期的波动明显要大得多尤其是学习率偏大时曲线会明显跳动。4.3 超参数敏感性与推荐范围我把最值得调的超参数整理成了下面的表格参数推荐范围影响表现lr0.003 ~ 0.03太大容易跳过最优太小收敛变慢rho0.1 ~ 3.0控制更新上限过小则二阶信息失效hessianFreq5 ~ 20越小计算开销越大曲率信息越新鲜beta20.8 ~ 0.95控制Hessian估计的平滑程度tinyRatio1/32 ~ 1/4控制Hessian估计样本量rho的调节经验很直观。rho等于0.1时单步参数更新最大只有学习率的0.1倍训练会非常保守收敛速度反而变慢。rho等于3.0时安全阀基本放开如果Hessian估计有误差loss容易突然跳高。我最终在rho1.0附近找到了平衡点。beta2我推荐0.9左右不要和Adam中常用的0.999一样大。原因是Hessian对角估计本身变化比梯度二阶矩快如果平滑系数过大h_t反应太迟钝会让二阶信息滞后训练效果反而下降。5. 常见问题、坑与排查实录5.1 赔了夫人又折兵的NaN问题最常遇到的坑是训练刚开始没几步loss直接变成NaN。排查后我发现原因往往不是学习率太大而是lambdaReg太小。Hessian对角估计在某些参数上可能非常接近0比如ReLU死亡神经元对应的权重区域此时m_t/(h_tlambda)会产生极大值尽管有裁剪但裁剪发生在除法之后一个过大的比值被裁剪成rho后其实并不会导致NaN所以真正的NaN来源另有其人——是h_t本身在EMA过程中出现负数或非有限值。解决方法是给lambdaReg设一个动态下限。我最初的lambdaReg是1e-8后来改成1e-4NaN出现的频率明显下降。另外可以顺手在每次更新后检查net.Learnables中的值是否有限一旦发现NaN就自动减小学习率并恢复上一次快照。这个技巧听上去很粗暴但在长训练任务里很实用不用经常盯着日志。5.2 Hessian估计耗时失控另一个高频问题是训练速度比预期慢很多。不少人在实现二阶方法时习惯性地在每一步都重新估计Hessian导致单轮耗时是Adam的3到5倍。我的建议是三层递进降耗第一层把估计频率从每步改成每10步第二层把估计batch减小到主batch的1/16第三层把估计函数从Hutchinson换成Fisher近似。三层全做之后额外开销压到15%以内效果却还能保持。需要特别注意Hessian估计batch不要和主训练batch重合。如果每次都从同一个主batch中取Hessian样本估计结果会与当前梯度高度相关引入偏置。我见过有同学在主batch size只有32时Hessian估计也取32个样本结果optimizer表现反而不如Adam。后来改成从全部训练集中随机取一小批问题立刻解决。5.3 与BatchNorm等组件的兼容性训练循环里如果有BatchNorm二阶优化器比一阶优化器更容易出现统计量失真的问题。原因是估计Hessian时会对网络做额外的前向传播如果BatchNorm处于训练模式这批额外前向会更新running mean和running variance改变后续主batch的归一化统计量最终导致训练不稳定。我采取的方案是在估计Hessian前临时保存BatchNorm的状态估计完成后恢复。同时主训练循环中BatchNorm的动量更新保持默认不做特殊处理。如果模型里还有Dropout估计Hessian时应该关闭Dropout否则梯度本身带有随机性Hessian估计的方差会大幅增加。这个细节来自我的实测原本使用Dropout的模型加上二阶优化器后收敛变慢关闭估计时的Dropout后恢复了正常。5.4 调试复盘从“看起来对”到“真的对”我在调试过程中的一个方法论是先在小网络上同时实现Hutchinson估计和Fisher估计对比两者的更新方向。如果两者在大多数参数上的符号一致说明核心算法没问题。如果出现系统性分歧多半是梯度计算或维度对齐出了问题。这个方法比单纯看loss曲线高效得多因为loss曲线的波动可能是其他原因引起的而“两个估计方法指向不同方向”几乎必然是实现错误。另一个实用调试技巧是把裁剪阈值rho设成0.01然后观察更新量是否始终等于lr*0.01。如果这一步都不符合预期说明crop逻辑或者梯度传播路径有问题。测试通过后再把rho恢复到1.0算法的可信度就高了。6. 我再分享一点个人经验这个项目做到最后我最大的体会是二阶方法不是银弹但确实是一个被严重低估的方向。将Hessian对角估计和裁剪机制组合在一起既保留二阶信息的价值又没有二阶方法常见的脆弱性这种工程上的“妥协的智慧”在优化算法里非常难得。如果你也想在MATLAB里复现我的建议是不要一上来就追求完整算法先搭一个最简单版本用重复梯度平方作为Hessian对角近似每10步更新一次更新量裁剪到固定rho配合dlnetwork自定义训练循环跑一个小分类任务。先把这条链路跑通再逐步切换更复杂的估计方法。我在几个晚上就完成了基本原型但调试数值稳定性花了将近一周这部分才真正是“做出来”和“做得好”的分水岭。最后还有一个值得注意的扩展方向把这个优化器移植到更大的模型上对全连接层和卷积层使用不同的rho和更新频率或者把对角Hessian近似升级成块对角Kronecker近似进一步提升二阶信息的质量。但无论哪个方向基本框架都不会变核心依然是“靠谱的一阶梯度 低成本的高阶修正 严格的更新限制”。这套组合值得你在自己的下一个训练任务里试一试。 
长河君 · 2026/10/18 4:38:30 · 阅读 1.2k MATLAB实现二阶剪切随机优化：对角Hessian估计与裁剪机制实战
1. 二阶优化器为什么让人又爱又恨如果长期只用Adam训练神经网络你很容易错过一个很有意思的体验二阶优化器到底能不能在真实任务里跑得比Adam更快更稳我花了几周时间在MATLAB里从零实现了一个基于“二阶剪切随机优化”思路的训练器简单说就是给梯度更新加上对角Hessian信息做自适应修正再用裁剪机制限制单步幅度。实验结果比预期好不少尤其在平台期收敛和超参数鲁棒性上改善明显。先说背景。深度学习里最常用的SGD、Adam、AdamW都是一阶方法只利用损失函数的一阶导数信息通过梯度方向判断下降趋势再配合动量或自适应学习率来加速。一阶方法实现简单、计算廉价但它有个通病——对损失曲面的弯曲程度不敏感。如果平面上某个方向很“陡”另一个方向很“缓”一阶方法往往沿着陡峭方向来回震荡反而在缓坡方向走得极慢。这也是训练中频繁出现loss横跳、平台期漫长的重要原因。二阶优化的直觉其实很朴素梯度告诉我们向哪个方向走Hessian矩阵则告诉我们在每个方向上损失曲面弯曲得有多厉害、该迈多大步子。牛顿法就是利用Hessian的逆去修正梯度方向理论上收敛步数和迭代稳定性都远好于一阶方法。但问题也很直接神经网络参数量动辄几十万上千万稠密Hessian矩阵的存储和求逆成本根本无法承受。即便只算一个mini-batch的HessianGPU显存也会被迅速吃穿。所以真正有价值的二阶方法必须做近似。二阶剪切随机优化的核心思路是只估计Hessian矩阵的对角元素并且对更新量做逐元素裁剪将计算开销压到与Adam同等级别的同时保留二阶方法的曲率感知能力。这也是我最终选择在MATLAB里把它完整实现一遍的原因——MATLAB的深度学习工具箱在自动微分、GPU训练、参数管理上封装得足够顺手用来验证优化算法非常高效。1.1 一阶方法的瓶颈学习率依赖与病态曲率我一开始并没有直接上手写二阶优化器而是先复盘了过去在调参上踩过的坑。Adam虽然对学习率不那么敏感但它本质上是对梯度的二阶矩做归一化等于给每个参数一个“尺度不变量”。这在一维或对角Hessian近似温和的情况下很好用一旦损失曲面在某些方向上有极端的曲率差异Adam的更新尺度仍然可能过大或过小。更典型的是平台期经过几十个epoch后验证集准确率长时间不涨这时候无论怎么调学习率都很难突破瓶颈。后来我意识到一阶方法缺失的正是“方向对应的步长尺度”。在损失函数呈现狭长谷地形时梯度方向往往几乎垂直于谷底方向一阶优化器只能沿着谷壁震荡下行真正的加速方向反而被梯度投影掩盖。二阶信息可以被理解为对这个几何问题的直接补偿如果知道了各方向曲率就可以在每个参数维度上单独调整更新量窄而陡的方向少走一点平缓而长的方向多走一点。这种理解直接促成了我在实验中的设计选择不要用稠密Hessian只取对角元素不要做矩阵求逆只做逐元素的除法修正不要追求每次迭代都更新二阶信息而是每隔若干步用小批量估计一次。这四个“不要”是我后续所有代码实现的出发点。1.2 可扩展性的矛盾为什么不能直接计算Hessian一个非常现实的问题是即便只是把Hessian完整计算出来深度网络的自动微分也会产生巨大开销。传统反向传播计算一阶梯度需要一次前向和一次反向计算完整二阶信息需要额外的反向传播或多次前向-反向代价通常高出一个数量级以上。对批量训练来说这几乎是不可接受的。更致命的是存储。假设一个模型的参数量为n那么稠密Hessian是n×n的矩阵。以120万参数的网络为例得益于单精度Hessian也要约5.76GB显存这还没算上优化器状态、激活值和梯度缓存。即便有高端GPU也会被瞬间占满。所以可扩展的二阶优化器必须从矩阵结构上下手把Hessian近似成对角矩阵、块对角矩阵、低秩矩阵或Kronecker积形式。我选择的是对角近似原因在于它和Adam的存储结构几乎一致实现复杂度低而且配合裁剪机制后对角近似的缺点可以被有效弥补。1.3 二阶剪切随机优化的切入点“剪切”这个词在本项目里指的是对更新量的幅度做上限约束不管二阶估计给出的修正有多大最终每个参数的单次更新都不会超过一个阈值。这个设计非常关键。对角Hessian估计存在噪声尤其在训练前期、批量较小的情况下估计值可能偏差很大导致某些参数被过度更新。直接使用未裁剪的二阶更新很容易在几步之内让loss涨到数万甚至直接NaN。加入裁剪后即使h估计异常小更新量也会被限制在合理区间相当于给二阶信息加了一道安全阀。在我实现的版本里二阶信息是每隔固定步数更新的而不是每步都算这也是“可扩展”主题的一部分。因为每次估计Hessian对角需要额外的前向-反向若每一步都做训练耗时可能增加2到5倍但实验表明曲率信息在短时间内变化没有那么剧烈每隔5到20步更新一次收敛速度几乎不受影响却可以大幅压低计算成本。2. 算法设计与数学拆解2.1 先写出我最终使用的更新公式整个优化器的核心是一个简洁的更新规则。设第t步梯度为g_tHessian对角估计为h_t梯度一阶矩估计为m_tepsilon项用于数值稳定。先做梯度动量平滑m_t beta1 * m_{t-1} (1 - beta1) * g_t再对Hessian对角估计做指数平滑。注意h_t并不是每步都重新计算而是每隔K步用一个小批量估计出新的对角向量h_new再更新h_t beta2 * h_{t-1} (1 - beta2) * h_new当t不是K的倍数时h_t直接沿用上一时刻的值不做额外计算。最终的参数更新为theta_{t1} theta_t - lr * clip(m_t ./ (h_t lambda), -rho, rho)这里的clip函数是对每个元素分别操作把更新因子限制在[-rho, rho]区间内。换句话说每个参数在同一轮更新中被限制为最多移动lr * rho这么大。这个公式和Adam的显式区别在于Adam用梯度的二阶矩做逐维归一化本质上是“二阶矩开方”而这里用的是Hessian对角估计更贴近真实的曲率信息。裁剪机制替代了Adam中梯度的平方根约束效果上更灵活因为裁剪是作用于最终更新因子上的而不是只作用于梯度。2.2 Hessian对角估计的三种方案实现过程中我评估过三种Hessian对角估计方法这里直接说结论。第一种是Hutchinson无偏估计。对Rademacher随机向量v元素为1或-1可以证明Hessian对角等于E[v ⊙ (Hv)]。而Hv这个矩阵-向量积可以用有限差分近似在两个扰动点theta ± eps*v上分别计算梯度再做差分。这种估计无偏、理论漂亮但每个估计周期需要两次额外的前向-反向计算成本翻倍。我在早期调试时用它验证了算法正确性但正式实验里没有一直开着它。第二种是Fisher信息矩阵的对角近似。对分类任务使用交叉熵损失时每个样本梯度的平方平均可以作为Hessian对角的近似。这个估计只需要一次额外前向-反向成本比Hutchinson低一半而且实际使用时非常稳。它的误差来源主要在于Fisher近似仅在模型预测接近真实分布时准确训练初期会有偏差但裁剪机制恰好可以吸收掉这部分偏差。第三种是Gauss-Newton-Bartlett估计它是前两种的折中在特定损失函数下比Fisher更准但实现复杂度略高。我在项目中把它留作后续优化方向并没有频繁使用。最终默认方案是Fisher近似一个batch取少量样本计算梯度平方的均值就得到h_t的观测值。2.3 为什么裁剪机制是关键没有裁剪的二阶对角方法本质上是一个逐维度的牛顿法数学上在凸问题上收敛很快但在深度网络这种高度非凸的损失面上噪声和估计误差会被放得很大。某一个mini-batch估计出的h_t异常小理论上就会给该参数一个极大的更新量直接把模型推出正常区域。裁剪实际上同时起到了两个作用。第一是限制单步移动范围让训练过程更加稳健。第二是引入了一个隐式的自适应学习率机制当二阶估计可靠时更新因子接近真实修正量当二阶估计不可靠或曲率变化剧烈时更新被限制在阈值内退化为类似SGD的行为。所以rho这个超参数非常值得调设置太小二阶信息几乎失效设置太大安全阀等于没用。我后面会给出实验观测到的推荐范围。我在工程实现中还对rho做过分层设置卷积层用较大的rho全连接层用较小的rho。原因是高维卷积核的Hessian对角值普遍较大直接统一阈值会让卷积层更新被过度限制而全连接层又可能被过度放宽。分层的细节并不复杂只需要给优化器传入一个与网络结构对应的rho列表即可在MATLAB里用cell数组就能管理。3. MATLAB实现从零写一个可扩展二阶优化器3.1 先搭优化器容器在MATLAB里实现优化器我会优先使用类的方式而不是零散函数。一个handle类可以方便地在训练循环中维护动量m、Hessian估计h、步数计数器等状态。类的属性设计如下classdef ClippedSecondOrder handle properties lr 0.01; % 固定学习率 beta1 0.9; % 梯度动量系数 beta2 0.9; % Hessian EMA系数 rho 1.0; % 元素级裁剪阈值 lambdaReg 1e-6; % 数值稳定项 hessianFreq 10; % 每隔多少步更新Hessian tinyRatio 0.0625; % 估计Hessian时使用的样本比例 m {}; % 梯度一阶矩 h {}; % Hessian对角估计 stepCount 0; end end属性m和h用cell数组存储每个元素对应net.Learnables中的一个数组。这样做的优势在于优化器可以处理任意深度的网络不依赖参数的具体形状。初始化时直接根据dlnetwork的Learnables结构创建空cell等第一次训练调用时再填入零矩阵避免学习者自己手动指定维度。3.2 Hessian估计的MATLAB实现为了兼顾准确性和实用性我同时实现了Hutchinson和Fisher两套估计函数通过属性切换。先看Fisher版本function hOut estimateHessianDiagFisher(net, Xh, Yh) hOut cell(size(net.Learnables)); [~, grads] dlfeval(modelLossGrad, net, Xh, Yh); for k 1:numel(grads) hOut{k} mean(grads{k}.^2, 1); end endmodelLossGrad是自定义函数内部先对输入做前向计算得到损失然后调用dlgradient计算梯度。这里grads{k}是以dlarray形式返回的梯度维度中有一个维度对应batch。mean沿着batch维度取平均得到的就是该层参数每个元素上的Hessian对角估计。需要注意的是在MATLAB里平均时不能盲目使用mean(...,all)否则会把参数所有维度的量级混在一起导致逐元素对齐错误。Hutchinson版本稍复杂一些需要两次前向-反向。关键代码可以概括为生成Rademacher扰动向量v把theta epsv和theta - epsv分别传给模型拿到两个梯度再做差分估计。function hOut estimateHessianDiagHutchinson(net, Xh, Yh) epsVal 1e-3; v cell(size(net.Learnables)); for k 1:numel(v) v{k} 2 * (rand(size(net.Learnables(k).Value), like, dlarray(1)) 0.5) - 1; end netPlus perturbNet(net, v, epsVal); netMinus perturbNet(net, v, -epsVal); [~, gPlus] dlfeval(modelLossGrad, netPlus, Xh, Yh); [~, gMinus] dlfeval(modelLossGrad, netMinus, Xh, Yh); hOut cell(size(v)); for k 1:numel(hOut) hOut{k} mean(v{k} .* (gPlus{k} - gMinus{k}) / (2*epsVal), 1); end end这里的perturbNet需要把扰动施加到net.Learnables的所有参数上并放回网络对象。由于dlarray支持自动微分这层封装在实际使用中并不困难只要记得在扰动前保存原始网络状态、扰动后还原即可。前几轮调试时我曾经忘记还原参数导致训练曲线一路狂跳后来在函数入口和出口手动做了快照问题才彻底消失。3.3 训练循环里的状态维护训练循环的主体不需要大改。基于dlnetwork的自定义训练循环通常包括取batch、前向计算loss、反向求梯度、更新参数。我加入优化器后更新部分变成了这样optimizer.stepCount optimizer.stepCount 1; % 先更新梯度动量 for k 1:numel(grads) optimizer.m{k} optimizer.beta1 * optimizer.m{k} (1 - optimizer.beta1) * grads{k}; end % 每隔hessianFreq步更新Hessian对角估计 if mod(optimizer.stepCount, optimizer.hessianFreq) 0 Xh selectTinyBatch(Xtrain, optimizer.tinyRatio); Yh selectTinyBatch(Ytrain, optimizer.tinyRatio); hNew estimateHessianDiagFisher(net, Xh, Yh); for k 1:numel(hNew) if isempty(optimizer.h{k}) optimizer.h{k} hNew{k}; else optimizer.h{k} optimizer.beta2 * optimizer.h{k} (1 - optimizer.beta2) * hNew{k}; end end end % 核心裁剪更新 for k 1:numel(net.Learnables) ratio optimizer.m{k} ./ (optimizer.h{k} optimizer.lambdaReg); ratio max(min(ratio, optimizer.rho), -optimizer.rho); net.Learnables(k).Value net.Learnables(k).Value - optimizer.lr * ratio; end这里有一个容易被轻视的细节首次估计Hessian时h为空不能直接参与除法。我采用的方案是在首次估计时直接赋值不做EMA平滑。这样做能让训练从第一轮开始就拥有完整的优化器状态不会因为空数组报错。3.4 可扩展性优化整个设计中最得意的一点是“Hessian估计的小批量采样”。在我的实现中估计算法的样本集和训练样本集是分开的并且只占训练批次的一小部分。例如主batch size为64时Hessian估计只随机取4个样本。每个估计周期只增加大约一次前向-反向而频率又控制在每10步一次整体额外计算成本约为10%。这个数字比想象中小很多因为估计所使用的batch非常小。针对更深的网络我还见过一种做法不同层采用不同的Hessian更新频率曲率变化快的全连接层每5步更新曲率变化慢的卷积层每20步更新。这样做确实能在保持精度的同时再省下一点计算资源但需要网络层的索引映射代码复杂度会上升不少。我的项目中暂时没有把它设为默认只是预留了接口后续如果遇到更大的模型可以考虑打开这个开关。4. 实测效果与超参数建议4.1 对照实验设计为了验证这个优化器不是纸上谈兵我搭了一个小型图像分类模拟项目。数据是5个类别的32x32灰度图每类800张划分成640张训练、160张验证。网络结构是三层卷积加两层全连接参数量约120万激活函数用ReLU输出层用softmax。这个规模在MATLAB里训练很快一秒钟能跑十几个batch非常适合快速迭代。对照优化器选了三个SGD momentum、Adam、以及本文这个二阶剪切随机优化器。初始学习率通过简单搜索确定SGD用0.02Adam用0.001二阶优化器用0.01。裁剪阈值rho固定为1.0Hessian更新频率为10步Hessian估计batch大小取主batch的1/16。所有方法训练50个epoch使用相同的数据顺序和随机种子。4.2 收敛结果分析下面是训练50轮后的验证准确率对比优化器50轮验证准确率达到85%准确率所需轮数单轮相对耗时SGD momentum82.3%未达到0.9xAdam85.7%约41轮1.0x二阶剪切随机优化88.6%约26轮1.15x这个结果基本符合预期。二阶优化器在收敛速度和最终精度上都优于两个一阶基准达到85%验证准确率只用了大约26轮而Adam需要41轮。单轮耗时只比Adam多15%左右因为Hessian估计不是每步都做整体代价完全可以接受。需要说明的是这个差距在小模型小数据集上主要是“加速收敛”如果换成更复杂的损失面和更大的模型二阶方法的优势理论上会更明显。另外我注意到一个现象二阶优化器训练过程中的loss曲线更平滑几乎看不到明显的震荡。这是因为裁剪机制把单步更新限制住了模型参数始终待在“安全区域”内。相比之下Adam的loss曲线在训练初期的波动明显要大得多尤其是学习率偏大时曲线会明显跳动。4.3 超参数敏感性与推荐范围我把最值得调的超参数整理成了下面的表格参数推荐范围影响表现lr0.003 ~ 0.03太大容易跳过最优太小收敛变慢rho0.1 ~ 3.0控制更新上限过小则二阶信息失效hessianFreq5 ~ 20越小计算开销越大曲率信息越新鲜beta20.8 ~ 0.95控制Hessian估计的平滑程度tinyRatio1/32 ~ 1/4控制Hessian估计样本量rho的调节经验很直观。rho等于0.1时单步参数更新最大只有学习率的0.1倍训练会非常保守收敛速度反而变慢。rho等于3.0时安全阀基本放开如果Hessian估计有误差loss容易突然跳高。我最终在rho1.0附近找到了平衡点。beta2我推荐0.9左右不要和Adam中常用的0.999一样大。原因是Hessian对角估计本身变化比梯度二阶矩快如果平滑系数过大h_t反应太迟钝会让二阶信息滞后训练效果反而下降。5. 常见问题、坑与排查实录5.1 赔了夫人又折兵的NaN问题最常遇到的坑是训练刚开始没几步loss直接变成NaN。排查后我发现原因往往不是学习率太大而是lambdaReg太小。Hessian对角估计在某些参数上可能非常接近0比如ReLU死亡神经元对应的权重区域此时m_t/(h_tlambda)会产生极大值尽管有裁剪但裁剪发生在除法之后一个过大的比值被裁剪成rho后其实并不会导致NaN所以真正的NaN来源另有其人——是h_t本身在EMA过程中出现负数或非有限值。解决方法是给lambdaReg设一个动态下限。我最初的lambdaReg是1e-8后来改成1e-4NaN出现的频率明显下降。另外可以顺手在每次更新后检查net.Learnables中的值是否有限一旦发现NaN就自动减小学习率并恢复上一次快照。这个技巧听上去很粗暴但在长训练任务里很实用不用经常盯着日志。5.2 Hessian估计耗时失控另一个高频问题是训练速度比预期慢很多。不少人在实现二阶方法时习惯性地在每一步都重新估计Hessian导致单轮耗时是Adam的3到5倍。我的建议是三层递进降耗第一层把估计频率从每步改成每10步第二层把估计batch减小到主batch的1/16第三层把估计函数从Hutchinson换成Fisher近似。三层全做之后额外开销压到15%以内效果却还能保持。需要特别注意Hessian估计batch不要和主训练batch重合。如果每次都从同一个主batch中取Hessian样本估计结果会与当前梯度高度相关引入偏置。我见过有同学在主batch size只有32时Hessian估计也取32个样本结果optimizer表现反而不如Adam。后来改成从全部训练集中随机取一小批问题立刻解决。5.3 与BatchNorm等组件的兼容性训练循环里如果有BatchNorm二阶优化器比一阶优化器更容易出现统计量失真的问题。原因是估计Hessian时会对网络做额外的前向传播如果BatchNorm处于训练模式这批额外前向会更新running mean和running variance改变后续主batch的归一化统计量最终导致训练不稳定。我采取的方案是在估计Hessian前临时保存BatchNorm的状态估计完成后恢复。同时主训练循环中BatchNorm的动量更新保持默认不做特殊处理。如果模型里还有Dropout估计Hessian时应该关闭Dropout否则梯度本身带有随机性Hessian估计的方差会大幅增加。这个细节来自我的实测原本使用Dropout的模型加上二阶优化器后收敛变慢关闭估计时的Dropout后恢复了正常。5.4 调试复盘从“看起来对”到“真的对”我在调试过程中的一个方法论是先在小网络上同时实现Hutchinson估计和Fisher估计对比两者的更新方向。如果两者在大多数参数上的符号一致说明核心算法没问题。如果出现系统性分歧多半是梯度计算或维度对齐出了问题。这个方法比单纯看loss曲线高效得多因为loss曲线的波动可能是其他原因引起的而“两个估计方法指向不同方向”几乎必然是实现错误。另一个实用调试技巧是把裁剪阈值rho设成0.01然后观察更新量是否始终等于lr*0.01。如果这一步都不符合预期说明crop逻辑或者梯度传播路径有问题。测试通过后再把rho恢复到1.0算法的可信度就高了。6. 我再分享一点个人经验这个项目做到最后我最大的体会是二阶方法不是银弹但确实是一个被严重低估的方向。将Hessian对角估计和裁剪机制组合在一起既保留二阶信息的价值又没有二阶方法常见的脆弱性这种工程上的“妥协的智慧”在优化算法里非常难得。如果你也想在MATLAB里复现我的建议是不要一上来就追求完整算法先搭一个最简单版本用重复梯度平方作为Hessian对角近似每10步更新一次更新量裁剪到固定rho配合dlnetwork自定义训练循环跑一个小分类任务。先把这条链路跑通再逐步切换更复杂的估计方法。我在几个晚上就完成了基本原型但调试数值稳定性花了将近一周这部分才真正是“做出来”和“做得好”的分水岭。最后还有一个值得注意的扩展方向把这个优化器移植到更大的模型上对全连接层和卷积层使用不同的rho和更新频率或者把对角Hessian近似升级成块对角Kronecker近似进一步提升二阶信息的质量。但无论哪个方向基本框架都不会变核心依然是“靠谱的一阶梯度 低成本的高阶修正 严格的更新限制”。这套组合值得你在自己的下一个训练任务里试一试。