Week 4 - Genetic Algorithm

Aim: Write a code in MATLAB to optimise the stalagmite function and find the global maxima of the function.

Theory behind GA (Genetic Algorithm):

Five phases are considered in a genetic algorithm.

1) Initial population

2) Fitness function

3) Selection

4) Crossover

5) Mutation

Genetic Algorithms are a method that work based on Darwin’s theory of natural selection: The fittest survives, and has a greater chance of breeding and passing forward its genetic information.Let’s say we have an equation unsolvable by analytical methods, or an optimisation problem that takes too long to solve. First, the program randomly generates a population of solutions. Then, these solutions are evaluated by a function adequate to the problem, determining how well they ‘fit’. The ones with the best fit will have a greater chance of being selected to breed and generate the next population. Finally, the selected solutions are recombined amongst themselves creating a new population. Random shifts, just like mutations, also happen in order to introduce more variability. Whereas the first generation was completely random, the second is biased towards solutions with better results. After several iterations, they tend to converge to the optimum value.

In essence, it is very close to what evolution proposes. A population is evaluated by the environment; the fittest have a greater chance to breed and pass their ‘good’ genes along, which will then be recombined and mixed. The next generation will now have more of the genetic characteristics that represent a good fit in the current environment.GA’s are especially attractive to rapidly find solutions to difficult high-dimension problems, with little information about the search space. As a downside, it does not always return the optimal value, as any heuristic. Instead, it returns a value as good as the problem and the stopping criteria permit.

Syntax in MATLAB

Note: ga function only give minimum value (However we want to calculate maxima, hence refer my fitness function for maxima)

[x,fval] = ga(function,no_of_variables,A,b,Aeq,beq,lb,ub,nonlcon,IntCon,options)

1) x=optimum solution/variables

2) no_of_variables=Number of variables, specified as a positive integer

3) A,b= Linear inequality constraints

For example, to specify

2x1 + 2x2 ≤ 13

x1 - x2 ≤ 20

x1 + 3x2 ≤ 30

give these constraints:

A = [2,2;1,-1;1,3]

b = [13;20;30]

4) Aeq,beq=Linear equality constraints

For example, to specify

2x1 + 2x2 = 13

x1 - x2 = 20

x1 + 3x2 = 30

give these constraints:

Aeq = [2,2;1,-1;1,3]

beq = [13;20;30]

5) lb,ub

lb=Lower bounds

ub=Upper bounds

For example, to specify

1 <= x(1) <= 6

-3 <= x(2) <= 8

give these constraints:

lb=[1 -3]

ub=[6 8]

6) nonlcon : Use for nonlinear constraints

7) IntCon : to get solution in integer.

8) options: Minimize with non default options (Refer below program)

options = optimoptions('function','Constraint',limit,);

% stalagmite function
function [F] = stalagmite(input_vector)
         x=input_vector(1) ;
         y=input_vector(2) ;
         
         f1_x = (sin(5.1*pi*x+0.5)).^6;
         f1_y = (sin(5.1*pi*y+0.5)).^6;
         f2_x = exp(-4*log(2)*((x-0.0667).^2)/0.64);
         f2_y = exp(-4*log(2)*((y-0.0667).^2)/0.64);
         [F] = (f1_x.*f1_y.*f2_x.*f2_y);
end

%Genetic_Algorithm_Case_1

% Step 1
% clearing workspace, command window & closing all plots
close all
clear all
clc

% Step 2
% Define our search space
x = linspace(0,0.6,150);
y = linspace(0,0.6,150);

% Step 3
% Calculating a 2Dimentionl mesh
[xx, yy] = meshgrid(x,y);

% Step 4
%Evaluating the stalagmite function
for i = 1:length(xx)
    for j = 1:length(yy)
        input_vector(1) = xx(i,j);
        input_vector(2) = yy(i,j);
        f(i,j) = stalagmite(input_vector);
    end
end

tic
% Step 5
% Study1- Statical behaviar
num_cases =  50;
for i=1:num_cases
    [inputs, fopt(i)] = ga(@stalagmite,30);
    xopt(i) = inputs(1);
    yopt(i) = inputs(2);
end
study1_time = toc
F=f;
figure(1)
subplot(2,1,1)
hold on
surfc(x, y, F)
shading interp
plot3(xopt,yopt,fopt,'marker','o','markersize',5,'markerfacecolor','r')
title('Unbounded inputs')
xlabel('x')
ylabel('y')
zlabel('f')
subplot(2,1,2)
plot(fopt)
xlabel('Iterations')
ylabel('Function maximum')

Step 1: Clearing workspace, command window & closing all plots

Step 2: Define our search space by using linspace command

Step 3: Calculating a 2 Dimentionl mesh using meshgrid command

Step 4: Evaluating the stalagmite function for each value of x & y

Step 5: using GA command with unbounded inputs means optimum value may be come from outside of search space, Refer attached graph.

Result:-

%Genetic_Algorithm_Case_2

% Step 1
% clearing workspace, command window & closing all plots
close all
clear all
clc

% Step 2
% Define our search space
x = linspace(0,0.6,150);
y = linspace(0,0.6,150);

% Step 3
% Calculating a 2Dimentionl mesh
[xx, yy] = meshgrid(x,y);

% Step 4
%Evaluating the stalagmite function
for i = 1:length(xx)
    for j = 1:length(yy)
        input_vector(1) = xx(i,j);
        input_vector(2) = yy(i,j);
        f(i,j) = stalagmite(input_vector);
    end
end

% Step 6
% Study2- Statical buhaviour with upper and lover bounds
options = optimoptions('ga','populationSize',50);
tic
num_cases =  50;
for i = 1:num_cases
    [inputs,fopt(i)] = ga(@stalagmite,2,[],[],[],[],[0;0],[0.6;0.6]);
    xopt(i) = inputs(1);
    yopt(i) = inputs(2);
end
Study2_time = toc
figure(2)
F=f;
subplot(2,1,1)
hold on
surfc(x, y, F)
shading interp
plot3(xopt,yopt,fopt,'marker','o','markersize',5,'markerfacecolor','r')
title('Bounded inputs')
xlabel('x')
ylabel('y')
zlabel('f')
subplot(2,1,2)
plot(fopt)
xlabel('Iterations')
ylabel('Function Maximum')

Step 6: using GA command with bounded inputs means optimum value will be come from inside of search space, Refer attached graph. However there will be huge variation in optimum value.

Result:-

%Genetic_Algorithm_Case_3

% Step 1
% clearing workspace, command window & closing all plots
close all
clear all
clc

% Step 2
% Define our search space
x = linspace(0,0.6,150);
y = linspace(0,0.6,150);

% Step 3
% Calculating a 2Dimentionl mesh
[xx, yy] = meshgrid(x,y);

% Step 4
%Evaluating the stalagmite function
for i = 1:length(xx)
    for j = 1:length(yy)
        input_vector(1) = xx(i,j);
        input_vector(2) = yy(i,j);
        f(i,j) = stalagmite(input_vector);
    end
end

% Step 7
% Study3- Increasing GA iteration
options = optimoptions('ga','populationSize',200);
tic
num_cases =  50;
for i = 1:num_cases
    [inputs, fopt(i)] = ga(@stalagmite,2,[],[],[],[],[0;0],[0.6;0.6],[],[],...
        options);
    xopt(i) = inputs(1);
    yopt(i) = inputs(2);
end
study3_time = toc
figure(3)
F=f;
subplot(2,1,1)
hold on
surfc(x, y, F)
shading interp
plot3(xopt,yopt,fopt,'marker','o','markersize',5,'markerfacecolor','r')
title('Bounded inputs with modified population size')
xlabel('x')
ylabel('y')
zlabel('f')
subplot(2,1,2)
plot(fopt)
axis([0 50 0 1.5]);
xlabel('Iterations')
ylabel('Function Maximum')
FMAX=-min(fopt);
X=mean(xopt);
Y=mean(yopt);

Step 7: using GA command with bounded inputs means optimum value will be come from inside of search space, Refer attached graph. However this time to avoid huge variation in optimum value, we have increased initial population size to 200. ( If we do not specify initial population in program, MATLAB take 50 as initial population)

Result:-