Genetic Algorithm: a highly robust inversion scheme for geophysical applications

An introduction to the basics of genetic algorithm along with a simple numerical example and solution of an earthquake location problem
Kumar, U. (2020, June 7). Genetic Algorithm explained in the context of earthquake location problem. https://doi.org/10.5281/zenodo.4679544

What is Optimization?

  • “The selection of a best element (with regard to some criteria) from some set of available alternatives” {Source: Wikipedia}
  • Picking the “best” option from several ways of accomplishing the same task.
  • We require a model that can give us the best result.
Optimum Route

Global vs. Local Optimization

Consider a hypothetical multi-modal function:

Hypothetical objective functions with global and local minima.

Genetic Algorithm

Genetic algorithm (GA), first proposed by John Holland in 1975 and described in Goldberg (1988), is based on analogies with biological evolution processes. The principle of GA is quite simple and mirrors what is perceived to occur in evolutionary mechanisms. The idea is to start with a given set of feasible trial solutions (either constrained or unconstrained) and iteratively evaluate objective (or fitness) function by keeping only those solutions that give the optimal value of the objective function and randomly mutate them to try and do even better and successive iterations. Thus, the beneficial mutations (that lead to optimal values) are kept while those that perform poorly are discarded, i.e., the survival of the fittest.

Consider the unconstrained optimization problem with the objective function, min f(\vec{x}) where \vec{x} is an n dimensional vector. Suppose m initial guesses are given for the values of \vec{x} so that

guess \quad j \quad is \quad x_j

Thus, m solutions are evaluated and compared with each other to see which of the solutions generate the smallest objective function since our goal (in general) is to minimize it. We can order the guesses such that the first (p < m) give the smallest values of f (\vec{x}). Arranging our data, we have:

keep \quad x_j \quad j = 1,2,...,p

discard \quad x_j = p+1,p+2,...,m

Since the first p solutions are the best, these are kept in the next generation. In addition, we now generate m-p new trial solutions that are randomly mutated from the p best solutions. This process is repeated through a finite number of iterations with the hope that convergence to the optimal solution is achieved.

Scheme

The steps involved in the Genetic Algorithm are:

  1. The initial population of model parameters is randomly chosen within the given search range.
  2. The simple GA undergoes a set of operations on the model population to produce the next generation: selection, coding, crossover, and mutation.
  3. In the selection scheme, the model parameters exhibiting the higher fitness value (lower objective function value relative to others) are selected and replicated with the given probability. The total population size remains constant; then, the population members are randomly paired among themselves.
  4. In the coding scheme, each population’s decimal values are converted to a binary system, forming a long bit string (analogous to a chromosome).
  5. In the crossover scheme, some part of the long bit string of binary model parameters is exchanged with their corresponding pair to produce a new population.
  6. In the mutation scheme, some randomly selected sites (with a given probability) of the new set of the binary model population are switched.
  7. These sets of operations will continue until some pre-defined termination criteria for the technique are satisfied.

The Optimum Size of a Coffee Mug for a Café Owner?

Here, the objective function is to minimize the loss of the Café owner. Let us make a very simple function for that.

y = 0.2 x^2 +50/x

where x is the size of the coffee mug. If the coffee mug’s size is too small, then the customers may not like it, and if the size is too large, then the owner will have to endure loss, hence the hypothetical function. We assumed that the price of the coffee is fixed.

The plot of the above function. It has a minimum of 15.0 units.

Solution

Let us solve this iteratively:

Generation 1

The average fitness in this generation is: 52.55.

Generation 2

The average fitness in this generation is 29.17. We can continue iterating as long as our termination criteria are satisfied. Some of the standard termination criteria are the tolerance value (the difference in the fitness of the generation compared to the previous generation), the maximum number of generations, etc.

Fitness with Generations

Fitness optimization with generations (lower is better in this case)

Earthquake Location Using Genetic Algorithm

In this series of earthquake location problems, we have used the Monte Carlo method for the earthquake location that gives us reliable results. The earthquake location problem was introduced in the previous post.

Now, I applied GA from the direct search toolbox of MATLAB to find the solution of the best earthquake location parameters for the earthquake problem defined in previous post for 30 stations. The GA does not improve the solutions but similar to the Monte Carlo method; it allows a large model space to be searched for solutions. I used the lower and upper limits for the search of best parameters as defined in Table 1. After 85 generations of a run and terminated based on the tolerance criterion, the best parameters found by the GA are [2.9750 2.5202 -2.5535 6.6696 -0.0839]. Since GA is a stochastic optimization scheme, it is best practice to make a reasonable number of acceptable solutions and consider its mean as the final solution. For the 50 runs of the GA for the earthquake location problem, the estimated mean and std of the model parameters are shown in Table 2. The difference from the actual model parameters arises because of random noise in the observed arrival time data.

Table 1: Lower and upper limits used for the earthquake location problem using genetic algorithm
The estimation of hypothetical earthquake location solution using Genetic Algorithm
(a) Estimated location of a hypothetical earthquake using Genetic algorithm. The yellow ellipsoid shows the 1𝜎 confidence interval in the hypocenter location. (b) Fitness value versus generations. The best score value (0.106) and mean score (0.106) versus generation for estimating earthquake location parameters using genetic algorithm scheme.
Table 2: Estimated model parameters for earthquake location using genetic algorithm after 50 runs
lower=[-3 -3 -3 5 -1];
upper=[3 3 0 7 1];
options = optimoptions('ga','PlotFcn', @gaplotbestf,'Display','iter'); 6. x=ga(@(pp)fit_arrival_times(pp),5,[],[],[],[],lower,upper,[],options)

function E=fit_arrival_times(pp)
filename = 'arrival_times.csv';
[dobs, dpre] = read_arrivaltimes(filename);

stationlocations = read_stnloc('station_locations.csv', 2, 31);
d_pre = zeros(length(stationlocations),1);

for i=1:length(stationlocations)
    dist = sqrt((pp(1) - stationlocations(i,1))^2 + (pp(2) - stationlocations( i,2))^2 + (pp(3) - stationlocations(i,3))^2);
    arr = dist/pp(4) + pp(5);
    d_pre(i) = arr;
end
E=sum((dobs-d_pre).^2);

Complete Codes

Complete codes can be downloaded from my GitHub repository.

Utpal Kumar
Utpal Kumar

Geophysicist | Geodesist | Seismologist | Open-source Developer
I am a geophysicist with a background in computational geophysics, currently working as a postdoctoral researcher at UC Berkeley. My research focuses on seismic data analysis, structural health monitoring, and understanding deep Earth structures. I have had the opportunity to work on diverse projects, from investigating building characteristics using smartphone data to developing 3D models of the Earth's mantle beneath the Yellowstone hotspot.

In addition to my research, I have experience in cloud computing, high-performance computing, and single-board computers, which I have applied in various projects. This includes working with platforms like AWS, Docker, and Kubernetes, as well as supercomputing environments such as STAMPEDE2, ANVIL, Savio and PERLMUTTER (and CORI). My work involves developing innovative solutions for structural health monitoring and advancing real-time seismic response analysis. I am committed to applying these skills to further research in computational seismology and structural health monitoring.

Articles: 32

Leave a Reply

Your email address will not be published. Required fields are marked *