Optimizing the stalagmite function and finding its global maxima

Genetic Algorithm

Genetic Algorithm (GA) is a search-based optimization technique based on the principles of Genetics and Natural Selection. They are commonly used to generate high-quality solutions for optimization problems and search problems. Genetic algorithms simulate the process of natural selection which means those species who can adapt to changes in their environment are able to survive and reproduce and go to next generation. This is in line with the Darwinian Theory of “Survival of the Fittest”.

The initial population of possible solutions to the problem or the initial generation is randomly selected. Once the initial generation is created, the algorithm evolves the generation and three types of children are formed for the next generation through three operators

1. Selection operator: Where individuals(or solutions) with best fitness values are passed onto the next generation.
2. Crossover operator: New children are produced from the crossover of two selected parents.
3. Mutation operator: A selected individula is mutated and passed to th enext generation

Objective function(fitness fucntion) is the function which is being optimized by the genetic algorithm and a fitness value is assigned to each individual(or solution) based on its obejective function value.

Program

The aim of this program is to optimize the stalagmite function and find its global maxima.

The stalagmite function is given below.

The syntax for the genetic algorithm command is

[x,fval]=ga(fun,nvars,A,b,Aeq,beq,lb,ub,nonclon,options) where

1. x: the local minimum.
2. fval: value of fitness function at x.
3. fun: the objective function.
4. nvars: number of design variables of fun.
5. A and b: when subject to the linear inequalities A*x<=b.
6. Aeq and beq: when subject to the linear inequalities Aeq*x=beq along with A*x<=b.
7. lb and ub: defines a set of lower and upper bounds on the design variables.
8. nonclon: for non linear constraints
9. options: the default optimization parameters replaced by values in options.

My program starts with the creation of two arrays using \"linspace\" command. A 2D mesh is created using these two arrays.

%defining search space
x=linspace(0,0.6,150);
y=linspace(0,0.6,150);

%creating 2D mesh
[xx,yy]=meshgrid(x,y);
num_cases=50;

The stalgamite function is then evaluated by passing \"input\" variable as its input. The \"x\" and \"y\" used in the stalagmite function are defined by \"input(1)\" and \"input(2)\".

%evaluating stalagmite function
for i=1:length(xx)
    for j=1:length(yy)
        input(1)=xx(i,j);
        input(2)=yy(i,j);
        f(i,j)=geneticalgo_stalagmite(input);
    end
end

The stalagmite function is defined as given above excpet with a slight modification. The final step, where f(x,y) is defined, is multiplied with -1. This is done because the \"ga\" function inherently minimizes the function and gives the global or local minima value. Since the stalagmite function when plotted produces stalagmites, multiplying it by -1 would give an inverted stalagmite for which the \"ga\" function can calculate the global minima which when multiplied by -1 again before plotting would give us the global or local maxima value and corresponding graph.

% modified stalagmite function
function result= geneticalgo_stalagmite(input)
x=input(1);
y=input(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);
result=(-1)*f1x*f1y*f2x*f2y; %multiplied by -1 to get inverted stalagmite
end

Now the statistical behaviour of GA is studied for three cases.

Case1: Unbounded inputs

The genetic algorithm is run 50 times and the minima values and their corresponding input values are stored in \"fopt\" and \"coord\" respectively. A 2x1 subplot is created and in the top half the stalagmite function is plotted(\"f\" is multiplied by -1 to get the original or erect stalagmite function) along with a 3D plot of the different input values used and the corresponding maxima values. In the bottom half the maxima vales are plotted against the number of iterations. Here \"fopt\" is multiplied by -1 to get the maxima values as the \"ga\" function gives only minima values.

%case 1: Statistical behaviour of GA
for i=1:num_cases
    [coord,fopt(i)]=ga(@geneticalgo_stalagmite,2);
    xopt(i)=coord(1);
    yopt(i)=coord(2);
end
figure(1)
subplot(2,1,1)
surfc(x,y,-f);% multipled by -1 to get erect stalagmite plot
xlabel(\"x value\");
ylabel(\"y value\");
hold on
shading interp;
plot3(xopt,yopt,-fopt,\"marker\",\"o\",\"markersize\",5); %multiplied by -1 to get maxima values
title(\"Unbounded inputs\");
hold off
subplot(2,1,2)
plot(-fopt);
xlabel(\"Iterations\");
ylabel(\"Function maximum\");

It can be seen from the plots below that the maxima values are highly varying and this is due to the random selection of input variables of the \"ga\" function.

Case2: Bounded Inputs

In the \"ga\" function the design variables or inputs are bounded by an upper and lower bound. Everything else is the same as case 1.

%case2: Statistical behaviour of GA with upper and lower bounds
for i=1:num_cases
    [coord,fopt(i)]=ga(@geneticalgo_stalagmite,2,[],[],[],[],[0;0],[1;1]); %upper and lower bound for the input variables
    xopt(i)=coord(1);
    yopt(i)=coord(2);
end
figure(2)
subplot(2,1,1)
surfc(x,y,-f);% multipled by -1 to get erect stalagmite plot
xlabel(\"x value\");
ylabel(\"y value\");
hold on
shading interp;
plot3(xopt,yopt,-fopt,\"marker\",\"o\",\"markersize\",5);%multiplied by -1 to get maxima values
title(\"Bounded inuts\");
hold off
subplot(2,1,2)
plot(-fopt);
xlabel(\"Iterations\");
ylabel(\"Function maximum\");

The variation in maxima values is much less as compared to case 1 as here the input variables have been bounded. The algorithm seems to be converging to a global maxima.

Case 3: Bound inputs along with increase population size

As we know the genetic algorithm requires an initial population which is randomly selected. By default this value for the \"ga\" function in MATLAB is 50. I have increased the population size to 250 so as to get the global minima instead of local minimas. This can be done using \"options \" in the \"ga\" function and options can be created using \"optimoptions\".

The \"rng\" command is used at the end to generate the same random values generated by MATLAB for the initial population. This way the program give sthe same output every time.

%case 3: Statistical behaviour of GA with upper and lower bounds and increased population size
options=optimoptions(\"ga\"); %gets default options for ga function
options=optimoptions(options,\"populationsize\",250); %returns options with altered values. Population size incresed to 250

for i=1:num_cases
    [coord,fopt(i)]=ga(@geneticalgo_stalagmite,2,[],[],[],[],[0;0],[1;1],[],[],options);
    xopt(i)=coord(1);
    yopt(i)=coord(2);
end
figure(3)
subplot(2,1,1)
surfc(x,y,-f);% multipled by -1 to get erect stalagmite plot
xlabel(\"x value\");
ylabel(\"y value\");
hold on
shading interp;
plot3(xopt,yopt,-fopt,\"marker\",\"o\",\"markersize\",5);%multiplied by -1 to get maxima values
title(\"Bounded inputs with modified population size\");
hold off
subplot(2,1,2)
plot(-fopt);
xlabel(\"Iterations\");
ylabel(\"Function maximum\");
rng(\"default\");

There is almost no variation in the maxima values generated by the \"ga\" function and it has converged on a global maxima value of 1. This is due to the bounding of the input variables along with the increased initial population size.

Code

clear all
close all
clc
%defining search space
x=linspace(0,0.6,150);
y=linspace(0,0.6,150);

%creating 2D mesh
[xx,yy]=meshgrid(x,y);
num_cases=50;

%evaluating stalagmite function
for i=1:length(xx)
    for j=1:length(yy)
        input(1)=xx(i,j);
        input(2)=yy(i,j);
        f(i,j)=geneticalgo_stalagmite(input);
    end
end

%case 1: Statistical behaviour of GA
for i=1:num_cases
    [coord,fopt(i)]=ga(@geneticalgo_stalagmite,2);
    xopt(i)=coord(1);
    yopt(i)=coord(2);
end
figure(1)
subplot(2,1,1)
surfc(x,y,-f);% multipled by -1 to get erect stalagmite plot
xlabel(\"x value\");
ylabel(\"y value\");
hold on
shading interp;
plot3(xopt,yopt,-fopt,\"marker\",\"o\",\"markersize\",5); %multiplied by -1 to get maxima values
title(\"Unbounded inputs\");
hold off
subplot(2,1,2)
plot(-fopt);
xlabel(\"Iterations\");
ylabel(\"Function maximum\");

%case2: Statistical behaviour of GA with upper and lower bounds
for i=1:num_cases
    [coord,fopt(i)]=ga(@geneticalgo_stalagmite,2,[],[],[],[],[0;0],[1;1]); %upper and lower bound for the input variables
    xopt(i)=coord(1);
    yopt(i)=coord(2);
end
figure(2)
subplot(2,1,1)
surfc(x,y,-f);% multipled by -1 to get erect stalagmite plot
xlabel(\"x value\");
ylabel(\"y value\");
hold on
shading interp;
plot3(xopt,yopt,-fopt,\"marker\",\"o\",\"markersize\",5);%multiplied by -1 to get maxima values
title(\"Bounded inuts\");
hold off
subplot(2,1,2)
plot(-fopt);
xlabel(\"Iterations\");
ylabel(\"Function maximum\");

%case 3: Statistical behaviour of GA with upper and lower bounds and increased population size
options=optimoptions(\"ga\"); %gets default options for ga function
options=optimoptions(options,\"populationsize\",250); %returns options with altered values. Population size incresed to 250

for i=1:num_cases
    [coord,fopt(i)]=ga(@geneticalgo_stalagmite,2,[],[],[],[],[0;0],[1;1],[],[],options);
    xopt(i)=coord(1);
    yopt(i)=coord(2);
end
figure(3)
subplot(2,1,1)
surfc(x,y,-f);% multipled by -1 to get erect stalagmite plot
xlabel(\"x value\");
ylabel(\"y value\");
hold on
shading interp;
plot3(xopt,yopt,-fopt,\"marker\",\"o\",\"markersize\",5);%multiplied by -1 to get maxima values
title(\"Bounded inputs with modified population size\");
hold off
subplot(2,1,2)
plot(-fopt);
xlabel(\"Iterations\");
ylabel(\"Function maximum\");
rng(\"default\");

% modified stalagmite function
function result= geneticalgo_stalagmite(input)
x=input(1);
y=input(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);
result=(-1)*f1x*f1y*f2x*f2y; %multiplied by -1 to get inverted stalagmite
end

References

1. https://in.mathworks.com/help/gads/how-the-genetic-algorithm-works.html
2. https://www.geeksforgeeks.org/genetic-algorithms/
3. https://www.tutorialspoint.com/genetic_algorithms/genetic_algorithms_introduction.htm
4. https://in.mathworks.com/help/gads/ga.html#d120e39545