Challenge

Genetic Algorithm using MATLAB

SSaurav Sankhe27 Jun 2020

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

Theory: A stalagmite is an upward-growing mound of mineral deposits that have precipitated from water dripping onto the floor of a cave. Most stalagmites have rounded or flattened tips.

Below figures shows how the actual stalagmite structure is modelled based on its equations in mathematics.

`f(x,y)=f_(1,x)∙f_(1,y)∙f_(2,x)∙f_(2,y)`

`f_(1,x)=[sin⁡(5.1πx+0.5) ]^6` `f_(1,y)=[sin⁡(5.1πy+0.5) ]^6`

`f_(2,x)=exp[-4ln⁡(2)∙(x-0.0667)^2/0.64]`

`f_(2,y)=exp[-4ln⁡(2)∙(y-0.0667)^2/0.64]`

Genetic Algorithm

A genetic algorithm (GA) is a method for solving both constrained and unconstrained optimization problems based on a natural selection process that mimics biological evolution. The algorithm repeatedly modifies a population of individual solutions. At each step, the genetic algorithm randomly selects individuals from the current population and uses them as parents to produce the children for the next generation. Over successive generations, the population "evolves" toward an optimal solution.

We can apply the genetic algorithm to solve problems that are not well suited for standard optimization algorithms, including problems in which the objective function is discontinuous, nondifferentiable, stochastic, or highly nonlinear.

The genetic algorithm differs from a classical, derivative-based, optimization algorithm in two main ways, as summarized in the following table. ref[1]

Classical Algorithm	Genetic Algorithm
Generates a single point at each iteration. The sequence of points approaches an optimal solution.	Generates a population of points at each iteration. The best point in the population approaches an optimal solution.
Selects the next point in the sequence by a deterministic computation.	Selects the next population by computation which uses random number generators.

How the Genetic Algorithm Works

Outline of the Algorithm ref[2]

The following outline summarizes 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.
The algorithm stops when one of the stopping criteria is met. Stopping Conditions for the Algorithm are:
1. Generations— The algorithm stops when the number of generations reaches the value of Generations.
2. Time limit— The algorithm stops after running for an amount of time in seconds equal to Time limit.
3. Fitness limit— The algorithm stops when the value of the fitness function for the best point in the current population is less than or equal to Fitness limit.
4. Stall generations— The algorithm stops when the average relative change in the fitness function value over Stall generations is less than Function tolerance.
5. Stall time limit— The algorithm stops if there is no improvement in the objective function during an interval of time in seconds equal to Stall time limit.
6. Stall test— The stall condition is either average change or geometric weighted. For geometric weighted, the weighting function is 1/2n, 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.
7. Function tolerance— The algorithm runs until the average relative change in the fitness function value over Stall generations is less than Function tolerance.
8. Constraint tolerance— The Constraint tolerance is not used as stopping criterion. It is used to determine the feasibility with respect to nonlinear constraints.

Genetic Algorithm Terminology ref[3]

Fitness Functions

The fitness function is the function we want to optimize. For standard optimization algorithms, this is known as the objective function. The toolbox software tries to find the minimum of the fitness function.

Write the fitness function as a file or anonymous function, and pass it as a function handle input argument to the main genetic algorithm function.

Individuals

An individual is any point to which we can apply the fitness function. The value of the fitness function for an individual is its score.

An individual is sometimes referred to as a genome and the vector entries of an individual as genes.

Populations and Generations

A population is an array of individuals. For example, if the size of the population is 100 and the number of variables in the fitness function is 3, we represent the population by a 100-by-3 matrix. The same individual can appear more than once in the population. For example, the individual (2, -3, 1) can appear in more than one row of the array.

At each iteration, the genetic algorithm performs a series of computations on the current population to produce a new population. Each successive population is called a new generation.

Diversity

Diversity refers to the average distance between individuals in a population. A population has high diversity if the average distance is large; otherwise it has low diversity. In the following figure, the population on the left has high diversity, while the population on the right has low diversity.

Diversity is essential to the genetic algorithm because it enables the algorithm to search a larger region of the space.

Fitness Values and Best Fitness Values

The fitness value of an individual is the value of the fitness function for that individual. Because the toolbox software finds the minimum of the fitness function, the best fitness value for a population is the smallest fitness value for any individual in the population.

Parents and Children To create the next generation, the genetic algorithm selects certain individuals in the current population, called parents, and uses them to create individuals in the next generation, called children. Typically, the algorithm is more likely to select parents that have better fitness values.

Syntax of Genetic Algorithm in MATLAB

x = ga(fun,nvars)

x = ga(fun,nvars,A,b,[],[],lb,ub,nonlcon,IntCon,options)

[x,fval] = ga(___)

where,

fun = function call for which the optimisation is to be carried out

nvars = number of variables present in the fitness function

A,b = Coefficients of linear inequalities, (Set A=[] and b=[] if no linear inequalities exist.) ga evaluates the matrix product Aeq*x as if x is transposed (Aeq*x').

lb,ub = lower and upper bounds

nonlcon = Coefficients of nonlinear inequalities (Set nonlcon = [] if no nonlinear constraints exist.) Create options using optimoptions.

[x,fval] = returns values of x and function value at which the fitness function is minimised. ref[4]

For more details regarding Genetic algorithm, please refer the following reference. ref[5]

Matlab Program to plot and generate Stalagmite Function

%% PLOTTING A STALAGMITE FUNCTION
clear all
close all
clc
%Defining our Earth Space
x = linspace(0, 0.6, 150);
y = linspace(0, 0.6, 150);
 
%creating 2D Mesh
[xx yy] = meshgrid(x,y);
 
%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
 
%Ploting the function
surfc(xx, yy, f)
shading interp
xlabel('x value')
ylabel('y value')
zlabel('z value')
title('Model of Stalagmite Function')

Stalagmite Function

function y = stalagmite(input_value)
f1_x = (sin(5.1*pi*input_value(1) + 0.5))^6;
f1_y = (sin(5.1*pi*input_value(2) + 0.5))^6;
f2_x = exp(-4*log(2)*(input_value(1)-0.0667)^2/0.64);
f2_y = exp(-4*log(2)*(input_value(2)-0.0667)^2/0.64);
y = (f1_x*f2_x*f1_y*f2_y);
end

Code Explanation

We start with defining x and y values which are varied using linspace This is also referred as an xy plane or earth space
Then we make a grid for plotting purpose using meshgrid
Once the surface of xy plane is generated through the previous step, function values are plotted for all the values of x and y present in the earth space through a for loop.
The function stalagmite is defined in a separate file that returns ye. output of stalagmite function for all the corresponding values of x and y defined in first step.
Finally, the graph is plotted using surfc command along with chart labels and titles.

Code output

Calculating maximum value of Stalagmite function using Genetic Algorithm

Using the code below along with the Stalagmite function code stated above, we get the following results:

[inputs, fval] = ga(@stalagmite,2)

As we can see, the minimum value (fval) of the function comes out to be zero.

We can cross-check the same results by using optimisation tool just by typing optimtool in command window and feeding necessary outputs as shown in the figure:

By setting up the highlighted inputs and starting the process of optimisation by clicking the start button, the function converges 51th iteration in my case. Note that Genetic Algorithm being a random sampling technique, we may get different results every time we try to optimise the function.

Also, from the plot we can see that the function minimises to fitness value close to zero, almost equal to zero.

As now we know that Genetic algorithm gives only the minimum value, so in order to calculate the maximum value of the function, we can approach to the problem by inverting the stalagmite function and calculating its minimum values using the same approach. Following plots shows the maximum value before inverting the plot.

Instead of converting the plot we will just introduce a negative sign while calculating the optimum values of the fitness function.

The code along with its Stepwise explanation is given below.

Matlab code for finding maximum value of fitness function

%% COMPLETE GENETIC ALGORITHM APPROACH
clear all
close all
clc

x = linspace(0, 0.6, 150);
y = linspace(0, 0.6, 150);
num_cases = 50;

%Creating 2D Mesh
[xx yy] = meshgrid(x,y);

% To make sure that the code gives the same output,
% even if it is made to run several times.Here rng refers to random number
% generator
rng default

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

%% Study1: Statistical behaviour
tic
for i = 1:num_cases
    [inputs1, fopt1(i)] = ga(@stalagmite,2);
    xopt1(i) = inputs1(1);
    yopt1(i) = inputs1(2);
end
study_time1 = toc; 


figure(1)
subplot(2,1,1)
hold on
surfc(xx,yy,-f)
% axis([-0.6,0.6],[-0.6,0.6])
xlabel('x value')
ylabel('y value')
zlabel('function value')
shading interp
plot3(xopt1,yopt1,-fopt1,'marker','o','markersize',5,'markerfacecolor','r')
title('Unbounded Statistical Behaviour')
subplot(2,1,2)
plot(-fopt1)
xlabel('Iterations')
ylabel('Function Maximum')


%% Study 2 Statistical Behaviour with Lower and Upper Bounds

tic
for i = 1:num_cases
    [inputs2, fopt2(i)] = ga(@stalagmite, 2,[],[],[],[],[0,0],[1;1]);
    xopt2(i) = inputs2(1);
    yopt2(i) = inputs2(2);
end
study_time2 = toc;

figure(2)
subplot(2,1,1)
hold on
surfc(xx,yy,-f)

xlabel('x value')
ylabel('y value')
zlabel('function value')
shading interp
plot3(xopt2,yopt2,-fopt2,'marker','o','markersize',5,'markerfacecolor','r')
title('Bounded Statistical Behaviour')
subplot(2,1,2)
plot(-fopt2)
xlabel('Iterations')
ylabel('Function Maximum')


%% Study 3 Increasing the population size

options = optimoptions('ga')
options = optimoptions(options,'populationsize',200);
tic
for i = 1:num_cases;
    [inputs3,fopt3(i)] = ga(@stalagmite,2,[],[],[],[],[0,0],[1,1],[],[],options);
    xopt3(i) = inputs3(1);
    yopt3(i)= inputs3(2);
end
study_time3 = toc;

figure(3)
subplot(2,1,1)
hold on
surfc(xx,yy,-f)

xlabel('x value')
ylabel('y value')
zlabel('function value')
shading interp
plot3(xopt3,yopt3,-fopt3,'marker','o','markersize',5,'markerfacecolor','r')
title('Bounded Statistical Behaviour with Increased Population Size')
subplot(2,1,2)
plot(-fopt3)
xlabel('Iterations')
ylabel('Function Maximum')

%% SUMMARY OF COMPLETE STATISTICS
clc
total_time = study_time1+study_time2+study_time3; 
fprintf('Total elapsed Time is %f seconds. \n',total_time)
fprintf('Following Table represents the Overall Summary. \n')

T = table([study_time1;study_time2;study_time3],[mean(xopt1);mean(xopt2);mean(xopt3)],...
          [mean(yopt1);mean(yopt2);mean(yopt3)],[mean(-fopt1);mean(-fopt2);mean(-fopt3)],...
          'VariableNames',{'Run_Time','x_value','y_value','Max_Function_value'},...
          'RowNames',{'Trial 1';'Trial 2';'Trial 3'})

Code Explanation

Initiated with the clearing all the previous commands
Initially setting the iterations to 50 for less run time.
Defining the Earth space as did earlier above.
Setting random number generator to default in order to get same results for multiple times of run.
Extracting all the optimum values for stalagmite function by genetic algorithm in a for loop till set 50 iterations. Note that this time the stalagmite function has been changed a little in order to get inverted values. The function for stalagmite is attached below this explanation.
Then two plots are plotted in the same figure using subplot command where first plot represents all the optimum points on the surface plot by using hold on command. The other subplot represents the variation of the optimum values with respect to the iterations carried out. As we can see that the scatter in optimum values are over the entire plot and thus this represents that every time the algorithm returns local minima but not the global minima.
In order to reduce these variations in spread, we set up bound in study 2 and then again carry out the previous calculation and compare the results. We find that now local minima are close between the global minima. But still it is difficult to guess the close value. Also, now the maximum function values are close to 1 as can be seen in second sub plot.
So, in study 3, we try to increase the population size from initially 50 to 170 and again try to compare the results with previous two studies. We see that now the function converges close towards the global minima with very less variations in the data. Thus, these reflects our close convergence to the minimum value. (Note that the algorithm calculates the minimum value only bot for our case this minimum value is actually the maximum value of the function, as we have inverted our stalagmite function.)
Note that by varying the converging potions we compromise towards the run time, the more the convergence options, longer will the code take to show the results, so just try to make sure that necessary options are set with their first priority.
We can try to set more convergence options as explained above based upon the need.
In the end, the values are compared with their studies.

Stalagmite function used

function y = stalagmite(input_value)
f1_x = (sin(5.1*pi*input_value(1) + 0.5))^6;
f1_y = (sin(5.1*pi*input_value(2) + 0.5))^6;
f2_x = exp(-4*log(2)*(input_value(1)-0.0667)^2/0.64);
f2_y = exp(-4*log(2)*(input_value(2)-0.0667)^2/0.64);
y = -(f1_x*f2_x*f1_y*f2_y);
end

Output Results

Results of Study 1

Results of Study 2

Results of Study 3

Rectifying the above results using Optimisation toolbox

Errors Encountered

Improper definition of the surface plot (function value of stalagmite function)

Later the function code was changed to get proper dimensioned matrix (here 150 by 150 vector)

Conclusion

From the above topic, I conclude the following takeaways:

Genetic Algorithm is very important in analysing complex optimisation problems.
Based upon the controls limits and many other factors, the algorithm run time depends upon user inputs.
Maximization problems of some cases can also be found out using genetic algorithm optimisation approach very easily.

References

Genetic Algorithm using MATLAB

Ship work like this
in your own course.

Talk to a career counsellor

Ship work like thisin your own course.

Talk to a career counsellor

Ship work like this
in your own course.