dtcms

支持向量机：鲸鱼算法优化支持向量机SVM预测原理与MATLAB代码

wuchangjian2021-11-13 07:47:05编程学习

目录

1. 鲸鱼优化算法

1.1 包围猎物

1.4 算法步骤

2 支持向量机基本理论

2.1 支持向量机简介

2.2 支持向量机回归SVR

3. 优化原理

3.1 优化模型构建

3.2 优化步骤

4 优化实例

4.1 实验数据

4.2 实验结果

4.3 MATLAB代码

1. 鲸鱼优化算法

    鲸鱼优化算法（Whale Optimization Algorithm，WOA），于2016 年由澳大利亚格里菲斯大学学者Mirjalili 等人提出，是一种模拟座头鲸特殊捕食行为的新型群智能优化算法。WOA 采用随机或最佳搜索代理来模拟捕猎行为，并使用螺旋来模拟座头鲸的泡泡网攻击机制。座头鲸喜欢捕食一群靠近水面的小鱼和小虾，遇到猎物后会先向下俯冲大约12 m，然后开始在猎物周围制造螺旋形气泡，最后游向水面捕食猎物。这个过程提取出三个数学模型即围捕猎物、螺旋气泡网捕食和寻找猎物。

1.1 包围猎物

    座头鲸在狩猎时会选择包围猎物，当p <0.5且| A | < 1时，执行包围猎物行为。p是[0，1]范围内产生的随机数。鲸鱼算法视当前最优候选解为猎物位置或者距离目标猎物最近的鲸鱼个体位置，其他候选解将向种群最优位置靠近。式（1）～（2）模拟了这个过程：

| D=\left|C X^{*}(t)-X(t)\right|       　　   （1）

  X(t)=X^{*}(t)-A D   　      　   （2）

    式中，A 和C 为系数；D是狩猎对象和当下鲸鱼个体之间的距离；t 为当前迭代次数，X*（t）为当前最佳的鲸鱼位置向量，X（t）为当前鲸鱼的位置向量。A和C 由式（3）～（5）计算得到：

A=2 a r_{1}-a　　　　     　　　 　（3）

　C=2 r_{2} \mathrm{~m}　　　　　　　　　　（4）

 a=2-2 t / T_{\max }         　  　            （5）

    式中，r1和r2 为（0，1）区间上的随机数，a 的值从2线性减小到0，t 为当前的迭代次数，Tmax为最大迭代次数。

1.2 狩猎行为

　　当p≥0.5，执行螺旋轨迹狩猎行为。在座头鲸的狩猎行为中，当座头鲸搜索猎物时，对目标猎物的位置以螺旋形式作为运动轨迹来捕猎食物，其数学模型如式（6）所示：

X(t+1)=X^{*}(t)+D_{p} \mathrm{e}^{b l} \cos (2 \pi l)　　 （6）

　　式中，D_{P}=\left|X^{*}(t)-X(t)\right|为鲸鱼个体与猎物之间的距离。b 为一个用来定义螺线形状的常数。

    鲸鱼以螺旋形路径游向猎物的同时还要收缩包围圈。假设鲸鱼有50% 的概率选择收缩包围机制和1% ～ 50% 的概率选择螺旋模型来更新位置，其数学模型如式（7）所示，P 为[0，1] 之间的随机数：

X(t+1)=\left\{\begin{array}{l} X^{*}(t)-A D \quad P<0.5 \\ D_{p} \mathrm{e}^{b l}(2 \pi l)+X^{*}(t) \quad P \geqslant 0.5 \end{array}\right.    （7）

1.3 搜索猎物

    当p <0.5且| A |>1时，执行随机猎物搜寻行为。此时随机选择一个鲸鱼个体更新下一个鲸鱼的位置，鲸鱼搜索猎物模型如式（8）～（9）所示：

　D=\left|C X_{\text {rand }}-X(t)\right|　　　　　　 　    　（8）

 X(t+1)=X_{\text {rand }}-A D 　　　　　 　 　　（9）

    式中，X_{rand} 为随机选择的鲸鱼位置。算法规定当A≥1时，在种群中随机选取出一个随机代理，根据随机选择的鲸鱼个体位置来更新其他个体的位置，从而使鲸鱼偏离当前目标，寻找一个更加合适的猎物，提高WOA 算法的全局搜

索能力。

1.4 算法步骤

鲸鱼优化算法的步骤为：

步骤1：初始化算法参数a、b，鲸鱼的种群规模和最大迭代次数；

步骤2：计算每只鲸鱼个体的适应度值，并进行比较，确定当前适应度值最优的鲸鱼个体，定义为X*（t）；

步骤3：进入算法主循环，局部搜索公式更新鲸鱼位置；

步骤4：再次对整个鲸鱼种群进行评价，找到全局最优的鲸鱼个体及其位置；

步骤5：若达到算法最大迭代次数，则循环结束；否则转到步骤2，继续进行算法迭代。

步骤6：输出全局最优个体X*（t）。

2 支持向量机基本理论

2.1 支持向量机简介

支持向量机（Support Vector Machine，SVM）是在20世纪90年代初提出的，它建立在统计学习理论的VC维概念和结构风险最小原理基础上，可以根据有限的样本信息实现计算量和计算能力之间的平衡性，是归属于统计学的新兴分支。其独特的对于数据的归类识别能力，使其能够有效地应用到机器学习的相关领域进行预测分析。

对于群体在线性可分情况下的分类情况示意如图1所示，由图可知，两类不同的样本被分类线H分离，且存在两条与H平行的直线H1和H2分别经过了两组样本最靠近分类线的样本。针对此结果，可得出最优分类线的定义：针对两样本存在一条分类线能将其正确分类，而且H1和H2的距离最大。

 图1. 最优分类线示意图

针对以上基本理论，将分类线H的方程设为wx+b=0，并将该方程归一化，有：

y_{i}\left(w \cdot x_{i}+b\right)-1 \geqslant 0,(1)

在方程（1）的样本集（xi ,yi）中有以下关系：xi∈Rd，yi={-1,+1}，i=1，2，…，n。据此可知分类间隔（margin）为2/ ||w||。若要在计算中实现分类间隔最大化，则通过约束实现||w||2最小即可。

据此，最优分类线在数学上的定义可描述为满足式（1）且使||w||2最小的方程表达式，而在H1和H2上进行训练的样本点则为对应的支持向量。

以上分析是基于二维平面概念而言的，那么在VC维空间里，则所分析的样本理论上分布在一个超球范围内，同样存在一超平面能够将超球内的样本进行有效分类，将超平面表示为 f（x,w,b）=sgn（wx+b）。若要实现有效分类，则该函数满足以下关系：

h \leqslant \min \left(A_{2} R_{2}, N\right)+1,(2)

其中，R和N分别为超球半径和空间的维数，且在超球中存在条件||w||≤A。与二维系统同理，当w最小时，VC维可取得最小值。

为了实现对最优分类面的求解，通常引用拉格朗日优化方法，利用该方法将此问题转化为对偶问题，可得：

Q(\alpha)=\sum_{i=1}^{n} \alpha_{i}-\frac{1}{2} \sum_{i, j=1}^{n} \alpha_{i} \alpha_{j} y_{i} y_{j}\left(x_{i} \cdot x_{j}\right),(3)

式（3）中的αi对应于每个数据样本中的拉格朗日乘子，且该因子有以下关系：

\sum_{i=1}^{n} \alpha_{i} y_{i}=0, \alpha_{i} \geq 0, i=1,2, \cdots, n,(4)

由于式（3）受到不等式的约束，所以能够求解出该方程的唯一解。同时，对于系列解的结果中存在少部分αi不为零的数据样本的合集便组成了支持向量。基于以上解的结果，可知分类曲面所对应的最优函数表达式为：

f(x)=\operatorname{sgn}(w \cdot x+b)=\operatorname{sgn}\left(\sum_{i=1}^{n} \alpha_{i}^{*} y_{i}\left(x_{i} \cdot x\right)+b^{*}\right),(5)

据上述分析可知，大部分的求解结果αi均为零，所以式（5）中的有效求和只是有效支持向量的解的和。对应的分类阈值b*可通过将式（5）代入任一支持向量求得。

2.2 支持向量机回归SVR

支持向量机回归（Support Vector Regression ， SVR）的主要思想是，引入核函数将样本空间中非线性回归映射为高维空间中的线性回归问题，即寻找一个最优分类超平面，使得所有训练样本离该最优超平面的误差最小。回归目的是得到一个能够尽量拟合训练集样本的模型，通常用的方法是构建一个样本标签与模型预测值的损失函数，使损失函数最小化，从而确定模型。

对于给定的样本数据集合𝑉 = (𝑥𝑖, 𝑦𝑖), 𝑖 = 1,2, … , 𝑛, 𝑥𝑖∈𝑅, 𝑦𝑖∈𝑅。

式中， 𝑥为输入样本值，𝑦为对应的输出值，𝑛为训练样本数。建立高维特征空间的回归模型：

y=\omega \cdot \phi(x)+b,(6)

参考文献：[1]罗默涵.基于优化支持向量机方法的风电场风速预测研究[J].机电信息,2021(26):6-7.

3. 优化原理

3.1 优化模型构建

优化参数选取为支持向量机的惩罚参数C和核惩罚参数g。在使用支持向量机预测SVR进行建模时，惩罚参数C与核参数g对模型的预测性能有着重要的影响。惩罚参数C用于权衡损失的权重，核参数g影响核函数的径向作用范围，决定训练样本数据的范围和分布特性。适应度值选用均方误差（Mean Square Error，MSE）来评估：

\text { Fitness }=\operatorname{argmin}\left(M S E_{\text {pridect }}\right),(1)

MSE误差越小表明预测的数据与原始数据重合度越高。

3.2 优化步骤

1）读取数据，进行数据预处理，并构建SVR模型。

2）初始化鲸鱼优化算法的种群参数，迭代次数，以及确定C、g的取值范围。

3）确定鲸鱼优化算法的适应度函数，计算适应度函数值，并且评估座头鲸的适应能力。依据鲸鱼优化算法的原理，寻找最优的适应度值，以及座头鲸所获得的最佳位置。

4）根据最佳鲸鱼优化个体位置，相应得到参数C、g的最优值。

5）将最优参数C、g赋给支持向量机，进行训练，得到优化的支持向量机预测模型。

6）输入测试样本，优化的SVR预测模型将输出测试样本的预测值，与优化前的支持向量机预测对比和分析误差。

4 优化实例

4.1 实验数据

利用优化的支持向量机算法，对2000组随机噪声的模拟数据进行非线性预测实验，部分数据如下表所示：

4.2 实验结果

 

 

4.3 MATLAB代码

获取全部源代码

WOA.m鲸鱼优化算法

% 鲸鱼优化算法
function WOA()

SearchAgents_no = 30;
Max_iter = 30;
lb=[0.01 0.0001];
ub=[1000 1000];
dim =2;

% initialize position vector and score for the leader，初始化位置向量和领导者得分
Leader_pos=zeros(1,dim);
Leader_score=inf; 

%Initialize the positions of search agents
Positions=initialization(SearchAgents_no,dim,ub,lb);%Positions，存放数个个体的多维位置。

Convergence_curve=zeros(1,Max_iter);%Convergence_curve收敛曲线

t=1;% Loop counter

% 开始循环
while t<Max_iter+1
        
    for i=1:size(Positions,1)%对每个个体一个一个检查是否越界
        
        % Return back the search agents that go beyond the boundaries of
        % the search space，返回超出搜索空间边界的搜索代理
        Flag4ub=Positions(i,:)>ub;
        Flag4lb=Positions(i,:)<lb;
        Positions(i,:)=(Positions(i,:).*(~(Flag4ub+Flag4lb)))+ub.*Flag4ub+lb.*Flag4lb;%超过最大值的设置成最大值，超过最小值的设置成最小值
        
        % Calculate objective function for each search agent，目标函数值的计算
        fitness=fobj(Positions(i,:));
        
        % Update the leader
        if fitness<Leader_score % Change this to > for maximization problem
            Leader_score=fitness; % Update alpha
            Leader_pos=Positions(i,:);
        end
        
    end
    
    a=2-t*((2)/Max_iter); % a decreases linearly fron 2 to 0 in Eq. (2.3)
    
    % a2 linearly dicreases from -1 to -2 to calculate t in Eq. (3.12)，有疑问？
    a2=-1+t*((-1)/Max_iter);
    
    % Update the Position of search agents，参数更新 
    for i=1:size(Positions,1)
        r1=rand(); % r1 is a random number in [0,1]
        r2=rand(); % r2 is a random number in [0,1]
        
        A=2*a*r1-a;  % Eq. (2.3) in the paper
        C=2*r2;      % Eq. (2.4) in the paper
                
        b=1;               %  parameters in Eq. (2.5)
        l=(a2-1)*rand+1;   %  parameters in Eq. (2.5)
        
        p = rand();        % p in Eq. (2.6)
        
        for j=1:size(Positions,2)%对每一个个体地多维度进行循环运算
            
            if p<0.5%收缩包围机制 
                if abs(A)>=1
                    rand_leader_index = floor(SearchAgents_no*rand()+1);%floor将 X 的每个元素四舍五入到小于或等于该元素的最接近整数
                    X_rand = Positions(rand_leader_index, :);
                    D_X_rand=abs(C*X_rand(j)-Positions(i,j)); % Eq. (2.7)
                    Positions(i,j)=X_rand(j)-A*D_X_rand;      % Eq. (2.8)
                    
                elseif abs(A)<1
                    D_Leader=abs(C*Leader_pos(j)-Positions(i,j)); % Eq. (2.1)
                    Positions(i,j)=Leader_pos(j)-A*D_Leader;      % Eq. (2.2)
                end
                
            elseif p>=0.5%螺旋更新位置
              
                distance2Leader=abs(Leader_pos(j)-Positions(i,j));
                % Eq. (2.5)
                Positions(i,j)=distance2Leader*exp(b.*l).*cos(l.*2*pi)+Leader_pos(j);
                
            end
            
        end
    end
    
    Convergence_curve(t)=Leader_score;
    t=t+1;

end

返回列表

上一篇：白夜追题之icpc一言南京版

下一篇：Math.random()方法确定区间[a,b]间的随机数

相关文章

详细讲解 Redis 的两种安装部署方式

详细讲解 Redis 的两种安装部署方式

Redis 是一款比较常用的 NoSQL 数据库，我们通常使用 Redis...

C语言 fgets读取中文出现乱码

...

计算机网络网络层协议之RIP协议、OSPF协议和BGP协议

...

费解的开关

...

统计了300课时的代码，得到学Python需要掌握的1024个英语单词

统计了300课时的代码，得到学Python需要掌握的1024个英语单词

对话框回复「导航」获取历史推送 与非学堂 专注青少年信息技术教育，探讨...

uni-app 80初始化总未读数

...

发表评论    

◎欢迎参与讨论，请在这里发表您的看法、交流您的观点。

Copyright Your WebSite.Some Rights Reserved.