function f = crowding_distance(x,problem)
% This function calculates the crowding distance
[N,M] = size(x);
switch problem
case 1
M = 2;
V = 6;
case 2
M = 3;
V = 12;
end
% Crowding distance for each front
for i = 1 : length(F(front).f)
y(i,:) = x(F(front).f(i),:);
end
for i = 1 : M
[sorted(i).individual,sorted(i).index] = sort(y(:,V + i));
distance(sorted(i).index(1)).individual = Inf;
distance(sorted(i).index(length(sorted(i).index))).individual = Inf;
end
[num,len] = size(y);
% Initialize all the distance of individuals as zero.
for i = 1 : M
for j = 2 : num - 1
distance(j).individual = 0;
end
objective(i).range = ...
sorted(i).individual(length(sorted(i).individual)) - ...
sorted(i).individual(1);
% Maximum and minimum objectives value for the ith objective
end
% Caluclate the crowding distance for front one.
for i = 1 : M
for j = 2 : num - 1
distance(j).individual = distance(j).individual + ...
(sorted(i).individual(j + 1) - sorted(i).individual(j - 1))/...
objective(i).range;
y(sorted(i).index(j),M + V + 2) = distance(j).individual;
end
end
</PRE
Published withMATLAB® 7.0
function f =genetic_operator(parent_chromosome,pro,mu,mum);
% This function is utilized to produce offsprings from parent chromosomes.
% The genetic operators corssover and mutation which are carried out with
% slight modifications from the original design. For more information read
% the document enclosed.
[N,M] = size(parent_chromosome);
switch pro
case 1
M = 2;
V = 6;
case 2
M = 3;
V = 12;
end
p = 1;
was_crossover = 0;
was_mutation = 0;
l_limit = 0;
u_limit = 1;
for i = 1 : N
if rand(1) < 0.9
child_1 = [];
child_2 = [];
parent_1 = round(N*rand(1));
if parent_1 < 1
parent_1 = 1;
end
parent_2 = round(N*rand(1));
if parent_2 < 1
parent_2 = 1;
end
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
parent_1 = parent_chromosome(parent_1,:);
parent_2 = parent_chromosome(parent_2,:);
for j = 1 : V
% SBX (Simulated Binary Crossover)
% 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
child_1(j) = ...
0.5*(((1 + bq(j))*parent_1(j)) + (1 - bq(j))*parent_2(j));
child_2(j) = ...
0.5*(((1 - bq(j))*parent_1(j)) + (1 + bq(j))*parent_2(j));
if child_1(j) > u_limit
child_1(j) = u_limit;
elseif child_1(j) < l_limit
child_1(j) = l_limit;
end
if child_2(j) > u_limit
child_2(j) = u_limit;
elseif child_2(j) < l_limit
child_2(j) = l_limit;
end
end
child_1(:,V + 1: M + V) = evaluate_objective(child_1,pro);
child_2(:,V + 1: M + V) = evaluate_objective(child_2,pro);
was_crossover = 1;
was_mutation = 0;
else
parent_3 = round(N*rand(1));
if parent_3 < 1
parent_3 = 1;
end
% Make sure that the mutation does not result in variables out of
% the search space. For both the MOP's the range for decision space
% is [0,1]. In case different variables have different decision
% space each variable can be assigned a range.
child_3 = parent_chromosome(parent_3,:);
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
child_3(j) = child_3(j) + delta(j);
if child_3(j) > u_limit
child_3(j) = u_limit;
elseif child_3(j) < l_limit
child_3(j) = l_limit;
end
end
child_3(:,V + 1: M + V) = evaluate_objective(child_3,pro);
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;
Published with MATLAB® 7.0
function f = initialize_variables(N,problem)
% function f = initialize_variables(N,problem)
% N - Population size
% problem - takes integer values 1 and 2 where,
% '1' for MOP1
% '2' for MOP2
%
% This function initializes the population with N individuals and each
% individual having M decision variables based on the selected problem.
% M = 6 for problem MOP1 and M = 12 for problem MOP2. The objective space
% for MOP1 is 2 dimensional while for MOP2 is 3 dimensional.
% Both the MOP's has 0 to 1 as its range for all the decision variables.
min = 0;
max = 1;
switch problem
case 1
M = 6;
K = 8;
case 2
M = 12;
K = 15;
end
for i = 1 : N
% Initialize the decision variables
for j = 1 : M
f(i,j) = rand(1); % i.e f(i,j) = min + (max - min)*rand(1);
end
% Evaluate the objective function
f(i,M + 1: K) = evaluate_objective(f(i,:),problem);
end
</PRE
Published withMATLAB® 7.0
Non-Donimation Sort
This function sort the current popultion based onnon-domination. All the individuals in the first front are given a rank of 1,the second front individuals are assigned rank 2 and so on. After assigning therank the crowding in each front is calculated.
[N,M] = size(x);
switch problem
case 1
M = 2;
V = 6;
case 2
M = 3;
V = 12;
end
front = 1;
% There is nothing to this assignment, used only to manipulate easily in
% MATLAB.
F(front).f = [];
individual = [];
for i = 1 : N
% Number of individuals that dominate this individual
individual(i).n = 0;
% Individuals which this individual dominate
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
% Find the subsequent fronts
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;
% Find the crowding distance for each individual in each front
for front = 1 : (length(F) - 1)
objective = [];
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;
% Sort each individual based on the objective
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();
Published with MATLAB® 7.0
Main Function
Main program to run the NSGA-II MOEA. Read thecorresponding documentation to learn more about multiobjective optimizationusing evolutionary algorithms. initialize_variables has two arguments; Firstbeing the population size and the second the problem number. '1' corresponds toMOP1 and '2' corresponds to MOP2.
Contents
- Initialize the variables
- Sort the initialized population
- Start the evolution process
- Result
- Visualize
Initialize the variables
Declare the variables and initialize their values pop -population gen - generations pro - problem number
pop = 200;
gen = 1;
pro = 1;
switch pro
case 1
% M is the number of objectives.
M = 2;
% V is the number of decision variables. In this case it is
% difficult to visualize the decision variables space while the
% objective space is just two dimensional.
V = 6;
case 2
M = 3;
V = 12;
end
% Initialize the population
chromosome = initialize_variables(pop,pro);
Sort the initialized population
Sort the population using non-domination-sort. This returnstwo columns for each individual which are the rank and the crowding distancecorresponding to their position in the front they belong.
chromosome = non_domination_sort_mod(chromosome,pro);
Start the evolution process
% The following are performed in each generation
% Select the parents
% Perfrom crossover and Mutation operator
% Perform Selection
for i = 1 : gen
% Select the parents
% Parents are selected for reproduction to generate offspring. The
% original NSGA-II uses a binary tournament selection based on the
% crowded-comparision operator. The arguments are
% pool - size of the mating pool. It is common to have this to be half the
% population size.
% tour - Tournament size. Original NSGA-II uses a binary tournament
% selection, but to see the effect of tournament size this is kept
% arbitary, to be choosen by the user.
pool = round(pop/2);
tour = 2;
parent_chromosome = tournament_selection(chromosome,pool,tour);
% Perfrom crossover and Mutation operator
% The original NSGA-II algorithm uses Simulated Binary Crossover (SBX) and
% Polynomial crossover. Crossover probability pc = 0.9 and mutation
% probability is pm = 1/n, where n is the number of decision variables.
% Both real-coded GA and binary-coded GA are implemented in the original
% algorithm, while in this program only the real-coded GA is considered.
% The distribution indeices for crossover and mutation operators as mu = 20
% and mum = 20 respectively.
mu = 20;
mum = 20;
offspring_chromosome = genetic_operator(parent_chromosome,pro,mu,mum);
% Intermediate population
% Intermediate population is the combined population of parents and
% offsprings of the current generation. The population size is almost 1 and
% half times the initial population.
[main_pop,temp] = size(chromosome);
[offspring_pop,temp] = size(offspring_chromosome);
intermediate_chromosome(1:main_pop,:) = chromosome;
intermediate_chromosome(main_pop + 1 : main_pop + offspring_pop,1 : M+V) = ...
offspring_chromosome;
% Non-domination-sort of intermediate population
% The intermediate population is sorted again based on non-domination sort
% before the replacement operator is performed on the intermediate
% population.
intermediate_chromosome = ...
non_domination_sort_mod(intermediate_chromosome,pro);
% Perform Selection
% Once the intermediate population is sorted only the best solution is
% selected based on it rank and crowding distance. Each front is filled in
% ascending order until the addition of population size is reached. The
% last front is included in the population based on the individuals with
% least crowding distance
chromosome = replace_chromosome(intermediate_chromosome,pro,pop);
if ~mod(i,10)
fprintf('%d\n',i);
end
end
Result
Save the result in ASCII text format.
save solution.txt chromosome -ASCII
Visualize
The following is used to visualize the result for the givenproblem.
switch pro
case 1
plot(chromosome(:,V + 1),chromosome(:,V + 2),'*');
title('MOP1 using NSGA-II');
xlabel('f(x_1)');
ylabel('f(x_2)');
case 2
plot3(chromosome(:,V + 1),chromosome(:,V + 2),chromosome(:,V + 3),'*');
title('MOP2 using NSGA-II');
xlabel('f(x_1)');
ylabel('f(x_2)');
zlabel('f(x_3)');
end
Published with MATLAB® 7.0
replace_chromosome(intermediate_chromosome,pro,pop)
This function replaces the chromosomes based on rank andcrowding distance. Initially until the population size is reached each front isadded one by one until addition of a complete front which results in exceedingthe population size. At this point the chromosomes in that front is addedsubsequently to the population based on crowding distance.
[N,V] = size(intermediate_chromosome);
switch pro
case 1
M = 2;
V = 6;
case 2
M = 3;
V = 12;
end
% Get the index for the population sort based on the rank
[temp,index] = sort(intermediate_chromosome(:,M + V + 1));
% 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
current_index = max(find(sorted_chromosome(:,M + V + 1) == i));
if current_index > pop
remaining = pop - previous_index;
temp_pop = ...
sorted_chromosome(previous_index + 1 : current_index, :);
[temp_sort,temp_sort_index] = ...
sort(temp_pop(:, M + V + 2),'descend');
for j = 1 : remaining
f(previous_index + j,:) = temp_pop(temp_sort_index(j),:);
end
return;
elseif current_index < pop
f(previous_index + 1 : current_index, :) = ...
sorted_chromosome(previous_index + 1 : current_index, :);
else
f(previous_index + 1 : current_index, :) = ...
sorted_chromosome(previous_index + 1 : current_index, :);
return;
end
previous_index = current_index;
end
Published with MATLAB® 7.0
function f =selection_individuals(chromosome,pool_size,tour_size)
% function selection_individuals(chromosome,pool_size,tour_size) is the
% selection policy for selecting the individuals for the mating pool. The
% selection is based on tournament selection. Argument 'chromosome' is the
% current generation population from which the individuals are selected to
% form a mating pool of size 'pool_size' after performing tournament
% selection, with size of the tournament being 'tour_size'. By varying the
% tournament size the selection pressure can be adjusted.
[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
Published with MATLAB® 7.0
function f = crowding_distance(x,problem)
% This function calculates the crowding distance
[N,M] = size(x);
switch problem
case 1
M = 2;
V = 6;
case 2
M = 3;
V = 12;
end
% Crowding distance for each front
for i = 1 : length(F(front).f)
y(i,:) = x(F(front).f(i),:);
end
for i = 1 : M
[sorted(i).individual,sorted(i).index] = sort(y(:,V + i));
distance(sorted(i).index(1)).individual = Inf;
distance(sorted(i).index(length(sorted(i).index))).individual = Inf;
end
[num,len] = size(y);
% Initialize all the distance of individuals as zero.
for i = 1 : M
for j = 2 : num - 1
distance(j).individual = 0;
end
objective(i).range = ...
sorted(i).individual(length(sorted(i).individual)) - ...
sorted(i).individual(1);
% Maximum and minimum objectives value for the ith objective
end
% Caluclate the crowding distance for front one.
for i = 1 : M
for j = 2 : num - 1
distance(j).individual = distance(j).individual + ...
(sorted(i).individual(j + 1) - sorted(i).individual(j - 1))/...
objective(i).range;
y(sorted(i).index(j),M + V + 2) = distance(j).individual;
end
end
</PRE
Published with MATLAB® 7.0
function f = genetic_operator(parent_chromosome,pro,mu,mum);
% This function is utilized to produce offsprings from parent chromosomes.
% The genetic operators corssover and mutation which are carried out with
% slight modifications from the original design. For more information read
% the document enclosed.
[N,M] = size(parent_chromosome);
switch pro
case 1
M = 2;
V = 6;
case 2
M = 3;
V = 12;
end
p = 1;
was_crossover = 0;
was_mutation = 0;
l_limit = 0;
u_limit = 1;
for i = 1 : N
if rand(1) < 0.9
child_1 = [];
child_2 = [];
parent_1 = round(N*rand(1));
if parent_1 < 1
parent_1 = 1;
end
parent_2 = round(N*rand(1));
if parent_2 < 1
parent_2 = 1;
end
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
parent_1 = parent_chromosome(parent_1,:);
parent_2 = parent_chromosome(parent_2,:);
for j = 1 : V
% SBX (Simulated Binary Crossover)
% 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
child_1(j) = ...
0.5*(((1 + bq(j))*parent_1(j)) + (1 - bq(j))*parent_2(j));
child_2(j) = ...
0.5*(((1 - bq(j))*parent_1(j)) + (1 + bq(j))*parent_2(j));
if child_1(j) > u_limit
child_1(j) = u_limit;
elseif child_1(j) < l_limit
child_1(j) = l_limit;
end
if child_2(j) > u_limit
child_2(j) = u_limit;
elseif child_2(j) < l_limit
child_2(j) = l_limit;
end
end
child_1(:,V + 1: M + V) = evaluate_objective(child_1,pro);
child_2(:,V + 1: M + V) = evaluate_objective(child_2,pro);
was_crossover = 1;
was_mutation = 0;
else
parent_3 = round(N*rand(1));
if parent_3 < 1
parent_3 = 1;
end
% Make sure that the mutation does not result in variables out of
% the search space. For both the MOP's the range for decision space
% is [0,1]. In case different variables have different decision
% space each variable can be assigned a range.
child_3 = parent_chromosome(parent_3,:);
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
child_3(j) = child_3(j) + delta(j);
if child_3(j) > u_limit
child_3(j) = u_limit;
elseif child_3(j) < l_limit
child_3(j) = l_limit;
end
end
child_3(:,V + 1: M + V) = evaluate_objective(child_3,pro);
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;
Published with MATLAB® 7.0
function f = initialize_variables(N,problem)
% function f = initialize_variables(N,problem)
% N - Population size
% problem - takes integer values 1 and 2 where,
% '1' for MOP1
% '2' for MOP2
%
% This function initializes the population with N individuals and each
% individual having M decision variables based on the selected problem.
% M = 6 for problem MOP1 and M = 12 for problem MOP2. The objective space
% for MOP1 is 2 dimensional while for MOP2 is 3 dimensional.
% Both the MOP's has 0 to 1 as its range for all the decision variables.
min = 0;
max = 1;
switch problem
case 1
M = 6;
K = 8;
case 2
M = 12;
K = 15;
end
for i = 1 : N
% Initialize the decision variables
for j = 1 : M
f(i,j) = rand(1); % i.e f(i,j) = min + (max - min)*rand(1);
end
% Evaluate the objective function
f(i,M + 1: K) = evaluate_objective(f(i,:),problem);
end
</PRE
Published with MATLAB® 7.0
Non-Donimation Sort
This function sort the current popultion based on non-domination. All the individuals in the first front are given a rank of 1, the second front individuals are assigned rank 2 and so on. After assigning the rank the crowding in each front is calculated.
[N,M] = size(x);
switch problem
case 1
M = 2;
V = 6;
case 2
M = 3;
V = 12;
end
front = 1;
% There is nothing to this assignment, used only to manipulate easily in
% MATLAB.
F(front).f = [];
individual = [];
for i = 1 : N
% Number of individuals that dominate this individual
individual(i).n = 0;
% Individuals which this individual dominate
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
% Find the subsequent fronts
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;
% Find the crowding distance for each individual in each front
for front = 1 : (length(F) - 1)
objective = [];
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;
% Sort each individual based on the objective
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();
Published with MATLAB® 7.0
Main Function
Main program to run the NSGA-II MOEA. Read the corresponding documentation to learn more about multiobjective optimization using evolutionary algorithms. initialize_variables has two arguments; First being the population size and the second the problem number. '1' corresponds to MOP1 and '2' corresponds to MOP2.
Contents
• Initialize the variables
• Sort the initialized population
• Start the evolution process
• Result
• Visualize
Initialize the variables
Declare the variables and initialize their values pop - population gen - generations pro - problem number
pop = 200;
gen = 1;
pro = 1;
switch pro
case 1
% M is the number of objectives.
M = 2;
% V is the number of decision variables. In this case it is
% difficult to visualize the decision variables space while the
% objective space is just two dimensional.
V = 6;
case 2
M = 3;
V = 12;
end
% Initialize the population
chromosome = initialize_variables(pop,pro);
Sort the initialized population
Sort the population using non-domination-sort. This returns two columns for each individual which are the rank and the crowding distance corresponding to their position in the front they belong.
chromosome = non_domination_sort_mod(chromosome,pro);
Start the evolution process
% The following are performed in each generation
% Select the parents
% Perfrom crossover and Mutation operator
% Perform Selection
for i = 1 : gen
% Select the parents
% Parents are selected for reproduction to generate offspring. The
% original NSGA-II uses a binary tournament selection based on the
% crowded-comparision operator. The arguments are
% pool - size of the mating pool. It is common to have this to be half the
% population size.
% tour - Tournament size. Original NSGA-II uses a binary tournament
% selection, but to see the effect of tournament size this is kept
% arbitary, to be choosen by the user.
pool = round(pop/2);
tour = 2;
parent_chromosome = tournament_selection(chromosome,pool,tour);
% Perfrom crossover and Mutation operator
% The original NSGA-II algorithm uses Simulated Binary Crossover (SBX) and
% Polynomial crossover. Crossover probability pc = 0.9 and mutation
% probability is pm = 1/n, where n is the number of decision variables.
% Both real-coded GA and binary-coded GA are implemented in the original
% algorithm, while in this program only the real-coded GA is considered.
% The distribution indeices for crossover and mutation operators as mu = 20
% and mum = 20 respectively.
mu = 20;
mum = 20;
offspring_chromosome = genetic_operator(parent_chromosome,pro,mu,mum);
% Intermediate population
% Intermediate population is the combined population of parents and
% offsprings of the current generation. The population size is almost 1 and
% half times the initial population.
[main_pop,temp] = size(chromosome);
[offspring_pop,temp] = size(offspring_chromosome);
intermediate_chromosome(1:main_pop,:) = chromosome;
intermediate_chromosome(main_pop + 1 : main_pop + offspring_pop,1 : M+V) = ...
offspring_chromosome;
% Non-domination-sort of intermediate population
% The intermediate population is sorted again based on non-domination sort
% before the replacement operator is performed on the intermediate
% population.
intermediate_chromosome = ...
non_domination_sort_mod(intermediate_chromosome,pro);
% Perform Selection
% Once the intermediate population is sorted only the best solution is
% selected based on it rank and crowding distance. Each front is filled in
% ascending order until the addition of population size is reached. The
% last front is included in the population based on the individuals with
% least crowding distance
chromosome = replace_chromosome(intermediate_chromosome,pro,pop);
if ~mod(i,10)
fprintf('%d\n',i);
end
end
Result
Save the result in ASCII text format.
save solution.txt chromosome -ASCII
Visualize
The following is used to visualize the result for the given problem.
switch pro
case 1
plot(chromosome(:,V + 1),chromosome(:,V + 2),'*');
title('MOP1 using NSGA-II');
xlabel('f(x_1)');
ylabel('f(x_2)');
case 2
plot3(chromosome(:,V + 1),chromosome(:,V + 2),chromosome(:,V + 3),'*');
title('MOP2 using NSGA-II');
xlabel('f(x_1)');
ylabel('f(x_2)');
zlabel('f(x_3)');
end
Published with MATLAB® 7.0
replace_chromosome(intermediate_chromosome,pro,pop)
This function replaces the chromosomes based on rank and crowding distance. Initially until the population size is reached each front is added one by one until addition of a complete front which results in exceeding the population size. At this point the chromosomes in that front is added subsequently to the population based on crowding distance.
[N,V] = size(intermediate_chromosome);
switch pro
case 1
M = 2;
V = 6;
case 2
M = 3;
V = 12;
end
% Get the index for the population sort based on the rank
[temp,index] = sort(intermediate_chromosome(:,M + V + 1));
% 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
current_index = max(find(sorted_chromosome(:,M + V + 1) == i));
if current_index > pop
remaining = pop - previous_index;
temp_pop = ...
sorted_chromosome(previous_index + 1 : current_index, :);
[temp_sort,temp_sort_index] = ...
sort(temp_pop(:, M + V + 2),'descend');
for j = 1 : remaining
f(previous_index + j,:) = temp_pop(temp_sort_index(j),:);
end
return;
elseif current_index < pop
f(previous_index + 1 : current_index, :) = ...
sorted_chromosome(previous_index + 1 : current_index, :);
else
f(previous_index + 1 : current_index, :) = ...
sorted_chromosome(previous_index + 1 : current_index, :);
return;
end
previous_index = current_index;
end
Published with MATLAB® 7.0
function f = selection_individuals(chromosome,pool_size,tour_size)
% function selection_individuals(chromosome,pool_size,tour_size) is the
% selection policy for selecting the individuals for the mating pool. The
% selection is based on tournament selection. Argument 'chromosome' is the
% current generation population from which the individuals are selected to
% form a mating pool of size 'pool_size' after performing tournament
% selection, with size of the tournament being 'tour_size'. By varying the
% tournament size the selection pressure can be adjusted.
[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
Published with MATLAB® 7.0
扫一扫
11万+
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
