智能优化算法
Matlab
下标从1开始！！
-
~=：相当于不等于!= -
eps表示的是一个数可以分辨的最小精度，返回1.0和下一个精度最高的双精度浮点数的差值， 即2^(-52)。
-
Inf和-Inf分别代表正无穷和负无穷
-
y =
length(x)函数计算指定向量或矩阵的长度y。如果参数变量x是向量，则返回其长度；如果参数变量是非空矩阵，则length(x)与max(size(x))等价
矩阵取值
- v(m , : )：取第m行
- v( : ，m)：取第m列
- v(x,y)：取第x行第y列
矩阵计算
- [SortFit1,Index]=
sort(Fit1)：对Fit进行排序，排序结果放入SortFit1矩阵，结果每位上的元素在原来列上的序号放入Index矩阵 - [aa,bb]=
size(A)：aa=行数，bb=列数 cumsum(x)：对矩阵x进行逐列累加，例：[a,b,c,d]=>[a,a+b,a+b+c,a+b+c+d]sum(x)：对矩阵x进行逐列求和，即每列之和
取整函数
- fix朝零方向取整，如fix(-1.3)=-1; fix(1.3)=1;
floor，顾名思义，就是地板，所以是取比它小的整数，即朝下取整，如floor(-1.3)=-2; floor(1.3)=1;floor(-1.8)=-2，floor(1.8)=1ceil，与floor相反，它的意思是天花板，也就是取比它大的最小整数，即朝上取整，如ceil(-1.3)=-1; ceil(1.3)=2;ceil(-1.8)=-1，ceil(1.8)=2round四舍五入到最近的整数，如round(-1.3)=-1;round(-1.52)=-2;round(1.3)=1;round(1.52)=2
生成随机数函数
- rand 生成均匀分布的伪随机数。分布在（0~1）之间
- rand(m,n)生成m行n列的均匀分布的伪随机数
- randn 生成标准正态分布的伪随机数（均值为0，方差为1）
- randi 生成均匀分布的伪随机整数
- randi(iMax)在闭区间[1，iMax]生成均匀分布的伪随机整数
- randi(iMax，m，n)在开区间[1，iMax]生成mXn型随机矩阵
- r = randi([iMin,iMax],m,n)在开区间[iMin，iMax]生成m*n型随机矩阵
- randperm 生成整数的随机排列
排序函数
B = sort(A) 按升序对 A 的元素进行排序。
- 如果
A是向量，则sort(A)对向量元素进行排序。 - 如果
A是矩阵，则sort(A)会将A的列视为向量并对每列进行排序，从上到下，从小到大 - 如果
A是多维数组，则sort(A)会沿大小不等于 1 的第一个数组维度计算，并将这些元素视为向量
复制函数
B = repmat(A,x,y)生成重复数矩阵
- 若A是一个数，则生成x*y的矩阵，全是A
- 若A是一个矩阵，将矩阵A复制2行3列
智能优化算法
代表人物汇总
| 算法名 | 英文名 | 代表人物 |
|---|---|---|
| 遗传算法GA | Genetic Algorithm | J.H.Holand |
| 差分进化算法DE | Differential evolution | Storn |
| 免疫算法IA | Immune algorithm | Burnet |
| 蚁群算法ACO | Ant Colony optimization | M.Dorigo，V.Maniezzo，A.Colorni |
| 粒子群算法PSO | Particle swarm optimization | James Kennedy，Rusell Eberhart |
| 模拟退火算法SA | Simulated annealing | Metropolis |
| 禁忌搜索算法TS | Tabu Search | Glover |
| 神经网络算法NN | Neural Newwork | McCulloch，Pitts，J.J.Hopfield |
算法特点汇总
| 算法名字 | 特点 | 优缺点 |
|---|---|---|
| 遗传算法 | 群体搜索策略和简单的遗传算子 | 全局搜索能力强，局部搜索能力较弱，早熟，算法收敛性无法保证 |
| 差分进化算法 | ||
| 免疫算法 | ||
| 蚁群算法 | ||
| 粒子群算法 | 收敛速度快但容易陷入局部最优解 | |
| 模拟退火算法 | ||
| 禁忌搜索算法 | ||
| 神经网络算法 |
遗传算法GA
遗传算法模拟生物在自然环境中的遗传和进化的过程，从而形成自适应全局优化搜索算法，它借用了生物遗传学的观点，通过自然选择，遗传，变异等机制，实现各个个体适应性的提高
名词解释
| 遗传学术语 | 遗传算法术语 |
|---|---|
| 群体 | 可行解集 |
| 个体 | 可行解 |
| 染色体 | 可行解的编码 |
| 基因 | 可行解编码的分量 |
| 基因形式 | 遗传编码 |
| 适应度 | 评价函数值 |
| 选择 | 选择操作 |
| 交叉 | 交叉操作 |
| 变异 | 变异操作 |
关键参数
- 群体规模Np：群体规模将影响遗传优化的最终结果以及遗传算法的执行效率。一般取10~200。
- 交叉概率Pc：交叉概率控制着交叉操作被使用的频度。一般取0.25~1.00。
- 变异概率Pm：变异的主要目的是保持群体的多样性，一般地频度的变异可防止群体中重要基因的丢失。一般取0.001~0.1。
- 遗传运算的终止进化代数G：终止遗传代数表示遗传算法运行到指定的进化代数之后就停止运行。一般取100~1000。
主要流程
选择
轮盘赌选择（Roulette Wheel Selection）
一种回放式随机采样方法。每个个体进入下一代的概率等于它的适应度值与整个种群中个体适应度值和的比例
轮盘赌选择法是最简单也是最常用的选择方法，在该方法中，各个个体的选择概率和其适应度值成比例，适应度越大，选中概率也越大。但实际在进行轮盘赌选择时个体的选择往往不是依据个体的选择概率，而是根据“累积概率”来进行选择。
例子：
- 假设群体 A：a, b, c, d, e
- 适应度 F1：0.45, 0.1, 0.1, 0.15, 0.2
- 对F1进行累加操作得到F2（累积概率）：0.45, 0.55, 0.65, 0.8, 1.0
- 生成排序随机队列MS：0.2, 0.3, 0.4 , 0.7 ,0.8 ,0.9
- 选择完成群体 B：a, a, a, d, e, e
合并新旧种群取最优
该选择策略，是在交配产生的新子代和父代中，通过合并子代和父代，并且按照个体的适应度进行排序，选择适应度最佳的前 NP 个个体作为新的种群，达到更好的收敛的效果。
例子：
- 父代个体parent : a,b,c
- 父代适应度parent_obj : 1 , 5 , 3
- 子代个体son : d , e , f
- 子代适应度son_obj : 2 6 4
- 合并个体total = a,b,c,d,e,f
- 合并适应度total_obj : 1,5,3,2,6,4
- 适应度排序order : e,b,f,c,d,a
- 前 NP=3个作为下一代：save = e,b,f
交叉
单点交叉
单点交叉通过选取两条染色体，在随机选择的位置点上进行分割并交换右侧的部分，从而得到两个不同的子染色体。
例子：
（空格处代表该处为随机生成的交叉点，往后的基因需要进行交换操作）
待交叉染色体： 11111 22222
33333 44444
交换后染色体： 11111 44444
33333 22222
君主交叉
君主交叉就是在种群内选择出最优个体，用最优个体作为君主染色体，并随机生成交叉位点，然后将君主染色体的特定位点上的基因遗传给子代染色体。
例子：
（中间空格分隔的基因代表随机生成的交叉位点，该位点上的基因需要遗传给子代染色体）
待交叉染色体： 101000 1 101
011001 0 011
交换后染色体： 101000 1 101
011001 1 011
多点奇偶交叉
- 先规定一个交叉率，然后循环奇数位
- 每次生成一个随机数，判断是否小于交叉率
- 如果是，则进行交叉
- 对二进制位上每一位取随机数，如果等于1，则与偶数位同位置的数进行交换
- 如果不是，则跳过交叉
- 如果是，则进行交叉
变异
通过变异率计算第i个个体的j个基因二进制位上会发生变异，即取反
域内随机值
描述：
对新种群内的所有染色体每个基因进行变异操作，先通过随机值和变异率去判断该基因是否会发生变异，再通过随机值在基因的取值范围内进行改变。
例子：
- 在所有染色体中随机挑选 NP * Pm 个染色体进行变异
- 假设挑选出染色体h1：11000011
- 随机挑选 L（染色体长度，基因个数） * Pm 个基因进行变异
- 挑选下标：1,2,4
- h2：00101011（下标1,2,4处基因取反）
代码
奇偶交叉GA1
%%%%%%%%%%%%%%%%%%%%标准遗传算法求函数极值%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%初始化参数%%%%%%%%%%%%%%%%%%%%%%%%%%
clear all; %清除所有变量
close all; %清图
clc; %清屏
NP=50; %种群数量
L=20; %二进制数串长度
Pc=0.8; %交叉率
Pm=0.1; %变异率
G=100; %最大遗传代数
Xs=10; %上限
Xx=0; %下限
f=randi([0,1],NP,L);
%%%%%%%%%%%%%%%%%%%%%%%%%遗传算法循环%%%%%%%%%%%%%%%%%%%%%%%%
for k=1:G
%%%%%%%%%%%%将二进制解码为定义域范围内十进制%%%%%%%%%%%%%%
for i=1:NP
U=f(i,:);
m=0;
for j=1:L
m=U(j)*2^(j-1)+m;
end
x(i)=Xx+m*(Xs-Xx)/(2^L-1);
Fit(i)= func1(x(i));
end
maxFit = max(Fit); %最大值
minFit = min(Fit); %最小值
rr = find(Fit==maxFit);
fBest = f(rr(1,1),:); %历代最优个体
xBest = x(rr(1,1));
Fit = (Fit-minFit)/(maxFit-minFit); %归一化适应度值
%%%%%%%%%%%%%%%%%%基于轮盘赌的复制操作%%%%%%%%%%%%%%%%%%%
sum_Fit=sum(Fit); %逐列求和
fitvalue=Fit./sum_Fit; %每个元素计算适应度百分比
fitvalue=cumsum(fitvalue); %逐列累加适应度百分比，求得累积概率
ms=sort(rand(NP,1)); %生成排序随机队列
fiti=1;
newi=1;
while newi<=NP %选择完成新群体，
if (ms(newi))<fitvalue(fiti) %当随机队列概率小于当前累积概率，则复制该个体为新群体上的新个体，新群体个体序号+1
nf(newi,:)=f(fiti,:);
newi=newi+1;
else %当随机队列概率大于等于当前累积概率，则累积概率数组序号+1
fiti=fiti+1;
end
end
%%%%%%%%%%%%%%%%%%%%%%基于概率的交叉操作%%%%%%%%%%%%%%%%%%
for i=1:2:NP %1，3，5...奇数位遍历
p=rand;
if p<Pc %如果小于交叉率
q=randi([0,1],1,L); %生成1到L的数组，每位上随机取0或1
for j=1:L
if q(j)==1; %若该位上是1，则与i+1上相同位上数进行swap
temp=nf(i+1,j);
nf(i+1,j)=nf(i,j);
nf(i,j)=temp;
end
end
end
end
%%%%%%%%%%%%%%%%%%%基于概率的变异操作%%%%%%%%%%%%%%%%%%%%%%%
i=1;
while i<=round(NP*Pm) %迭代数=变异率x总数=变异染色体数量
h=randi([1,NP],1,1); %随机选取一个需要变异的染色体
for j=1:round(L*Pm) %在一个染色体上随机变异一定数量的基因
g=randi([1,L],1,1); %计算随机需要变异的基因在第几位
nf(h,g)=~nf(h,g); %将其取反=即是发生变异
end
i=i+1;
end
f=nf;
f(1,:)=fBest; %保留最优个体在新种群中
trace(k)=maxFit; %历代最优适应度
end
xBest; %最优个体
君主交叉GA2
%%%%%%%%%%%%%%%%%%%%实值遗传算法求函数极值%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%初始化%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear all; %清除所有变量
close all; %清图
clc; %清屏
D=10; %基因数目
NP=100; %染色体数目
Xs=20; %上限
Xx=-20; %下限
G=100; %最大遗传代数
f=zeros(D,NP); %初始种群赋空间
nf=zeros(D,NP); %子种群赋空间
Pc=0.8; %交叉概率
Pm=0.1; %变异概率
f=rand(D,NP)*(Xs-Xx)+Xx; %随机获得初始种群
%%%%%%%%%%%%%%%%%%%%%%按适应度升序排列%%%%%%%%%%%%%%%%%%%%%%%
for np=1:NP
FIT(np)=func2(f(:,np));
end
[SortFIT,Index]=sort(FIT);
Sortf=f(:,Index);
%%%%%%%%%%%%%%%%%%%%%%%遗传算法循环%%%%%%%%%%%%%%%%%%%%%%%%%%
for gen=1:G
%%%%%%%%%%%%%%采用君主方案进行选择交叉操作%%%%%%%%%%%%%%%%
Emper=Sortf(:,1); %君主染色体
NoPoint=round(D*Pc); %每次交叉基因的个数=round(基因数x交叉率)
PoPoint=randi([1 D],NoPoint,NP/2); %交叉基因的位置，NoPoint x NP/2的随机矩阵
nf=Sortf;
for i=1:NP/2
nf(:,2*i-1)=Emper; %奇数位=君主
nf(:,2*i)=Sortf(:,2*i); %偶数位=原基因
for k=1:NoPoint
nf(PoPoint(k,i),2*i-1)=nf(PoPoint(k,i),2*i);
nf(PoPoint(k,i),2*i)=Emper(PoPoint(k,i));
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%变异操作%%%%%%%%%%%%%%%%%%%%%%%%%
for m=1:NP
for n=1:D
r=rand(1,1); %随机生成1个数在(0,1]之间
if r<Pm %若随机概率小于变异率则发生变异
nf(n,m)=rand(1,1)*(Xs-Xx)+Xx; %即用概率乘以区间差+下限重新生成新的值
end
end
end
%%%%%%%%%%%%%%%%%%%%%子种群按适应度升序排列%%%%%%%%%%%%%%%%%%
for np=1:NP
NFIT(np)=func2(nf(:,np)); %计算子种群适应度
end
[NSortFIT,Index]=sort(NFIT); %子种群排序
NSortf=nf(:,Index); %子种群
%%%%%%%%%%%%%%%%%%%%%%%%%产生新种群%%%%%%%%%%%%%%%%%%%%%%%%%%
f1=[Sortf,NSortf]; %子代和父代合并
FIT1=[SortFIT,NSortFIT]; %子代和父代的适应度值合并
[SortFIT1,Index]=sort(FIT1); %适应度按升序排列
Sortf1=f1(:,Index); %按适应度排列个体
SortFIT=SortFIT1(1:NP); %取前NP个适应度值
Sortf=Sortf1(:,1:NP); %取前NP个个体
trace(gen)=SortFIT(1); %历代最优适应度值
end
Bestf=Sortf(:,1); %最优个体
trace(end) %最优值
GA算法求解旅行商问题
%%%%%%%%%%%%%%%%%%%%%%%%%遗传算法解决TSP问题%%%%%%%%%%%%%%%%%%%%%%%
clear all; %清除所有变量
close all; %清图
clc; %清屏
C=[1304 2312;3639 1315;4177 2244;3712 1399;3488 1535;3326 1556;...
3238 1229;4196 1044;4312 790;4386 570;3007 1970;2562 1756;...
2788 1491;2381 1676;1332 695;3715 1678;3918 2179;4061 2370;...
3780 2212;3676 2578;4029 2838;4263 2931;3429 1908;3507 2376;...
3394 2643;3439 3201;2935 3240;3140 3550;2545 2357;2778 2826;...
2370 2975]; %31个省会城市坐标
N=size(C,1); %TSP问题的规模,即城市数目，也就是基因数
D=zeros(N); %任意两个城市距离间隔矩阵
%%%%%%%%%%%%%%%%%%%%%求任意两个城市距离间隔矩阵%%%%%%%%%%%%%%%%%%%%%
for i=1:N
for j=1:N
D(i,j)=((C(i,1)-C(j,1))^2+(C(i,2)-C(j,2))^2)^0.5;%计算距离
end
end
NP=200; %种群规模
G=1000; %最大遗传代数
f=zeros(NP,N); %用于存储种群
F=[]; %种群更新中间存储
for i=1:NP
f(i,:)=randperm(N); %随机生成初始种群
end
R=f(1,:); %存储最优种群
len=zeros(NP,1); %存储路径长度
fitness=zeros(NP,1); %存储归一化适应值
gen=0;
%%%%%%%%%%%%%%%%%%%%%%%%%遗传算法循环%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
while gen<G
%%%%%%%%%%%%%%%%%%%%%计算路径长度%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for i=1:NP
len(i,1)=D(f(i,N),f(i,1)); %终点到起点的距离
for j=1:(N-1)
len(i,1)=len(i,1)+D(f(i,j),f(i,j+1)); %路径长度累加
end
end
maxlen=max(len); %最长路径
minlen=min(len); %最短路径
%%%%%%%%%%%%%%%%%%%%%%%%%更新最短路径%%%%%%%%%%%%%%%%%%%%%%%%%%
rr=find(len==minlen); %最短路径位置数组
R=f(rr(1,1),:); %最短路径
%%%%%%%%%%%%%%%%%%%%%计算归一化适应值%%%%%%%%%%%%%%%%%%%%%%%%%%
for i=1:length(len)
fitness(i,1)=(1-((len(i,1)-minlen)/(maxlen-minlen+0.001)));
end
%%%%%%%%%%%%%%%%%%%%%%%%%%选择操作%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
nn=0;
for i=1:NP
if fitness(i,1)>=rand %适应度越高越有机会被选中
nn=nn+1; %被选中个体数加1
F(nn,:)=f(i,:); %用F记录被选中的个体
end
end
[aa,bb]=size(F);
while aa<NP
nnper=randperm(nn); %对选择的群体里随机挑2个个体
A=F(nnper(1),:); %随机序列第一位对应的个体
B=F(nnper(2),:); %随机序列第二位对应的个体
%%%%%%%%%%%%%%%%%%%%%%%交叉操作%%%%%%%%%%%%%%%%%%%%%%%%%%%%
W=ceil(N/10); %交叉点个数31/10 = 4个
p=unidrnd(N-W+1); %随机选择交叉范围，从p到p+W
for i=1:W %4个交叉点
%先保存交叉点所在值在对方的位置，然后交换交叉点值，再交换失去的值在对方的位置，即保证31个城市不会因为交叉导致没走完
x=find(A==B(p+i-1)); %B中交叉点所在值在A的位置
y=find(B==A(p+i-1));
temp=A(p+i-1);
A(p+i-1)=B(p+i-1);
B(p+i-1)=temp;
temp=A(x);
A(x)=B(y);
B(y)=temp;
end
%%%%%%%%%%%%%%%%%%%%%%%%%%变异操作%%%%%%%%%%%%%%%%%%%%%%%%%
p1=floor(1+N*rand());
p2=floor(1+N*rand());
while p1==p2 %生成两个不同城市位置序号
p1=floor(1+N*rand());
p2=floor(1+N*rand());
end
tmp=A(p1); %在A与B中进行该位置上城市的交换，即发生变异
A(p1)=A(p2);
A(p2)=tmp;
tmp=B(p1);
B(p1)=B(p2);
B(p2)=tmp;
F=[F;A;B]; %合并取出来并经过交叉，变异的个体到原选择数组上
[aa,bb]=size(F);
end
if aa>NP
F=F(1:NP,:); %保持种群规模为NP
end
f=F; %更新种群
f(1,:)=R; %保留每代最优个体
clear F;
gen=gen+1
Rlength(gen)=minlen;
end
figure
for i=1:N-1
plot([C(R(i),1),C(R(i+1),1)],[C(R(i),2),C(R(i+1),2)],'bo-');
hold on;
end
plot([C(R(N),1),C(R(1),1)],[C(R(N),2),C(R(1),2)],'ro-');
title(['优化最短距离:',num2str(minlen)]);
figure
plot(Rlength)
xlabel('迭代次数')
ylabel('目标函数值')
title('适应度进化曲线')
差分进化算法DE
差分演化算法是一种基于群体差异的启发式随机搜索算法，将问题的求解表示成"染色体"的适者生存过程，通过"染色体"群的一代代不断进化，包括复制、交叉和变异等操作，最终收敛到"最适应环境"的个体，从而求得问题的最优解或满意解。
关键参数
• 种群数量Np：一般来说，种群数量越大，种群的多样性也就越大，寻优能力也就越强，但也会增加计算的难度。为了算法具有足够的不同的变异向量，Np >= 4。一般取5D~10D。（D为变量维数）
• 变异算子F：变异算子 F∈[0, 2] 决定偏差向量的缩放比例。F过小，可能造成算法“早熟”，容易陷入局部最优；F过大，算法很难快速收敛到最优值。F = 0.5通常是一个较好的初始选择。如果种群过早收敛，那么 F 或 Np 应该增大。
• 交叉算子CR：交叉算子 CR∈[0, 1]控制着一个试验向量参数来自随机选择的变异向量而不是原来的向量的概率。CR越大，发生交叉的概率越大。CR的一个较好选择是 0.1，但较大的 CR通常会加速收敛。为了尝试获得一个快速解，可以先尝试 CR = 0.9 或 CR = 0.1。
主要流程
变异操作
DE/rand/1/bin
描述：
- 对需要变异的染色体 m，通过随机值选择 3 条互不相同且与染色体 m 也不相同的其他染色体 r1, r2, r3
- 利用变异算子对 r2, r3 的差异进行缩小，
- r2, r3与 r1 进行合并后替换到原染色体 m 上，完成变异过程。
DE/rand/1/bin + 自适应变异
描述：
在 DE/rand/1/bin 的基础上，针对变异算子进行自适应计算。此自适应计算的方法可以基于变异代数，对于变异代数越高的变异操作，可以考虑减小变异算子，从而做到更好的收敛效果。
自适应算子计算
将变异算子中随机选择的三个个体Xb,Xm,Xw进行计算得变异种群数组
交叉
DE 交叉
描述：
同时采用随机选择和概率遗传两种方法，随机选择通过随机选择一位基因进行变异替换，概率遗传通过概率值决定是否采用变异后的基因进行交叉。这样可以确保子代至少有一个基因来自变异。
二项式交叉
选择
贪婪准则：对原种群和变异后种群每个个体一对一进行判断，取适应度更高的作为下一代种群
代码
DE/rand/1/bin + 自适应变异
%%%%%%%%%%%%%%%%%差分进化算法求函数极值%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%初始化%%%%%%%%%%%%%%%%%%%%%%%%%
clear all; %清除所有变量
close all; %清图
clc; %清屏
NP=50; %个体数目
D=10; %变量的维数
G=200; %最大进化代数
F0=0.4; %初始变异算子
CR=0.1; %交叉算子
Xs=20; %上限
Xx=-20; %下限
yz=10^-6; %阈值
MAXRUN = 10; %独立测试次数
%%%%%%%%%%%%%%%%%%%%%%%%%赋初值%%%%%%%%%%%%%%%%%%%%%%%%
x=zeros(D,NP); %初始种群
v=zeros(D,NP); %变异种群
u=zeros(D,NP); %选择种群
for run = 1:MAXRUN
%%%%%%%%%%%%%%%%一次独立测试%%%%%%%%%%%%%%
x=rand(D,NP)*(Xs-Xx)+Xx; %赋初值
%%%%%%%%%%%%%%%%%%%%计算目标函数%%%%%%%%%%%%%%%%%%%%
for m=1:NP
Ob(m)=func1(x(:,m)); %计算初始种群每位上的目标函数值
end
trace(1)=min(Ob); %获取最优值
%%%%%%%%%%%%%%%%%%%%%%%差分进化循环%%%%%%%%%%%%%%%%%%%%%
for gen=1:G
%%%%%%%%%%%%%%%%%%%%%%变异操作%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%自适应变异算子%%%%%%%%%%%%%%%%%%%
lamda=exp(1-G/(G+1-gen));
F=F0*2^(lamda); %自适应变异算子公式
%%%%%%%%%%%%%%%%%r1,r2,r3和m互不相同%%%%%%%%%%%%%%%%
for m=1:NP
r1=randi([1,NP],1,1);
while (r1==m) %防止r1==m
r1=randi([1,NP],1,1);
end
r2=randi([1,NP],1,1);
while (r2==m)|(r2==r1) %防止r2==r1,r2==m
r2=randi([1,NP],1,1);
end
r3=randi([1,NP],1,1);
while (r3==m)|(r3==r1)|(r3==r2) %防止r3==r1,r3==r2,r3==m
r3=randi([1,NP],1,1);
end
v(:,m)=x(:,r1)+F*(x(:,r2)-x(:,r3)); %变异向量公式3.3，通过对r1，r2，r3列操作后放入变异种群数组的第m列
end
%%%%%%%%%%%%%%%%%%%%%%交叉操作%%%%%%%%%%%%%%%%%%%%%%%
r=randi([1,D],1,1); %整个交叉操作的随机值，确保交叉操作肯定会选择变异种群其中一个作为新种群，随机选择法
for n=1:D
cr=rand(1); %每个变量上的随机值，概率遗传法
if (cr<=CR)|(n==r) %u为新选择种群，v为变异种群，x为初始种群
%如果每个变量上的随机值小于变异算子，或n==整个操作的随机值，则进行交叉，即用当前序号变异种群作为选择种群
u(n,:)=v(n,:);
else %否则，则不进行交叉，即用当前序号初始种群作为选择种群
u(n,:)=x(n,:);
end
end
%%%%%%%%%%%%%%%%%%%边界条件的处理%%%%%%%%%%%%%%%%%%%%%
for n=1:D
for m=1:NP
if (u(n,m)<Xx)|(u(n,m)>Xs) %超过边界的重新生成值
u(n,m)=rand*(Xs-Xx)+Xx;
end
end
end
%%%%%%%%%%%%%%%%%%%%%%选择操作%%%%%%%%%%%%%%%%%%%%%%%
for m=1:NP
Ob1(m)=func1(u(:,m)); %重新计算新选择种群每位上的目标函数值
end
for m=1:NP
if Ob1(m)<Ob(m) %若小于之前的初始种群的该位上的目标函数值(更优)，初始种群上该位被选择种群替换
x(:,m)=u(:,m); %即是适者生存
end
end
for m=1:NP
Ob(m)=func1(x(:,m)); %重新计算现在的初始种群每位上的目标函数值
end
trace(gen+1)=min(Ob); %取最优值放入轨迹
if min(Ob(m))<yz %如果最优值小于阈值，则跳出
break
end
end
[SortOb,Index]=sort(Ob);
x=x(:,Index);
X=x(:,1); %最优变量
Y=min(Ob) %最优值
end
%%%%%%%%%%%%%%%%%%%%%%%%%画图%%%%%%%%%%%%%%%%%%%%%%%%%%
figure
plot(trace);
xlabel('迭代次数')
ylabel('目标函数值')
title('适应度进化曲线')
免疫算法IA
免疫算法是模仿生物免疫机制，结合基因的进化机理，人工构造出的一种新型智能优化算法。采用群体搜索策略，通过迭代计算，最终以较大的概率得到问题的最优解。相比于其他算法，免疫算法利用自身产生多样性和维持机制的特点，保证了群体多样性，克服了‘早熟’问题，可以得到全局最优解。具有自适应性，随机性并行性，全局收敛性，种群多样性等特点。
| 生物免疫系统 | 免疫算法 |
|---|---|
| 抗原 | 优化问题 |
| 抗体（B细胞） | 优化问题的可行解 |
| 亲和度 | 可行解的质量 |
| 细胞活化 | 免疫选择 |
| 细胞分化 | 个体克隆 |
| 亲和度成熟 | 变异 |
| 克隆抑制 | 克隆抑制 |
| 动态维持平衡 | 种群刷新 |
关键参数
- 抗体种群大小Np：保留了免疫细胞的多样性，种群越大，算法全局搜索能力越好，但算法每代的计算量也就相应增大。一般取10~100，且200以内。
-
- 免疫选择比例：免疫选择的抗体数量越多，产生抗体越多，搜索能力越强，但增加每代计算量。一般取0.1Np
- 0.5Np
- 抗体克隆扩增的倍数：决定了克隆扩增的细胞的数量，从而决定算法的搜索能力。数值越大，局部搜索能力越强，全局搜索能力也有一定提高，但是计算量增大。一般取5~10。
- 种群刷新比例：细胞的淘汰和更新是产生抗体多样性的重要机制，从而对免疫算法的全局搜索能力产生重要影响。一般不超过0.5Np。
- 算子
- 亲和度评价算子：相当于遗传算法中的适应度，与具体问题相关
- 抗体浓度评价算子：表征抗体种群的多样性好坏，浓度过高意味非常类似个体大量存在；不利于全局优化
- 激励度算子：对抗体质量的最终评价结果，通常亲和度大，浓度低的抗体会得到较大激励度
- 克隆抑制算子：用于对经过变异后的克隆体进行再选择，抑制亲和度低的抗体(重新生成随机新生种群-种群刷新)，保留亲和度高的抗体作为免疫种群，最后免疫种群与新生种群进行进行合并
主要流程
免疫选择
前 NP / 2 个激励度高的个体指针对 NP/2 个激励度高的个体进行交叉和变异，目的是试图产生激励度更高的个体。随机生成指针对另一半的个体产生采用随机生成的方式，保证了种群中不断有新个体的加入。
抗体间亲和度计算
对于实数编码，亲和度通常可以使用抗体向量间的欧氏距离来计算
抗体浓度计算
抗体浓度公式
其中N为种群规模，S(abi, abj)表示抗体间的相似度，具体表示如下
相似度检查
若小于相似度阈值，则设为1
激励度计算
激励度算子就是抗体的综合评价，计算方式如下，二者取其一，其本质目的是为了共同考虑亲和度和浓度，从而筛选下一代的抗体
变异操作
变异源邻域
描述：
在需要变异的个体染色体中，通过添加一个扰动，从而使其偏离原来的位置落入原个体相邻的另外一个位置中，从而实现变异源领域的搜索。
实数编码算法变异计算
交叉操作
克隆抑制
描述：
对于克隆产生的所有个体进行亲和度的计算，然后保留亲和度最高的个体，从而实现交叉操作。
例子：
- 原个体亲和度为 10
- 克隆后新个体亲和度分别为 [12，6，34，43，2，15，3，8，19，22]
- 将10个新个体和原个体合并后取亲和度最高的个体，即亲和度为 43 的个体
代码
%%%%%%%%%%%%%%%%%免疫算法求函数极值%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%初始化%%%%%%%%%%%%%%%%%%%%%%%%%
clear all; %清除所有变量
close all; %清图
clc; %清屏
D=10; %免疫个体维数
NP=100; %免疫个体数目
Xs=20; %取值上限
Xx=-20; %取值下限
G=500; %最大免疫代数
pm=0.7; %变异概率
alfa=1; %激励度系数
belta=1; %激励度系数
detas=0.2; %相似度阈值
gen=0; %免疫代数
Ncl=10; %克隆个数
deta0=1*Xs; %邻域范围初值
%%%%%%%%%%%%%%%%%%%%%%%初始种群%%%%%%%%%%%%%%%%%%%%%%%%
f=rand(D,NP)*(Xs-Xx)+Xx; %随机生成初始值数组10行100列
for np=1:NP
FIT(np)=func1(f(:,np));
end
%%%%%%%%%%%%%%%%%计算个体浓度和激励度%%%%%%%%%%%%%%%%%%%
for np=1:NP
for j=1:NP
nd(j)=sqrt(sum((f(:,np)-f(:,j)).^2));%基于欧式距离的抗体间亲和度计算
%nd一开始存储抗体间亲和度，后来计算出抗体间相似度，从而取平均值获得抗体浓度，最后通过
if nd(j)<detas %如果小于相似度阈值，则说明相似，设为1
nd(j)=1;
else
nd(j)=0; %如果大于等于相似度阈值，则说明不相似，设为0
end
end
ND(np)=sum(nd)/NP; %取平均值得抗体浓度
end
FIT = alfa*FIT- belta*ND; %4.6 抗体abi激励度=激励度系数ax初始激励度-激励度系数bx抗体浓度
%%%%%%%%%%%%%%%%%%%激励度按升序排列%%%%%%%%%%%%%%%%%%%%%%
[SortFIT,Index]=sort(FIT);
Sortf=f(:,Index);
%%%%%%%%%%%%%%%%%%%%%%%%免疫循环%%%%%%%%%%%%%%%%%%%%%%%%
while gen<G
for i=1:NP/2
%%%%%%%%选激励度前NP/2个体进行免疫操作%%%%%%%%%%%
a=Sortf(:,i); %取第i列
Na=repmat(a,1,Ncl); %生成重复第i列的1xNcl的矩阵
deta=deta0/gen; %邻域范围初值/当前循环代数=扰动因子，代数越大，扰动越低，减小误差
for j=1:Ncl %迭代克隆个体数
for ii=1:D
%%%%%%%%%%%%%%%%%变异%%%%%%%%%%%%%%%%%%%
if rand<pm %如果小于变异率，则发生变异
Na(ii,j)=Na(ii,j)+(rand-0.5)*deta; %变异算子
end
%%%%%%%%%%%%%%边界条件处理%%%%%%%%%%%%%%%
if (Na(ii,j)>Xs) | (Na(ii,j)<Xx)
Na(ii,j)=rand * (Xs-Xx)+Xx; %如果超出边界则重新生成
end
end
end
Na(:,1)=Sortf(:,i); %保留克隆源个体
%%%%%%%%%%克隆抑制，保留亲和度最高的个体%%%%%%%%%%
for j=1:Ncl
NaFIT(j)=func1(Na(:,j)); %亲和度计算，得数组
end
[NaSortFIT,Index]=sort(NaFIT); %亲和度排序
aFIT(i)=NaSortFIT(1); %取函数最小值，即亲和度最高值
NaSortf=Na(:,Index);
af(:,i)=NaSortf(:,1); %亲和度最高的个体
end
%%%%%%%%%%%%%%%%%%%%免疫种群激励度%%%%%%%%%%%%%%%%%%%
for np=1:NP/2
for j=1:NP/2
nda(j)=sqrt(sum((af(:,np)-af(:,j)).^2));
if nda(j)<detas
nda(j)=1;
else
nda(j)=0;
end
end
aND(np)=sum(nda)/NP/2;
end
aFIT = alfa*aFIT- belta*aND;
%%%%%%%%%%%%%%%%%%%%%%%种群刷新%%%%%%%%%%%%%%%%%%%%%%%
bf=rand(D,NP/2)*(Xs-Xx)+Xx; %生成D行NP/2列在取值范围内的新种群
for np=1:NP/2
bFIT(np)=func1(bf(:,np)); %计算新种群函数评价
end
%%%%%%%%%%%%%%%%%%%新生成种群激励度%%%%%%%%%%%%%%%%%%%%
for np=1:NP/2
for j=1:NP/2
ndc(j)=sqrt(sum((bf(:,np)-bf(:,j)).^2));
if ndc(j)<detas
ndc(j)=1;
else
ndc(j)=0;
end
end
bND(np)=sum(ndc)/NP/2;
end
bFIT = alfa*bFIT- belta*bND;
%%%%%%%%%%%%%%免疫种群与新生种群合并%%%%%%%%%%%%%%%%%%%
f1=[af,bf]; %新旧种群合并
FIT1=[aFIT,bFIT]; %新旧种群激励度数组合并
[SortFIT,Index]=sort(FIT1); %新旧种群激励度合并数组排序
Sortf=f1(:,Index); %取激励度最高的一列，即这一代最优解
gen=gen+1; %免疫代数+1
trace(gen)=func1(Sortf(:,1)); %加入最优解轨迹数组
end
%%%%%%%%%%%%%%%%%%%%%%%输出优化结果%%%%%%%%%%%%%%%%%%%%%%%%
Bestf=Sortf(:,1); %最优变量
trace(end); %最优值
figure,plot(trace)
xlabel('迭代次数')
ylabel('目标函数值')
title('亲和度进化曲线')
蚁群算法ACO
蚁群算法是一种仿生学算法，是由自然界中蚂蚁觅食的行为而启发的。
在自然界中，蚂蚁觅食过程中，蚁群总能够按照寻找到一条从蚁巢和食物源的最优路径。蚂蚁在寻找食物的过程中往往是随机选择路径的，但它们能感知当前地面上的信息素浓度， 并倾向于往信息素浓度高的方向行进。信息素由蚂蚁自身释放，是实现蚁群内间接通信的物质。
由于较短路径上蚂蚁的往返时间比较短，单位时间内经过该路径的蚂蚁多，所以信息素的积累速度比较长路径快。因此，当后续蚂蚁在路口时，就能感知先前蚂蚁留下的信息，并 倾向于选择一条较短的路径前行。
这种正反馈机制使得越来越多的蚂蚁在巢穴与食物之间的 最短路径上行进。由于其他路径上的信息素会随着时间蒸发，最终所有的蚂蚁都在最优路径上行进。
- 构建解
- 状态转移概率
- 启发式信息
- 启发因子：城市间距离的倒数
- 信息素
- 启发式信息
- 状态转移概率
- 信息素
- 蒸发
- 增量
- 最优解
- 走过部分解
关键参数
- 信息素启发式因子α：代表信息量对是否选择当前路径的影响程度，反应蚂蚁在运动过程中所积累的信息素在知道蚁群搜索中的相对重要程度。一般取α ∈[1, 4]。
- 期望启发因子β：表示在搜索时路径上的信息素在指导蚂蚁选择路径时的向导性，反映蚁群在搜索最优路径的过程中的先验性和确定性因素的作用强度。一般取β ∈[3, 5]。β越大，蚂蚁在某个局部点上选择局部最短路径可能性越大，算法收敛速度越快，但随机性越弱，容易陷入局部最优解
- 信息蒸发系数p：p ∈[0,1]，表示信息素的蒸发程度，反映了蚂蚁群体中个体之间相互影响的强弱。p过小时，影响算法的随机性能和全局搜索能力；p过大是，降低算法的收敛速度。
- 信息素强度Q：代表计算信息素增量的一个取值
- 蚂蚁数目m：蚂蚁数量增多，提高算法的全局搜索能力和稳定性，但正反馈作用不明显，收敛速度减慢；蚂蚁数量减少，收敛速度加快，但算法全局性能降低，稳定性差，容易出现过早停滞现象。一般取10~50。
主要流程
选择
采用的是轮盘赌的方式，其概率的计算是基于信息素和距离权重的，公式如下（ 为信息素， 为距离权重）
信息素更新
信息素增量计算公式，其中 Q 信息素强度系数， 为第 k 只蚂蚁在本轮周游中所走过的路径的长度
信息素蒸发量计算公式，其中  表示路径上信息素的蒸发系数，表示信息素的持久性系数
代码
蚁群算法解决旅行商问题
%%%%%%%%%%%%%%%%%%%%蚁群算法解决TSP问题%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%初始化%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear all; %清除所有变量
close all; %清图
clc; %清屏
m=50; %蚂蚁个数
Alpha=1; %信息素重要程度参数
Beta=5; %启发式因子重要程度参数
Rho=0.1; %信息素蒸发系数
G_max=200; %最大迭代次数
Q=100; %信息素增加强度系数
C=[1304 2312;3639 1315;4177 2244;3712 1399;3488 1535;3326 1556;...
3238 1229;4196 1044;4312 790;4386 570;3007 1970;2562 1756;...
2788 1491;2381 1676;1332 695;3715 1678;3918 2179;4061 2370;...
3780 2212;3676 2578;4029 2838;4263 2931;3429 1908;3507 2376;...
3394 2643;3439 3201;2935 3240;3140 3550;2545 2357;2778 2826;...
2370 2975]; %31个省会城市坐标
%%%%%%%%%%%%%%%%%%%%%%%%第一步：变量初始化%%%%%%%%%%%%%%%%%%%%%%%%
n=size(C,1); %n表示问题的规模（城市个数）
D=zeros(n,n); %D表示两个城市距离间隔矩阵
for i=1:n
for j=1:n
if i~=j %如果不是同一个城市
D(i,j)=((C(i,1)-C(j,1))^2+(C(i,2)-C(j,2))^2)^0.5; %计算两个城市距离
else
D(i,j)=eps; %设为一个数可以分辨的最小精度
end
D(j,i)=D(i,j); %更新两个城市距离
end
end
Eta=1./D; %Eta为启发因子，这里设为距离的倒数
Tau=ones(n,n); %Tau为信息素矩阵
Tabu=zeros(m,n); %存储并记录路径的生成
NC=1; %迭代计数器
R_best=zeros(G_max,n); %各代最佳路线
L_best=inf.*ones(G_max,1); %各代最佳路线的长度
figure(1);%优化解
while NC<=G_max
%%%%%%%%%%%%%%%%%%第二步：将m只蚂蚁放到n个城市上%%%%%%%%%%%%%%%%
Randpos=[];
for i=1:(ceil(m/n)) %求得每个城市放m/n只蚂蚁，总共需要几次迭代
Randpos=[Randpos,randperm(n)]; %每次迭代取n只蚂蚁随机去往n个城市
end
Tabu(:,1)=(Randpos(1,1:m))'; %最后取存储m只蚂蚁去往的城市序号的数组
%%%%%第三步：m只蚂蚁按概率函数选择下一座城市，完成各自的周游%%%%%%
for j=2:n
for i=1:m
visited=Tabu(i,1:(j-1)); %已访问的城市数组
J=zeros(1,(n-j+1)); %待访问的城市数组，未初始化，全为0
P=J; %待访问城市的选择概率分布
Jc=1; %待访问的城市数组的遍历变量
for k=1:n
if length(find(visited==k))==0 %find查找已访问城市为k的，生成序号数组，并计算有几个(length)
J(Jc)=k; %若为0，则待访问，将城市序号用来初始化待访问的城市数组
Jc=Jc+1; %遍历变量++
end
end
%%%%%%%%%%%%%%%%%%计算待选城市的概率分布%%%%%%%%%%%%%%%%
for k=1:length(J)
P(k)=(Tau(visited(end),J(k))^Alpha)...
*(Eta(visited(end),J(k))^Beta);
end
P=P/(sum(P)); %状态转移概率公式计算5.1
%%%%%%%%%%%%%%%%按概率原则选取下一个城市%%%%%%%%%%%%%%%%
Pcum=cumsum(P); %轮盘赌概率累加
Select=find(Pcum>=rand); %随机生成一个概率，查找概率累加数组里大于等于该随机概率的数，生成数组Select
to_visit=J(Select(1)); %取Select数组第一个数作为下一个城市
Tabu(i,j)=to_visit;
end
end
if NC>=2
Tabu(1,:)=R_best(NC-1,:);
end
%%%%%%%%%%%%%%%%%%%第四步：记录本次迭代最佳路线%%%%%%%%%%%%%%%%%%
L=zeros(m,1);
for i=1:m
R=Tabu(i,:); %第i行路线放入R
for j=1:(n-1)
L(i)=L(i)+D(R(j),R(j+1)); %遍历n个城市计算第i行路线长度
end
L(i)=L(i)+D(R(1),R(n)); %最后补上终点到起点的距离
end
L_best(NC)=min(L); %这一代最佳路径的长度记录
pos=find(L==L_best(NC)); %查找哪一行(路线)长度==最佳路径长度
R_best(NC,:)=Tabu(pos(1),:); %这一代最佳路径的记录
%%%%%%%%%%%%%%%%%%%%%%%%%第五步：更新信息素%%%%%%%%%%%%%%%%%%%%%%
Delta_Tau=zeros(n,n);
for i=1:m
for j=1:(n-1)
Delta_Tau(Tabu(i,j),Tabu(i,j+1))=... %ant-cycle模型
Delta_Tau(Tabu(i,j),Tabu(i,j+1))+Q/L(i);%Q/L(i)=信息素增加强度系数/本次周游所走路径长度=迭代新增信息素
end
Delta_Tau(Tabu(i,n),Tabu(i,1))=...
Delta_Tau(Tabu(i,n),Tabu(i,1))+Q/L(i);
end
Tau=(1-Rho).*Tau+Delta_Tau; %5.2信息素公式=保留下来的信息素+迭代新增信息素，Rho代表蒸发系数，1-Rho代表持久系数
%%%%%%%%%%%%%%%%%%%%%%%第六步：禁忌表清零%%%%%%%%%%%%%%%%%%%%%%
Tabu=zeros(m,n);
%%%%%%%%%%%%%%%%%%%%%%%%%历代最优路线%%%%%%%%%%%%%%%%%%%%%%%%%%
for i=1:n-1
plot([ C(R_best(NC,i),1), C(R_best(NC,i+1),1)],...
[C(R_best(NC,i),2), C(R_best(NC,i+1),2)],'bo-'); %画路线
hold on;
end
plot([C(R_best(NC,n),1), C(R_best(NC,1),1)],...
[C(R_best(NC,n),2), C(R_best(NC,1),2)],'ro-');
title(['优化最短距离:',num2str(L_best(NC))]);
hold off;
pause(0.005);
NC=NC+1;
end
%%%%%%%%%%%%%%%%%%%%%%%%%%第七步：输出结果%%%%%%%%%%%%%%%%%%%%%%%%%%
Pos=find(L_best==min(L_best));
Shortest_Route=R_best(Pos(1),:); %最佳路线
Shortest_Length=L_best(Pos(1)); %最佳路线长度
figure(2),
plot(L_best)
xlabel('迭代次数')
ylabel('目标函数值')
title('适应度进化曲线')
粒子群算法PSO
粒子群算法模拟鸟群的捕食行为，每个优化问题的解都是搜索空间中的一只鸟。我们称之为“粒子”。所有的粒子都有一个由被优化的函数决定的适应值(fitnessvalue)，每个粒子还有一个速度决定他们飞翔的方向和距离。然后粒子们就追随当前的最优粒子在解空间中搜索。
关键参数
- 粒子种群规模N：粒子数目越大，算法搜索的空间范围就越大，也就更容易发现全局最优解，但算法运行时间也就越长。一般10个粒子已经可以取得比较好的结构。对于比较难的问题或特定的问题，粒子的数量可以取到100或200。
- 惯性权重w：代表对粒子当前速度继承的多少，用来控制算法的开发和探索能力。权重较大时，全局寻优能力较强，局部寻优能力较弱；权重较小时，全局寻优能力较弱，局部寻优能力较强。一般取0.8~1.2。
- 加速常数c1和c2：分别调节向 pbest 和 gbest方向飞行的最大步长，他们分别决定粒子个体经验和群体经验对粒子运行轨迹的影响。通常可以取c1 = c2 = 1.5。
主要流程
PSO初始化为一群随机粒子(随机解)，然后通过迭代找到最优解，在每一次迭代中，粒子通过跟踪两个“极值”来更新自己。
第一个就是粒子本身所找到的最优解，这个解叫做个体极值pBest，
另一个极值是整个种群找到的最优解，这个极值是全局极值gBest。另外也可以不用整个种群而只是用其中一部分最优粒子的邻居，那么在所有邻居中的极值就是局部极值。
- 粒子的速度及位置更新的方式如下：
- 其中
- 粒子的速度更新公式由三部分组成：粒子先前速度+个体最优+全局最优
是一个非负数，称为惯性因子，对算法的收敛起到很大的作用，其值越大，粒子飞跃的范围就越广，更容易找到全局最优，但是也会错失局部搜寻的能力。
分别为局部和全局最优位置。
- 加速常数
也是非负常数，是调整局部最优值和全局最优值权重的参数，如果前者为0说明搜寻过程中没有自身经验只有社会经验，容易陷入局部最优解；若后者为0，即只有社会经验，没有自身经验，常常会陷入局部最优解中，不能飞越该局部最优区域。
是[0,1]范围之内的随机数，
是约束因子，目的是控制速度的权重。
代码
%%%%%%%%%%%%%%%%%粒子群算法求函数极值%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%初始化%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear all; %清除所有变量
close all; %清图
clc; %清屏
global fes; %fes统计函数计算次数
N = 100; %群体粒子个数
D = 10; %粒子维数
c1 = 1.5; %加速常数1
c2 = 1.5; %加速常数2
w = 0.8; %惯性权重
a = 1.0; %约束因子
Xmax = 20; %位置最大值
Xmin = -20; %位置最小值
Vmax = 10; %速度最大值
Vmin = -10; %速度最小值
MAXRUN = 20; %独立测试次数
MAXFES = 10^5; %最大函数迭代次数
error = 10^-6; %误差阈值
success = 0;
for run = 1:MAXRUN
fes = 0;
%%%%%%%%%%%%%%%%初始化种群个体（限定位置和速度）%%%%%%%%%%%%%%%%
x = rand(N, D) * (Xmax - Xmin) + Xmin;
v = rand(N, D) * (Vmax - Vmin) + Vmin;
%%%%%%%%%%%%%%%%%%初始化个体最优位置和最优值%%%%%%%%%%%%%%%%%%%
p = x;
pbest = ones(N,1);
for i = 1:N
pbest(i) = func(x(i, :));
end
%%%%%%%%%%%%%%%%%%%初始化全局最优位置和最优值%%%%%%%%%%%%%%%%%%
g = ones(1, D);
gbest = inf; %初始设为正无穷
for i = 1:N
if(pbest(i) < gbest) %如果个体最优小于全局最优
g = p(i, :); %更新全局最优数组
gbest = pbest(i); %更新全局最优值
end
end
gen = 1;
%%%%%%%%%%%按照公式依次迭代直到满足精度或者迭代次数%%%%%%%%%%%%%
while fes <= MAXFES
for j = 1:N
%%%%%%%%%%%%%%更新个体最优位置和最优值%%%%%%%%%%%%%%%%%
temp = func(x(j, :));
if (temp < pbest(j))
p(j, :) = x(j, :);
pbest(j) = temp;
end
%%%%%%%%%%%%%%%%更新全局最优位置和最优值%%%%%%%%%%%%%%%
if(pbest(j) < gbest)
g = p(j, :);
gbest = pbest(j);
end
%%%%%%%%%%%%%%%%%更新位置和速度值%%%%%%%%%%%%%%%%%%%%%
v(j, :) = w * v(j, :) + c1 * rand * (p(j, :) - x(j, :))... %粒子的速度更新公式
+ c2 * rand * (g - x(j, :));
x(j, :) = x(j, :) + a*v(j, :); %粒子的位置更新公式
%%%%%%%%%%%%%%%%%%%%边界条件处理%%%%%%%%%%%%%%%%%%%%%%
for ii = 1:D
if (v(j, ii) > Vmax) | (v(j, ii) < Vmin) %对粒子的速度，位置进行边界处理
v(j, ii) = rand * (Vmax - Vmin) + Vmin;
end
if (x(j, ii) > Xmax) | (x(j,ii) < Xmin)
x(j, ii) = rand * (Xmax - Xmin) + Xmin;
end
end
end
%%%%%%%%%%%%%%%%%%%%记录历代全局最优值%%%%%%%%%%%%%%%%%%%%%
trace(gen) = gbest;
gen = gen + 1;
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%画图%%%%%%%%%%%%%%%%%%%%%%%%%%
figure
plot(trace)
xlabel('迭代次数');
ylabel('适应度值');
title('适应度进化曲线')
模拟退火算法
模拟退火算法来源于固体退火原理，将固体加温至充分高，再让其徐徐冷却，加温时，固体内部粒子随温升变为无序状，内能增大，而徐徐冷却时粒子渐趋有序，在每个温度都达到平衡态，最后在常温时达到基态，内能减为最小，即所求值
能量即目标函数，能量最低态即最优解，
温度是Metropolis算法的一个重要控制参数，开始时T大，可以接受较差恶化解；随着T减小，只能接受较好的恶化解了
根据Metropolis准则，粒子在温度T时趋于平衡的概率为 exp(-ΔE/(kT)) ，其中 E 为温度 T 时的内能，ΔE为其改变数, k 为 Boltzmann 常数。 Metropolis 准则常表示为
主要流程
- 系统从一个能量状态变化到另一个能量状态时(根据当前温度产生新解：邻域随机选取+边界处理)，相应的能量从Eold变化到Enew
- 若新状态是全局最优，则更新全局最优解，保留上一个最优解
- Metropolis过程
- 若状态向下(局部最优)，则接受新解
- 若状态向上(全局搜索)，则以一定概率接受(
exp(-ΔE/(kT))>rand)
- 在系统每次变化时，T也在冷却，
T(n+1)=KxT(n)。当T趋向于0(设置一个小值，如0.001)时，求得问题整体最优解
禁忌搜索算法
基本思想：
采用邻域选优的搜索方法，为了逃离局部最优解，算法必须能够接受劣解，也就是每一次得到的解不一定优于原来的解。
但是，一旦接受了劣解，算法迭代即可能陷入循环。为了避免循环，算法将最近接受的一些移动放在禁忌表中，在以后的迭代中加以禁止。即只有不再禁忌表中的较好解（可能比当前解差）才能接受作为下一代迭代的初始解。
随着迭代的进行，禁忌表不断更新，经过一定的迭代次数后，最早进入禁忌表的移动就从禁忌表中解禁退出。
主要流程
算法记忆汇总
- 遗传算法
- 最早是由美国的 John Holland于20世纪70年代提出
- 能有效求解NP问题(所有的非确定性多项式时间可解的判定问题构成NP类问题)，以及非线性，多峰函数优化和多目标优化问题
- 差分进化算法
- 由Storn和Price于1995年首次提出
- 免疫算法
- 1958年澳大利亚学者Burnet率先提出克隆选择原理
- 1973年Jerne提出免疫系统的数学框架
- 蚁群算法
- 由Marco Dorigo于1992年在他的博士论文中提出
- 蚁群算法是一种用来寻找优化路径的概率型算法，灵感来源于蚂蚁在寻找食物过程中发现路径的行为
- 粒子群算法
- 由Eberhart博士和kennedy博士发明
- 模拟退火算法
- 最早的思想是由Metropolis等人于1953年提出。1983 年,Kirkpatrick 等成功地将退火思想引入到组合优化领域。
- 模拟退火算法从某一较高初温出发，伴随温度参数的不断下降,结合概率突跳特性在解空间中随机寻找目标函数目标函数)的全局最优解，即在局部最优解能概率性地跳出并最终趋于全局最优。
- 禁忌搜索算法
- 由Glover教授于1986年在一篇论文中首次提出
- 从一个初始可行解出发，选择一系列的特定搜索方向（移动）作为试探，选择实现让特定的目标函数值变化最多的移动。为了避免陷入局部最优解，TS搜索中采用了一种灵活的“记忆”技术，对已经进行的优化过程进行记录和选择，指导下一步的搜索方向，这就是Tabu表的建立。
- 神经网络算法
- 最早由McCulloch与Pitts提出形式神经元数学模型
- 由Hopfield提出具有联想记忆功能的Hopfield神经网络，标志神经网络取得重大进展
- 用大量的简单计算单元（神经元）构成非线性系统，在一定程度上模拟了大脑的信息处理、存储和检索等功能。
- 常见网络结构
- 前馈神经网络
- 反馈神经网络
- BP网络的误差反向后传学习算法，是最常用的神经网络算法。它利用输出后的误差来估计输出层的直接前导层误差，再利用这个误差更新前一层的误差，如此反传下去获得所有其他各层的误差估计。
- 自组织网络
测试输出汇总
-
tall：运行计时
-
mean：平均
-
gen：代数
-
fes：函数评价次数
-
trace：最优值
注意输出信息最终应该包括：
-
测试的函数名Fname算法参数取值，如D, Pc, Pm, NP,MAXG, MAXRUN等
-
表3最优函数值bestfun, 最差函数值worstfun，平均值meanfun,标准差stdfun，算法总体平均用时meantotaltimes(s)
-
表4首次获得至今最优解值对应的平均代数meanbestGen，对应的函数评价次数meanbestFEs，对应的用时meanbesttime(s)，统计得到的成功次数success%
//单次测试F1run
//每次运行记录F1totoal
//汇总测试total求最优值，平均值，平均用时
global fes %fes统计函数计算次数
MAXRUN = 10; %独立测试次数
funfile = fopen(‘F1total.txt’,‘a’); %w 改为a表示添加（w表示覆盖原有内容）
tallrecord = []; %用于记录每次独立测试的运行用时
bestfrecord = []; %用于记录每次独立测试得到的最优函数值
allfile = fopen(‘total.txt’,‘a’); %w 改为a表示添加（w表示覆盖原有内容）
for run = 1:MAXRUN
%....忽略一段代码…..
fprintf('第%d次运行\t第%d代\t%f\n',run,0,trace(1));
onerunfile = fopen([‘F1_run’ num2str(run) ‘.txt’],‘w’); %建立第run次运行的记录文档
fprintf(onerunfile,‘%d\t%d\t%g\r\n’,0,fes,trace(1)); %写入第0代种群结
t0 = clock; %记录当前时刻
tall = 0; %统计计算总用时（排除写文件用时）
for gen=1:G %差分进化循环
%....忽略好多段代码…..
tall = tall + etime(clock,t0); %累加计算用时（单位秒，排除写文件用时）
fprintf(onerunfile,'%d\t%d\t%g\r\n',gen,fes,trace(gen+1)); %写入第gen代结果
t0 = clock; %记录当前时刻
end
fclose(onerunfile);
tallrecord = [tallrecord,tall]; %添加第run次独立测试的运行用时
bestfrecord = [bestfrecord, trace(gen+1)]; %添加第run次独立测试得到的最优函数值
fprintf(funfile,‘%d\t%g\t%g\r\n’,run,trace(gen+1),tall); %输出第几次运行最优值和用时
End
fclose(funfile);
fprintf(allfile,'%s\t%g\t%g\t%g\r\n','F1',min(bestfrecord),mean(bestfrecord),mean(tallrecord));
%输出汇总信息：【函数名，多次运行最优值，平均值，平均用时】
