|
|
EDA365欢迎您登录!
您需要 登录 才可以下载或查看,没有帐号?注册
x
+ S% X- C& Z& ]. cNSGA2的算法不具有普遍性,下面参考课国外的课题小组的代码重新修改了内部冗余内容,使之能够自定义优化函数。7 b$ I* p8 g2 K, k" J: z; ?
8 N% Q3 K/ [# w1 q
- C- h$ e' ]# c# }
( A6 t: b8 D9 q. y; l
NSGA2的过程为:3 P* W/ n: A! S/ m
8 W0 [: ^/ A! ~5 T- t5 o
1、随机产生一个初始父代Po,在此基础上采用二元锦标赛选择、交叉和变异操作产生子代Qo, Po 和Qo群体规模均为N$ Z& [* m: U4 v5 E6 B1 D6 u7 r& Q
( z: o) D& s1 m6 S; ?+ q# O
2、将Pt和Qt并入到Rt中(初始时t=0),对Rt进行快速非支配解排序,构造其所有不同等级的非支配解集F1、F2……..& i5 Y' K9 f- n2 u( x
$ C" M5 U) J: u# b8 Q4 G n! F, J
3、按照需要计算Fi中所有个体的拥挤距离,并根据拥挤比较运算符构造Pt+1,直至Pt+1规模为N,图中的Fi为F31 F6 N* W4 J1 _# Y% A
% m: [3 r5 u7 r8 R6 N4 |
, |7 K; n# B0 o6 N下面是完整版的代码:
: R/ s7 d0 \2 y+ M* v% i
; H6 P0 ^. _5 i3 U! ~) Y①nsga2-optimization.m
2 L9 h7 B, k4 d3 b5 ?; V7 Q/ {2 A
1 u! t6 S/ I4 q5 U9 M- function nsga_2_optimization
- %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
- %此处可以更改
- %更多机器学习内容请访问omegaxyz.com
- pop = 500; %种群数量
- gen = 500; %迭代次数
- M = 2; %目标数量
- V = 30; %维度
- min_range = zeros(1, V); %下界
- max_range = ones(1,V); %上界
- %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
- chromosome = initialize_variables(pop, M, V, min_range, max_range);
- chromosome = non_domination_sort_mod(chromosome, M, V);
- for i = 1 : gen
- pool = round(pop/2);
- tour = 2;
- parent_chromosome = tournament_selection(chromosome, pool, tour);
- mu = 20;
- mum = 20;
- offspring_chromosome = genetic_operator(parent_chromosome,M, V, mu, mum, min_range, max_range);
- [main_pop,~] = size(chromosome);
- [offspring_pop,~] = size(offspring_chromosome);
- clear temp
- intermediate_chromosome(1:main_pop,:) = chromosome;
- intermediate_chromosome(main_pop + 1 : main_pop + offspring_pop,1 : M+V) = offspring_chromosome;
- intermediate_chromosome = non_domination_sort_mod(intermediate_chromosome, M, V);
- chromosome = replace_chromosome(intermediate_chromosome, M, V, pop);
- if ~mod(i,100)
- clc;
- fprintf('%d generations completed\n',i);
- end
- end
- if M == 2
- plot(chromosome(:,V + 1),chromosome(:,V + 2),'*');
- xlabel('f_1'); ylabel('f_2');
- title('Pareto Optimal Front');
- elseif M == 3
- plot3(chromosome(:,V + 1),chromosome(:,V + 2),chromosome(:,V + 3),'*');
- xlabel('f_1'); ylabel('f_2'); zlabel('f_3');
- title('Pareto Optimal SuRFace');
- end, \4 W* o/ S! M+ p7 A! {" L
8 d6 o* @/ H3 [) }
: B! l' {, ]) \1 a* L# }
②initialize_variables.m
9 g( X* H/ e' y$ {, q# }: y4 a* F: I1 A. j6 ]$ m
- function f = initialize_variables(N, M, V, min_range, max_range)
- min = min_range;
- max = max_range;
- K = M + V;
- for i = 1 : N
- for j = 1 : V
- f(i,j) = min(j) + (max(j) - min(j))*rand(1);
- end
- f(i,V + 1: K) = evaluate_objective(f(i,:), M, V);
- end
" J5 c7 L# A6 D1 u" q, h9 d/ r; J
# d1 `$ X9 t, t9 l( m! q7 K6 e/ k; P9 x8 f. S+ X1 I
③non_domination_sort_mod.m6 k5 r& ~5 J4 {; V3 `, R& M7 P
* C7 r5 z$ t% z$ K- `5 _3 Q- function f = non_domination_sort_mod(x, M, V)
- [N, ~] = size(x);
- clear m
- front = 1;
- F(front).f = [];
- individual = [];
- for i = 1 : N
- individual(i).n = 0;
- individual(i).p = [];
- for j = 1 : N
- dom_less = 0;
- dom_equal = 0;
- dom_more = 0;
- for k = 1 : M
- if (x(i,V + k) < x(j,V + k))
- dom_less = dom_less + 1;
- elseif (x(i,V + k) == x(j,V + k))
- dom_equal = dom_equal + 1;
- else
- dom_more = dom_more + 1;
- end
- end
- if dom_less == 0 && dom_equal ~= M
- individual(i).n = individual(i).n + 1;
- elseif dom_more == 0 && dom_equal ~= M
- individual(i).p = [individual(i).p j];
- end
- end
- if individual(i).n == 0
- x(i,M + V + 1) = 1;
- F(front).f = [F(front).f i];
- end
- end
- while ~isempty(F(front).f)
- Q = [];
- for i = 1 : length(F(front).f)
- if ~isempty(individual(F(front).f(i)).p)
- for j = 1 : length(individual(F(front).f(i)).p)
- individual(individual(F(front).f(i)).p(j)).n = ...
- individual(individual(F(front).f(i)).p(j)).n - 1;
- if individual(individual(F(front).f(i)).p(j)).n == 0
- x(individual(F(front).f(i)).p(j),M + V + 1) = ...
- front + 1;
- Q = [Q individual(F(front).f(i)).p(j)];
- end
- end
- end
- end
- front = front + 1;
- F(front).f = Q;
- end
- [temp,index_of_fronts] = sort(x(:,M + V + 1));
- for i = 1 : length(index_of_fronts)
- sorted_based_on_front(i,:) = x(index_of_fronts(i),:);
- end
- current_index = 0;
- %% Crowding distance
- for front = 1 : (length(F) - 1)
- distance = 0;
- y = [];
- previous_index = current_index + 1;
- for i = 1 : length(F(front).f)
- y(i,:) = sorted_based_on_front(current_index + i,:);
- end
- current_index = current_index + i;
- sorted_based_on_objective = [];
- for i = 1 : M
- [sorted_based_on_objective, index_of_objectives] = ...
- sort(y(:,V + i));
- sorted_based_on_objective = [];
- for j = 1 : length(index_of_objectives)
- sorted_based_on_objective(j,:) = y(index_of_objectives(j),:);
- end
- f_max = ...
- sorted_based_on_objective(length(index_of_objectives), V + i);
- f_min = sorted_based_on_objective(1, V + i);
- y(index_of_objectives(length(index_of_objectives)),M + V + 1 + i)...
- = Inf;
- y(index_of_objectives(1),M + V + 1 + i) = Inf;
- for j = 2 : length(index_of_objectives) - 1
- next_obj = sorted_based_on_objective(j + 1,V + i);
- previous_obj = sorted_based_on_objective(j - 1,V + i);
- if (f_max - f_min == 0)
- y(index_of_objectives(j),M + V + 1 + i) = Inf;
- else
- y(index_of_objectives(j),M + V + 1 + i) = ...
- (next_obj - previous_obj)/(f_max - f_min);
- end
- end
- end
- distance = [];
- distance(:,1) = zeros(length(F(front).f),1);
- for i = 1 : M
- distance(:,1) = distance(:,1) + y(:,M + V + 1 + i);
- end
- y(:,M + V + 2) = distance;
- y = y(:,1 : M + V + 2);
- z(previous_index:current_index,:) = y;
- end
- f = z();
1 I. K+ m# Q6 i2 r, i$ A H
& _# K0 j( x4 C* F& z8 p: m$ s6 k4 k8 K" Y! ]( F3 q1 N; I! v; j% c
④tournament_selection.m& a* d$ N& v8 e6 w& l4 Z/ T, ]
# d8 J8 ?0 _9 R4 Z7 o
- function f = tournament_selection(chromosome, pool_size, tour_size)
- [pop, variables] = size(chromosome);
- rank = variables - 1;
- distance = variables;
- for i = 1 : pool_size
- for j = 1 : tour_size
- candidate(j) = round(pop*rand(1));
- if candidate(j) == 0
- candidate(j) = 1;
- end
- if j > 1
- while ~isempty(find(candidate(1 : j - 1) == candidate(j)))
- candidate(j) = round(pop*rand(1));
- if candidate(j) == 0
- candidate(j) = 1;
- end
- end
- end
- end
- for j = 1 : tour_size
- c_obj_rank(j) = chromosome(candidate(j),rank);
- c_obj_distance(j) = chromosome(candidate(j),distance);
- end
- min_candidate = ...
- find(c_obj_rank == min(c_obj_rank));
- if length(min_candidate) ~= 1
- max_candidate = ...
- find(c_obj_distance(min_candidate) == max(c_obj_distance(min_candidate)));
- if length(max_candidate) ~= 1
- max_candidate = max_candidate(1);
- end
- f(i,:) = chromosome(candidate(min_candidate(max_candidate)),:);
- else
- f(i,:) = chromosome(candidate(min_candidate(1)),:);
- end
- end
( o S: d2 C5 E7 V* r0 z5 S 6 k& z P2 @- M% {0 P& P
6 y. B) i9 C, _, B( s⑤genetic_operator.m' z* c4 j& V2 c7 o$ T
8 T5 G( p9 H" Z9 f9 J# R7 R9 y- function f = genetic_operator(parent_chromosome, M, V, mu, mum, l_limit, u_limit)
- [N,m] = size(parent_chromosome);
- clear m
- p = 1;
- was_crossover = 0;
- was_mutation = 0;
- for i = 1 : N
- % With 90 % probability perform crossover
- if rand(1) < 0.9
- % Initialize the children to be null vector.
- child_1 = [];
- child_2 = [];
- % Select the first parent
- parent_1 = round(N*rand(1));
- if parent_1 < 1
- parent_1 = 1;
- end
- % Select the second parent
- parent_2 = round(N*rand(1));
- if parent_2 < 1
- parent_2 = 1;
- end
- % Make sure both the parents are not the same.
- while isequal(parent_chromosome(parent_1,:),parent_chromosome(parent_2,:))
- parent_2 = round(N*rand(1));
- if parent_2 < 1
- parent_2 = 1;
- end
- end
- % Get the chromosome information for each randomnly selected
- % parents
- parent_1 = parent_chromosome(parent_1,:);
- parent_2 = parent_chromosome(parent_2,:);
- % Perform corssover for each decision variable in the chromosome.
- for j = 1 : V
- % SBX (Simulated Binary Crossover).
- % For more information about SBX refer the enclosed pdf file.
- % Generate a random number
- u(j) = rand(1);
- if u(j) <= 0.5
- bq(j) = (2*u(j))^(1/(mu+1));
- else
- bq(j) = (1/(2*(1 - u(j))))^(1/(mu+1));
- end
- % Generate the jth element of first child
- child_1(j) = ...
- 0.5*(((1 + bq(j))*parent_1(j)) + (1 - bq(j))*parent_2(j));
- % Generate the jth element of second child
- child_2(j) = ...
- 0.5*(((1 - bq(j))*parent_1(j)) + (1 + bq(j))*parent_2(j));
- % Make sure that the generated element is within the specified
- % decision space else set it to the appropriate extrema.
- if child_1(j) > u_limit(j)
- child_1(j) = u_limit(j);
- elseif child_1(j) < l_limit(j)
- child_1(j) = l_limit(j);
- end
- if child_2(j) > u_limit(j)
- child_2(j) = u_limit(j);
- elseif child_2(j) < l_limit(j)
- child_2(j) = l_limit(j);
- end
- end
- child_1(:,V + 1: M + V) = evaluate_objective(child_1, M, V);
- child_2(:,V + 1: M + V) = evaluate_objective(child_2, M, V);
- was_crossover = 1;
- was_mutation = 0;
- % With 10 % probability perform mutation. Mutation is based on
- % polynomial mutation.
- else
- % Select at random the parent.
- parent_3 = round(N*rand(1));
- if parent_3 < 1
- parent_3 = 1;
- end
- % Get the chromosome information for the randomnly selected parent.
- child_3 = parent_chromosome(parent_3,:);
- % Perform mutation on eact element of the selected parent.
- for j = 1 : V
- r(j) = rand(1);
- if r(j) < 0.5
- delta(j) = (2*r(j))^(1/(mum+1)) - 1;
- else
- delta(j) = 1 - (2*(1 - r(j)))^(1/(mum+1));
- end
- % Generate the corresponding child element.
- child_3(j) = child_3(j) + delta(j);
- % Make sure that the generated element is within the decision
- % space.
- if child_3(j) > u_limit(j)
- child_3(j) = u_limit(j);
- elseif child_3(j) < l_limit(j)
- child_3(j) = l_limit(j);
- end
- end
- child_3(:,V + 1: M + V) = evaluate_objective(child_3, M, V);
- % Set the mutation flag
- was_mutation = 1;
- was_crossover = 0;
- end
- if was_crossover
- child(p,:) = child_1;
- child(p+1,:) = child_2;
- was_cossover = 0;
- p = p + 2;
- elseif was_mutation
- child(p,:) = child_3(1,1 : M + V);
- was_mutation = 0;
- p = p + 1;
- end
- end
- f = child;1 t* e" x: P l2 \9 A( O
K% C" }; n" Y) F⑥replace_chromosome.m
* U' f( g! f: Y/ j$ m* q8 O
, y) F" x u3 b- function f = replace_chromosome(intermediate_chromosome, M, V,pop)
- [N, m] = size(intermediate_chromosome);
- % Get the index for the population sort based on the rank
- [temp,index] = sort(intermediate_chromosome(:,M + V + 1));
- clear temp m
- % Now sort the individuals based on the index
- for i = 1 : N
- sorted_chromosome(i,:) = intermediate_chromosome(index(i),:);
- end
- % Find the maximum rank in the current population
- max_rank = max(intermediate_chromosome(:,M + V + 1));
- % Start adding each front based on rank and crowing distance until the
- % whole population is filled.
- previous_index = 0;
- for i = 1 : max_rank
- % Get the index for current rank i.e the last the last element in the
- % sorted_chromosome with rank i.
- current_index = max(find(sorted_chromosome(:,M + V + 1) == i));
- % Check to see if the population is filled if all the individuals with
- % rank i is added to the population.
- if current_index > pop
- % If so then find the number of individuals with in with current
- % rank i.
- remaining = pop - previous_index;
- % Get information about the individuals in the current rank i.
- temp_pop = ...
- sorted_chromosome(previous_index + 1 : current_index, :);
- % Sort the individuals with rank i in the descending order based on
- % the crowding distance.
- [temp_sort,temp_sort_index] = ...
- sort(temp_pop(:, M + V + 2),'descend');
- % Start filling individuals into the population in descending order
- % until the population is filled.
- for j = 1 : remaining
- f(previous_index + j,:) = temp_pop(temp_sort_index(j),:);
- end
- return;
- elseif current_index < pop
- % Add all the individuals with rank i into the population.
- f(previous_index + 1 : current_index, :) = ...
- sorted_chromosome(previous_index + 1 : current_index, :);
- else
- % Add all the individuals with rank i into the population.
- f(previous_index + 1 : current_index, :) = ...
- sorted_chromosome(previous_index + 1 : current_index, :);
- return;
- end
- % Get the index for the last added individual.
- previous_index = current_index;
- end
) d$ ~2 a2 |* Q% M2 S/ p- e , j F4 j$ K) a$ p* @
" r; W0 G( R% @- k0 a
⑦自定义评价函数(我选用的ZDT1函数)$ \0 }3 r |8 o, T: Z1 K" }
% X, s' y! } e
- function f = evaluate_objective(x, M, V)
- f = [];
- f(1) = x(1);
- g = 1;
- sum = 0;
- for i = 2:V
- sum = sum + x(i);
- end
- g = g + 9*sum / (V-1));
- f(2) = g * (1 - sqrt(x(1) / g));
- end# r2 p& N) s+ m- H
$ ?; ]9 |; ?2 l2 s2 K y6 u) {3 J
, {# _5 a% {: {# D u500个种群运行500代的结果:& S% U2 Z6 y. o4 ~
! J8 C- K' ]! ~, D2 d! M7 O
4 Z+ e* ]4 p2 o* V% t" d3 A, u
- y; x; O: R, g
$ f9 X$ ?5 K& V* t3 a( S
$ O$ t$ ]9 j. H3 N: I# A
( E/ o0 O7 J- Q, u3 f
+ n" o* O0 D6 U- u. h" U
( [' a+ Y) B( d9 v; w! I A
# F/ {0 r# ?4 b. ], ]8 S
$ B, }0 r' F8 x% k |
|