抱歉，您的浏览器无法访问本站
本页面需要浏览器支持（启用）JavaScript
了解详情 >

NingLi's Blog

Stay foolish, Stay hungry

主页

分类

差分进化(DE)算法

research

发布于：2024年3月22日

差分进化算法(differential evolution, DE)是一类基于群体的自适应全局优化算法，有较强的鲁棒性和全局寻优能力，主要用于求解实数优化问题。它从数学角度看是一种随机搜索算法，从工程角度看是一种自适应的迭代寻优过程。差分进化算法凭借其独特的优势在数据挖掘、模式识别、电磁学等领域被广泛应用。

算法流程

差分进化算法采用实数编码，基于差分的简单变异操作和“一对一”的竞争生存策略，具体步骤如下

  1. 随机初始化种群规模N，进化代数t=1;
  2. 计算初始种群中的每个个体的目标函数值；
  3. 判断是否达到终止条件或达到最大进化代数。若满足，则进化停止，将此时的最佳个体作为解输出；否则执行步骤4；
  4. 进行变异操作和交叉操作，并对边界条件进行处理，得到临时种群；
  5. 对临时种群进行评价，计算临时种群中每个个体的适应度值；
  6. 进行选择操作，得到新的种群；
  7. 进化代数t=t+1，并用新的种群代替旧种群，返回步骤3.

DE算法流程图

伪代码

1
2
3
4
5
6
7
8
9
10
11
12
13
14
初始化N个个体
初始化参数D、N_iter
i=1
while(t<=N_iter)
	for i=1:N
		for j=1:D
			进行变异操作和交叉操作，得到临时种群
		end for j
		计算临时种群中每个个体的适应度值
		对临时种群进行选择操作，得到新种群
	end for i
	t = t+1
end while
输出结果

MATLAB仿真实例

1.计算函数\(f(x)= \sum_{i=1}^n x_i^2(-20\leq x_i\leq20)\)的最小值，其中个体\(x\)的维数\(n=10\)。这是一个简单的平方和函数，只有一个极小点\(x=(0,0,\cdots,0)\)，理论最小值\(f(0,0,\cdots,0)=0\)

仿真过程如下：

（1）初始化个体数目为\(NP=50\)，变量维数为\(D=10\)，最大进化代数为\(G=200\)，初始变异算子\(F_0=0.4\)，交叉算子\(CR=0.1\)，阈值\(yz=10^{-6}\)。

（2）产生初始种群，计算个体目标函数；进行变异操作、交叉操作、边界条件处理，产生临时种群，其中变异操作采用自适应变异算子，边界条件处理采用在可行域中随机产生参数向量的方式。

（3）计算临时种群个体目标函数，与原种群对应个体进行“一对一”选择操作，产生新种群。

（4）判断是否满足终止条件：若满足，则结束搜索过程，输出优化值；若不满足，则继续进行迭代优化。

优化结束后，DE目标函数曲线如下图所示

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
%% 差分进化算法求极值
clc; clear; close all;
NP = 50;  % 种群数量
D = 10;  % 变量的维数
G = 200;  % 最大的进化代数
F0  = 0.4; % 初始变异算子
CR = 0.1;  % 交叉算子
Xs = 20;   % 上限
Xx = -20;  % 下限
yz = 10^-6;  % 阈值
% 赋初值
x = zeros(D, NP);  % 初始种群
v = zeros(D, NP);  % 变异种群
u = zeros(D, NP);  % 选择种群
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
    %%%%%%%%%%%%%%%%%%%%%%%%% 变异操作%%%%%%%%%%%%%%%%
    %%%%%%%%%%%%%%%%%%%%% 自适应变异算子%%%%%%%%%%%%%%%
    lambda = exp(1-G/(G+1-gen));
    F = F0*2^(1-lambda);
    %%%%%%%%%%%%%%%% r1, r2, r3和m互不相同%%%%%%%%%%%%%%%%
    for m = 1:NP
        r1 = randi([1, NP], 1, 1);
        while(r1==m)
            r1 = randi([1, NP], 1, 1);
        end
        r2 = randi([1, NP], 1, 1);
        while(r2==m) || (r2==r1)
            r2 = randi([1, NP], 1, 1);
        end
        r3 = randi([1, NP], 1, 1);
        while (r3==m) || (r3==r1) || (r3==r2)
            r3 = randi([1, NP], 1, 1);
        end
        v(:,m) = x(:,r1) + F*(x(:,r2)-x(:,r3));
    end
    %%%%%%%%%%%%%%%%%%%%%%%交叉操作%%%%%%%%%%%%%%%%%%%%%%%%%%%
    r = randi([1, D], 1, 1);
    for n = 1:D
        cr = rand(1);
        if(cr<=CR) || (n==r)
            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);  % 最优值
%%%%%%%%%%%%%%%%%%画图%%%%%%%%%%%%%%%%%%
figure
plot(trace);
xlabel('迭代次数');
ylabel('目标函数值');
title('DE目标函数曲线')
%%%%%%%%%%%%适应度函数%%%%%%%%%%%%%%%%%%%
function result = func1(x)
summ = sum(x.^2);
result  = summ;
end
  1. 求函数\(f(x,y)=3\cos(xy)+x+y\)的最小值，其中x的取值范围为［-4，4］，y的取值范围为［-4，4］。这是一个有多个局部极值的函数，其函数值图形如下图所示，其MATLAB实现程序如下：
1
2
3
4
5
6
7
8
9
10
11
12
13
clc; clear; close all;
%%%%%%%%%%%%%%%f(x,y) =3cos(xy)+x+y%%%%%%%%%%%%%%%%
x = -4:0.02:4;
y = -4:0.02:4;
N = size(x, 2);
for i = 1:N
    for j = 1:N
       z(i,j) = 3*cos(x(i)*y(j))+x(i)+y(j);
    end
end
mesh(x, y, z);
xlabel('x');
ylabel('y');

仿真过程如下：

（1）初始化个体数目为NP=20，变量维数为D=2，最大进化代数为G=100，变异算子F=0.5，交叉算子CR=0.1；

（2）产生初始种群，计算个体目标函数；进行变异操作、交叉操作、边界条件处理，产生临时种群，其中边界条件处理采用边界吸收方式；

（3）计算临时种群个体目标函数，与原种群对应个体进行“一对一”选择操作，产生新种群；

（4）判断是否满足终止条件：若满足，则结束搜索过程，输出优化值；若不满足，则继续进行迭代优化。、

优化结束后，DE目标函数曲线如下图所示

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
clc; clear; close all;
%%%%%%%%%%%%%%%f(x,y) =3cos(xy)+x+y%%%%%%%%%%%%%%%%
x = -4:0.02:4;
y = -4:0.02:4;
N = size(x, 2);
for i = 1:N
    for j = 1:N
       z(i,j) = 3*cos(x(i)*y(j))+x(i)+y(j);
    end
end
mesh(x, y, z);
xlabel('x');
ylabel('y');

%% 差分进化算法求极值
NP = 20;  % 种群数量
D = 2;  % 变量的维数
G = 100;  % 最大的进化代数
F  = 0.5; % 初始变异算子
CR = 0.1;  % 交叉算子
Xs = 4;   % 上限
Xx = -4;  % 下限
yz = 10^-6;  % 阈值
% 赋初值
x = zeros(D, NP);  % 初始种群
v = zeros(D, NP);  % 变异种群
u = zeros(D, NP);  % 选择种群
x = rand(D, NP) * (Xs - Xx) + Xx;  % 赋初值

%%%%%%%%%%%%%%计算目标函数%%%%%%%%%%%%%%%%
for m=1:NP
    Ob(m) = func2(x(:, m));
end
trace(1) = min(Ob);
%%%%%%%%%%%%%%差分进化循环%%%%%%%%%%%%%%%%
for gen = 1:G
    %%%%%%%%%%%%%%%%%%%%%%%%% 变异操作%%%%%%%%%%%%%%%%
    %%%%%%%%%%%%%%%% r1, r2, r3和m互不相同%%%%%%%%%%%%%%%%
    for m = 1:NP
        r1 = randi([1, NP], 1, 1);
        while(r1==m)
            r1 = randi([1, NP], 1, 1);
        end
        r2 = randi([1, NP], 1, 1);
        while(r2==m) || (r2==r1)
            r2 = randi([1, NP], 1, 1);
        end
        r3 = randi([1, NP], 1, 1);
        while (r3==m) || (r3==r1) || (r3==r2)
            r3 = randi([1, NP], 1, 1);
        end
        v(:,m) = x(:,r1) + F*(x(:,r2)-x(:,r3));
    end
    %%%%%%%%%%%%%%%%%%%%%%%交叉操作%%%%%%%%%%%%%%%%%%%%%%%%%%%
    r = randi([1, D], 1, 1);
    for n = 1:D
        cr = rand(1);
        if(cr<=CR) || (n==r)
            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) = Xx;
            end
            if(u(n,m)>Xs)
                u(n,m) = Xs;
            end
        end
    end
    %%%%%%%%%%%%%%%%%%%%%选择操作%%%%%%%%%%%%%%%%%%%%%%%%%%
    for m = 1:NP
        Ob1(m) = func2(u(:, m));
    end
    for m = 1:NP
        if Ob1(m) < Ob(m)
            x(:, m) = u(:, m);
        end
    end
    for m = 1:NP
        Ob(m) = func2(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);  % 最优值
%%%%%%%%%%%%%%%%%%画图%%%%%%%%%%%%%%%%%%
figure
plot(trace);
xlabel('迭代次数');
ylabel('目标函数值');
title('DE目标函数曲线')
%%%%%%%%%%%%适应度函数%%%%%%%%%%%%%%%%%%%
function result = func2(x)
result  = 3*cos(x(1)*x(2))+x(1)+x(2);
end
  1. 用离散差分进化算法求函数\(f(x,y)=-((x^2+y-1)^2+(x+y^2-7)^2)/200+10\)的最大值，其中\(x\)的取值为-100～100之间的整数，\(y\)的取值为-100～100之间的整数，其函数值图形如下图所示

1
2
3
4
5
6
7
8
9
10
11
12
13
%%%%F(x,y)=-((x~2+y-1).^2+(x+y~2-7)^2)/200+10%%%%
clc; clear; close all;
x=-100:1:100;
y=-100:1:100;
N=size(x,2);
for i=1:N
    for j=1:N
        z(i,j)=-((x(i)^2+y(j)-1).^2+(x(i) +y(j)^2-7)^2)/200+10;
    end
end
mesh(x,y,z);
xlabel('x');
ylabel('y');

仿真过程如下：

（1）初始化个体数目为NP=20，变量维数为D=2，最大进化代数为G=100，变异算子F=0.5，交叉算子CR=0.1

（2）产生数值为［-100，100］内整数的初始种群，计算个体目标函数；进行变异操作，对变异后的种群内数值进行取整操作，然后进行交叉操作、边界条件处理操作，产生临时种群，其中边界条件处理采用边界吸收方式

（3）计算临时种群个体目标函数，与原种群对应个体进行“一对一”选择操作，产生新种群

（4）判断是否满足终止条件：若满足，则结束搜索过程，输出优化值；若不满足，则继续进行迭代优化

优化结束后，DE目标函数曲线如下图所示，优化后的结果为：x=-2，y=-3，函数f（x，y）的最大值为10。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
%%%%F(x,y)=-((x~2+y-1).^2+(x+y~2-7)^2)/200+10%%%%
clc; clear; close all;
x=-100:1:100;
y=-100:1:100;
N=size(x,2);
for i=1:N
    for j=1:N
        z(i,j)=-((x(i)^2+y(j)-1).^2+(x(i) +y(j)^2-7)^2)/200+10;
    end
end
mesh(x,y,z);
xlabel('x');
ylabel('y');

%% 差分进化算法求极值
NP = 20;  % 种群数量
D = 2;  % 变量的维数
G = 100;  % 最大的进化代数
F  = 0.5; % 初始变异算子
CR = 0.1;  % 交叉算子
Xs = 100;   % 上限
Xx = -100;  % 下限
yz = 10^-6;  % 阈值
% 赋初值
x = zeros(D, NP);  % 初始种群
v = zeros(D, NP);  % 变异种群
u = zeros(D, NP);  % 选择种群
x = randi([Xx, Xs], D, NP);  % 赋初值

%%%%%%%%%%%%%%计算目标函数%%%%%%%%%%%%%%%%
for m=1:NP
    Ob(m) = func3(x(:, m));
end
trace(1) = max(Ob);
%%%%%%%%%%%%%%差分进化循环%%%%%%%%%%%%%%%%
for gen = 1:G
    %%%%%%%%%%%%%%%%%%%%%%%%% 变异操作%%%%%%%%%%%%%%%%
    %%%%%%%%%%%%%%%% r1, r2, r3和m互不相同%%%%%%%%%%%%%%%%
    for m = 1:NP
        r1 = randi([1, NP], 1, 1);
        while(r1==m)
            r1 = randi([1, NP], 1, 1);
        end
        r2 = randi([1, NP], 1, 1);
        while(r2==m) || (r2==r1)
            r2 = randi([1, NP], 1, 1);
        end
        r3 = randi([1, NP], 1, 1);
        while (r3==m) || (r3==r1) || (r3==r2)
            r3 = randi([1, NP], 1, 1);
        end
        v(:,m) = floor(x(:,r1) + F*(x(:,r2)-x(:,r3)));
    end
    %%%%%%%%%%%%%%%%%%%%%%%交叉操作%%%%%%%%%%%%%%%%%%%%%%%%%%%
    r = randi([1, D], 1, 1);
    for n = 1:D
        cr = rand(1);
        if(cr<=CR) || (n==r)
            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) = Xx;
            end
            if(u(n,m)>Xs)
                u(n,m) = Xs;
            end
        end
    end
    %%%%%%%%%%%%%%%%%%%%%选择操作%%%%%%%%%%%%%%%%%%%%%%%%%%
    for m = 1:NP
        Ob1(m) = func3(u(:, m));
    end
    for m = 1:NP
        if Ob1(m) > Ob(m)
            x(:, m) = u(:, m);
        end
    end
    for m = 1:NP
        Ob(m) = func3(x(:, m));
    end
    trace(gen+1) = max(Ob);
%     if min(Ob(m)) <yz
%         break
%     end
end
[SortOb, Index] = sort(Ob);
x = x(:, Index);
X = x(:, end);  % 最优变量
Y = max(Ob);  % 最优值
%%%%%%%%%%%%%%%%%%画图%%%%%%%%%%%%%%%%%%
figure
plot(trace);
xlabel('迭代次数');
ylabel('目标函数值');
title('DE目标函数曲线')
%%%%%%%%%%%%适应度函数%%%%%%%%%%%%%%%%%%%
function result = func3(x)
result  = -((x(1)^2+x(2)-1).^2+(x(1) +x(2)^2-7)^2)/200+10;
end

参考

[1] 包子阳等. 智能优化算法及其MATLAB实例(第3版)[M]. 北京: 电子工业出版社, 2020.

更新于：2024年3月27日

评论