MATLAB Program to optimize the stalagmite function and find global maxima using Genetic algorithm

Aim:- To optimize the stalagmite function and find the global maxima of the function.

Theory:-

Genetic Algorithm (GA) is a search-based optimization technique based on the principles of Genetics and Natural Selection. It is used in optimization problems wherein we have to maximize or minimize a given objective function value under a given set of constraints.

How the genetic algorithm works:

• The algorithm begins by creating a random initial population.

• The algorithm then creates a sequence of new populations. At each step, the algorithm uses the individuals in the current generation to create the next population. To create the new population, the algorithm performs the following steps:

  1. Scores each member of the current population by computing its fitness value. These values are called the raw fitness scores.

  2. Scales the raw fitness scores to convert them into a more usable range of values. These scaled values are called expectation values.

  3. Selects members, called parents, based on their expectation.

  4. Some of the individuals in the current population that have lower fitness are chosen as elite. These elite individuals are passed to the next population.

  5. Produces children from the parents. Children are produced either by making random changes to a single parent—mutation—or by combining the vector entries of a pair of parents—crossover.

  6. Replaces the current population with the children to form the next generation.

Syntax of ga

Main code:-

clear all
close all
clc

%defining our search space
x = linspace(0,.6,150);
y = linspace(0,.6,150);

[X Y]=meshgrid(x,y);     %creating 2d mesh
num_cases=50;            %no. of iterations

f=get_stalagmite;        %calling stalagmite function

% Study-1 Statistical Behavior
tic                               %calculation of time
for i = 1:num_cases               %using loop for finding global maxima of the stalagmite function
    [inputs, fopt(i)] = ga(@stalagmite, 2);
    xopt(i) = inputs(1);
    yopt(i) = inputs(2);
end
study_time_1 = toc

% Figure 1 is plotted from results of study 1
figure(1)
subplot(2 ,1 ,1)
hold on
grid on
surfc(X,Y,-f)    %creating 3D surface plot of function
view(17,10)      %views plot by horizontal rotation of 17deg & vertical elevation of 10deg.
shading interp
plot3(xopt,yopt,-fopt,'marker','d','markersize',5,'markerfacecolor','r') %marking points of max. value.
title('Study 1: Unbounded inputs')
xlabel('X-axis')
ylabel('Y-axis')
subplot(2,1,2)
plot(-fopt)   
grid on
xlabel('Iterations')
ylabel('Function maximum')


% Study-2 Statistical behavior with upper and lower bound
tic
for i = 1:num_cases
    [inputs, fopt(i)] = ga(@stalagmite,2,[],[],[],[],[0;0],[1;1]);
    xopt(i) = inputs(1);
    yopt(i) = inputs(2);
end
study_time_2 = toc

% Figure 2 is plotted from results of study 2
figure(2)
subplot(2 ,1 ,1)
hold on
grid on
surfc(X,Y,-f)    %creating 3D surface plot of function
view(17,10)
shading interp
plot3(xopt,yopt,-fopt,'marker','d','markersize',5,'markerfacecolor','r') %marking points of max. value.
title('Study 2: With upper and lower bound inputs')
xlabel('X-axis')
ylabel('Y-axis')
subplot(2,1,2)
plot(-fopt)
grid on
xlabel('Iterations')
ylabel('Function maximum')


% Study-3 Statistical behaviour with upper and lower bound and increasing the ga iterations 
options = optimoptions('ga');
options = optimoptions(options,'populationsize',180);

tic
for  i = 1:num_cases
     [inputs , fopt(i)] = ga(@stalagmite,2,[],[],[],[],[0;0],[1;1],[],options);
     xopt(i) = inputs(1);
     yopt(i) = inputs(2);
end
study_time_3 = toc

% Figure 3 is plotted from results of study 3
figure(3)
subplot(2 ,1 ,1)
hold on
grid on
surfc(X,Y,-f)    %creating 3D surface plot of function
view(17,10)
shading interp
plot3(xopt,yopt,-fopt,'marker','d','markersize',5,'markerfacecolor','r') %marking points of max. value.
title('Study 3: By incresing number of iterations')
xlabel('X-axis')
ylabel('Y-axis')
subplot(2,1,2)
plot(-fopt)
grid on
xlabel('Iterations')
ylabel('Function maximum')

• First we define the search space(x,y) between 0 to 0.6 with 150 intervals using linspace() and generate 2D array using meshgrid() then define no. of iterations.
• Then we call the stalagmite function and store its value in 'f'.

Study 1

• tic-toc command is used to find the time taken to run the case_study.
• num_cases = 50, that will run 50 iterations by using for loop.
• Output values of variables are stored in xopt(i) & yopt(i) and maximum value in fopt(i).
• The ga syntax used is 'x = ga(fun,nvars)'.
• Then by using subplot() we stalagmite function by marking the solution points using plot3(), marker and the fopt values.

Study 2

• tic-toc command is used to find the time taken to run the case_study.
• num_cases = 50, that will run 50 iterations by using for loop and using upper and lower bounds.
• Output values of variables are stored in xopt(i) & yopt(i) and maximum value in fopt(i).
• The ga syntax used is 'x = ga(fun,nvars,A,b,Aeq,beq,lb,ub)' where lb, ub ar lower and upper bound respetively.
• Then by using subplot() we stalagmite function by marking the solution points using plot3(), marker and the fopt values.

Study 3

• Using optimoptions we increase the population size to 180 which is assigned to options.
• tic-toc command is used to find the time taken to run the case_study.
• num_cases = 50, that will run 50 iterations by using for loop and we also use upper and lower bounds.
• Output values of variables are stored in xopt(i) & yopt(i) and maximum value in fopt(i).
• The ga syntax used is 'x = ga(fun,nvars,A,b,[],[],lb,ub,nonlcon,IntCon,options)'.
• Then by using subplot() we stalagmite function by marking the solution points using plot3(), marker and the fopt values.

Function codes:-

• Stalagmite function

%stalagmite function

function [f]=stalagmite(input_vector)

x=input_vector(1);
y=input_vector(2);
f1x=(sin(5.1*pi*x+0.5))^6;
f1y=(sin(5.1*pi*y+0.5))^6;
f2x=exp(-4*log(2)*((x-0.0667)^2)/0.64);
f2y=exp(-4*log(2)*((y-0.0667)^2)/0.64);

f = -(f1x*f2x*f1y*f2y);
end

• Function named stalagmite was defined with input_vector as argument that accepts the vector passed from function call.
• input_vector conatin the x and y values respectively.
• The function is split into 4 terms and later all are multiplied and this is given as output.
• Here we use -ve sign for output because ga only minimizes the function but we need to find maxima.

• Calling the stalagmite function

%Function for calling stalagmite function

function f = get_stalagmite
x = linspace(0,0.6,150);
y = linspace(0,0.6,150);
[X Y] = meshgrid(x,y);

% Evaluating the stalagmite function
for i = 1:length(X)
    for j = 1:length(Y)
        input_vector(1)=X(i,j);
        input_vector(2)=Y(i,j);
        f(i,j)= stalagmite(input_vector);
    end
end

• This function is used to call stalagmite function.
• First, we define the search space(x,y) between 0 to 0.6 with 150 intervals using linspace() and generate 2D array using meshgrid().
• Then at each point the stalagmite function is evaluated and output is stored in 'f'.

Output:-

Study 1

We can see that the maximum value is continuiusly varying between 0 to 1 for each iteration. This is overcomed by using bounds.

Study 2

In this also we can see that the maximum value is varying but only between 0.5 to 1. This is because of using bounds. We can still optimize it by increasing the population size.

Study 3

We can that the maximum value is almost constant and equal to '1' and it has changed to '0.85' only thrice in 50 iterations. Hence we can say that the value of global maxima is '1'.

Workspace

We can see that the study time in study 2 (inputs with upper and lower bounds) increased slightly compared to study 1 (unbounded inputs). But it increased by alomost 5 times in study 3 (increased the population size).

References:-

• https://in.mathworks.com/help/gads/fitness-function-forms.html#gafitness-3
• https://in.mathworks.com/help/gads/ga.html
• https://in.mathworks.com/help/gads/how-the-genetic-algorithm-works.html