Explore Genetic Algorithm using MATLAB

Aim: Explore Genetic Alorithim and calculate the golbal maxima of stalagmite function using Genetic Algorithm in MATLAB.

OBJECTIVE:

1. Write code for stalagmite function
2. Using Genetic algorithm builtin function write matlab code to search for global maxima
3. plot graph based on the different options used in Genetic Algorithm.

Working principle of Genetic Algorithm

1. The algorithm begines by creating a random initial population
2. The algorithm then creates a sequence of new population. 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:

a.Scores each member ofthe current population by computing its fitness value.These values are called the raw fitness scores.
b.scales the raw fitness scores to convert them into a more usable range of values. These scaled values are called expectation values.
c.select members, called parents, based on their expectation.
d.some of the individuals in the current population that have lower fitness are choosen as elite. These elite indivisulas are passed to the next population.
e. Produces children from 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.
f.Repalces the current population with the children to form the next generation
3. The algorithm stops when one of the stopping criteria is meet.

Genetic Algorithm uses the following conditions to determine when to stop

Generations - The algorithm stops when the number of generations reaches the value of Generations.

Time limit - The algorithm stops after running for an amount of time in seconds equal to Time limit

Fitness limit - The algorithm stops when the value ofthe fitness function for the best point in the current population is less than or equal to Fitness limit

Stall generations - The algorithm stops when the average relative change in the fitness function value over Stall generations is less than Funtion tolerance.

Stall time limit - The algorithms stops ifthere is no improvement in the objective function during an interval of time in seconds equal to Stall time limit.

Stall test- The stall condition is either average change or geometric weighted. For geometric weighted, the weighting function is (1/2)^n, where "n" is the number of generations prior to the current. Both Stall conditions apply to the relative change in the fitness function over Stall generations.

Function tolerance - The algorithm runs until the average relative change in the fitness function value over Stall generations is less than Function tolorance.

Constraint tolorance - The constraint tolorance is not used as stopping criterion. It is used to determine the feasibility with respect to nonlinear constraints.Also, max(sqrt(eps),ConstraintTolerance) determines feasibility with respect to linear constraints.

The algorithm stops as soon as any one of these conditions is meet.

Stalagmite function

Stalagmite is a type of rock formation that rises from the floor of the cave due to the accumulation of material deposited on the floor from ceiling drippings.

Stalagmite function is matematical replica with such geometrical similarities.

stalagmite function is a mathematical model which helps to generate local minima and maxima.
It is a three dimensional function and the way the function is designed is that it has a lot of local maxima or minima and one global maxima or minima.

MATLAB CODE

clear all
close all
clc

% Defining search space

x=linspace(0,0.6,150); % this command generates 150 data points from 0 to 0.6
y=linspace(0,0.6,150); % this command generates 150 data points from 0 to 0.6

[X,Y]=meshgrid(x,y);

% [X,Y]=meshgrid(x,y); returns 2D grid coordinates based on the coordinates contained in vector x and y. X is a matrix where each row is a copy of x, and Y is a matrix where each column is a copy of y. The grid %represented by the coordinates X and Y has length(y) rows and length(x) columns.

num_cases=50; % Number of iterations for which the function will evaluate the fitness

% creating the surface plot of stalagemite 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

% ------------------------------------------------------

tic %calculating the time elapsed for searching
% Study case-1 Statistical behaviour with unbounded inputs

for i=1:num_cases

[inputs,fopt(i)]=ga(@stalagmite,2); % stalagmite is a sub-routine written based on provided functions.
xopt(i)=inputs(1);
yopt(i)=inputs(2);
end

study1_time=toc;

figure(1) % figure creates a new figure window using default property values.

subplot(2,1,1);

% subplot(m,n,p) divides the current figure into an m-by-n grid & creates axes in the position specified by p.

hold on;

surf(X,Y,f);

shading interp;

plot3(xopt,yopt,fopt,'marker','o','markersize',5,'markerfacecolor','r') % 3d plot

title('unbound inputs')
subplot(2,1,2); % figure(1) is used to display two plots using subplot (2,1,1) and (2,1,2)
plot(fopt);
xlabel('iteration');
ylabel('funtion minimum');

% plot for above set of initial code lines.

% observe searching is +/- data availability range. So not effective method to obtain single solution.

% Further tuning code for better search with bounded range.

% Study case-2 Statistical behaviour with bounded inputs

for i=1:num_cases

[inputs,fopt(i)]=ga(@stalagmite,2,[],[],[],[],[0;0],[1;1]);

% limiting seaarch b/w (0,1) to (0,1) ->(x1,x2) & (y1,y2)

xopt(i)=inputs(1);
yopt(i)=inputs(2);
end
study1_time=toc;
figure(1)
subplot(2,1,1);
hold on;
surf(X,Y,f);

shading interp;

plot3(xopt,yopt,fopt,'marker','o','markersize',5,'markerfacecolor','r')

title('bound inputs')
subplot(2,1,2);
plot(fopt);
xlabel('iteration');
ylabel('funtion minimum');

% plot for above set of updated code lines.

% It is well observed that in modified code search lies with data range.

% however still single value of global maximum is not obtained. need to further tune the code.

% In next step, option is used with increased population data set from 150 to 300.

% Study case-3 Statistical behaviour with bounded inputs and incresed sample cases.

options=optimoptions('ga');
options=optimoptions(options,'populationsize',300);

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
study3_time=toc;
figure(3)
subplot(2,1,1);
hold on;
surf(X,Y,f);

shading interp;

plot3(xopt,yopt,fopt,'marker','o','markersize',5,'markerfacecolor','r')

title('bound inputs')
subplot(2,1,2);
plot(fopt);
xlabel('iteration');
ylabel('funtion minimum');

% with modified set of code now global maximum value has been searched which was aim of this Challenge.

% Also in subplot2, observe function minimum is close to 0.9999. variation is in 7th decimel point onwards.

	Time for search, sec
Study-1	2.3456
Study-2	3.1462
Study-3	12.2953

% With change in code line to get Global max value, time for searching is also increased.

---------------------------------------------------------------------------------------------------

% Sub routine program for Stalagmite

function[f]=stalagmite(input_vector)

x=input_vector(1);
y=input_vector(2);

ff_1x=(sin(5.1*pi*x+0.5))^6;
ff_1y=(sin(5.1*pi*y+0.5))^6;

ff_2x=exp((-4*log(2))*(((x-0.0667)^2)/(0.64)));
ff_2y=exp((-4*log(2))*(((y-0.0667)^2)/(0.64)));

f=-(ff_1x*ff_1y*ff_2x*ff_2y);

% -ve sign is included to ensure code finds the maximum value.

% by default GA is built to get minimum value.

--------------------------------------------------------------------------------------------------

MATLAB FIlE LINK :

https://drive.google.com/drive/folders/1sQDa4AVK5xpbBHxFWt3ZrjT5THW8i5Cx?usp=sharing