简介:鲸鱼优化算法WOA优化BP神经网络回归预测MATLAB代码,面向需要进行非线性回归预测的MATLAB开发者与算法学习者。以座头鲸捕食策略为灵感,通过WOA对BP神经网络的权重和阈值进行全局寻优,可有效缓解传统BP网络易陷入局部最优的问题,提升预测精度与稳定度。压缩包共5个文件,包含3个m脚本、1个mat数据文件和1个xlsx数据集,代码覆盖数据预处理、网络构建、WOA优化、训练测试与误差计算等环节,用户可直接运行主程序并替换Excel数据完成自定义回归任务。目前已有4203人学习下载,整体约197KB,结构简洁清晰。借助该代码,读者可同时理解生物启发式算法的优化机制与BP网络的应用流程,也可在此基础上调整适应度函数或网络层数,适配更多实际场景。
1. 为什么鲸鱼优化算法能让BP回归预测不再"看运气"
做过BP神经网络回归的人基本都遇到过同一个问题:同样的数据、同样的网络结构,多跑几次结果却天差地别。原因在于BP网络本质上是梯度下降驱动的,初始权重和阈值一旦落在误差曲面的平坦区或者局部极小点附近,训练过程就很难跳出来,预测精度基本靠"初始化运气"。鲸鱼优化算法(WOA)解决的就是这件事——它把BP的权重和阈值当作一组待寻优的参数,先用座头鲸的包围、气泡网攻击和随机搜索三种策略在解空间里做全局探索,找到一组较优的初值后再交给BP做局部精修。这套MATLAB代码把WOA和BP封装成了可直接运行的工程结构,只要准备好Excel格式的数据,main.m跑完就能拿到优化前后的误差对比。适合正在做回归预测课题的学生、需要快速验证算法效果的工程师,也适合想搞懂元启发式算法怎么和神经网络结合的入门研究者。
2. WOA寻优机制与BP网络参数空间的映射关系
2.1 座头鲸捕食策略的数学抽象
WOA的核心是三个阶段的数学模型:包围猎物、气泡网攻击、随机搜索。包围阶段用当前最优解作为参考位置,其它鲸鱼个体向它收缩;气泡网攻击阶段通过螺旋更新位置模拟座头鲸向上吐气泡圈的过程;随机搜索则让部分鲸鱼在全局范围内游荡,避免种群过早收敛。三种策略在迭代中按下式切换:
% 核心参数:a从2线性递减到0,r1、r2为[0,1]随机数 a = 2 - it * (2 / MaxIt); % 收敛因子,控制搜索步长 A = 2 * a * r1 - a; % 包围步长系数 C = 2 * r2; % 随机权重,增强探索性 p = rand; % 决定走螺旋更新还是包围收缩当A的绝对值大于1时,鲸鱼个体远离当前最优解,对应全局探索;小于1时向最优解靠拢,对应局部开发。这个A值的动态变化决定了WOA不需要像遗传算法那样单独设置交叉变异概率,算法本身的自适应机制就是从探索到开发的自然过渡。和粒子群算法相比,WOA的螺旋更新路径让种群多样性维持得更久,在处理高维参数优化时不容易早熟。
2.2 待优化参数与网络结构的关系
BP网络需要优化的参数分为两部分:输入层到隐藏层的权重矩阵W1、阈值b1,隐藏层到输出层的权重矩阵W2、阈值b2。这些参数的总数由网络结构唯一确定,计算公式为:
N_params = (n_in * n_hidden + n_hidden) + (n_hidden * n_out + n_out)假设输入特征为8维,隐藏层神经元10个,输出1维,参数总数就是(8×10+10)+(10×1+1)=101。这101个值按顺序展开成一个行向量,就是WOA中每条鲸鱼的位置向量。fitness.m做的就是把行向量重新拆回W1、b1、W2、b2,赋值给BP网络,然后用训练集计算均方误差作为适应度值。
2.2.1 适应度函数设计的关键
适应度函数必须用训练误差而不是测试误差,否则就相当于拿答案去考试,模型在测试集上的表现会失真。常见做法是把训练集的均方误差(MSE)作为适应度,WOA迭代结束后挑出适应度最小的鲸鱼个体,把它对应的参数作为BP网络的初始权重和阈值。这样BP要做的就不是从随机起点开始摸索,而是在WOA给出的优质解附近做局部精修,收敛速度和精度都会明显提升。
2.3 为什么要"WOA初值 + BP精修"两步走
WOA的全局搜索能力强,但精细搜索能力弱,迭代后期在最优解附近的微调不够细腻。BP的梯度下降恰恰擅长局部精修,但依赖初值。两者结合的逻辑很简单:先用WOA找出误差曲面上"大范围看已经很低"的区域,再用BP在这个区域内顺着梯度滑到谷底。这套思路和"先用全局优化器找好初始点,再用局部优化器精修"的通用框架一致,也符合MATLAB优化工具箱中GlobalSearch和fmincon配合使用的思想。
3. 代码结构与核心文件的执行逻辑
3.1 文件清单与调用关系
这套代码的文件结构非常清晰,核心文件就是main.m、fitness.m、calc_error.m,外加数据文件数据.xlsx和data1.mat。main.m是总入口,负责数据读取、WOA参数设置、种群初始化和主循环调用;fitness.m是适应度函数,被WOA主循环反复调用;calc_error.m负责把解码后的权重阈值赋给BP网络并计算误差。
%% main.m 主干流程 data = xlsread('数据.xlsx'); % 读取Excel数据,最后一列为输出 X = data(:, 1:end-1); Y = data(:, end); % 数据归一化到[0,1] X = mapminmax(X', 0, 1)'; Y = mapminmax(Y', 0, 1)'; % 划分训练集和测试集 ratio = 0.8; n_train = floor(size(X, 1) * ratio); X_train = X(1:n_train, :); Y_train = Y(1:n_train, :); X_test = X(n_train+1:end, :); Y_test = Y(n_train+1:end, :); % WOA参数 SearchAgents_no = 30; % 种群规模 MaxIt = 100; % 最大迭代次数 dim = 101; % 参数维度,由网络结构决定 lb = -5 * ones(1, dim); % 参数下界 ub = 5 * ones(1, dim); % 参数上界归一化这一步不能省,因为BP的激活函数(通常是tansig和purelin)对输入量级敏感,特征数值差异过大会让梯度计算不稳定。mapminmax把数据压到[0,1]区间,反向映射的时候用mapminmax('reverse', ...)还原预测值。lb和ub设置成[-5,5]是经验值,在这个范围内BP网络的权重初始化不会让神经元过早饱和。
3.2 WOA主循环的MATLAB实现
%% WOA主循环(简化版) Positions = lb + rand(SearchAgents_no, dim) .* (ub - lb); for it = 1:MaxIt for i = 1:SearchAgents_no % 边界处理 Positions(i, :) = max(Positions(i, :), lb); Positions(i, :) = min(Positions(i, :), ub); % 计算适应度 fitness(i) = fitness(Positions(i, :), X_train, Y_train, hiddennum); end [best_fit, idx] = min(fitness); best_pos = Positions(idx, :); % 更新每头鲸鱼位置(按A值选择包围或螺旋) for i = 1:SearchAgents_no if p < 0.5 if abs(A) < 1 Positions(i, :) = best_pos - A .* abs(C .* best_pos - Positions(i, :)); else rand_idx = randi(SearchAgents_no); Positions(i, :) = Positions(rand_idx, :) - A .* abs(C .* Positions(rand_idx, :) - Positions(i, :)); end else D = abs(best_pos - Positions(i, :)); Positions(i, :) = D .* exp(b) .* cos(2*pi*l) + best_pos; end end end每次迭代需要计算整个种群的适应度,而每次适应度调用都要建一次BP网络并做前向传播。如果数据量大、隐藏层神经元多,这个计算开销会非常可观。常见做法是把fitness函数写成不依赖神经网络工具箱的纯矩阵运算版本,直接计算隐藏层输出和最终误差,这样会比net = newff(...) + train(...)的方式快3到5倍。
3.3 fitness.m和calc_error.m的分工
fitness.m接收一个参数向量,内部调用calc_error.m完成误差计算。calc_error.m要做的是把参数向量按结构拆解,并完成BP的前向传播:
function err = calc_error(params, X_train, Y_train, hiddennum) % 拆解参数向量 n_in = size(X_train, 2); n_out = size(Y_train, 2); w1 = reshape(params(1:n_in*hiddennum), hiddennum, n_in); b1 = params(n_in*hiddennum+1 : n_in*hiddennum+hiddennum)'; w2 = reshape(params(n_in*hiddennum+hiddennum+1 : end-n_out), n_out, hiddennum); b2 = params(end-n_out+1 : end)'; % 前向传播 h_in = X_train * w1' + repmat(b1, size(X_train, 1), 1); h_out = 2 ./ (1 + exp(-2*h_in)) - 1; % tansig激活 y_pred = h_out * w2' + repmat(b2, size(X_train, 1), 1); err = mean(sum((y_pred - Y_train).^2, 2)); end这里用tansig公式替代了MATLAB的tansig函数调用,避免在fitness循环里频繁创建网络对象。repmat展开偏置项是MATLAB矩阵运算的常见技巧,如果直接用广播机制要注意版本兼容性,R2016b之前的版本对隐式扩展支持不好。
4. 数据适配、网络结构选择与参数调节策略
4.1 Excel数据格式与读取注意事项
数据.xlsx的默认格式是最后一列为输出Y,前面所有列是输入特征X。这里有一个经常踩的坑:如果Excel数据里包含文本表头,xlsread默认会把文本读成NaN,导致后续计算全部出错。解决办法是读取时指定数据范围,或者在Excel里就把表头删掉。另一个坑是缺失值,MATLAB的xlsread遇到空单元格会填NaN,WOA的适应度计算遇到NaN会连锁出错,所以数据预处理阶段必须检查是否有缺失值,用isnan(sum(X,2))找出含NaN的行并删掉或插值。
data1.mat是已经保存好的MAT格式数据,这和Excel数据是二选一的用法。如果你有现成的.mat文件,直接load('data1.mat')就能拿到变量,省去xlsread那一步。动态读取数据可以这样写:
if exist('数据.xlsx', 'file') data = xlsread('数据.xlsx'); else load('data1.mat'); % 变量名以实际文件中为准 end这样写的好处是别人拿到代码后不管用Excel还是MAT数据都能直接跑通。
4.2 隐藏层神经元数怎么定
隐藏层节点数没有标准公式,但有一个工程上常用的经验范围:输入维度加输出维度取中间值,再上下浮动。更科学一点用以下公式估算:
n_hidden = floor((n_in + n_out) / 2) + sqrt(n_train) / 10实际调试时可以从小到大依次试,记录每个隐藏层节点数对应的测试集MSE,选MSE最小的那个点。但注意隐藏层节点数不是越多越好,节点过多会带来过拟合问题,让模型在训练集上表现好而测试集上急剧退化。另外隐藏层节点数直接决定了WOA优化参数的维度——节点越多,dim越大,WOA的搜索空间呈线性增长,需要的种群规模和迭代次数也相应增加。
4.2.1 种群规模和迭代次数的平衡
SearchAgents_no和MaxIt是WOA最敏感的两个参数。种群规模太小,全局探索能力不足;规模太大,每次迭代的计算成本成倍增加。经验来看,参数维度在100左右时,种群规模30到50、迭代次数100到200是比较稳妥的起点。如果数据量不大(几百条样本),每次fitness计算只要几毫秒,200次迭代加30个种群也就是几十秒的耗时。但数据量过万时,每次前向传播就要遍历全部样本,计算量会陡增。
4.3 参数边界lb和ub的敏感性分析
WOA的位置向量在[lb, ub]范围内初始化并更新,这个边界设置对寻优结果有直接影响。BP的权重通常不需要很大的绝对值,tansig激活函数在输入绝对值大于3之后就接近饱和,梯度趋近于零。所以权重边界设置在[-5, 5]或[-3, 3]是比较合理的。边界设得过大会让WOA在无效区域浪费算力,边界过小又会限制搜索空间,可能漏掉最优解。
%% 动态计算维度,避免网络结构调整后忘记同步dim dim = (size(X_train, 2) + 1) * hiddennum + (hiddennum + 1) * size(Y_train, 2);用这段代码可以自动根据输入输出维度和隐藏层节点数计算参数维度,改网络结构时不用再手动更新dim,能省掉不少低级错误。
5. 收敛性分析、评估指标与调试技巧
5.1 怎么判断WOA优化是否有效
跑完main.m之后,代码会画出WOA的收敛曲线并对优化前后的预测结果做对比。收敛曲线应该是单调下降的,在迭代早期下降速度快,后期趋于平缓。如果你的收敛曲线是一条接近水平的直线,大概率是lb和ub设置不当或者fitness函数写错了。此时先用一个简单的测试函数验证WOA代码本身是否正确,比如用Rastrigin函数或者Sphere函数,能收敛到已知最优解,说明WOA主循环没问题,再回头查BP的前向传播和参数解码逻辑。
5.2 回归预测的常用评估指标
只用MSE衡量模型好坏不够全面,工程上一般同时看三个指标:相关系数R²衡量预测值与真实值的趋势一致性,均方根误差RMSE衡量绝对误差大小,平均绝对百分比误差MAPE衡量相对误差。R²越接近1越好,RMSE和MAPE越小越好。在calc_error.m里只用了MSE作为WOA的适应度函数,但最终评估时可以在main.m里补充计算这三个指标,代码只加几行:
R2 = 1 - sum((Y_test_true - Y_test_pred).^2) / sum((Y_test_true - mean(Y_test_true)).^2); RMSE = sqrt(mean((Y_test_true - Y_test_pred).^2)); MAPE = mean(abs((Y_test_true - Y_test_pred) ./ Y_test_true)) * 100;注意MAPE在真实值接近0时会变得非常大,数据包含接近0的样本时建议改用SMAPE或者直接忽略MAPE。
5.3 把fitness函数矩阵化的性能优化技巧
如果样本量上万,fitness.m里的循环会成为性能瓶颈。优化思路是把整个训练集的前向传播改成矩阵乘法,一次算出所有样本的预测值。关键步骤是把权重矩阵和输入矩阵的乘法用矩阵运算完成,避免for循环逐条样本计算。另一种常见的优化尝试是使用parfor替代for循环,但需要额外做数据切片和变量管理,收益不一定明显。
5.3.1 消除MapMinMax的隐藏坑
很多人在代码跑完后发现预测效果极差,排查了半天发现是归一化和反归一化的映射参数没有保存。mapminmax的默认模式会生成一个结构体记录训练集的归一化参数,测试集必须用同一组参数做映射,而不是用测试集自己的统计量重新归一化。正确做法是先获取归一化参数,再分别映射训练集和测试集:
[X_norm, ps] = mapminmax(X'); % ps保存归一化参数 X_test_norm = mapminmax('apply', X_test', ps); Y_pred = mapminmax('reverse', Y_pred_norm, ps_out);5.4 用不同激活函数的对比实验方法
Tansig是隐藏层的默认激活函数,它在[-1,1]区间内近似线性,两端饱和平滑,适合大多数回归场景。如果你想对比logsig和purelin的效果,只需在calc_error.m里改动激活函数那两行。logsig的输出范围是[0,1],如果归一化时把数据压到[-1,1],logsig和输出层的衔接会不匹配。实际做对比实验时,先确定归一化区间,再选择激活函数,两者必须配对使用。
5.4.1 换激活函数后需要同步调整的细节
激活函数改了,参数边界lb和ub也要对应调整。Logsig的输出范围是[0,1],中间层权重不需要很大的绝对值,边界可以缩小到[-3,3];纯线性输出层配合tansig隐藏层时,输出层偏置可能比较大,边界可以放宽到[-8,8]。这些细微调整对最终精度的影响有时比调迭代次数还明显,值得做一次完整的参数扫描实验。
5.5 从WOA优化结果反推网络设计问题
如果WOA优化后的BP预测精度依然不理想,优先检查数据而非检查算法。常见的三种情况:一是特征维度太高但样本量太少,模型学不到有效映射,这种情况减少特征或者加入正则项;二是输出值本身存在极端离群点,归一化后仍会拉偏误差曲面,建议提前做离群点剔除;三是输入特征中存在和输出完全无关的噪声维度,WOA虽然能降低权重但无法完全消除影响,用相关性分析删掉无效特征通常比加大迭代次数更有效。
本文还有配套的精品资源,点击获取