DEGWO-BP混合算法优化神经网络:回归预测稳定复现MATLAB实现
做回归预测的MATLAB用户应该都遇到过这种窘境同样的代码同样的数据集BP神经网络连续训练几次结果一次一个样。运气好几十轮迭代就收敛到不错的误差运气不好直接卡在局部极值里出不来预测曲线跟一条水平直线差不多。我第一次用BP做多输入单输出预测时为了复现一个能看的实验结果足足重复跑了二十多遍从那以后我就一直在找能稳定复现的改进方案。后来我把灰狼优化GWO和差分进化DE组合起来用这个混合算法去替代BP随机初始化的权值和阈值做成了这套基于DEGWO-BP的回归预测MATLAB代码。主程序直接读取EXCEL表格数据一次跑完会输出优化前后的对比曲线和各项误差指标不需要你手动去凑初始权值也不用反复重跑碰运气。这篇文章我想把整套代码的架构、关键函数、运行方式和调参思路完整讲一遍。适合正在做回归预测仿真、写论文需要加算法对比、或者想给BP神经网络加上智能优化算法的MATLAB用户参考。如果你手里已经有自己的数据集改一下Excel路径和输入输出列就能直接用。先说一下我的结论DEGWO-BP不是那种“花架子”优化它解决的是BP最让人头疼的稳定性问题在中等规模数据集上效果非常明显。1. DEGWO-BP到底在解决什么问题1.1 传统BP的三个老毛病你中了几个BP神经网络全称是反向传播神经网络核心机制是误差反向传播加梯度下降。它在函数拟合、回归预测、模式识别里用得非常多但用过的都知道它有几个根深蒂固的问题。第一个是初始权值和阈值随机生成导致的不可复现。BP在训练之前要随机初始化一组连接权值和偏置这个随机性直接决定了最终收敛到什么位置。哪怕你数据完全一样代码完全一样先后跑两次结果都可能差出一大截。做论文的时候审稿人最怕看到的就是这种“不可复现”的实验结果。第二个是梯度下降法本质上是一个局部搜索方法。误差曲面在参数空间里是高低起伏的BP的更新规则是沿着负梯度方向走也就是说它每一步都在找“当前点上下降最快的方向”。可问题是下降最快的方向只保证当前这一步在下山不保证最后能走到整个山脚。稍微复杂一点的误差曲面就会有一堆局部极小值坑BP很容易掉进其中一个坑里出不来。第三个是对超参数敏感。学习率设大了训练过程震荡误差不降反升设小了收敛慢得让人抓狂。隐含层节点数设置不合理要么网络容量不够拟合不了数据要么过拟合训练集误差很低但测试集一塌糊涂。这三个问题叠加在一起导致传统BP在小样本回归预测任务中经常表现得像个“盲盒”。我做DEGWO-BP这套东西最直接的目的就是解决第一个和第二个问题。用群智能优化算法去搜索一组更好的初始权值阈值让BP一开始就站在一个“比较好的山坡”上从源头减少随机初始化带来的不稳定性。1.2 灰狼优化和差分进化为什么天生适合组队灰狼优化GWO是2014年提出的一种群智能优化算法模拟的是灰狼群体的社会等级和捕猎行为。种群里有alpha、beta、delta、omega四个等级alpha对应当前最优解beta和delta对应次优解剩下的个体都是omega。每次迭代所有omega都会根据alpha、beta、delta三只“头狼”的位置来更新自己的位置相当于整个狼群同时朝几个较优方向逼近。GWO的优点很突出参数少实现简单收敛速度快。它只需要调两个系数A和C其中A由控制参数a线性递减刚开始a比较大算法倾向全局搜索后期a变小算法集中在局部精细搜索。这个机制让它比粒子群、遗传算法更容易上手。但GWO也有一个毛病后期所有个体都在向头狼靠拢种群多样性快速下降比较容易早熟收敛。简单说就是跑到后面大家挤在一起谁也跳不出当前区域一旦这个区域不是全局最优算法就只能困在那里。差分进化DE恰恰能补上这个短板。DE的变异算子是从种群中随机抽取三个互不相同的个体用它们的加权差生成变异向量然后再做交叉和选择。这种“随机差分引导”的方式能让种群始终保持探索活力全局搜索能力很强不容易陷入局部极值但纯DE的收敛速度相对GWO要慢一些。把两者组合起来逻辑就非常顺了GWO负责快速逼近较优区域DE负责在每轮迭代中制造扰动、维持多样性最后用贪心选择保证新个体不差于旧个体。我实际跑下来的感受是DEGWO的收敛曲线比单GWO更平滑测试集上的误差也更稳定。这个“稳定”在回归预测里非常重要因为它意味着你换一批数据、换个随机种子结果不会天差地别。1.3 这套代码从数据到结果一共走了几步整套主程序的工作流程其实并不复杂一共八步。第一步用xlsread读取Excel文件里的数据数据区域是A2到F2009也就是2008条样本、6列数据。第二步用mapminmax把所有特征归一化到[0,1]区间消除量纲影响。第三步按行划分训练集和测试集前1500条做训练后面508条做预测验证。第四步初始化DEGWO参数随机生成N个个体每个个体是一组BP权值和阈值编码。第五步进入迭代循环每轮先做GWO位置更新再做DE变异、交叉、选择同时用fobj函数计算每个个体的适应度也就是训练集误差。第六步迭代结束后取出历史最优个体Alpha_pos用这组权值阈值去构建BP神经网络完成训练和测试集预测。第七步同时再跑一个随机初始化的普通BP作为“优化前”对照组。第八步输出收敛曲线、回归拟合图和误差指标对比。这套流程最大的好处是通用性强。无论你的数据是什么领域——电力负荷预测、房价预测、气象温度预测、材料性能回归——只要整理成Excel表格最后一列是输出、前面几列是输入就可以直接套用。影响范围很广这也是我当初愿意把代码和思路完整分享出来的原因。2. 核心代码逐段拆解不是背代码是看懂代码2.1 数据读取与归一化xlsread和mapminmax的搭配细节主程序最开始是清空环境然后读取Excel数据代码长这样%% 清空环境变量 warning off % 关闭报警信息 close all % 关闭开启的图窗 clear % 清空变量 %% 读取数据 filename 数据.xlsx; sheet 1; xlRange A2:F2009; data xlsread(filename, sheet, xlRange);xlsread的第三个参数是单元格区域这里从A2开始是有讲究的。如果Excel第一行是表头也就是“日期、特征1、特征2……结果”这样的文字说明那读取时必须从第二行开始否则表头文字会被转成NaN后面归一化直接报错。读取完数据之后是归一化%% 数据归一化 [data_m, data_ps] mapminmax(data, 0, 1); data data_m;这里有个特别容易踩的坑mapminmax是按行处理数据的每一行是一个样本每一列是一个特征。但咱们从Excel读出来的data恰好是每一行一个样本所以要先把data转置让每个特征变成一行归一化完成后再转置回来。如果忘了转置归一化就会按样本维度去算结果完全错误。归一化到[0,1]是目前回归预测里最常用的选择。像Sigmoid这类激活函数输入在0附近时梯度最大学习效率高所以把数据压到[0,1]对BP训练非常友好。当然你也可以归一化到[-1,1]这个不影响主流程改一下参数就行。接着是训练集和测试集的划分%% 训练集和测试集划分 train_x data(1:1500, 1:end-1); train_y data(1:1500, end); test_x data(1501:end, 1:end-1); test_y data(1501:end, end);代码默认最后一列是输出前面所有列是输入。前1500条训练后面的都做测试。这种按时间或顺序的划分方式在回归预测里很常见但不适合做随机划分。如果你的数据集是乱序的建议改成randperm随机打乱后再划分避免某一段数据分布异常影响模型评估。2.2 适应度函数优化算法拿什么当“好吃”的标准fobj是整个DEGWO算法里的核心评价函数。DEGWO每生成一个新个体都要调用fobj去算这个个体的适应度适应度越小代表这组权值阈值越好。完整代码如下function error fobj(Positions) inputnum 13; % 输入层节点数 hiddennum 4; % 隐含层节点数 outputnum 1; % 输出层节点数 % 从个体中提取神经网络的权值和阈值 w1 Positions(1 : inputnum * hiddennum); B1 Positions(inputnum * hiddennum 1 : inputnum * hiddennum hiddennum); w2 Positions(inputnum * hiddennum hiddennum 1 : ... inputnum * hiddennum hiddennum hiddennum * outputnum); B2 Positions(inputnum * hiddennum hiddennum hiddennum * outputnum 1 : end); % 创建BP网络 net newff(train_x, train_y, hiddennum, {tansig, purelin}, trainlm); net.trainParam.epochs 1000; net.trainParam.goal 1e-5; net.trainParam.lr 0.01; % 网络权值阈值赋值 net.iw{1,1} reshape(w1, hiddennum, inputnum); net.lw{2,1} reshape(w2, outputnum, hiddennum); net.b{1} reshape(B1, hiddennum, 1); net.b{2} reshape(B2, outputnum, 1); % 网络训练 net train(net, train_x, train_y); % 仿真预测并计算误差 y_lj sim(net, train_x); error mean(abs(y_lj - train_y)); end这个函数看起来长其实核心逻辑就三块。第一块是从一个向量Position里按顺序截取出输入层到隐含层的权值w1、隐含层阈值B1、隐含层到输出层权值w2、输出层阈值B2。注意这里的截取顺序必须和主程序里个体编码的顺序完全一致否则切出来的值位置错乱网络直接报废。第二块是用newff创建BP网络tansig是隐含层激活函数purelin是输出层激活函数trainlm是Levenberg-Marquardt训练算法。第三块是把截取出的权值阈值填进网络用训练集训练然后对训练集做预测返回平均绝对误差mean(abs(y_lj - train_y))作为适应度。这里有一个关于维度的关键点必须提醒Dim的值必须等于BP网络全部权值和阈值的总数。按inputnum13、hiddennum4、outputnum1来计算总数是13×4 4×1 4 1也就是61。但你从网上下载的这份代码里写的是Dim26这说明原作者很可能按一个更小规模的网络结构留下的示例值。你复制代码后必须按自己的输入层节点数、隐含层节点数、输出层节点数重新算一遍Dim否则主程序初始化个体时维度不对fobj里的Position切片就会对不上运行必然报错。另一个我自己特别在意的设计细节是适应度用的是训练集误差而不是测试集误差。很多人为了追求“好看的结果”会把测试集放进优化过程里当反馈这其实是严重的方法论错误。优化算法看到测试集误差之后就相当于把测试集的信息泄露给了模型选择过程最后得到的测试集精度是虚高的论文里这么写非常容易被审稿人质疑。正确的做法是优化全程只看训练集反馈测试集留到最优参数确定之后一次性评估。2.3 DEGWO主循环GWO骨架 DE扰动 贪心选择主循环是整段代码的核心我拆成三部分来讲种群初始化、GWO位置更新、DE差分进化操作。初始化部分是这样N 30; % 种群规模 Dim 26; % 空间维度必须按实际网络结构重新计算 Max_iter 50; % 最大迭代次数 ub 1; lb -1; % 初始化种群位置和适应度 for i 1:N for j 1:Dim Positions(i,j) rand(1) * (ub - lb) lb; end Fitness(i) fobj(Positions(i,:)); end这一步的物理含义是生成30个候选解每个候选解都是一组BP权值阈值取值在[-1,1]范围内。随后把每只狼的适应度算出来并选出初始的Alpha_pos、Beta_pos、Delta_pos分别对应适应度最小的三个个体。进入迭代循环后第一大步是边界处理和等级更新然后计算衰减系数aa 2 - iter * ((2) / Max_iter);a从2线性递减到0这是GWO里最经典的操作。a越大A的取值范围越大狼的移动步长越大算法倾向全局探索a越小步长越小算法倾向在局部精细搜索。这个简单的线性衰减相当于给算法做了一个“先广后精”的节奏控制。第二大步是每个个体按灰狼追捕公式更新位置for i 1:N for j 1:Dim r1 rand(); r2 rand(); A1 2 * a * r1 - a; C1 2 * r2; D_alpha abs(C1 * Alpha_pos(j) - Positions(i,j)); X1 Alpha_pos(j) - A1 * D_alpha; % 同理计算X2、X3 Positions(i,j) (X1 X2 X3) / 3; end end简单说就是每个个体分别朝alpha、beta、delta这三个当前较优解的方向走一步然后取三步的平均位置。C系数是用来给距离加上随机权重避免算法太机械地贴近头狼。第三大步是DE的变异、交叉与选择这部分是DEGWO和纯GWO最大的区别。DE的变异公式是mutant Positions(r1,:) F * (Positions(r2,:) - Positions(r3,:));其中F是缩放因子r1、r2、r3是三个随机选出的、和当前个体i互不相同的个体编号。这相当于从种群中随机抽三个样本用它们之间的差分生成一个新方向这个方向完全是种群内部的信息交换和GWO的头狼引导方向无关所以能有效对抗GWO后期的多样性下降。生成变异向量之后做二项式交叉jrand randi([1, Dim]); for j 1:Dim if rand() CR || j jrand trial(j) mutant(j); else trial(j) Positions(i,j); end endCR是交叉概率jrand是保证至少有一个维度来自变异向量的保险机制。交叉完之后要做越界处理把超出[lb,ub]的值拉回边界最后用贪心选择trial_fit fobj(trial); if trial_fit Fitness(i) Positions(i,:) trial; Fitness(i) trial_fit; if trial_fit Alpha_score Alpha_score trial_fit; Alpha_pos trial; end end贪心选择的意思是如果试验个体trial的适应度比当前个体好就替换掉当前个体。这样每一轮迭代都能保证种群质量不下降。这一步非常关键因为它同时保证了DE带来的探索不会把算法带偏——好就收不好就丢。从我的使用经验来看F取0.5、CR取0.9是相对稳妥的默认值。F太大变异方向容易跳得太远收敛曲线会出现剧烈震荡F太小DE的扰动效果会变得微弱起不到增加多样性的作用。CR取0.9能让交叉后的试验向量更多地保留变异成分适合回归预测这种连续参数优化问题。2.4 优化完成后的网络重建与结果输出迭代结束之后Alpha_pos就是DEGWO找到的最优权值阈值组合。接下来要做的是重建网络、正式训练、预测、反归一化然后跑一个随机初始化的BP做对比。反归一化这一块是新手最容易漏的。前面把原始数据压到了[0,1]预测出来的结果也就在[0,1]范围内必须用mapminmax的反变换才能还原成原始量纲DEGWOpredict mapminmax(reverse, DEGWOpredict, output_ps);这里output_ps是归一化时保存的映射参数。注意如果你是在划分数据之前对全量data做了归一化那output_ps对应的是最后一列输出变量的映射如果你只对train_y做了归一化就得单独保存那一份output_ps。两种做法都可以但一定不能搞混。评价指标建议至少算四个均方根误差RMSE、平均绝对误差MAE、平均绝对百分比误差MAPE、决定系数R2。R2 1 - sum((test_y - predict).^2) / sum((test_y - mean(test_y)).^2); RMSE sqrt(mean((test_y - predict).^2)); MAE mean(abs(test_y - predict)); MAPE mean(abs((test_y - predict) ./ test_y)) * 100;R2越接近1说明模型解释能力越强RMSE和MAE越小说明预测误差越低。优化前后对比就是看你用随机初始化的BP跑出来的指标和用DEGWO优化后跑出来的指标差距。一般来说DEGWO-BP的R2会更高RMSE会更低而且多次运行时结果波动明显更小。绘图部分按个人喜好来我习惯画四张子图收敛曲线、训练集拟合图、测试集拟合图、误差对比柱状图。这样一屏就能看完全部信息。3. 实操运行与调参把代码变成你自己的工具3.1 跑代码前必须检查的5个位置第一文件路径。Excel数据文件和工作脚本必须在同一个MATLAB当前工作目录下或者你在代码里把filename写成完整路径。中文路径在部分MATLAB版本里偶尔会出问题尽量用英文目录别在桌面直接双击运行。第二数据区域。先打开Excel看一眼数据到底有多少行、多少列。如果数据是2008行、6列第一行是表头那xlRange就该写成A2:F2009。如果数据没有表头从第一行就是数值那就改成A1:F2008。这一行报错是新手最常见的问题。第三输入输出列。代码默认最后一列是y前面所有列是x。如果你的Excel布局不是这样要么改读取区域要么调整train_x和train_y的提取方式。第四网络结构参数和Dim。inputnum、hiddennum、outputnum三个值先在fobj里定好然后在主程序里按权值阈值总数重新计算Dim。一遍写不对就拆开算w1加B1加w2加B2的总个数就是Dim。第五样本量是否满足网络训练要求。输入维度很高但样本只有几百条BP很容易学不进去。这种情况下要么增加数据要么先做特征筛选减少inputnum。3.2 关键参数怎么调直接可抄的经验值参数调节是整个流程里最玄学但也最有规律的部分。我把常用参数的经验值整理成表格方便你对照调整参数常见范围我推荐的经验值参数影响N 种群规模20到5030太小多样性差太大训练时间成倍增加Max_iter 最大迭代数50到20080回归预测问题50代以后收益明显变小lb、ub 解空间边界-1到1或-3到3-1、1权值初始值范围过大会导致网络训练发散F 缩放因子0.4到0.90.5F过大容易震荡F过小DE扰动变弱CR 交叉概率0.7到0.90.9越大试验向量保留变异成分越多hiddennum 隐含层节点数经验公式计算试凑后确定直接影响网络拟合能力和过拟合风险隐含层节点数有一个常见的经验公式hiddennum sqrt(inputnum outputnum) a其中a取1到10之间的整数。另一个常用做法是hiddennum 2 * inputnum 1。你要是不确定就从小往大试先在训练集上看到误差明显下降再留意测试集误差是否回升。测试集误差开始回升就是过拟合的典型信号这时候要减少节点数或者增大训练样本量。学习率lr和训练目标goal这两个参数我一般保持lr0.01、goal1e-5不动。如果训练过程不收敛把学习率降到0.005看看如果收敛太快导致精度不够把goal从1e-5收紧到1e-6同时增加训练轮数epochs。3.3 优化前后对比图到底应该怎么读很多人跑完代码看到两张拟合图和一个误差表不知道从哪里判断DEGWO到底有没有起作用。我的建议是分三步看。第一步看收敛曲线。DEGWO的收敛曲线应当是平滑下降的前期快速下降后期趋于水平。如果曲线在后期还在剧烈上下跳说明F设置偏大或者种群规模N偏小这种情况下优化结果不稳定多跑几次结果会差很多。如果曲线下降很慢迭代到50代还在明显下降那就加Max_iter让它收敛到平台期再做对比。第二步看测试集拟合图。理想情况下散点应该紧贴yx对角线分布。优化后的图至少应该比优化前的图“收拢”一些尤其是远离对角线的异常点如果变少了说明DEGWO确实帮BP找到了更合理的初始位置。第三步看误差指标。RMSE、MAE、MAPE三个指标优化后应当同时下降R2应当上升。这里要特别提醒一下如果你的数据本身就比较简单线性关系很强原始BP加上多次重跑也可能跑出和DEGWO-BP差不多的结果这不代表代码有问题只说明这个数据集对初始值不敏感。这时候建议用rng固定随机种子多对比几次运行结果看DEGWO-BP的指标波动范围是否明显小于原始BP。稳定性本身就是优势。4. 常见报错与排查思路实录4.1 读取Excel数据时最容易踩的坑第一个坑是“文件不存在”报错。大概率是文件名写错、扩展名写漏或者工作目录没有切换。建议在MATLAB里用cd命令切到数据所在目录再执行脚本。第二个坑是数据区域不对。第一行有表头却写成A1:F2009或者数据实际只有1500行却写了F2009都会导致读取出的数据矩阵尺寸不符合预期。后面的train_x、train_y划分一旦越界立刻报“索引超出矩阵维度”。第三个坑是Excel单元格里有文本、空值或公式生成的异常值。xlsread会把这些转成NaN而NaN进入mapminmax之后整个数据列都可能变成NaN后面BP训练直接崩溃。排查方法是在读取数据之后加一句检查。我喜欢加一行data data(~any(isnan(data), 2), :);直接删掉包含NaN的行简单粗暴但有效。当然更稳妥的做法是把Excel里的空行、备注文字全部清掉。第四个坑是新版MATLAB对xlsread的警告。R2019a之后的版本官方更推荐readmatrix、readtable函数xlsread虽然还能用但会提示Warning。这个Warning不影响运行介意的话可以换成readmatrix返回的是数值矩阵用法基本一致。4.2 newff、train、mapminmax的版本和维度问题newff是经典老函数在较新版本的MATLAB里依然兼容但如果你用的是新版且同时安装了深度学习工具箱newff有时会弹出版本提醒。遇到“输入参数的数目不足”或“不能打开网络”这类报错先排查是不是手滑把hiddennum、输入输出层数写成了0或者负数。trainlm是默认训练函数但它需要计算雅可比矩阵内存占用比较大。如果你的训练样本量很大比如超过几万条trainlm可能直接报内存不足。这种情况下建议换成trainbr贝叶斯正则化或trainscg比例共轭梯度法收敛速度可能会稍慢但内存占用低很多泛化效果往往更好。reshape维度不匹配是我自己实操时最常遇到的报错。fobj函数里用的是reshape(w1, hiddennum, inputnum)如果w1的长度不等于hiddennum×inputnumMATLAB立刻报错。这个错误的根源几乎都是主程序的Dim算错了。记住一句话Dim必须和fobj里切分Positions的方式严格对应多一个少一个都不行。mapminmax的维度问题也值得一提。mapminmax的默认处理方式是“每个特征一行”如果直接对行样本做归一化等于把不同特征的取值范围混在一起算结果失去意义。代码里先做data转置归一化之后再转置回来这套顺序不要打乱。反归一化时同理预测值必须先转置再reverse得到的结果才是和原始标签同量纲的值。4.3 预测结果异常时的排查线路结果异常比报错更让人头大因为程序明明跑完了但输出根本没法用。我总结了四条常见线路。第一条预测值几乎是一个常数。这说明BP没有学到有效特征最可能的原因是输入数据没有归一化或归一化出错另一个可能是学习率太大导致训练过程直接震荡发散。检查思路是先把学习率调小到0.001再把归一化逻辑重新检查一遍。第二条训练集误差很低但测试集误差高得离谱。这是典型的过拟合。解决办法是减少hiddennum或者在BP训练时开设验证集做早停或者增加训练数据量。如果数据量确实没法增加考虑加正则化项。第三条收敛曲线下降很快但最终精度不高。这说明算法早熟了大概率是种群规模N太小或者初始解空间范围设得不够合理。可以尝试把N从30提到50把Max_iter从50提到100让算法有更多机会跳出早期的小坑。第四条多跑几次结果差异巨大。先用rng(1)固定随机种子排除随机性影响。如果固定种子后差异消失说明算法本身没问题只是你之前对比的时候没有控制变量。如果固定种子后测试集精度依然时好时坏说明优化过程没有充分收敛按第三条的思路加大N和Max_iter。4.4 报错与现象速查表日常被问到最多的问题我整理成一个表放在这里。遇到问题先对照这个表排查一轮80%的情况都能自己解决。报错或现象可能原因解决办法未定义变量或函数xlsread数据文件路径错误或MATLAB版本兼容性问题检查工作目录和文件名或改用readmatrix索引超出矩阵维度Dim与网络参数不匹配或数据行数不够划分范围重新计算Dim检查train_x/test_x划分范围输出结果全是NaN数据中含NaN或网络训练发散清洗数据降低学习率检查归一化reshape维度不匹配w1等切片长度与目标矩阵元素数不一致核对Dim计算公式和Position切片逻辑trainlm内存不足样本量或网络规模过大改用trainbr或trainscg训练算法R2为负数模型比直接取均值还差检查数据划分、网络结构、归一化是否正确收敛曲线剧烈震荡F太大或种群规模太小F降到0.3到0.4之间N提到40以上多次运行结果差异大随机种子未固定或优化未收敛用rng固定随机种子加大Max_iter训练效果好但测试效果差过拟合减少隐含层节点增加训练样本考虑早停最后再分享一个我自己的习惯拿到这类优化代码先别急着换自己的数据跑。先原样跑通确认每一步都不报错再把自己的Excel数据复制到项目目录下改文件名、改数据区域、改输入输出列、改Dim一步一步替换。我见过太多人一上来就把数据换成自己的报错之后根本分不清是数据问题还是代码问题。调试成本最低的方式永远是“小步改动逐步替换”。这套DEGWO-BP说到底只是把BP的初始解选得更合理它改变不了数据本身的质量。先把基础模型跑通再做优化才是效率最高的路径。