简介MATLAB实现的多核最小二乘支持向量机LSSVM优化代码包面向机器学习和数据挖掘方向的研究者与学生用于处理复杂非线性数据的分类与回归建模。代码引入粒子群优化PSO对多核LSSVM的正则化参数和核函数参数进行自动寻优从而提升模型泛化能力与预测精度。压缩包共7个文件以4个mat数据文件和3个m脚本文件组成整体仅10KBmat文件存放训练/测试样本及对应标签m文件则涵盖PSO参数优化、多核LSSVM模型构建以及自适应改进函数等核心模块。目前已有877人学习下载。资源在原单核LSSVM基础上扩展到多核版本并结合PSO自动调参既适合入门者理解多核学习与智能优化算法结合的原理也可为进阶研究提供可运行的轻量级代码参考便于进行模型对比与二次开发是一份实用价值较高的算法学习资料。1. 多核LSSVM与PSO-LSSVM到底在解决什么问题做回归预测建模时经常遇到一个很尴尬的情况LSSVM的单一RBF核在局部区间拟合得很好换一个区间就出现系统性偏差把核换成多项式整体趋势对了局部细节又丢掉。与其来回换核函数不如把多个核组合起来。多核LSSVM让核矩阵变成多个基核的加权和但代价是待调参数从两三个变成五六个网格搜索直接爆炸。PSO-LSSVM的做法就是用粒子群算法去搜这一组参数包括LSSVM的正则化参数γ、基核宽度σ和各核权重w。下面这套代码可以直接在MATLAB里运行不依赖matlab优化工具箱手写PSO主循环用3折交叉验证做适应度评估最后输出预测结果和最优核权重。适合需要在回归、软测量这类任务里快速搭一套可解释调参模型的人。2. 多核LSSVM的建模原理与待优化参数2.1 从SVM到LSSVM二次规划换成线性方程组LSSVMLeast Squares Support Vector Machine和标准SVM的核心差别在于误差项的写法。标准SVM用ε-不敏感带加松弛变量求解的是带不等式约束的二次规划LSSVM直接把误差写成等式约束min J(w, e) 1/2 * ||w||^2 γ/2 * Σ e_i^2 s.t. y_i w^T * φ(x_i) b e_i, i 1,...,n每个等式约束配一个拉格朗日乘子α_iKKT条件整理后得到一个 (n1) 维线性方程组[ 0 1^T ] * [ b ] [ 0 ] [ 1 Ω γ^{-1} I] [ α ] [ y ]其中 Ω(i,j) K(x_i, x_j)1 是 n 维全 1 向量。这个方程组直接用 MATLAB 的反斜杠\求解即可不需要任何优化迭代。这里要注意 γ 的物理含义γ 出现在核矩阵对角线的偏移项1/γ上γ 越大对角线偏移越小模型对训练样本的拟合越激进γ 越小模型越平滑。这也是后面对 γ 做对数编码的原因——它实际起作用的是数量级。LSSVM 没有稀疏解几乎每个训练样本的 α_i 都不为零这是它和标准 SVM 的最大差异。对于回归预测这类不追求支持向量稀疏性的任务这个代价完全可接受换来的是训练速度快一个量级。预测公式为y_new Σ α_i * K(x_i, x_new) b在 MATLAB 里这等价于测试核矩阵乘上 α 再加 b。理解了这条主线多核 LSSVM 就只剩一个改动把 K 从单一核换成多个核的加权和。2.2 多核LSSVM的核组合机理多核 LSSVM 的做法是把核矩阵写成若干基核的加权和K(i,j) Σ w_k * K_k(x_i, x_j), w_k ≥ 0, Σ w_k 1权重非负且和为 1 不是随便定的。核函数要满足 Mercer 条件才能保证对应的核矩阵半正定而非负系数的线性组合能保持这一性质。如果你把某个权重设成负数核矩阵可能出现负特征值LSSVM 的线性方程组求解就不稳定了。选择什么样的基核跟数据分布相关。我一般会固定三到五个差异足够大的基核而不是堆十几个。差异大意味着它们的特征映射不在同一个方向上组合起来才有意义。常见组合是 RBF 核 多项式核 线性核RBF 负责局部细节多项式负责全局趋势线性核在特征维度不高时作为保底项。如果三个核本身高度相似多核组合的收益就很有限。多核真正带来的问题是参数维度上升。RBF 核有宽度 σ多项式核有次数和偏移每个核还有权重 w。再加上 LSSVM 本身的 γ参数空间很容易到 5 维以上。网格搜索在这种维度下不可行所以需要 PSO 这种不依赖梯度、对目标函数解析形式没有要求的实数优化算法。2.3 PSO-LSSVM的待优化参数与选型把整套 PSO-LSSVM 的待优化参数列出来方便后面写代码时一一对应参数作用典型范围编码方式γLSSVM 正则化控制拟合与平滑10^-3 ~ 10^3log10(γ)σ_RBFRBF 核宽度控制局部作用半径0.1 ~ 5可视特征尺度调整w_1, w_2, w_3各基核权重≥ 0和为 1softmax 映射nPPSO 种群大小20 ~ 50直接设定nIterPSO 迭代代数20 ~ 100直接设定PSO 的选型理由很直接目标函数是交叉验证 RMSE它关于 γ 和 σ 不可导而且存在大量局部极小。粒子群的全局搜索能力不依赖初始点实现也就几十行代码。相比 matlab 自带的优化工具箱手写 PSO 能把每一次迭代里核矩阵的计算复用做到极致——这一点在下一章会体现出来。3. 多核LSSVM的MATLAB核心代码与核矩阵预计算3.1 数据准备与基核矩阵预计算先造一个带非线性、非平稳特性的回归数据集方便验证多核的效果。这里用 150 个样本输入是一维的输出包含周期项、趋势项和噪声clear; clc; rng(42); N 150; X linspace(-4, 4, N); Y 0.8 * sin(2*X) ./ (0.6 abs(X)) 0.25*X 0.1*randn(N, 1); % 训练/测试 60/40 划分 idx_train randperm(N, round(0.6*N)); idx_test setdiff(1:N, idx_train); Xtr X(idx_train); Ytr Y(idx_train); Xte X(idx_test); Yte Y(idx_test);数据生成里sin(2*X)./(0.6abs(X))保证了 X0 处不会出现除零同时周期和趋势叠加的结构会让单核模型“偏科”RBF 核能抓住局部周期但趋势外推弱多项式核趋势好但局部拟合不细腻。这正是多核发挥作用的场景。接下来预计算三个基核矩阵。这一步是整个方案里最关键的提速点PSO 每次评估适应度时只需要对三个固定好的矩阵做加权求和不需要重新算核函数。% 训练集核矩阵 K_rbf exp(-pdist2(Xtr, Xtr).^2 / (2*1.0^2)); % RBF核宽度固定为1.0 K_poly (Xtr * Xtr 1).^2; % 多项式核次数2偏移1 K_lin Xtr * Xtr; % 线性核 % 测试集对训练集的核矩阵用于最终预测 K_rbf_te exp(-pdist2(Xte, Xtr).^2 / (2*1.0^2)); K_poly_te (Xte * Xtr 1).^2; K_lin_te Xte * Xtr; K_cell {K_rbf, K_poly, K_lin}; K_cell_te {K_rbf_te, K_poly_te, K_lin_te};这里把 RBF 宽度固定为 1.0 是有意的。第一版先让 PSO 只搜 γ 和核权重核宽度后面再放开。pdist2 是 MATLAB 基础函数不需要任何工具箱。三个核矩阵的尺寸都是 90×90训练集约 90 个样本测试核矩阵是 60×90内存开销完全可以忽略。3.2 LSSVM训练与预测的MATLAB函数核心求解器写在两个局部函数里。MATLAB R2016b 及以后版本允许把函数直接放在脚本末尾不需要单独建文件function [alpha, b] lssvm_train(K, y, gamma) n size(K, 1); A [0, ones(1, n); ones(n, 1), K eye(n)/gamma]; sol A \ [0; y(:)]; b sol(1); alpha sol(2:end); end function yp lssvm_predict(K_pred, alpha, b) yp K_pred * alpha b; endA \ [0; y]就是求解那个 (n1) 维线性方程组。矩阵 A 的第一行第一列是 0对应约束 Σα_i 0右下角块K eye(n)/gamma就是前面公式里的Ω γ^{-1}I。gamma 越大eye(n)/gamma越小核矩阵越接近原始 K模型越容易过拟合。预测函数里K_pred是测试样本与训练样本的核矩阵尺寸 m×n乘上 α 再加上偏置 b 得到 m×1 预测值。注意y(:)先把标签拉成列向量避免行向量输入导致广播维度出错。这是 MATLAB 里很常见的坑——如果你从表格或元胞数组里取数据y 可能是行向量。3.3 核权重的softmax解码避免边界约束直接优化 w1、w2、w3 有个麻烦PSO 是无约束优化粒子的位置更新后很可能跑出 [0,1] 范围导致权重和不为 1。虽然可以用“钳制归一化”处理但边界上的梯度信息很别扭。我一般用 softmax 把无约束的实数 t 映射成有效权重function w softmax_weights(t) z [t(:), 0]; % 第3个核的t固定为0作为基准 ez exp(z - max(z)); % 减max防溢出 w ez / sum(ez); end这个函数接收一个 2 维向量 t输出 3 维权重且恒满足 w_k ≥ 0、Σw_k 1。PSO 搜索的其实是 t1、t2而不是权重本身。理由是 softmax 的输入落在实数轴任意位置都能映射到合法权重粒子在任何区域的解都可用。3.4 待搜索参数与范围设置这一版 PSO 搜索 3 个维度每个维度的意义和边界如下维度含义搜索范围解码方式x(1)log10(γ)[-3, 3]γ 10^x(1)x(2)t1控制 RBF 核权重[-10, 10]softmax 映射x(3)t2控制多项式核权重[-10, 10]softmax 映射t 的范围给 [-10, 10] 看似很宽但 softmax 输出对输入差值不敏感±10 已经能把权重压缩到接近 0 或接近 1不需要精确搜索边界。log10(γ) 范围对应 γ ∈ [0.001, 1000]基本覆盖了回归任务里“欠平滑到过拟合”的完整区间。4. PSO寻优多核LSSVM的MATLAB主程序与交叉验证4.1 粒子编码与PSO参数范围粒子向量为pos [log10(γ), t1, t2]维度 dim 3。PSO 的固定参数设置如下nP 30; nIter 30; dim 3; lb [-3, -10, -10]; ub [3, 10, 10]; pos rand(nP, dim) .* (ub - lb) lb; % 均匀初始化 vel zeros(nP, dim); % 初始速度为0惯性权重在这里需要随迭代线性衰减从 0.9 降到 0.4。前期大惯性权重保证全局探索后期小权重收敛到局部精细搜索。c1 和 c2 都取 1.5这是 PSO 经验参数里最常用的默认值对 3~6 维问题表现稳定。在进行 PSO 迭代之前要先固定交叉验证的折划分。否则每次调用适应度函数都随机重划分同一个粒子两次评估的 RMSE 会因数据划分不同而抖动PSO 会误把噪声当作适应度变化收敛明显变差。4.2 适应度函数3折交叉验证的RMSE我用 3 折交叉验证计算 RMSE 作为适应度。折数太少评估噪声大折数太多每次要训练 5 次 LSSVM30 个粒子迭代 30 代就是 4500 次训练没必要。3 折是效率与稳定的平衡点。fold mod(randperm(length(Ytr)), 3) 1; % 固定划分保证可复现 function rmse_cv fitness_A(x, Ytr, K_cell, fold) gamma 10^x(1); w softmax_weights(x(2:3)); % 核矩阵加权组合 K zeros(size(K_cell{1})); for k 1:3 K K w(k) * K_cell{k}; end sse 0; for f 1:3 ite (fold f); % 测试折 itr ~ite; % 训练折 [alpha, b] lssvm_train(K(itr, itr), Ytr(itr), gamma); yp K(ite, itr) * alpha b; sse sse sum((Ytr(ite) - yp).^2); end rmse_cv sqrt(sse / length(Ytr)); end核矩阵组合只做三次矩阵数乘再加起来共用size(K_cell{1})的零矩阵。三个基核预计算后适应度评估里最重的运算是那个 60×60 线性方程组的求解速度很快。归一化 RMSE 用所有样本的总误差除以总样本数而不是把每折的 RMSE 做算术平均这样折大小不一致时评估更稳。提示如果适应度出现 NaN先检查 gamma 是否解析为 Inf 或 0。log10 编码一般不会但粒子速度过大时位置可能溢出边界钳制要加上。4.3 PSO主循环与边界钳制完整的 PSO 主循环如下。每代先评估所有粒子的适应度更新个体最优 pbest 和全局最优 gbest再用标准速度更新公式推进粒子pbest pos; pbest_fit inf(nP, 1); gbest pos(1, :); gbest_fit inf; for it 1:nIter w_i 0.9 - (0.9 - 0.4) * (it - 1) / (nIter - 1); % 惯性权重线性衰减 for i 1:nP fit fitness_A(pos(i, :), Ytr, K_cell, fold); if fit pbest_fit(i) pbest_fit(i) fit; pbest(i, :) pos(i, :); end if fit gbest_fit gbest_fit fit; gbest pos(i, :); end end r1 rand(nP, dim); r2 rand(nP, dim); vel w_i * vel 1.5 * r1 .* (pbest - pos) 1.5 * r2 .* (gbest - pos); pos pos vel; pos min(max(pos, lb), ub); % 边界钳制 fprintf(iter %d, best RMSE_cv %.5f\n, it, gbest_fit); end这里速度更新里1.5 * r1 .* (pbest - pos)是认知项1.5 * r2 .* (gbest - pos)是社会项。两者分别把粒子拉向自身历史最优和全局最优。注意速度vel没有显式上限所以位置钳制到[lb, ub]是必要的否则粒子可能直接飞到 log10(γ) 20 这种无意义区域softmax 输入过大时 exp 溢出。迭代 30 代后gbest 里保存的就是最优的[log10(γ), t1, t2]。可以在迭代过程中把gbest_fit画出来观察收敛情况正常情况下前 5 代 RMSE 快速下降后 20 代缓慢波动波动幅度逐步收窄。4.4 用gbest重建多核LSSVM并对比单核模型PSO 得到的最优参数还需要重建模型在留出的测试集上做一次最终评估。同时用相同 γ 跑一个纯 RBF 核的 LSSVM 作为 baseline判断多核到底带来多少收益gamma_best 10^gbest(1); w_best softmax_weights(gbest(2:3)); % 重建最优多核矩阵 K_best zeros(size(K_cell{1})); K_best_te zeros(size(K_cell_te{1})); for k 1:3 K_best K_best w_best(k) * K_cell{k}; K_best_te K_best_te w_best(k) * K_cell_te{k}; end % 多核模型 [alpha, b] lssvm_train(K_best, Ytr, gamma_best); Yp lssvm_predict(K_best_te, alpha, b); rmse_mk sqrt(mean((Yte - Yp).^2)); % 单RBF核baseline [alpha_rbf, b_rbf] lssvm_train(K_cell{1}, Ytr, gamma_best); Yp_rbf lssvm_predict(K_cell_te{1}, alpha_rbf, b_rbf); rmse_rbf sqrt(mean((Yte - Yp_rbf).^2)); fprintf(多核LSSVM test RMSE %.5f, weights [%.3f %.3f %.3f]\n, ... rmse_mk, w_best); fprintf(单RBF核LSSVM test RMSE %.5f\n, rmse_rbf);这里有个细节很容易被忽略baseline 用的 γ 和 PSO 搜出来的最优 γ 相同。这样对比的是“多核权重”带来的增益而不是“调参带来的增益”。如果你发现单 RBF 核的 RMSE 反而更低先不要急着否定多核可能是三个基核的宽度设置不合理RBF 核宽度 1.0 在这个数据上恰好太合适了。这种场景下真正有价值的做法是把 σ 也放进粒子搜索下一章展开。5. PSO-LSSVM的验证与进阶核宽度搜索与稳健求解5.1 先做两组对比验证核权重起了作用多核权重容易给人“参数越多越好”的错觉。要验证权重真正有用建议做两组对照第一组是“等权重 最优 γ”第二组是“PSO 搜出的权重 最优 γ”。如果两者测试 RMSE 差距很小说明三个基核本身冗余多核收益有限如果差距明显说明 PSO 确实找到了更合理的核组合方向。我在实际项目中更关注权重的分布形态某个核权重被压到 0.05 以下基本可以判断这个核贡献不足下一步应该换别的基核而不是继续加权重。5.2 核矩阵数值问题与稳健求解LSSVM 的线性方程组在 gamma 偏大时容易病态尤其当样本存在重复或高度相关特征时。A\b内部会做列主元消去大多数情况没问题但如果求解结果出现符号抖动或预测值异常大可以把核矩阵显式对称化H K eye(n) / gamma; H (H H) / 2; % 消除数值不对称 sol [0, ones(1,n); ones(n,1), H] \ [0; y(:)];还可以观察核矩阵的条件数。条件数超过 1e12 时建议先对输入特征做标准化再考虑调小 gamma。PSO 里对 gamma 做对数编码本身就能缓解这个问题因为 log10 范围限制了粒子不可能直接取到 1e8 这种危险值。5.3 让PSO连同核宽度一起搜索固定核宽度只是第一版方案。当核宽度本身需要调优时把粒子扩展到 5 维[log10(γ), log10(σ1), log10(σ2), t1, t2]。σ1 是 RBF 核宽度σ2 是多项式核的偏移系数。这时核矩阵每次都要按 σ 重新计算之前预计算的方案失效需要用缓存避免重复计算cache containers.Map(KeyType, char, ValueType, any); key sprintf(s%.4f, sigma); if cache.isKey(key) K_rbf cache(key); else K_rbf exp(-pdist2(Xtr, Xtr).^2 / (2*sigma^2)); cache(key) K_rbf; end核矩阵缓存让同一个 σ 在多次适应度评估中只算一次。粒子数 30、迭代 30 代时缓存命中率通常能到 60% 以上。σ 的搜索范围我一般设为 log10 域 [-1, 1]对应 σ ∈ [0.1, 10]覆盖从“几乎线性”到“高度局部”的完整过渡这个范围再放宽后你在测试集上看到的 RMSE 下降会比单纯调核权重更大一些。本文还有配套的精品资源点击获取