使用灰狼优化算法(GWO)优化最小二乘支持向量回归(LSSVR)参数的MATLAB实现
classdef GWO_LSSVR% 灰狼算法优化的LSSVR模型% 使用灰狼优化算法寻找LSSVR的最佳参数(C, gamma)propertiesC; % 正则化参数gamma; % RBF核参数X_train; % 训练数据特征y_train; % 训练数据标签X_test; % 测试数据特征y_test; % 测试数据标签model; % 训练好的LSSVR模型fitnessHistory; % 适应度历史记录endmethodsfunction obj = GWO_LSSVR(X_train, y_train, X_test, y_test)% 构造函数obj.X_train = X_train;obj.y_train = y_train;obj.X_test = X_test;obj.y_test = y_test;obj.fitnessHistory = [];endfunction [C_opt, gamma_opt] = train(obj, numWolves, maxIter)% 使用灰狼算法优化LSSVR参数% 参数设置dim = 2; % 优化参数维度 (C, gamma)lb = [0.1, 0.01]; % 参数下界 [C_min, gamma_min]ub = [1000, 100]; % 参数上界 [C_max, gamma_max]% 初始化灰狼种群wolves = initializeWolves(numWolves, dim, lb, ub);alpha_pos = zeros(1, dim); % α狼位置alpha_score = inf; % α狼适应度beta_pos = zeros(1, dim); % β狼位置beta_score = inf; % β狼适应度delta_pos = zeros(1, dim); % δ狼位置delta_score = inf; % δ狼适应度% 迭代优化for iter = 1:maxItera = 2 - iter * (2 / maxIter); % 线性递减参数afor i = 1:numWolves% 边界检查wolves(i, :) = boundConstraint(wolves(i, :), lb, ub);% 计算适应度fitness = obj.evaluateFitness(wolves(i, 1), wolves(i, 2));% 更新α, β, δ狼if fitness < alpha_scorealpha_score = fitness;alpha_pos = wolves(i, :);endif fitness > alpha_score && fitness < beta_scorebeta_score = fitness;beta_pos = wolves(i, :);endif fitness > alpha_score && fitness > beta_score && fitness < delta_scoredelta_score = fitness;delta_pos = wolves(i, :);endend% 更新每只狼的位置for i = 1:numWolvesfor j = 1:dim% 更新α狼影响r1 = rand();r2 = rand();A1 = 2*a*r1 - a;C1 = 2*r2;D_alpha = abs(C1*alpha_pos(j) - wolves(i, j));X1 = alpha_pos(j) - A1*D_alpha;% 更新β狼影响r1 = rand();r2 = rand();A2 = 2*a*r1 - a;C2 = 2*r2;D_beta = abs(C2*beta_pos(j) - wolves(i, j));X2 = beta_pos(j) - A2*D_beta;% 更新δ狼影响r1 = rand();r2 = rand();A3 = 2*a*r1 - a;C3 = 2*r2;D_delta = abs(C3*delta_pos(j) - wolves(i, j));X3 = delta_pos(j) - A3*D_delta;% 位置更新wolves(i, j) = (X1 + X2 + X3) / 3;endend% 记录最佳适应度obj.fitnessHistory(iter) = alpha_score;% 显示迭代信息if mod(iter, 10) == 0fprintf('Iteration %d: Best Fitness = %.4f, C = %.4f, gamma = %.4f\n', ...iter, alpha_score, alpha_pos(1), alpha_pos(2));endend% 保存最优参数C_opt = alpha_pos(1);gamma_opt = alpha_pos(2);obj.C = C_opt;obj.gamma = gamma_opt;% 使用最优参数训练最终模型obj.model = obj.trainLSSVR(C_opt, gamma_opt);endfunction fitness = evaluateFitness(obj, C, gamma)% 评估适应度函数 (MSE)model = obj.trainLSSVR(C, gamma);y_pred = obj.predict(model, obj.X_test);mse = mean((obj.y_test - y_pred).^2);fitness = mse; % 最小化MSEendfunction model = trainLSSVR(obj, C, gamma)% 训练LSSVR模型% 使用RBF核函数sigma = 1/sqrt(2*gamma); % RBF核参数转换% 计算核矩阵n = size(obj.X_train, 1);K = zeros(n, n);for i = 1:nfor j = 1:ndiff = obj.X_train(i, :) - obj.X_train(j, :);K(i, j) = exp(-norm(diff)^2/(2*sigma^2));endend% 构建线性方程组Omega = [K + eye(n)/C, ones(n, 1); ones(n, 1)', 0];Y = [obj.y_train; 0];% 求解线性方程组sol = Omega \ Y;alpha = sol(1:n);b = sol(n+1);% 保存模型model.alpha = alpha;model.b = b;model.X_train = obj.X_train;model.gamma = gamma;model.C = C;model.sigma = sigma;endfunction y_pred = predict(obj, model, X)% 使用LSSVR模型进行预测n_train = size(model.X_train, 1);n_test = size(X, 1);y_pred = zeros(n_test, 1);for i = 1:n_tests = 0;for j = 1:n_traindiff = X(i, :) - model.X_train(j, :);k_val = exp(-norm(diff)^2/(2*model.sigma^2));s = s + model.alpha(j)*k_val;endy_pred(i) = s + model.b;endendfunction plotResults(obj)% 绘制优化过程和预测结果figure;% 适应度曲线subplot(2, 2, 1);plot(obj.fitnessHistory, 'LineWidth', 2);title('灰狼优化适应度曲线');xlabel('迭代次数');ylabel('均方误差(MSE)');grid on;% 预测结果与真实值对比y_pred = obj.predict(obj.model, obj.X_test);subplot(2, 2, 2);plot(obj.y_test, 'b-o', 'LineWidth', 1.5, 'DisplayName', '真实值');hold on;plot(y_pred, 'r-*', 'LineWidth', 1.5, 'DisplayName', '预测值');title('预测结果与真实值对比');xlabel('样本索引');ylabel('目标值');legend;grid on;% 误差分布subplot(2, 2, 3);errors = obj.y_test - y_pred;histogram(errors, 20);title('预测误差分布');xlabel('预测误差');ylabel('频数');grid on;% 参数优化结果subplot(2, 2, 4);text(0.1, 0.8, sprintf('最优参数 C = %.4f', obj.C), 'FontSize', 12);text(0.1, 0.6, sprintf('最优参数 gamma = %.4f', obj.gamma), 'FontSize', 12);text(0.1, 0.4, sprintf('最终MSE = %.4f', mean((obj.y_test - y_pred).^2)), 'FontSize', 12);text(0.1, 0.2, sprintf('R² = %.4f', 1 - sum((obj.y_test - y_pred).^2)/sum((obj.y_test - mean(obj.y_test)).^2)), 'FontSize', 12);axis off;title('优化结果总结');endfunction metrics = evaluateModel(obj)% 评估模型性能y_pred = obj.predict(obj.model, obj.X_test);% 计算性能指标mse = mean((obj.y_test - y_pred).^2);rmse = sqrt(mse);mae = mean(abs(obj.y_test - y_pred));r2 = 1 - sum((obj.y_test - y_pred).^2)/sum((obj.y_test - mean(obj.y_test)).^2);% 显示结果fprintf('\n===== 模型性能评估 =====\n');fprintf('均方误差(MSE): %.4f\n', mse);fprintf('均方根误差(RMSE): %.4f\n', rmse);fprintf('平均绝对误差(MAE): %.4f\n', mae);fprintf('决定系数(R²): %.4f\n', r2);% 返回指标metrics = struct('MSE', mse, 'RMSE', rmse, 'MAE', mae, 'R2', r2);endend
end%% 辅助函数
function wolves = initializeWolves(numWolves, dim, lb, ub)% 初始化灰狼种群wolves = zeros(numWolves, dim);for i = 1:dimwolves(:, i) = lb(i) + (ub(i) - lb(i)) * rand(numWolves, 1);end
endfunction x = boundConstraint(x, lb, ub)% 边界约束处理for i = 1:length(x)if x(i) < lb(i)x(i) = lb(i);elseif x(i) > ub(i)x(i) = ub(i);endend
end%% 数据生成和演示脚本
function demoGWO_LSSVR()% 演示灰狼优化LSSVRclear; close all; clc;% 生成示例数据rng(42); % 设置随机种子nSamples = 200;nFeatures = 5;% 生成非线性关系数据X = rand(nSamples, nFeatures) * 10;y = sin(X(:,1)) + 0.5*cos(2*X(:,2)) + 0.3*X(:,3).^2 + 0.1*X(:,4).*X(:,5) + randn(nSamples, 1)*0.1;% 划分训练集和测试集splitRatio = 0.8;nTrain = round(splitRatio * nSamples);idx = randperm(nSamples);trainIdx = idx(1:nTrain);testIdx = idx(nTrain+1:end);X_train = X(trainIdx, :);y_train = y(trainIdx);X_test = X(testIdx, :);y_test = y(testIdx);% 创建并训练GWO-LSSVR模型gwo_lssvr = GWO_LSSVR(X_train, y_train, X_test, y_test);[C_opt, gamma_opt] = gwo_lssvr.train(20, 50); % 20只狼,50次迭代% 评估模型metrics = gwo_lssvr.evaluateModel();% 绘制结果gwo_lssvr.plotResults();% 与网格搜索比较compareWithGridSearch(X_train, y_train, X_test, y_test);
endfunction compareWithGridSearch(X_train, y_train, X_test, y_test)% 与网格搜索调参比较fprintf('\n===== 网格搜索调参比较 =====\n');% 定义参数范围C_list = logspace(-1, 3, 10); % 0.1到1000gamma_list = logspace(-2, 2, 10); % 0.01到100best_mse = inf;best_C = 0;best_gamma = 0;% 网格搜索for i = 1:length(C_list)for j = 1:length(gamma_list)C = C_list(i);gamma = gamma_list(j);% 训练临时模型temp_model = trainLSSVR(X_train, y_train, C, gamma);y_pred = predictLSSVR(temp_model, X_test);mse = mean((y_test - y_pred).^2);if mse < best_msebest_mse = mse;best_C = C;best_gamma = gamma;endendendfprintf('网格搜索最佳参数: C = %.4f, gamma = %.4f\n', best_C, best_gamma);fprintf('网格搜索最佳MSE: %.4f\n', best_mse);% 使用最佳参数训练最终模型final_model = trainLSSVR(X_train, y_train, best_C, best_gamma);y_pred_grid = predictLSSVR(final_model, X_test);r2_grid = 1 - sum((y_test - y_pred_grid).^2)/sum((y_test - mean(y_test)).^2);fprintf('网格搜索R²: %.4f\n', r2_grid);
endfunction model = trainLSSVR(X_train, y_train, C, gamma)% 训练LSSVR模型 (简化版)sigma = 1/sqrt(2*gamma);n = size(X_train, 1);K = zeros(n, n);for i = 1:nfor j = 1:ndiff = X_train(i, :) - X_train(j, :);K(i, j) = exp(-norm(diff)^2/(2*sigma^2));endendOmega = [K + eye(n)/C, ones(n, 1); ones(n, 1)', 0];Y = [y_train; 0];sol = Omega \ Y;model.alpha = sol(1:n);model.b = sol(n+1);model.X_train = X_train;model.sigma = sigma;
endfunction y_pred = predictLSSVR(model, X)% 预测函数 (简化版)n_train = size(model.X_train, 1);n_test = size(X, 1);y_pred = zeros(n_test, 1);for i = 1:n_tests = 0;for j = 1:n_traindiff = X(i, :) - model.X_train(j, :);k_val = exp(-norm(diff)^2/(2*model.sigma^2));s = s + model.alpha(j)*k_val;endy_pred(i) = s + model.b;end
end
程序功能说明
1. 核心类:GWO_LSSVR
-
属性:
-
C,gamma:LSSVR的正则化参数和核参数 -
X_train,y_train:训练数据 -
X_test,y_test:测试数据 -
model:训练好的LSSVR模型 -
fitnessHistory:优化过程中的适应度历史
-
-
方法:
-
train():使用灰狼算法优化LSSVR参数 -
evaluateFitness():评估给定参数的适应度(MSE) -
trainLSSVR():训练LSSVR模型 -
predict():使用训练好的模型进行预测 -
plotResults():可视化优化过程和结果 -
evaluateModel():计算模型性能指标
2. 灰狼优化算法实现
-
初始化:在参数空间内随机生成狼群位置
-
迭代优化:
-
更新参数
a线性递减 -
根据α、β、δ狼的位置更新其他狼的位置
-
边界约束处理
-
-
适应度评估:使用测试集MSE作为适应度函数
3. LSSVR实现
-
核函数:径向基函数(RBF)
-
模型训练:求解线性方程组得到拉格朗日乘子和偏置
-
预测:使用核函数和学到的参数进行预测
4. 性能评估
-
计算多种回归指标:MSE、RMSE、MAE、R²
-
可视化预测结果与真实值对比
-
误差分布直方图
-
优化过程收敛曲线
5. 演示脚本
-
生成非线性回归数据集
-
划分训练集和测试集
-
训练GWO-LSSVR模型
-
与网格搜索调参方法进行性能比较
使用示例
% 运行演示程序
demoGWO_LSSVR();% 自定义数据使用
% 1. 准备数据
load('your_dataset.mat'); % 加载数据
X = features; % 特征矩阵
y = targets; % 目标变量% 2. 划分数据集
splitRatio = 0.8;
nSamples = size(X, 1);
idx = randperm(nSamples);
nTrain = round(splitRatio * nSamples);
X_train = X(idx(1:nTrain), :);
y_train = y(idx(1:nTrain), :);
X_test = X(idx(nTrain+1:end), :);
y_test = y(idx(nTrain+1:end), :);% 3. 创建并训练模型
gwo_lssvr = GWO_LSSVR(X_train, y_train, X_test, y_test);
[C_opt, gamma_opt] = gwo_lssvr.train(30, 100); % 30只狼,100次迭代% 4. 评估模型
metrics = gwo_lssvr.evaluateModel();% 5. 可视化结果
gwo_lssvr.plotResults();
参数调整建议
-
灰狼算法参数:
-
numWolves:狼群数量(通常10-50) -
maxIter:最大迭代次数(通常50-500) -
增加狼群数量和迭代次数可提高优化精度,但会增加计算时间
-
-
LSSVR参数范围:
-
lb和ub:参数搜索范围 -
根据具体问题调整范围:
-
对于简单问题:C∈[0.1, 100], gamma∈[0.01, 10]
-
对于复杂问题:C∈[1, 1000], gamma∈[0.1, 100]
-
-
-
核函数选择:
-
当前实现使用RBF核
-
可扩展其他核函数(多项式核、线性核等)
参考代码 使用灰狼算法优化的LSSVR程序 www.youwenfan.com/contentcns/100255.html
扩展功能建议
1. 多核学习支持
function model = trainMultiKernelLSSVR(obj, C, gamma1, gamma2, weight)% 训练多核LSSVR% 结合RBF核和多项式核n = size(obj.X_train, 1);K1 = zeros(n, n); % RBF核K2 = zeros(n, n); % 多项式核for i = 1:nfor j = 1:ndiff = obj.X_train(i, :) - obj.X_train(j, :);K1(i, j) = exp(-norm(diff)^2/(2*(1/sqrt(2*gamma1))^2));K2(i, j) = (1 + obj.X_train(i, :) * obj.X_train(j, :)') ^ gamma2;endend% 加权组合核矩阵K = weight * K1 + (1-weight) * K2;% 其余部分与单核相同...
end
2. 并行计算加速
function fitness = evaluateFitnessParallel(obj, wolves, C_list, gamma_list)% 并行评估适应度nWolves = size(wolves, 1);fitness = zeros(nWolves, 1);parfor i = 1:nWolvesC = wolves(i, 1);gamma = wolves(i, 2);model = obj.trainLSSVR(C, gamma);y_pred = obj.predict(model, obj.X_test);fitness(i) = mean((obj.y_test - y_pred).^2);end
end
3. 混合优化策略
function [C_opt, gamma_opt] = hybridOptimization(obj, numWolves, maxIter)% 混合灰狼算法和局部搜索% 第一阶段:灰狼全局优化[C_gwo, gamma_gwo] = obj.train(numWolves, maxIter/2);% 第二阶段:局部搜索精调C_range = linspace(C_gwo*0.5, C_gwo*1.5, 10);gamma_range = linspace(gamma_gwo*0.5, gamma_gwo*1.5, 10);best_fitness = inf;for i = 1:length(C_range)for j = 1:length(gamma_range)fitness = obj.evaluateFitness(C_range(i), gamma_range(j));if fitness < best_fitnessbest_fitness = fitness;C_opt = C_range(i);gamma_opt = gamma_range(j);endendend
end
4. 时间序列预测应用
function forecast = timeSeriesForecast(obj, data, forecastHorizon)% 时间序列预测% data: 历史时间序列数据% forecastHorizon: 预测步长% 构造特征和目标X = [];y = [];for i = 1:(length(data)-obj.nFeatures)X(i, :) = data(i:i+obj.nFeatures-1);y(i) = data(i+obj.nFeatures);end% 训练模型obj.X_train = X;obj.y_train = y;obj.train(20, 50);% 递归预测lastValues = data(end-obj.nFeatures+1:end);forecasts = zeros(forecastHorizon, 1);for i = 1:forecastHorizonnextPred = obj.predict(obj.model, lastValues);forecasts(i) = nextPred;lastValues = [lastValues(2:end); nextPred];endforecast = forecasts;
end
实际应用建议
-
数据预处理:
-
标准化/归一化特征
-
处理缺失值和异常值
-
特征选择和降维(PCA等)
-
-
参数调优:
-
根据问题复杂度调整搜索范围
-
使用敏感性分析确定关键参数
-
考虑参数间的耦合关系
-
-
模型集成:
-
结合多个优化算法的优点
-
集成多个LSSVR模型的预测
-
使用bagging或boosting方法
-
-
实时应用:
-
增量学习更新模型
-
滑动窗口处理新数据
-
模型压缩减少计算开销
