Genetic Programming Book PDF
Genetic Programming Book PDF
Alavi
Conor Ryan Editors
Handbook
of Genetic
Programming
Applications
Handbook of Genetic Programming Applications
Amir H. Gandomi • Amir H. Alavi • Conor Ryan
Editors
Handbook of Genetic
Programming Applications
123
Editors
Amir H. Gandomi Amir H. Alavi
BEACON Center for the Study Department of Civil and Environmental
of Evolution in Action Engineering
Michigan State University Michigan State University
East Lansing, MI, USA East Lansing, MI, USA
Conor Ryan
Department of Computer Science
and Information Systems
University of Limerick
Limerick, Ireland
vii
Contents
ix
x Contents
Part IV Tools
22 GPTIPS 2: An Open-Source Software Platform
for Symbolic Data Mining . . . . . . . . . . . . . . . . . . . . . . . . . . . . .. . . . . . . . . . . . . . . . . . . . 551
Dominic P. Searson
23 eCrash : a Genetic Programming-Based Testing Tool
for Object-Oriented Software . . . . . . . . . . . . . . . . . . . . . . . . .. . . . . . . . . . . . . . . . . . . . 575
José Carlos Bregieiro Ribeiro, Ana Filipa Nogueira, Francisco
Fernández de Vega, and Mário Alberto Zenha-Rela
Part I
Overview of Genetic Programming
Applications
Chapter 1
Graph-Based Evolutionary Art
1.1 Introduction
Electronic supplementary material The online version of this chapter (doi: 10.1007/978-3-319-
20883-1_1) contains supplementary material, which is available to authorized users.
P. Machado () • J. Correia • F. Assunção
CISUC, Department of Informatics Engineering, University of Coimbra,
3030 Coimbra, Portugal
e-mail: [email protected]; [email protected]; [email protected]
Although there are noteworthy expression-based evolutionary art systems (e.g. Sims
(1991); World (1996); Unemi (1999); Machado and Cardoso (2002); Hart (2007)),
systems that allow the evolution of images that are composed of a set of distinct
graphic elements such as lines, shapes, colors and textures are extremely scarce.
Among the exceptions to the norm, we can cite the work of: Baker and Seltzer
(1994), who uses a Genetic Algorithm (GA) operating on strings of variable size
to evolve line drawings; den Heijer and Eiben (2011) who evolve Scalable Vector
Graphics (SVG), manipulating directly SVG files through a set of specifically
designed mutation and recombination operators. Unlike GP approaches, where
the representation is procedural, the representations adopted in these works are,
essentially, descriptive—in the sense that the genotypes describe the elements of
the images in a relatively directed way instead of describing a procedure, i.e.
program, that once executed or interpreted produces the image as output.
In addition to our early work on this topic (Machado et al. 2010; Machado
and Nunes 2010), there are two examples of the use of CFDG for evolutionary
1 Graph-Based Evolutionary Art 5
1.3 Representation
Fig. 1.2 Examples of images produced by the CFDG depicted in Fig. 1.1
In this Section we describe the genetic operators designed to manipulate the graph-
based representation of CFDGs, namely: initialization, mutation and crossover.
The creation of the initial population for the current evolutionary engine is of huge
importance, being responsible for generating the first genetic material that will be
evolved through time. In our previous works on the evolution of CFDGs the initial
population was supplied to the evolutionary engine: the first population was either
composed of human-created grammars (Machado and Nunes 2010) or of a single
minimal grammar (Machado et al. 2010). Although both those options have merit,
the lack of an initialization procedure for the creation of a random population of
CFDGs was a limitation of the approach.
In simple terms, the procedure for creating a random CFDG can be described as
follows: we begin by randomly determining the number of non-terminal symbols
and the number of production rules for each of the symbols (i.e. the number of
different options for its expansion). Since this defines the nodes of the graph,
the next step is the random creation of connections among nodes and calls to
non-terminal symbols. The parameters associated with the calls to terminal and non-
terminal symbols are also established randomly. Finally, once all productions have
been created, we randomly select a starting node and background color. Algorithm 1
details this process, which is repeated until the desired number of individuals is
reached. Figure 1.3 depicts a sample of a random initial population created using
this method.
The crossover operator used for the experiments described in this Chapter is similar
to the one used in our previous work on the same topic (Machado et al. 2010;
Machado and Nunes 2010). The rational was to develop a crossover operator that
would promote the meaningful exchange of genetic material between individuals.
Given the nature of the representation, this implied the development of a graph-
based crossover operator that is aware of the structure of the graphs being
manipulated. The proposed operator can be seen as an extension of the one presented
by Pereira et al. (1999). In simple terms, this operator allows the exchange of
subgraphs between individuals.
The crossover of the genetic code of two individuals, a and b, implies: (1)
selecting one subgraph from each parent; (2) swapping the nodes and internal edges
8 P. Machado et al.
of the subgraphs, i.e., edges that connect two subgraph nodes; (3) establishing a
correspondence between nodes; (4) restoring the outgoing and incoming edges, i.e.,
respectively, edges from nodes of the subgraph to non-subgraph nodes and edges
from non-subgraph nodes to nodes of the subgraph.
Subgraph selection—randomly selects for each parent, a and b, one crossover
node, va and vb , and a subgraph radius, ra and rb . Subgraph sra is composed
of all the nodes, and edges among them, that can be reached in a maximum of ra
steps starting from node va . Subgraph srb is defined analogously. Two methods
were tested for choosing va and vb , one assuring that both va and vb are in the
connected part of the graph and one without restrictions. The radius ra and rb
were randomly chose being the maximum allowed value the maximum depth of
the graph.
Swapping the subgraphs—swapping sra and srb consists in replacing sra by srb
(and vice-versa). After this operation the outgoing and the incoming edges are
destroyed. Establishing a correspondence between nodes repairs these connec-
tions.
Correspondence of Nodes—let sraC1 and srbC1 be the subgraphs that would be
obtained by considering a subgraph radius of ra C 1 and rb C 1 while performing
the subgraph selection. Let msta and mstb be the minimum spanning trees (MSTs)
with root nodes va and vb connecting all sraC1 and srbC1 nodes, respectively.
For determining the MSTs all edges are considered to have unitary cost. When
several MSTs exist, the first one found is the one considered. The correspondence
between the nodes of sraC1 and srbC1 is established by transversing msta and mstb ,
starting from their roots, as described in Algorithm 2.
Restoring outgoing and incoming edges—the edges from a … sra to sra are
replaced by edges from a … srb to srb using the correspondence between the nodes
established in the previous step (e.g. the incoming edges to va are redirected to
vb , and so on). Considering a radius of ra C 1 and rb C 1 instead of ra and rb in
the previous step allows the restoration of the outgoing edges. By definition, all
outgoing edges from sa and sb link to nodes that are at a minimum distance of
ra C 1 and rb C 1, respectively. This allows us to redirect the edges from sb to b
… sb to a … sa using the correspondence list.
The mutation operators were designed to attend two basic goals: allowing the
introduction of new genetic material in the population and ensuring that the search
space is fully connected, i.e., that all of its points are reachable from any starting
point through the successive application of mutation operators. This resulted in the
use of a total of ten operators, which are succinctly described on the following
paragraphs.
10 P. Machado et al.
Fitness assignment implies interpreting and rendering the CFDG. This is accom-
plished by calling the Context Free (Horigan and Lentczner 2009) application.
Grammars with infinite recursive loops are quite common. As such, it was necessary
to establish an upper bound to the number of steps that a CFDG is allowed to make
before its expansion is considered complete. The original version of Context Free
only allows the definition of an upper bound for the number of drawn shapes. This
is insufficient for our goals, because it allows endless loops, provided that no shapes
are drawn. As such, it was necessary to introduce several changes to the source
code of Context Free (which is open source) to accommodate our needs. When
calling Context Free we give as input (1) the CFDG to be interpreted and rendered,
1 Graph-Based Evolutionary Art 11
(2) the rendering size, (3) the maximum number of steps (4) the rendering seed. We
receive as output an image file. The maximum number of steps was set to 100;000
for all the experiments described in this Chapter. The “rendering seed” defines the
seed of the random number generator used by Context Free during the expansion
of the CFDGs. The rendering of the same CFDG using different rendering seeds
can, and often does, result in different images (see Sect. 1.3). We performed tests
using fixed and randomly generated rendering seeds. The results of those tests will
be described in Sect. 1.6.
We use six different hardwired fitness functions based on evolutionary art
literature and conduct tests using each of these functions to guide evolution. In a
second stage, we perform runs using a combination of these measures to assign
fitness. In the reminder of this Section we describe each of the functions and the
procedure used to combine them.
The image returned by Context Free is encoded in JPEG format using the maximum
quality settings. The size of the JPEG file becomes the fitness of the individual.
The rationale is that complex images, with abrupt transitions of color are harder to
compress and hence result in larger file sizes, whereas simple images will result in
small file sizes (Machado and Cardoso 2002; Machado et al. 2007). Although this
assignment scheme is rather simplistic, it has the virtue of being straightforward to
implement and yield results that are easily interpretable. As such, it was used to
assess the ability of the evolutionary engine to complexify and to establish adequate
experimental settings.
1
fitness D (1.1)
1 C jtargetvalue observedvalue j
We use the target values of 1:3 and 0:90 for fractal dimension and lacunarity,
respectively. These values were established empirically by calculating the fractal
dimension and lacunarity of images that we find to have desirable aesthetic qualities.
Fig. 1.4 Example of the transformation from the input color image (left image) to the back-
ground/foreground image (right image) used for the Fractal Dimension and Lacunarity estimates
1 Graph-Based Evolutionary Art 13
1.5.4 Complexity
This fitness function, based on the work of Machado and Cardoso (2002); Machado
et al. (2007, 2005), assesses several characteristics of the image related with
complexity. In simple terms, the rationale is valuing images that constitute a
complex visual stimulus but that are, nevertheless, easy to process. A thorough
discussion of the virtues and limitations of this approach is beyond the scope of
this Chapter, as such, we focus on practical issues pertaining its implementation.
The approach relies on the notion of compression complexity, which is defined as
calculated using the following formula:
s.scheme.i//
C.i; scheme/ D RMSE.i; scheme.i// (1.2)
s.i/
where i is the image being analysed, scheme is a lossy image compression scheme,
RMSE stands for the root mean square error, and s is the file size function.
To estimate the complexity of the visual stimulus (IC.i/) they calculate the
complexity of the JPEG encoding of the image (i.e. IC.i/ D C.i; JPEG/). The
processing complexity (PC.i/) is estimated using a fractal (quadratic tree based)
encoding of the image (Fisher 1995). Considering that as time passes the level
of detail in the perception of the image increases, the processing complexity is
estimated for different moments in time (PC.t0; i/, PC.t1; i/) by using fractal image
compression with different levels of detail. In addition to valuing images with
high visual complexity and low processing complexity, the approach also values
images where PC is stable for different levels of detail. In other words, according
to this approach, an increase in description length should be accompanied by an
increase in image fidelity. Taking all of these factors into consideration, Machado
and Cardoso (2002); Machado et al. (2007, 2005) propose the following formula for
fitness assignment:
IC.i/a
(1.3)
.PC.t0; i/ PC.t1; i//b . PC.t1;i/PC.t0;i/
PC.t1;i/ /c
where ˛, and ı operate as target values for IC.i/, .PC.t0; i/ PC.t1; i/ and
PC.t1; i/ PC.t0; i/, which were set to 6, 24 and 1:1, respectively. These values
were determined empirically through the analysis of images that we find to be
14 P. Machado et al.
desirable. Due to the limitations of the adopted fractal image compression scheme
this approach only deals with greyscale images. Therefore, all images are converted
to greyscale before being processed.
1.5.5 Bell
This fitness function is based on the work of Ross et al. (2006) and relies on the
observation that many fine-art works exhibit a normal distribution of color gradients.
Following Ross et al. (2006) the gradients of each color channel are calculated, one
by one, in the following manner:
where ri;j is the image pixel intensity values for position .i; j/ and d is a scaling
factor that allows to compare images of different size; this value was set to 0:1 % of
half the diagonal of the input image (based on Ross et al. (2006)). Then the overall
gradient Si;j is computed as follows:
q
Si;j D jrri;j j2 C jrgi;j j2 C jrbi;j j2 (1.6)
Si;j
Ri;j D log (1.7)
S0
P
2 i;j Ri;j .Ri;j /2
D P (1.9)
i;j Ri;j
At this step we introduce a subtle but important change to (Ross et al. 2006)
work: we consider a lower bound for the 2 , which was empirically set to 0.7. This
prevents the evolutionary engine to converge to monochromatic images that, due
to the use of a small number of colors, trivially match a normal distribution. This
change has a profound impact in the experimental results, promoting the evolution
of colorful images that match a normal distribution of gradients.
1 Graph-Based Evolutionary Art 15
Using , 2 and the values of Ri;j a frequency histogram with a bin size of
=100 is created, which allows calculating the deviation from normality (DFN). The
DFN is computed using qi , which is the observed probability and pi , the expected
probability considering a normal distribution. Ross et al. (2006) uses:
X pi
DFN D 1000 pi log (1.10)
qi
Which measures the squares of the differences between expected and observed
probabilities. Therefore, in the experiments described in this Chapter Bell fitness is
assigned according to the following formula: 1=.1 C DFNs /.
In addition to the tests where the fitness functions described above were used
to guide evolution, we conducted several experiments where the goal was to
simultaneously maximize several of these functions. This implied producing a
fitness score from multiple functions, which was accomplished using the following
formula:
Y
combinedfitness .i/ D log .1 C fj .i// (1.12)
j
where i is the image being assessed and fj refers to the functions being considered.
Thus, to assign fitness based on the Complexity and Bell functions we compute:
log.1 C Complexity.i// log.1 C Bell.i//. By adopting logarithmic scaling and
a multiplicative fitness function we wish to promote the discovery of images that
maximize all the measures being considered in the experiment.
The evolutionary engine has several novel characteristics that differentiate it from
conventional GP approaches. Therefore, it was necessary to conduct a series of tests
to assess the adequacy of the engine for the evolution of CFDGs and to determine a
reasonable set of configuration parameters. These tests were conducted using JPEG
Size as fitness function and allowed us to establish the experimental parameters
16 P. Machado et al.
Table 1.1 Parameters used for the experiments described in this chapter
Initialization (see Algorithm 1) Values
min, max number of symbols (1,3)
min, max number of rules (1,3)
min, max calls per production rule (1,2)
Evolutionary Engine Values
Number of runs 30
Number of generations 100
Population size 100
Crossover probability 0.6
Mutation probability 0.1
Tournament size 10
Elite size Top 2 % of the population
CFDG Parameters Values
Maximum number expansion steps 100,000
Limits of the geometric transformations rotate 2 [0,359], size 2 [-5,5]
x 2 [-5,5], y 2 [-5,5], z 2 [-5,5]
flip 2 [-5,5], skew 2 [-5,5]
Limits of the color transformations hue 2 [0,359], saturation 2 [-1,1]
brightness 2 [-1,1], alpha 2 [-1,1]
Terminal symbols SQUARE, CIRCLE, TRIANGLE
summarized in Table 1.1, which are used throughout all the experiments described
herein. In general, the results show that the engine is not overly sensitive to the
configuration parameters, depicting an adequate behavior for a wide set of parameter
configurations. Although the optimal parameters settings are likely to depend on the
fitness function, a detailed parametric study is beyond the scope of this Chapter.
Therefore, we did not attempt to find an optimal combination of parameters.
The use of a graph-based representation and genetic operators is one of the
novel aspects of our approach. The use of such operators may introduce changes
to the graph that may make some of the nodes (i.e. some production variables)
unreachable from the starting node. For instance, a mutation of the node Edera of
Fig. 1.1 may remove the call to node Ciglio making most of the graph unreachable.
Although, unreachable nodes have no impact on the phenotype, their existence may
influence the evolutionary process. On one hand they may provide space for neutral
variations and promote evolvability (unreachable nodes may become reattached by
subsequent genetic operators), on the other they may induce bloat since they allow
protection from destructive crossover. To study the impact of unreachable nodes in
the evolutionary process we considered three variations of the algorithm:
Unrestricted—the crossover points are chosen randomly;
Restricted—the crossover points are chosen randomly from the list of reachable
nodes of each parent;
1 Graph-Based Evolutionary Art 17
25K
20K
← FITNESS →
15K
10K
5K
0
0 10 20 30 40 50 60 70 80 90 100
← GENERATION →
Fig. 1.5 Best and average fitness values for different implementations of the genetic operators
using JPEG Size as fitness function. The results are averages of 30 independent runs
125
75
50
25
0
0 10 20 30 40 50 60 70 80 90 100
← GENERATION →
Fig. 1.6 Evolution of the average number of reachable and unreachable nodes across populations
for different implementations of the genetic operators using JPEG Size as fitness function. The
results are averages of 30 independent runs
passed to Context Free. We considered two scenarios: using a fixed rendering seed
for all individuals; randomly generating the rendering seed whenever genotype to
phenotype occurs. The second option implies that the fitness of a genotype may, and
often does, vary from one evolution to the other, since the phenotype may change.
Figure 1.7 summarizes the results of these tests in terms of the evolution of fitness
through time. As expected, using a fixed rendering seed yields better fitness, but
the differences between the approaches are surprisingly small and decrease as the
number of generations increases. To better understand this result we focused on the
analysis of the characteristics of the CFGDs being evolved. Figure 1.8 depicts box
plots of fitness values of the fittest individuals of each of the 30 evolutionary runs
using different setups:
Fixed—individuals evolved and evaluated using fixed rendering seeds; Random—
individuals evolved using random rendering seeds and evaluated using the same
seeds as the ones picked randomly during evolution;
Fixed Random—individuals evolved using fixed rendering seeds and evaluated with
30 random seeds each;
Random Random—individuals evolved using random rendering seeds and evaluated
with 30 random seeds each.
In other words, we take the genotypes evolved in a controlled static environment
(fixed random seed) and place them in different environments, proceeding in
the same way for the ones evolved in a changing environment. The analysis of
the box plots shows that, in the considered experimental settings, the fitness of
1 Graph-Based Evolutionary Art 19
25K
20K
← FITNESS →
15K
10K
5K
0
0 10 20 30 40 50 60 70 80 90 100
← GENERATION →
Avg Fixed Max Fixed Avg Random Max Random
Fig. 1.7 Evolution of the best and average fitness across populations when using fixed and random
rendering seeds using JPEG Size as the fitness function. The results are averages of 30 independent
runs
30K
25K
20K
← FITNESS →
15K
10K
5K
0
Fixed Random Fixed Random
Random Random
Fig. 1.8 Box plots of fitness values of the fittest individuals of each of the 30 evolutionary runs
using different rendering seed setups
the individuals evolved in a fixed environment may change dramatically when the
environmental conditions are different. Conversely, using a dynamic environment
promotes the discovery of robust individuals that perform well under different
conditions. Although this result is not unexpected, it was surprising to notice how
fast the evolutionary algorithm was able to adapt to the changing conditions and
find robust individuals. In future tests we wish to explore, and exploit, this ability.
20 P. Machado et al.
Nevertheless, for the purposes of this Chapter, and considering that the use of a
fixed rendering seed makes the analysis and reproduction of the experimental results
easier, we adopt a fixed rendering seed in all further tests presented in this Chapter.
After establishing the experimental conditions for the evolutionary runs we con-
ducted a series of tests using each of the fitness functions described in Sect. 1.5
to guide evolution. In a second step, based on the results obtained, we combined
several of these measures performing further tests. The results of using each of the
measures individually are presented in Sect. 1.7.1 while those resulting from the
combination of several are presented in Sect. 1.7.2.
1.0 1.0
0.8 0.8
¬ FITNESS ®
¬ FITNESS ®
0.6 0.6
0.4 0.4
0.2 0.2
0 0
0 20 40 60 80 100 120 140 160 180 200 0 20 40 60 80 100 120 140 160 180 200
¬ GENERATION ® ¬ GENERATION ®
0.8 0.8
¬ FITNESS ®
¬ FITNESS ®
0.6 0.6
0.4 0.4
0.2 0.2
0 0
0 20 40 60 80 100 120 140 160 180 200 0 20 40 60 80 100 120 140 160 180 200
¬ GENERATION ® ¬ GENERATION ®
Complexity Bell
1.0 1.0
0.8 0.8
¬ FITNESS ®
¬ FITNESS ®
0.6 0.6
0.4 0.4
0.2 0.2
0 0
0 20 40 60 80 100 120 140 160 180 200 0 20 40 60 80 100 120 140 160 180 200
¬ GENERATION ® ¬ GENERATION ®
Fig. 1.9 Evolution of the fitness of the best individual across populations. The fitness function
used to guide evolution is depicted in the title of each chart. The other values are presented for
reference. The results are averages of 30 independent runs for each chart
during the early stages of the runs, but the number of Contrasting Colors tends to
decrease as the number of generations progresses. The Complexity fitness function
operates on a greyscale version of the images, as such it is not sensitive to changes
of color. Furthermore, abrupt changes from black to white create artifacts that are
hard to encode using JPEG compression, resulting in high IC estimates. Fractal
image compression, which is used to estimate PC, is less sensitive to these abrupt
changes. Therefore, since the approach values images with high IC and low PC,
and since it does not take color information into consideration, the convergence
to images using a reduced palette of contrasting colors is expected. Like for the
other measures, Complexity and Bell appear to be unrelated. Finally, maximizing
Bell promotes an increase of JPEG Size, Contrasting Colors and Complexity during
22 P. Machado et al.
the first generations. It is important to notice that this behavior was only observed
after enforcing a lower bound for 2 (see Sect. 1.5). Without this limit, maximizing
Bell results in the early convergence to simplistic monochromatic images (typically
a single black square on a white background). The adoption of a quadratic DFN
estimate (DFNs ) also contributed to the improvement of the visual results.
Figures 1.10, 1.11, 1.12, 1.13, 1.14, and 1.15 depict the best individual of each
evolutionary run using the different fitness functions individually. A degree of
subjectivity in the analysis of the visual results is unavoidable. Nevertheless, we
Fig. 1.10 Best individual of each of the 30 runs using JPEG Size as fitness function
1 Graph-Based Evolutionary Art 23
Fig. 1.11 Best individual of each of the 30 runs using Contrasting Colors as fitness function
believe that most of the findings tend to be consensual. When using JPEG Size
to guide evolution, the evolutionary engine tended to converge to colorful circular
patterns, with high contrasts of color (see Fig. 1.10). The tendency to converge to
circular patterns, which is observed in several runs, is related with the recursive
nature of the CFDGs and the particularities of the Context Free rendering engine.
For instance, repeatedly drawing and rotating a square while changing its color
will generate images that are hard to encode. Furthermore, the rendering engine
automatically “zooms in” the shapes drawn cropping the empty regions of the
24 P. Machado et al.
Fig. 1.12 Best individual of each of the 30 runs using Fractal Dimension as fitness function
canvas. As such, rotating about a fixed point in space tends to result in images that
fill the entire canvas, maximizing the opportunities for introducing abrupt changes
and, therefore, maximizing file size. Additionally, these CFDGs tend to be relatively
stable and robust, which further promotes the convergence to this type of image.
Unsurprisingly, the results obtained when using Contrasting Colors are char-
acterized by the convergence to images that are extremely colorful. Although
some exceptions exist, most runs converged to amorphous unstructured shapes,
which contrasts with circular patterns found when using JPEG Size. In our opinion
1 Graph-Based Evolutionary Art 25
Fig. 1.13 Best individual of each of the 30 runs using Lacunarity as fitness function
this jeopardizes the aesthetic appeal of the images, that tend to have a random
appearance, both in terms of shape and color.
As anticipated by the data pertaining the evolution of fitness, the visual results
obtained using Fractal Dimension and Lacunarity (Figs. 1.12 and 1.13 are disap-
pointing. None of the runs converged to images of fractal nature. These results
reinforce earlier findings using expression based evolutionary art systems, indicating
that these measures are not suitable for aesthetically driven evolution (Machado
et al. 2007).
26 P. Machado et al.
Fig. 1.14 Best individual of each of the 30 runs using Complexity as fitness function
Fig. 1.15 Best individual of each of the 30 runs using Bell as fitness function
easily attainable and provides the conditions for reaching a natural distribution of
color gradients. Although this is not visible in Fig. 1.15 the individuals reaching the
highest fitness values tend to use a large color palette.
1.0 1.0
0.8 0.8
¬ FITNESS ®
¬ FITNESS ®
0.6 0.6
0.4 0.4
0.2 0.2
0 0
0 20 40 60 80 100 120 140 160 180 200 0 20 40 60 80 100 120 140 160 180 200
¬ GENERATION ® ¬ GENERATION ®
0.8 0.8
¬ FITNESS ®
¬ FITNESS ®
0.6 0.6
0.4 0.4
0.2 0.2
0 0
0 20 40 60 80 100 120 140 160 180 200 0 20 40 60 80 100 120 140 160 180 200
¬ GENERATION ® ¬ GENERATION ®
JPEG Size Contrasting Colors Fractal Dimension Lacunarity Complexity Bell Combined
Fig. 1.16 Evolution of the fitness of the best individual across populations using a combination of
measures. The combination used to guide evolution is depicted in the title of each chart. The other
values are presented for reference, but have no influence in the evolutionary process. The results are
averages of 30 independent runs for each chart and have been normalized to improve readability
followed by a slow, but steady, increase of both measures throughout the runs.
The combination of the three measures further establishes Bell as the measure that
is most difficult to address, since the improvements of fitness are mostly due to
increases in the other two measures. Significantly longer runs would be necessary
to obtain noteworthy improvements in Bell.
Figure 1.17 depicts the best individual of each evolutionary run using as fitness
a combination of the Contrasting Colors and Complexity measures. As it can
be observed, in most cases, the neat structures that characterize the runs using
Complexity (see Fig. 1.14) continue to emerge. However, due to the influence of
the Contrasting Colors measure, they tend to be colorful instead of monochromatic.
Thus, the visual results appear to depict a good combination of both measures. The
same can be stated for the images resulting from using Contrasting Colors and Bell.
As can be observed in Fig. 1.18, they are more colorful than those evolved using Bell
(see Fig. 1.15) but retain a natural distribution of color gradients, deviating from the
“random” coloring schemes that characterize the images evolved using Contrasting
Colors (see Fig. 1.11).
The images obtained when using Complexity and Bell simultaneously (Fig. 1.19)
are less colorful than expected. Visually, the impact of Complexity appears to
overshadow the impact of Bell. Nevertheless, a comparison between these images
and those obtained using Complexity alone (Fig. 1.14) reveals the influence of Bell
30 P. Machado et al.
Fig. 1.17 Best individual of each of the 30 runs using the combination of Contrasting Colors with
Complexity as fitness function
in the course of the runs: the monochromatic images are replaced by ones with a
wider number of color gradients, and these color changes tend to be subtler.
Finally, as expected, the images obtained in the runs using the three measures
(Fig. 1.20) often depict, simultaneously, the features associated with each of them.
As previously, the influence of the Bell measure is less obvious than the others,
but a comparison with the results depicted in Fig. 1.17 highlights the influence of
this measure. Likewise, the structures that emerge from runs using Complexity and
1 Graph-Based Evolutionary Art 31
Fig. 1.18 Best individual of each of the 30 runs using the combination of Contrasting Colors with
Bell as fitness function
the colorful images that characterize runs using Contrasting Colors are also less
often. Thus, although the influence of each measure is observable, we consider that
significantly longer runs would be necessary to enhance their visibility.
32 P. Machado et al.
Fig. 1.19 Best individual of each of the 30 runs using the combination of Complexity with Bell as
fitness function
1.8 Conclusions
We have presented a graph-based approach for the evolution of Context Free Design
Grammars. This approach contrasts with the mainstream evolutionary art practices
by abandoning expression-based evolution of images and embracing the evolution
of images created through the combination of basic shapes. Nevertheless, the
procedural nature of the representation, which characterizes Genetic Programming
1 Graph-Based Evolutionary Art 33
Fig. 1.20 Best individual of each of the 30 runs using the combination of Contrasting Colors,
Complexity and Bell as fitness function
and the influence of the environment in the robustness of the individuals. In the
considered experimental settings, we find that restricting crossover to the portions of
the genome that are expressed and cleaning unexpressed code is advantageous, and
that dynamic environmental conditions promote the evolution of robust individuals.
In a second step, we conducted runs using each of the six fitness functions
individually. The results show that Fractal Dimension and Lacunarity are ill-suited
for aesthetic evolution. The results obtained with the remaining fitness functions are
satisfactory and correspond to our expectations. Finally, we conducted runs using
a combination of the previously described measures to assign fitness. Globally, the
experimental results illustrate the ability of the system to simultaneously address
the different components taken into consideration for fitness assignment. They also
show that some components are harder to optimize than others, and that runs using
several fitness components tend to require a higher number of generations to reach
good results.
One of the most prominent features of the representation adopted herein is its
non-deterministic nature. Namely, the fact that a genotype may be mapped into a
multitude of phenotypes, i.e. images, produced from different expansions of the
same set of rules. As such, each genotype represents a family of shapes that, by
virtue of being generated using the same set of rules, tend to be aesthetically and
stylistically similar. The ability of the system to generate multiple phenotypes from
one genotype was not explored in this Chapter, and will be addressed in future work.
Currently we are conducting experiments where the fitness of a genotype depends
on a set of phenotypes generated from it. The approach values genotypes which are
able to consistently produce fit and diverse individuals, promoting the discovery of
image families that are simultaneously coherent and diverse.
Acknowledgements This research is partially funded by: Fundação para a Ciência e Tecnologia
(FCT), Portugal, under the grant SFRH/BD/90968/2012; project ConCreTe. The project ConCreTe
acknowledges the financial support of the Future and Emerging Technologies (FET) programme
within the Seventh Framework Programme for Research of the European Commission, under
FET grant number 611733. We acknowledge and thank the contribution of Manuel Levi who
implemented the Contrasting Colors fitness function.
References
E. Baker and M. Seltzer. Evolving line drawings. In Proceedings of the Fifth International
Conference on Genetic Algorithms, pages 91–100. Morgan Kaufmann Publishers, 1994.
J. Bird, P. Husbands, M. Perris, B. Bigge, and P. Brown. Implicit fitness functions for evolving a
drawing robot. In Applications of Evolutionary Computing, pages 473–478. Springer, 2008.
J. Bird and D. Stokes. Minimal creativity, evaluation and fractal pattern discrimination. Programme
Committee and Reviewers, page 121, 2007.
A. Borrell. CFDG Mutate. http://www.wickedbean.co.uk/cfdg/index.html, last accessed in Novem-
ber 2014.
C. Coyne. Context Free Design Grammar. http://www.chriscoyne.com/cfdg/, last accessed in
November 2014.
1 Graph-Based Evolutionary Art 35
E. den Heijer and A. E. Eiben. Evolving art with scalable vector graphics. In N. Krasnogor and
P. L. Lanzi, editors, GECCO, pages 427–434. ACM, 2011.
Y. Fisher, editor. Fractal Image Compression: Theory and Application. Springer, London, 1995.
D. A. Hart. Toward greater artistic control for interactive evolution of images and animation.
In Proceedings of the 2007 EvoWorkshops 2007 on EvoCoMnet, EvoFIN, EvoIASP, EvoIN-
TERACTION, EvoMUSART, EvoSTOC and EvoTransLog, pages 527–536, Berlin, Heidelberg,
2007. Springer-Verlag.
J. Horigan and M. Lentczner. Context Free. http://www.contextfreeart.org/, last accessed in
September 2009.
C. Incorporated. GIF Graphics Interchange Format: A standard defining a mechanism for the
storage and transmission of bitmap-based graphics information. Columbus, OH, USA, 1987.
A. Karperien. Fraclac for imagej. In http:// rsb.info.nih.gov/ ij/ plugins/ fraclac/ FLHelp/
Introduction.htm, 1999–2013.
P. Machado and A. Cardoso. All the truth about NEvAr. Applied Intelligence, Special Issue on
Creative Systems, 16(2):101–119, 2002.
P. Machado and H. Nunes. A step towards the evolution of visual languages. In First International
Conference on Computational Creativity, Lisbon, Portugal, 2010.
P. Machado, H. Nunes, and J. Romero. Graph-based evolution of visual languages. In C. D. Chio,
A. Brabazon, G. A. D. Caro, M. Ebner, M. Farooq, A. Fink, J. Grahl, G. Greenfield, P. Machado,
M. ONeill, E. Tarantino, and N. Urquhart, editors, Applications of Evolutionary Computation,
EvoApplications 2010: EvoCOMNET, EvoENVIRONMENT, EvoFIN, EvoMUSART, and Evo-
TRANSLOG, Istanbul, Turkey, April 7–9, 2010, Proceedings, Part II, volume 6025 of Lecture
Notes in Computer Science, pages 271–280. Springer, 2010.
P. Machado, J. Romero, A. Cardoso, and A. Santos. Partially interactive evolutionary artists.
New Generation Computing – Special Issue on Interactive Evolutionary Computation,
23(42):143–155, 2005.
P. Machado, J. Romero, and B. Manaris. Experiments in computational aesthetics: an iterative
approach to stylistic change in evolutionary art. In J. Romero and P. Machado, editors, The Art
of Artificial Evolution: A Handbook on Evolutionary Art and Music, pages 381–415. Springer
Berlin Heidelberg, 2007.
J. McCormack. Facing the future: Evolutionary possibilities for human-machine creativity. In
J. Romero and P. Machado, editors, The Art of Artificial Evolution: A Handbook on Evolu-
tionary Art and Music, pages 417–451. Springer Berlin Heidelberg, 2007.
T. Mori, Y. Endou, and A. Nakayama. Fractal analysis and aesthetic evaluation of geometrically
overlapping patterns. Textile research journal, 66(9):581–586, 1996.
M. O’Neill, J. McDermott, J. M. Swafford, J. Byrne, E. Hemberg, A. Brabazon, E. Shotton,
C. McNally, and M. Hemberg. Evolutionary design using grammatical evolution and shape
grammars: Designing a shelter. International Journal of Design Engineering, 3(1):4–24, 2010.
M. O’Neill and C. Ryan. Grammatical evolution: evolutionary automatic programming in an
arbitrary language, volume 4. Springer, 2003.
M. O’Neill, J. M. Swafford, J. McDermott, J. Byrne, A. Brabazon, E. Shotton, C. McNally,
and M. Hemberg. Shape grammars and grammatical evolution for evolutionary design. In
Proceedings of the 11th Annual Conference on Genetic and Evolutionary Computation,
GECCO ’09, pages 1035–1042, New York, NY, USA, 2009. ACM.
F. B. Pereira, P. Machado, E. Costa, and A. Cardoso. Graph based crossover – a case study with
the busy beaver problem. In Proceedings of the 1999 Genetic and Evolutionary Computation
Conference, 1999.
B. J. Ross, W. Ralph, and Z. Hai. Evolutionary image synthesis using a model of aesthetics. In G. G.
Yen, S. M. Lucas, G. Fogel, G. Kendall, R. Salomon, B.-T. Zhang, C. A. C. Coello, and T. P.
Runarsson, editors, Proceedings of the 2006 IEEE Congress on Evolutionary Computation,
pages 1087–1094, Vancouver, BC, Canada, 16–21 July 2006. IEEE Press.
R. Saunders and K. Grace. Teaching evolutionary design systems by extending “Context Free”. In
EvoWorkshops ’09: Proceedings of the EvoWorkshops 2009 on Applications of Evolutionary
Computing, pages 591–596. Springer-Verlag, 2009.
36 P. Machado et al.
K. Sims. Artificial evolution for computer graphics. ACM Computer Graphics, 25:319–328, 1991.
B. Spehar, C. W. G. Clifford, N. Newell, and R. P. Taylor. Universal aesthetic of fractals. Computers
and Graphics, 27(5):813–820, Oct. 2003.
G. Stiny and J. Gips. Shape grammars and the generative specification of paintings and sculpture. In
C. V. Freiman, editor, Information Processing 71, pages 1460–1465, Amsterdam, 1971. North
Holland Publishing Co.
T. Unemi. SBART2.4: Breeding 2D CG images and movies, and creating a type of collage. In
The Third International Conference on Knowledge-based Intelligent Information Engineering
Systems, pages 288–291, Adelaide, Australia, 1999.
L. World. Aesthetic selection: The evolutionary art of Steven Rooke. IEEE Computer Graphics
and Applications, 16(1), 1996.
Chapter 2
Genetic Programming for Modelling
of Geotechnical Engineering Systems
Mohamed A. Shahin
2.1 Introduction
Geotechnical engineering deals with materials (e.g., soil and rock) that, by their
very nature, exhibit varied and uncertain behaviour due to the imprecise physical
processes associated with the formation of these materials. Modelling the behaviour
of such materials is complex and usually beyond the ability of most traditional forms
of physically-based engineering methods (e.g., analytical formulations and limit
equilibrium methods). Artificial intelligence (AI) is becoming more popular and
particularly amenable to modelling the complex behaviour of most geotechnical
engineering materials as it has demonstrated superior predictive ability when
compared to traditional methods. AI is a computational method that attempts to
mimic, in a very simplistic way, the human cognition capability to solve engineering
problems that have defied solution using conventional computational techniques
(Flood 2008). The essence of AI techniques in solving any engineering problem
is to learn by examples of data inputs and outputs presented to them so that the
subtle functional relationships among the data are captured, even if the underlying
relationships are unknown or the physical meaning is difficult to explain. Thus,
AI models are data-driven approaches that rely on the data alone to determine
the structure and parameters that govern a phenomenon (or system), without the
need for making any assumptions about the physical behavior of the system. This
is in contrast to most physically-based models that use the first principles (e.g.,
physical laws) to derive the underlying relationships of the system, which usually
justifiably simplified with many assumptions and require prior knowledge about
Electronic supplementary material The online version of this chapter (doi: 10.1007/978-3-319-
20883-1_2) contains supplementary material, which is available to authorized users.
M.A. Shahin ()
Department of Civil Engineering, Curtin University, Perth, WA 6845, Australia
e-mail: [email protected]
the nature of the relationships among the data. This is one of the main benefits of
AI techniques when compared to most physically-based empirical and statistical
methods. Examples of the available AI techniques are artificial neural networks
(ANNs), genetic programming (GP), support vector machines (SVM), M5 model
trees, and k-nearest neighbors (Elshorbagy et al. 2010). Of these, ANNs are by far
the most commonly used AI technique in geotechnical engineering and interested
readers are referred to Shahin et al. (2001), where the pre-2001 ANN applications in
geotechnical engineering are reviewed in some detail, and Shahin et al. (2009) and
Shahin (2013), where the post-2001 papers of ANN applications in geotechnical
engineering are briefly examined. More recently, GP has been frequently used
in geotechnical engineering and has proved to be successful. The use of GP in
geotechnical engineering is the main focus of this book chapter.
Despite the success of ANNs in the analysis and simulation of many geotech-
nical engineering applications, they have some drawbacks such as the lack of
transparency and knowledge extraction, leading this technique to be criticised as
being black boxes (Ahangar-Asr et al. 2011). Model transparency and knowledge
extraction are the feasibility of interpreting AI models in a way that provides
insights into how model inputs affect outputs. Figure 2.1 shows a representation of
the classification of modelling techniques based on colours (Giustolisi et al. 2007)
in which the higher the physical knowledge used during model development, the
better the physical interpretation of the phenomenon that the model provides to
the user. It can be seen that the colour coding of mathematical modelling can be
classified into white-, black-, and grey-box models, each of which can be explained
as follows (Giustolisi et al. 2007). White-box models are systems that are based
on first principles (e.g., physical laws) where model variables and parameters are
known and have physical meaning by which the underlying physical relationships
of the system can be explained. Black-box models are data-driven or regressive
systems in which the functional form of relationships between model variables
are unknown and need to be estimated. Black-box models rely on data to map
the relationships between model inputs and corresponding outputs rather than to
B
Grey Box
Black Box
Functional nodes
+
–
4 x1 x2 x3 Terminal nodes
Fig. 2.2 Typical example of genetic programming (GP) tree representation for the function:
Œ.4 x1 / = .x2 C x3 /2
2 Genetic Programming for Modelling of Geotechnical Engineering Systems 41
models in the existing population, and the models are selected according to their
fitness. Reproduction is copying a computer model from an existing population
into the new population without alteration. Cross-over is genetically recombining
(swapping) randomly chosen parts of two computer models. Mutation is replacing
a randomly selected functional or terminal node with another node from the same
function or terminal set, provided that a functional node replaces a functional node
and a terminal node replaces a terminal node. The evolutionary process of evaluating
the fitness of an existing population and producing new population is continued
until a termination criterion is met, which can be either a particular acceptable
error or a certain maximum number of generations. The best computer model
that appears in any generation designates the result of the GP process. There are
currently three variants of GP available in the literature including the linear genetic
programming (LGP), gene expression programming (GEP), and multi-expression
programming (MEP) (Alavi and Gandomi 2011). More recently, the multi-stage
genetic programming (MSGP) (Gandomi and Alavi 2011) and multi-gene genetic
programming (MGGP) (Gandomi and Alavi 2012) are also introduced. However,
GEP is the most commonly used GP method in geotechnical engineering and is
thus described in some detail below.
Gene expression programming was developed by Ferreira (2001) and utilises
evolution of mathematical equations that are encoded linearly in chromosomes of
fixed length and expressed non-linearly in the form of expression trees (ETs) of
different sizes and shapes. The chromosomes are composed of multiple genes, each
gene is encoded a smaller sub-program or sub-expression tree (Sub-ET). Every
gene has a constant length and consists of a head and a tail. The head can contain
functions and terminals (variables and constants) required to code any expression,
whereas the tail solely contains terminals. The genetic code represents a one-to-
one relationship between the symbols of the chromosome and the function or
terminal. The process of information decoding from chromosomes to expression
trees is called translation, which is based on sets of rules that determine the spatial
organisation of the functions and terminals in the ETs and the type of interaction
(link) between the Sub-ETs (Ferreira 2001). The main strength of GEP is that the
creation of genetic diversity is extremely simplified as the genetic operators work at
the chromosome level. Another strength is regarding the unique multi-genetic nature
of GEP, which allows the evolution of more powerful models/programs composed
of several sub-programs (Ferreira 2001).
The major steps in the GEP procedure are schematically represented in Fig. 2.3.
The process begins with choosing sets of functions F and terminals T to randomly
create an initial population of chromosomes of mathematical equations. One
could choose, for example, the four basic arithmetic operators to form the set of
functions, i.e., F D fC, , , /g, and the set of terminals will obviously consist
of the independent variables of a particular problem, for example, for a problem
that has two independent variables, x1 and x2 would be T D fx1 , x2 g. Choosing
the chromosomal architecture, i.e., the number and length of genes and linking
functions (e.g., addition, subtraction, multiplication, and division), is also part of
this step. The chromosomes are then expressed as expression trees of different sizes
42 M.A. Shahin
X
Ct
ˇ ˇ
fi D M ˇC.i;j/ Tj ˇ (2.1)
jD1
where M is the range of selection, C(i.j) is the value returned by the individual
chromosome i for fitness case j (out of Ct fitness cases), and Tj is the target value
for the fitness case j. There are, of course, other fitness functions available that
can be appropriate for different problems. If the desired results (according to the
measured errors) are satisfactory, the GEP process is stopped, otherwise, some
chromosomes are selected and mutated to reproduce new chromosomes, and the
process is repeated for a certain number of generation or until the desired fitness
score is obtained.
Figure 2.4 shows a typical example of a chromosome with one gene, and its ET
and corresponding mathematical equation. It can be seen that, while the p head of a
gene contains arithmetic and trigonometric functions (e.g., C, , , /, , sin, cos),
the tail includes constants and independent variables (e.g., 1, a, b, c). The ET is
codified reading the ET from left to right in the top line of the tree and from top to
bottom.
More recently, a genetic programming based technique called evolutionary
polynomial regression (EPR) was developed and used in geotechnical engineering.
2 Genetic Programming for Modelling of Geotechnical Engineering Systems 43
Fig. 2.4 Schematic representation of a chromosome with one gene and its expression tree (ET)
and corresponding mathematical equation (Kayadelen 2011)
EPR is a hybrid regression technique that was developed by Giustolisi and Savic
(2006). It constructs symbolic models by integrating the soundest features of
numerical regression, with genetic programming and symbolic regression (Koza
1992). The following two steps roughly describe the underlying features of the
EPR technique, aimed to search for polynomial structures representing a system. In
the first step, the selection of exponents for polynomial expressions is carried out,
employing an evolutionary searching strategy by means of GA (Goldberg 1989). In
the second step, numerical regression using the least square method is conducted,
aiming to compute the coefficients of the previously selected polynomial terms. The
general form of expression in EPR can be presented as follows (Giustolisi and Savic
2006):
X
m
yD F X; f .X/; aj C ao (2.2)
jDi
where y is the estimated vector of output of the process, m is the number of terms
of the target expression, F is a function constructed by the process, X is the matrix
of input variables, f is a function defined by the user, and aj is a constant. A typical
example of EPR pseudo-polynomial expression that belongs to the class of Eq. (2.2)
is as follows (Giustolisi and Savic 2006):
X
m h i
b
Y D ao C aj : .X1 /ES.j;1/ : : : .Xk /ES.j;k/ :f .X1 /ES.j;kC1/ : : : .Xk /ES.j;2k/ (2.3)
jDi
44 M.A. Shahin
where Ŷ is the vector of target values, m is the length of the expression, aj is the
value of the constants, Xi is the vector(s) of the k candidate inputs, ES is the matrix
of exponents, and f is a function selected by the user.
EPR is suitable for modelling physical phenomena, based on two features (Savic
et al. 2006): (1) the introduction of prior knowledge about the physical system/pro-
cess, to be modelled at three different times, namely before, during, and after EPR
modelling calibration; and (2) the production of symbolic formulas, enabling data
mining to discover patterns that describe the desired parameters. In the first EPR
feature (1) above, before the construction of the EPR model, the modeller selects
the relevant inputs and arranges them in a suitable format according to their physical
meaning. During the EPR model construction, model structures are determined by
following user-defined settings such as general polynomial structure, user-defined
function types (e.g., natural logarithms, exponentials, tangential hyperbolics), and
searching strategy parameters. The EPR starts from true polynomials and also
allows for the development of non-polynomial expressions containing user-defined
functions (e.g., natural logarithms). After EPR model calibration, an optimum
model can be selected from among the series of models returned. The optimum
model is selected based on the modeller’s judgement, in addition to statistical
performance indicators such as the coefficient of determination. A typical flow
diagram of the EPR procedure is shown in Fig. 2.5, and a detailed description of
the technique can be found in Giustolisi and Savic (2006).
Fig. 2.5 Typical flow diagram of the evolutionary polynomial regression (EPR) procedure
(Rezania et al. 2011)
cohesionless soils utilizing EPR technique. Using a large database that contains 187
data records of field measurements of settlement of shallow foundations as well
as the corresponding information regarding the footings and soil, Shahin (2015)
developed an EPR model that was found to outperform the most commonly used
traditional methods. The data were obtained from the literature and cover a wide
range of variation in footing dimensions and cohesionless soil types and properties.
Details of the references from which the data were obtained can be found in Shahin
et al. (2002a). The model was trained using five inputs representing the footing
width, net applied footing pressure, average blow count obtained from the standard
penetration test (SPT) over the depth of influence of the foundations as a measure
of soil compressibility, footing length, and footing embedment depth. The single
model output was the foundation settlement. The EPR returned several different
models and the one selected to be optimal is as follows (Shahin 2015):
p p
q q B q B qDf
SpEPR D 8:327 2 C 8:849 2 C 2:993 0:651 C 2:883 (2.4)
N L N N N
where Sp (mm) is the predicted settlement, B (m) is the footing width, q (kPa) is the
net applied footing pressure, N is the average SPT blow count, L (m) is the footing
length, and Df (m) is the footing embedment depth.
The results between the predicted and measured settlements obtained by utilising
GP model were compared with those obtained from an artificial neural networks
(ANN) model previously developed by the author (Shahin et al. 2002b), and three
traditional methods, namely, Meyerhof (1965), Schultze and Sherif (1973), and
Schmertmann (1978). Comparisons of the results obtained using the GP model and
the methods used for comparison in the validation set are given in Table 2.1. It can be
seen that the EPR model performs better than the other methods, including the ANN
model, in all performance measures used including the coefficient of correlation, r,
coefficient of determination, R2 , root mean squared error, RMSE, mean absolute
error, MAE, and ratio of average measured to predicted outputs, .
Using the same database of Shahin et al. (2002a) and similar model inputs
and outputs used above, Rezania and Javadi (2007) and Shahnazari et al. (2014)
Table 2.1 Comparison of EPR model and other methods in the validation set for settlement of
shallow foundations on cohesionless soils (Shahin 2015)
Method
ANNs Schultze
EPR (Shahin et al. Meyerhof and Sherif Schmertmann
Performance measure (Shahin 2015) 2002a, b) (1965) (1973) (1978)
r 0.923 0:905 0:440 0:729 0:838
R2 0.844 0:803 0:014 0:185 0:153
RMSE (mm) 9.83 11:04 24:71 22:48 22:91
MAE (mm) 6.99 8:78 16:91 11:29 16:23
1.03 1:10 0:91 1:73 0:79
48 M.A. Shahin
The above GP models represented by Eqs. (2.5) and (2.7) were compared with the
traditional methods and found to outperform most available methods.
Dqctip p
u.steeldriven/ D 2:277 q
QEPR C 0:096DL C 1:714 104 D2 qctip L
qcshaft f sshaft
q
6:279 109 D2 L2 qctip f stip C 243:39 (2.7)
2 Genetic Programming for Modelling of Geotechnical Engineering Systems 49
Dqctip
Qu.concretedriven/ D 2:277 q C 0:096DL
qcshaft f sshaft
p
C 1:714 104 D2 qctip L 6:279
q
109 D2 L2 qctip f stip C 486:78 (2.8)
where D (mm) is the pile perimeter/ (for driven piles) or pile stem diameter (for
drilled shafts), L (m) is the pile embedment length, B (mm) is the drilled shaft base
diameter, qctip (MPa) is the weighted average cone point resistance over pile tip
failure zone, f stip (kPa) is the weighted average cone sleeve friction over pile tip
failure zone, qcshaft (MPa) is the weighted average cone point resistance over pile
embedment length, and f sshahft (kPa) is the weighted average cone sleeve friction
over pile embedment length.
The performance of the above EPR models, represented by Eqs. (2.7)(2.9),
was compared with four other models in the validation set and the results are given
in Table 2.2. For driven piles, the methods considered for comparison include an
ANN model developed by Shahin (2010), the European method (de Ruiter and
Beringen 1979), LCPC method (Bustamante and Gianeselli 1982), and Eslami and
Fellenius (1997) method. For drilled shafts, the methods considered for comparison
include an ANN model (Shahin 2010), Schmertmann (1978) method, LCPC method
(Bustamante and Gianeselli 1982), and Alsamman (1995) method. It can been seen
from Table 2.2 that the performance of the EPR models is as good as the ANN
model, or better, and outperforms the other available methods with the possible
exception of Alsamman (1995).
The application of GP in estimating the capacity of pile foundations was carried
out by Alkroosh and Nikraz (2011, 2012). Correlation models for predicting the
relationship between pile axial capacity and CPT data using gene expression
programming (GEP) technique were developed. The GEP models were developed
for bored piles as well as driven piles (a model for each of concrete and steel piles).
The performance of the GEP models was evaluated by comparing their results with
experimental data as well as the results of a number of currently used CPT-based
methods. The results indicated the potential ability of GEP models in predicting the
bearing capacity of pile foundations and outperformance of the developed models
over existing methods. More recently, Alkroosh and Nikraz (2014) developed GEP
50 M.A. Shahin
Table 2.2 Comparison of EPR model and other methods in the validation set for bearing
capacity of pile foundations (Shahin 2015)
Methods for driven piles
de Ruiter
EPR ANNs and Eslami and
Performance (Shahin (Shahin Beringen Bustamante and Fellenius
measure 2015) 2010) (1979) Gianeselli (1982) (1997)
r 0.848 0.837 0.799 0.809 0.907
R2 0.745 0.753 0.219 0.722 0.681
RMSE (kN) 249.0 244.0 435.0 260.0 278.0
MAE (kN) 185.0 203.0 382.0 219.0 186.0
1.00 0.97 1.36 1.11 0.94
Methods for drilled shafts
EPR ANNs
Performance (Shahin (Shahin Schmertmann Bustamante and Alsamman
measure 2015) 2010) (1978) Gianeselli (1982) (1995)
r 0.990 0.970 0.901 0.951 0.984
R2 0.944 0.939 0.578 0.901 0.939
RMSE (kN) 511.0 533.0 1404.0 681.0 534.0
MAE (kN) 347.0 374.0 702.0 426.0 312.0
1.03 1.02 1.33 0.97 1.03
model that correlates the pile capacity with the dynamic input and SPT data. The
performance of the model was assessed by comparing its predictions with those
calculated using two commonly used traditional methods and an ANN model. It was
found that the GEP model performed well with a coefficient of determination of 0.94
and 0.96 in the training and testing sets, respectively. The results of comparison with
other available methods showed that the GEP model predicted the pile capacity more
accurately than existing traditional methods and ANN model. Another successful
application of genetic programming in pile capacity prediction was carried out by
Gandomi and Alavi (2012), who used a multi-gene genetic programming (MGGP)
method for the assessment of the undrained lateral load capacity of driven piles and
undrained side resistance alpha factor of drilled shafts.
necessary tool for design of civil engineering structures located on active zones
of earthquakes. Hence, many researchers used GP to develop new models for
prediction of liquefaction potential of soils and induced deformation.
Alavi and Gandomi (2011) and Gandomi and Alavi (2012) developed gener-
alized GP models including LGP, GEP, MEP, and MGGP for the classification
of several liquefied and non-liquefied case records. Soil and seismic parameters
governing the soil liquefaction potential were used for model development including
the CPT cone tip resistance, qc (MPa), sleeve friction ratio, Rf (%), effective stress
0
at the depth of interest, v (kPa), total stress at the same depth, v (kPa), maximum
horizontal ground surface acceleration, amax (g), and earthquake moment magnitude,
Mw . the existence of the liquefaction (LC) was represented by binary variables,
non-liquefied and liquefied cases were represented by 0 and 1, respectively. The
GP models were developed based on CPT database that contains 226 case records,
with 133 liquefied cases and 93 non-liquefied cases. Out of the available data,
170 case records were used for model training and 56 case records were used for
model validation. The LGP, GEP, MEP, and MGGP models used to classify the non-
liquefied and liquefied cases, LC, are given as follows (Alavi and Gandomi 2011;
Gandomi and Alavi 2012):
1 02 0 0 0
LCLGP D a max v 4q c v 9R f v C 54v C 9v 54M w 378 (2.10)
0 2v
1 Rf qc .Mw qc / Rf 4 Rf
LCGEP D amax 0 5amax C C (2.11)
v amax 2qc v 3 qc C 2
1 9 9Rf 0
LCMEP D amax C 4amax C 4Mw 4qc C 0 C v
v0 4v 4 4
!
2.qc amax Mw /2 C Mw
4Rf (2.12)
v0
When the return of Eqs (2.10)–(2.13) is greater than or equal to 0.5, the case is
marked as “liquefied”, otherwise, it is marked as “non-liquefied”. The accuracy
of the GP models in the training and validation sets were, respectively, LGP
(training D 90 % and validation D 94.64 %), GEP (training D 88.82 % and Vali-
dation D 92.86 %), MEP (training D 86.47 % and validation D 85.71 %), MGGP
(training D 90 % and validation D 96.4 %). These results clearly indicate that the GP
52 M.A. Shahin
models are efficiently capable of classifying the liquefied and non-liquefied cases.
The best classification results are obtained by both LGP and MGGP models, which
yielded similar performance, followed by the GEP and MEP models.
Javadi et al. (2006) introduced GP models for determination of liquefaction
induced lateral spreading. The models were trained and validated using SPT-based
case histories. Separate models were presented to estimate the lateral displacements
for free face as well as gently sloping ground conditions. It was shown that
the GP models are capable of learning the complex relationship between lateral
displacement and its contributing factors, in the form of a high accuracy prediction
function. It was also shown that the attained function can be used to generalize the
learning to predict liquefaction induced lateral spreading for new cases that have not
been used in model calibration. The results of the developed GP models were also
compared with one of the most commonly used available methods in the literature,
i.e., multi linear regression (MLR) model (Youd et al. 2002), and the advantages of
the proposed GP models were highlighted. It was shown that the GP models offer
an improved determination of the lateral spreading over the most commonly used
MLR method.
Another successful application of genetic programming in soil liquefaction
potential was carried out by Rezania et al. (2010), who used CPT results and
EPR method for determination of liquefaction potential in sands. Furthermore,
Kayadelen (2011) used GEP method to forecast the safety factor of soil liquefaction
using standard penetration test (SPT) results. Both of the above GP models were
found to provide more accurate results compared to the conventional available
methods. More recently, Gandomi and Alavi (2013) developed a robust GP model,
coupled with orthogonal east squares, for predicting the soil capacity energy
required to trigger soil liquefaction, and Gandomi (2014) presented a short review
for use of soft computing, including GP, in earthquake engineering.
results by presenting new training examples as new data become available. These
factors combine to make AI techniques a powerful modelling tool in geotechnical
engineering.
In contrast to most AI techniques, GP does not suffer from the problem of lack of
transparency and knowledge extraction. GP has the ability to generate transparent,
compact, optimum and well-structured mathematical formulations of the system
being studied, directly from raw experimental or field data. Furthermore, prior
knowledge about the underlying physical process based on engineering judgement
or human expertise can also be incorporated into the learning formulation, which
greatly enhances the usefulness of GP over other AI techniques. It was evident
from the review presented in this chapter that GP has been applied successfully to
several applications in geotechnical engineering. Based on the results of the studies
reviewed, it can be concluded that genetic programming models provide high level
of prediction capability and outperform most traditional methods.
References
ADARSH, S., DHANYA, R., KRISHNA, G., MERLIN, R. & TINA, J. 2012. Prediction of
ultimate bearing capacity of cohesionless soils using soft computing techniques. ISRN Artificial
Intelligence, 2012, 10pp.
ADARSH, S. A. & JANGAREDDY, M. 2010. Slope stability modeling using genetic program-
ming. International Journal of Earth Sciences and Engineering, 3, 1–8.
ADELI, H. 2001. Neural networks in civil engineering: 1989-2000. Computer-Aided Civil and
Infrastructure Engineering, 16, 126–142.
AHANGAR-ASR, A., FARAMARZI, A. & JAVADI, A. 2010. A new approach for prediction
of the stability of soil and rock slopes. Engineering Computations: International Journal of
Computer-Aided Engineering and Software, 27, 878–893.
AHANGAR-ASR, A., FARAMARZI, A., MOTTAGHIFARD, N. & JAVADI, A. A. 2011. Mod-
eling of permeability and compaction characteristics of soils using evolutionary polynomial
regression. Computers and Geosciences, 37, 1860–1869.
ALAVI, A. H. & GANDOMI, A. H. 2011. A robust data mining approach for formulation
of geotechnical engineering systems. Engineering Computations: International Journal of
Computer-Aided Engineering and Software, 28, 242–274.
ALAVI, A. H. & GANDOMI, A. H. 2012. Energy-based models for assessment of soil liquefaction.
Geoscience Frontiers.
ALAVI, A. H., GANDOMI, A. H., NEJAD, H. C., MOLLAHASANI, A. & RASHED, A. 2013.
Design equations for prediction of pressuremeter soil deformation moduli utilizing expression
programming systems. Neural Computing and Applications, 23, 1771–1786.
ALAVI, A. H., GANDOMI, A. H., SAHAB, M. G. & GANDOMI, M. 2010. Multi expression pro-
gramming: a new approach to formulation of soil classification. Engineering with Computers,
26, 111–118.
ALAVI, A. H., MOLLAHASANI, A., GANDOMI, A. H. & BAZA, J. B. 2012a. Formulation of
secant and reloading soil deformation moduli using multi expression programming. Engineer-
ing Computations: International Journal of Computer-Aided Engineering and Software, 29,
173–197.
ALAVI, A. M., GANDOMI, A. H., BOLURY, J. & MOLLAHASANI, A. 2012b. Linear and
tree-based genetic programming for solving geotechnical engineering problems. In: YANG,
54 M.A. Shahin
GANDOMI, A. H. & ALAVI, A. H. 2013. Hybridizing genetic programming with orthogonal least
squares for modeling of soil liquefaction. International Journal of Earthquake Engineering and
Hazard Mitigation, 1, 2–8.
GANDOMI, A. H., ALAVI, A. H. & GUN, J. Y. 2011. Formulation of uplift capacity of
suction caissons using multi expression programming. KSCE Journal of Civil Engineering, 15,
363–373.
GARG, A., GARG, A., TAI, K. & SREEDEEP, S. 2014a. Estimation of pore water pressure of soil
using genetic programming. Geotechnical and Geological Engineering, 32, 765–772.
GARG, A., GARG, A., TAI, K. & SREEDEEP, S. 2014b. An integrated SRP-multi-gene genetic
programming approach for prediction of factor of safety of 3-D soil nailed slopes. Engineering
Applications of Artificial Intelligence, 30, 30–40.
GIUSTOLISI, O., DOGLIONI, A., SAVIC, D. A. & WEBB, B. W. 2007. A multi-model approach
to analysis of environmental phenomena. Environmental Modelling and Software, 22, 674–682.
GIUSTOLISI, O. & SAVIC, D. A. 2006. A symbolic data-driven technique based on evolutionary
polynomial regression. Journal of Hydroinformatics, 8, 207–222.
GOLDBERG, D. E. 1989. Genetic Algorithms in Search Optimization and Machine Learning,
Mass, Addison - Wesley.
HOLLAND JH. 1975 Published. Adaptation in natural and artificial systems. 1975 University of
Michigan
JAVADI, A., AHANGAR-ASR, A., JOHARI, A., FARAMARZI, A. & TOLL, D. 2012a. Mod-
elling stress-strain and volume change behaviour of unsaturated soils using an evolutionary
based data mining technique, and incremental approach. Engineering Applications of Artificial
Intelligence, 25, 926–933.
JAVADI, A., FARAMARZI, A. & AHANGAR-ASR, A. 2012b. Analysis of behaviour of soils
under cyclic loading using EPR-based finite element method. Finite Elements in Analysis and
Design, 58, 53–65.
JAVADI, A. & REZANIA, M. 2009. Intelligent finite element method: An evolutionary approach
to constitutive modelling. Advanced Engineering Informatics, 23, 442–451.
JAVADI, A., REZANIA, M. & MOUSAVI, N. M. 2006. Evaluation of liquefaction induced lateral
displacements using genetic programming. Computers and Geotechnics, 33, 222–233.
JOHARI, A., HABIBAGAHI, G. & GHAHRAMANI, A. 2006. Prediction of soil-water char-
acteristic curve using genetic programming. Journal of Geotechnical and Geoenvironmental
Engineering, 132, 661–665.
KAYADELEN, C. 2011. Soil liquefaction modeling by genetic expression programming and
neuro-fuzzy. Expert Systems with Applications, 38, 4080–4087.
KOZA, J. R. 1992. Genetic programming: on the programming of computers by natural selection,
Cambridge (MA), MIT Press.
MEYERHOF, G. G. 1965. Shallow foundations. Journal of Soil Mechanics & Foundation
Engineering Division, 91, 21–31
MOLLAHASANI, A., ALAVI, A. H. & GANDOMI, A. H. 2011. Empirical modeling of plate
load test moduli of soil via gene expression programming. Computers and Geotechnics, 38,
281–286.
MOUSAVI, S. M., ALAVI, A. H., MOLLAHASANI, A., GANDOMI, A. H. & ESMAEILI, M. A.
2013. Formulation of soil angle of resistance using a hybrid GP and OLS method. Engineering
with Computers, 29, 37–53.
MUDULI, P. K. & DAS, S. K. 2013. SPT-based probabilistic method for evaluation of liquefaction
potential of soil using multi-gene genetic programming. International Journal of Geotechnical
Earthquake Engineering, 4, 42–60.
MUDULI, P. K. & DAS, S. K. 2014. CPT-based seismic liquefaction potential evaluation using
multi-gene genetic programming approach. Indian Geotechnical Journal, 44, 86–93.
NADERI, N., ROSHANI, P., SAMANI, M. Z. & TUTUNCHIAN, M. A. 2012. Application of
genetic programming for estimation of soil compaction parameters. Applied Mechanics and
Materials, 147, 70–74.
56 M.A. Shahin
3.1 Introduction
In the real world, there are several natural and artificial phenomenons which follow
some rules. These rules can model in a mathematical and/or logical form consider-
ing simple or complex equation/s may be difficult in some systems. Moreover, it is
sometimes necessary to model just some parts of system without considering whole
system information. Data-driven models are a programming paradigm that employs
a sequence of steps to achieve best connection between data sets.
The contributions from artificial intelligence, data mining, knowledge discovery
in databases, computational intelligence, machine learning, intelligent data analysis,
soft computing, and pattern recognition are main cores of data-driven models with
a large overlap in the disciplines mentioned.
GP is a data-driven tool which applies computational programming to achieve the
best relation in a system. This tool can set in the inner or outer of system modeling
which makes it more flexible to adapt different system states.
Electronic supplementary material The online version of this chapter (doi: 10.1007/978-3-319-
20883-1_3) contains supplementary material, which is available to authorized users.
E. Fallah-Mehdipour
Department of Irrigation & Reclamation Engineering, Faculty of Agricultural Engineering
& Technology, College of Agriculture & Natural Resources, University of Tehran, Karaj,
Tehran, Iran
National Elites Foundation, Tehran, Tehran, Iran
e-mail: [email protected]
O. Bozorg Haddad ()
Department of Irrigation & Reclamation Engineering, Faculty of Agricultural Engineering
& Technology, College of Agriculture & Natural Resources, University of Tehran, Karaj,
Tehran, Iran
e-mail: [email protected]
Real (observed)
Input data set/s Mathematical model
output data set
Minimization
GP
this difference
a b
/
+
x 47 y
y x
47 exp(y)
f (x, y) = sin(x) + f (x,y) =
y cos(x)
Fig. 3.2 Two GP expressions in the tree structure
To generate the next tree set, trees with the better fitness values are selected using
techniques such as roulette wheel, tournament, or ranking methods (Orouji et al.
2014). In following, crossover and mutation as the two genetic operators as same
as GA operators create new trees using the selected trees. In the crossover operator,
two trees are selected and sub-tree crossover randomly (and independently) selects
a crossover point (a node) in each parent tree. Then, two new trees are produced
by replacing the sub-tree rooted at the crossover point in a copy of the first parent
with a copy of the sub-tree rooted at the crossover point in the second parent, as
illustrated in Fig. 3.3 (Fallah-Mehdipour et al. 2012).
In the mutation operator, point mutation is applied on a per node basis. That is,
some node/s are randomly selected, it is exchanged by another random terminal or
function, as it is presented in Fig. 3.4. The produced trees using genetic operators are
the input trees for the next iteration and the GP process continues up to a maximum
number of iterations or minimum of error performance.
+ /
sin /
+ exp cos
x 47 y
y x
⇓
+ /
x y 47 y x
/
/
/ cos ⇒ / sin
47 y x
47 x x
A watershed is a hydrologic unit in which surface water from rain, melting snow
and/or ice converges to a single point at a lower elevation, usually the exit of the
basin. Commonly, water that moves to external point and join another water body,
such as river, lake or sea. Figure 3.6 presents schematic of a watershed.
When rain falls on watershed, water that called runoff, flows on it. A rainfall-
runoff model is a mathematical model describing relations between rainfall and
runoff for a watershed. In this case, conceptual models are usually used to
obtain both short- and long-term forecasts of runoff. These models are applied
several variables such as climate parameters, topography and land use variables to
determine runoff volume. Thus, that volume depends directly on the accuracy of
each aforementioned variable estimation. On the other hand, some global circulation
model (GCM) that is used for runoff calculation apply for large scale and runoff
volume for smaller scale should be extracted by extra processes.
Input/s
Although conceptual models can calculate runoff for a watershed, their processes
are long and expensive. Therefore, to overcome these problems, Savic et al. (1999)
applied GP to estimate runoff volume for Kirkton catchment in Scotland.
Rainfall on the Kirkton catchment is estimated using a network of 11 period
gauges and 3 automatic weather stations at different altitudes. The daily average
rainfall is calculated from weighted domain areas for each gauge. Stream flow is
measured by a weir for which the rating has been adjusted after intensive current
metering (Savic et al. 1999). They compared obtained results with HYRROM,
one conceptual model by Eeles (1994) that applied 9 and 35 parameters for
runoff estimation considering different land use variables. Moreover, GP employed
different combinations rainfall, runoff and evaporation for one, two and three
previous periods and rainfall at current period as the input data to estimate runoff of
current period as the output data. Results showed that GP can present better solution
even by fewer input data sets than other conceptual models by Eeles (1994).
When rain falls, extra surface water and runoff moves under earth and forms
groundwater. In groundwater, soil pore spaces and fractures of rock formations fill
from water and called an aquifer. The depth at which soil pores and/or fractures
become completely saturated with water is water table or groundwater level.
Groundwater contained in aquifer systems is affected by various processes, such
as precipitation, evaporation, recharge, and discharge. Groundwater level is typically
measured as the elevation that the water rises in, for example, a test well.
66 E. Fallah-Mehdipour and O. Bozorg Haddad
Rt D F1 .St ; Qt / (3.2)
in which, Rt , St and Qt are release, storage and inflow at tth period. Moreover, F1
is linear or nonlinear function for transferring storage volume and inflow to the
released water from the reservoir at each period.
The common pattern of aforementioned decision rule which is a linear decision
rule that a, b and c are the decision variables (e.g., Mousavi et al. 2007; Bolouri-
Yazdeli et al. 2014):
R t D a Qt C b S t C c (3.3)
GP equation as a
Generate random decision rule Release Evaluate objective
trees Rt=f(St,Qt) function
SMin<St<SMax
Continuity equation
St+1=St+Qt-Rt
process by penalty. This penalty is added and subtracted in the minimization and
maximization objective for each violation unit from feasible bound. The other GP
process (selection, crossover and mutation) are continues to satisfy stopping criteria.
Acknowledgement Authors thank Iran’s National Elites Foundation for financial support of this
research.
3 Application of Genetic Programming in Hydrology 69
References
Batista Celeste, A., Suzuki, K., and Kadota, A. (2008). “Integrating long-and short-term reservoir
operation models via stochastic and deterministic optimization: case study in Japan.” Journal
of Water Resources Planning and Management (ASCE), 134(5), 440–448.
Bolouri-Yazdeli, Y., Bozorg Haddad, O., Fallah-Mehdipour, E., and Mariño, M.A. (2014).
“Evaluation of real-time operation rules in reservoir systems operation.” Water Resources
Management, 28(3), 715–729.
Bozorg Haddad, O., Rezapour Tabari, M.M., Fallah-Mehdipour, E., and Mariño, M.A. (2013).
“Groundwater model calibration by meta-heuristic algorithms.” Water Resources Management,
27(7), 2515–2529.
Cramer, N.L. (1985). “A representation for the adaptive generation of simple sequential programs.”
In Proceedings of an International Conference on Genetic Algorithms and the Applications,
Grefenstette, John J. (ed.), Carnegie Mellon University. 24-26 July, 183–187.
Eakins, B.W. and Sharman, G.F. (2010). Volumes of the World’s Oceans from ETOPO1, NOAA
National Geophysical Data Center, Boulder, CO, 2010.
Eeles CWO: (1994) Parameter optimization of conceptual hydrological models, PhD Thesis, Open
University, Milton Keynes, U.K.
Fallah-Mehdipour, E., Bozorg Haddad, O., and Mariño, M. A. (2012). “Real-time operation of
reservoir system by genetic programming.” Water Resources Management, 26(14), 4091–4103.
Fallah-Mehdipour, E., Bozorg Haddad, O., and Mariño, M. A. (2013a). “Prediction and Simulation
of Monthly Groundwater level by Genetic Programming.” Journal of Hydro-environment
Research, 7(4), 253–260.
Fallah-Mehdipour, E., Bozorg Haddad, O., and Mariño, M. A. (2013b). “Developing reservoir
operational decision rule by genetic programming.” Journal of Hydroinformatics, 15(1),
103–119.
Gandomi, A.H., Yang, X.S., Talatahari, S., and Alavi, A. H. (2013). “Metaheuristic Applications
in Structures and Infrastructures” Elsevier. 568 pages.
Guven, A., and Kisi, O. (2011). “Daily pan evaporation modeling using linear genetic program-
ming technique.” Irrigation Science, 29(2), 135–145.
Izadifar, Z., and Elshorbagy, A. (2010). “Prediction of hourly actual evapotranspiration using
neural network, genetic programming, and statistical models.” Hydrological Processes, 24(23),
3413–3425.
Koza, J. R. (1992). Genetic programming: on the programming of computers by means of natural
selection. MIT Press, Cambridge, MA.
Koza, J. R. (1994). Genetic Programming II: Automatic Discovery of Reusable Programs. MIT
293 Press. Cambridge, MA.
Li, L., Liu, P., Rheinheimer, D.E., Deng, C., and Zhou, Y. (2014). “Identifying explicit formulation
of operating rules for multi-reservoir systems using genetic programming.” Water Resources
Management, 28(6), 1545–1565.
Mousavi, S. J., Ponnambalam, K., and Karray, F. (2007). “Inferring operating rules for reservoir
operations using fuzzy regression and ANFIS.” Fuzzy Sets and Systems, 158(10), 1064–1082.
Nasseri, M., Moeini, A., and Tabesh, M. (2011). “Forecasting monthly urban water demand using
extended Kalman filter and genetic programming.” Expert Systems with Applications, 38(6),
7387–7395.
Orouji, H., Bozorg Haddad, O., Fallah-Mehdipour, E., and Mariño, M.A. (2014). “Flood routing
in branched river by genetic programming.” Proceedings of the Institution of Civil Engineers:
Water Management, 167(2), 115–123.
Savic, D. A., Walters, G. A., and Davidson, J. W. (1999). “A genetic programming approach to
rainfall-runoff modeling.” Water Resources Management, 13(3), 219–231.
Sivapragasam, C., Vasudevan, G., Maran, J., Bose, C., Kaza, S., and Ganesh, N. (2009). “Modeling
evaporation-seepage losses for reservoir water balance in semi-arid regions.” Water Resources
Management, 23(5), 853–867.
70 E. Fallah-Mehdipour and O. Bozorg Haddad
Traore, S., and Guven, A. (2012). “Regional-specific numerical models of evapotranspiration using
gene-expression programming interface in Sahel.” Water Resources Management, 26(15),
4367–4380.
Traore, S., and Guven, A. (2013). “New algebraic formulations of evapotranspiration extracted
from gene-expression programming in the tropical seasonally dry regions of West Africa.”
Irrigation Science, 31(1), 1–.10.
Yang, X.S., Gandomi, A.H., Talatahari, S., and Alavi, A. H. (2013a). “Metaheuristis in Water,
Geotechnical and Transportation Engineering” Elsevier. 496 pages
Yang, X.S., Cui, Z., Xiao, R., Gandomi, A.H., and Karamanoglu, M. (2013b). “Swarm Intelligence
and Bio-Inspired Computation: Theory and Applications”, Elsevier. 450 pages.
Chapter 4
Application of Gene-Expression Programming
in Hydraulic Engineering
4.1 Introduction
Electronic supplementary material The online version of this chapter (doi: 10.1007/978-3-319-
20883-1_4) contains supplementary material, which is available to authorized users.
A. Zahiri () • A.A. Dehghani
Department of Water Engineering, Gorgan University of Agricultural Sciences
and Natural Resources, Gorgan, Iran
e-mail: [email protected]; [email protected]; [email protected]
H.Md. Azamathulla
Associate Professor of Civil Engineering, Faculty of Engineering, University of Tabuk, Tabuk,
Saudi Arabia
e-mail: [email protected]
In which q is unit discharge over the spillway, H1 is total head, R is radius of the
bucket, is lip angle of the bucket, dw is tail water depth, d50 is median sediment size
and g is acceleration due to gravity. Results of these equations were compared with
the regression equation formulae and neural network approach. The comparison
revealed that the GEP models (Eqs. 4.1–4.3) have higher accuracy.
Mujahid et al. (2012) used GEP for estimation of bridge pier scour. The following
explicit relation was obtained and its performance was compared with artificial
neural networks (ANNs) and conventional regression-based techniques. The results
showed that GEP gives more accurate results than the other models.
4 Application of Gene-Expression Programming in Hydraulic Engineering 73
!
ds b d50 b 2 d50
D 0:595 Fr C Fr Fr C 0:063
Y Y Y Y Y
0 0 0 . 111
2 3:24d50 Y
d50 b
. 1/ C Fr @ @Fr @ AAA (4.4)
Y Y .Fr 1/
In the above equation, Y is approach flow depth, Fr is Froude number, b is pier width
and ¢ is standard deviation of particle grain size distribution.
Wang et al. (2013) using GEP model, presented Eq. (4.5) for estimating pier
scours depth based on available experimental data from various researches. Four
main dimensionless parameters such as pier width (D/d50 ), approaching flow depth
(Y/d50 ), threshold flow velocity ((V2 –Vc 2 )/(4gd50 )), and pier scour depth (ds /D)
were used as independent variables in Eq. (4.5).
0 , 1
, !! B D V 2 Vc 2
2
log@ A
ds V Vc 2 D Y d50
gd50
D ˛
D
gd d50 d50
0 0 , 1 10:156
, ! , B D C
B V 2 V 2 2:05 D log@ 1:44A C
B c D d50 C
B 4:94 ˛ C (4.5)
@
gd50 d50 d50 A
The results showed that GEP as an effective tool can be used for estimating
of daily discharge data in flood events. For developing flow rating curves, also
Guven and Aytek (2009) used GEP technique in two stations of Schuylkill River
(Pennsylvania). The performance of GEP approach was compared with more
conventional methods, common stage rating curve (SRC) and multiple linear
regression (MLR) techniques. The results showed that GEP gives more accuracy
74 A. Zahiri et al.
results than SRC and MLR models. Equations (4.7) and (4.8) were obtained for
discharge as a function of stage for upstream (Berne) and downstream stations
(Philadelphia), respectively:
Zakaria et al. (2010) used GEP model to predict total load transport in three rivers
(i.e., Kurau, Langat, and Muda). The explicit formulation of GEP for total bed
material load is presented in Eq. (4.9):
p
Qs D 0:39Rh Y 0 S0 = .0:72 CS0 / C Rh C eSin.QVRh /
p 3
C Tan1 .0:16Rh B/ C Rh Q C .d50 3:39/ d50 S0 C Tan1 .V/ eV
log .6:93 Y0 / ..!B/ = .2:075//
(4.9)
(Sharifi 2009). The branch of Genetic based algorithms states that the survival of
an organism is affected by the rule of “survival of the strongest species”. This
family of algorithms is often viewed as function optimizers, although the range of
problems to which genetic algorithms have been applied is quite broad (Whitley
1991). Among the different methods belong to this family, the Gene expression
programming (GEP) is the up-to-date technique. GEP was invented by Ferreira
in 1999, and is the natural development of EAs. The great insight of GEP was
the invention of chromosomes capable of representing any expression tree; GEP
surpasses the genetic programming (GP) system by a factor of 60,000 (Ferreira
2001). In GEP, complex relations are encoded in simpler, linear structures of a fixed
length called chromosomes. The chromosomes consist of a linear symbolic string
of a fixed length composed of one or more genes (Sattar 2014). To express the
genetic information encoded in the gene, Ferreira (2001) used expression tree (ET)
representations. Due to the simple rules that determine the structure of the ET, it
is possible to infer the gene composition given the ET and vice versa using the
unequivocal Karva language. The Karva language represents genes in a sequence
that begins with a start codon, continues with amino acid codons, and ends with a
termination codon (Sattar 2014). Consider, for example, the following mathematical
expression:
sin.x/
cos.x/C
zDe cos.y/ (4.11)
* %
X X Y
76 A. Zahiri et al.
(2001) and Sattar (2014). GEP fitting computations for experimental or field data is
can be performed using some commercial nonlinear data-mining softwares such as
GeneXProTools (www.gepsoft.com). The fitness function of a program fi in GEP is:
1
fi D 1000 (4.12)
1 C RRSEi
The fitness function ranges from 0 to 1000, with 1000 corresponding to a perfect fit.
In the above equation, RRSE is root relative square error of an individual program i
(i-th offspring) and is defined by the following equation:
v
u n
uX ˇ ˇ
u ˇYij Xj ˇ2
u
u jD1
RRSEi D u uX (4.13)
u
n
ˇ ˇ2
t ˇX X j ˇ
jD1
where X and Y, respectively, are the actual and predicted targets, Yij is the value
predicted by the program i for fitness case j, Xj is the target value for fitness case j,
X is the average of the measured outputs (X) and n is the number of samples (Sattar
2014).
The following procedure has been used to develop the final GEP equation (Sattar
2014);
1. Choosing an initial set of control variables as terminals for GEP.
2. Defining the chromosome architecture (number of genes, head size, functions)
and mutation rates for the initial work environment of GEP.
3. Producing several first-generation offspring by GEP through randomly formu-
lation of the parent program’s chromosomes and implementation of genetic
operators.
4. Selection of the fittest offspring by using the fitness criteria (Eq. 4.12). This
offspring represents the solution to the problem in the first generation.
5. Producing several second-generation offspring by GEP using the fittest offspring
as new parent the then implementing genetic operators.
6. Repeating steps 3–5 until the required program fitness is met. Unfortunately,
there is no specific range for GEP fitness indicator fi , however, a domain of
600–800 has been suggested by Ferreira (2001) for suitable model predictions.
The final GEP-model (the fittest offspring of generation i) is scored on a set of
performance indicators.
7. Repeating steps 1–7 with a different set of control variables to produce another
GEP-model.
4 Application of Gene-Expression Programming in Hydraulic Engineering 77
To evaluate the accuracy of the final explicit GEP equation in both training and test-
ing phases, some common statistical measures including correlation coefficient (R),
the mean squared error (MSE), the mean absolute error (MAE) and performance
index (), are used in this study as follows:
X
xy
R D qX X (4.14)
x2 y2
X
.X Y/2
MSE D (4.15)
N
X jX Yj
MAE D X (4.16)
N
p
MSE 1
D (4.17)
X 1CR
where xD(X X), yD(Y Y), X is the mean of X (measured outputs), Y is the
mean of Y (predicted outputs) and N is the data point’s number for GEP evaluation
(experimental data). The last statistical parameter () is a new criterion proposed by
Gandomi and Roke (2013) which combines both correlation and error functions.
4.3 Applications
Rivers are vital carriers of water and sediments. At extreme discharge conditions
floods may occur that could damage nearby infrastructure and also cause casualties
(Huthoff 2007). Over more than three decades, hydraulics of compound channels
has been extensively investigated by many researchers.
Flow hydraulic characteristics are completely different in main channel and
floodplains, and hence, it’s necessary to treat the channel into subsections for any
analysis and computation (Lambert and Myers 1998). For dividing of compound
78 A. Zahiri et al.
Fig. 4.2 A typical compound channel (a) with common dividing channel methods (b) (Bousmar
2002)
u2 u1
channels, there are three common approaches, vertical, horizontal and diagonal
planes. For hydraulic modeling of river compound channels, the first approach
has the most applications in one-dimensional commercial mathematical packages
such as MIKE11, HEC-RAS, ISIS and SOBEC (Huthoff et al. 2008). However,
this method has great over-prediction error for discharge estimation in field and
laboratory compound sections.
In Fig. 4.2a, a typical compound channel with associated important parameters
is shown. Also, in Fig. 4.2b, the division ways for vertical, horizontal and diagonal
dividing channel methods have been illustrated. As mentioned by many researchers,
it’s assumed, in all these simple methods, that there is no shear stress and momentum
transfer at the division lines between main channel and floodplains. Results of
experimental works carried out in compound channels, have revealed that this
assumption isn’t correct and therefore, these methods are maybe very erroneous
(Martin and Myers 1991; Ackers 1992).
Unreliability of dividing channel methods’ assumption was demonstrated
through Van Prooijen et al. (2005) experimental work. It’s seen from Fig. 4.3, that
due to lateral shear stress, a fully turbulent flow with high momentum transfer
is induced at the main channel/floodplain interface. This shear stress is maybe
comparable, in magnitude, with the bed shear stress. Furthermore, it’s interesting
to note that the flow velocity in the floodplain is considerably less than the main
channel and it’s insufficient to cause major movement of suspended sediment.
In accordance to large error of traditional divided channel methods, several
modified approaches have been provided by many researchers (Wormleaton and
4 Application of Gene-Expression Programming in Hydraulic Engineering 79
Merrett 1990; Ackers 1992; Bousmar and Zech 1999; Atabay 2001; Huthoff et al.
2008; Yang et al. 2014). Even though these methods perform well, but in most of
these studies, the main aim is precise prediction of total flow discharge or average
velocity. However, in many practical situations, distribution of flow rate in the main
channel and the flood plains is also important. Ackers (1992) states that any suitable
method should has such ability to predict flow discharge in subsections, especially
in the main channel, with sufficient accuracy.
In overbank flows the river system not only behaves as a conveyance but also as a
storage or pond. It is recognized that for sediment transport, only the flow discharge
in main channel is effective and floodplain’s discharge is nearly negligible (Ackers
1992). In fact the floodplains, due to their high capacities, play an important role in
flood water level reduction, water retention and sediment deposition. These features
are essential for wetlands restoration and preserve of river ecology as well as for
success of flood mitigation works. The main channel flow discharge determination
also covers the main input data for several hydraulic and morphologic computations
such as pollutant dispersion, sediment transport and bed shear stress distribution in
river compound channels.
The bank-full level is defined as the level at which the water has its maximum
power to move sediment. In flood event and when the water rises above the bank-
full level, flow spills onto the floodplain, which the average flow velocity and
consequently stream power dramatically reduce. As the stream power is reduced, so
too is its capacity sediment transport. Thus, for better monitoring of river behavior
during flood events, accurate computation of flow velocity and hence sediment
transport capacity of both main channel and floodplains are needed.
It should be noted that for computation of sediment transport capacity in flooded
rivers, one initially should separate the main channel flow discharge from the total
flow rate and then put it into a suitable empirical sediment transport equation.
Towards the finding suitable methods for prediction of main channel and
floodplain discharges, first, the traditional methods are reviewed.
In Fig. 4.2b, main channel and floodplains sections separated by three dividing
methods (e.g., horizontal, vertical, and diagonal) are shown. Total flow discharge is
the sum of discharges calculated separately in each subsection using an appropriate
conventional friction formula, for example, Manning’s equation (Chow 1959):
3
X 3
X 2=3 1=2
Ai R S 0
QDCM D Qi D i
(4.18)
iD1 iD1
ni
where QDCM is total flow discharge in compound channel, A is area. In this equation,
i refers to each subsection (main channel or floodplains).
80 A. Zahiri et al.
where P is the wetted perimeter and * denotes the ratio of floodplain to main
channel’s value.
Table 4.1 Overview of data sets used for development and assess-
ment of GEP model
Variable definition Variable range Mean
Bank-full height, h (m) 0.05–0.2 0.103
Flow depth, H (m) 0.058–0.32 0.1482
Main channel width, bc (m) 0.05–1.6 0.89
Floodplain width, bf (m) 0.16–6 1.49
Bank-full discharge, Qb (m3 /s) 0.0023–0.2162 0.096
Total flow discharge, Qt (m3 /s) 0.003–1.1142 0.2145
Main channel flow discharge, Qmc (m3 /s) 0.00233–0.6271 0.1499
Floodplain flow discharge, Qf (m3 /s) 0.00046–0.6340 0.064
Bed slope 0.00099–0.013 0.0021
4 Application of Gene-Expression Programming in Hydraulic Engineering 81
where Qmc and Qfp are flow discharges in main channel and floodplain, respectively.
In Fig. 4.4a, b, results of three methods of dividing compound channels are shown
for flow discharges in main channel and floodplains. As can be seen, for main
channel discharge (Fig. 4.4a), vertical and horizontal approaches have still over and
under predictions, respectively. Errors of these methods are growing with increasing
flow discharges, especially for vertical method. Among these approaches, the
diagonal planes, has a suitable result, although for large main channel’s discharges,
the errors are increasing. For floodplains (Fig. 4.4b), both vertical and diagonal
dividing planes produce considerably better predictions than the horizontal case.
It is interesting to note, that even for large discharges, these two methods have good
results. Furthermore, for floodplains, the vertical divided method under-predicts the
flow discharge, opposite to the main channels case. This is due to the interaction
effect that causes the actual discharge to decrease in the main channel and increase in
the floodplains. This flow exchange isn’t considered in the vertical divided method,
as well as for both horizontal and diagonal methods.
The formulations of GEP model for main channel and floodplains flow dis-
charges, as a function of Dr, COH and vertical divided discharges, were obtained as
following:
a VDCM
.
DDCM
Computed Qmc (m/s) HDCM
. Perfect Line
.
.
. . . .
Observed Qmc (m/s)
b
.
Computed Qfp (m/s)
.
.
VDCM
. DDCM
HDCM
Perfect Line
. . . .
Observed Qfp (m/s)
Fig. 4.4 Flow discharge calculation of traditional divided channel methods for main channels
(a) and floodplains (b)
The performance of the GEP model was compared with the traditional vertical
divided method. Figure 4.5a, b show the observed and estimated main channel and
floodplains flow discharges of the all used data. As can be seen, the GEP model
produced much enhanced results, especially for floodplaindischarges.
Table 4.2 presentsa comparison of R, MSE, MAE and for predicted main
channel flow discharges obtained from different models. It can be concluded that
according to the error functions, especially the mean absolute error, the GEP model
gives much better results than the other approaches. The GEP model produces the
least errors (MSE D 0.0003 and MAE D 2.1 %). Among traditional methods, the
4 Application of Gene-Expression Programming in Hydraulic Engineering 83
a .
.
. VDCM
GEP
Perfect Line
. . . .
Observed Qmc (m/s)
b
.
Computed Qfp (m/s)
.
.
. VDCM
GEP
Perfect Line
. . . .
Observed Qfp (m/s)
Fig. 4.5 Comparison of flow discharge calculation of traditional vertical divided channel method
and GEP model for main channels (a) and floodplains (b)
vertical approach, which is currently used in many engineering packages, with mean
absolute error of 19 %, has the lowest accuracy. On the other hands, the diagonal
approach is the best one.
In rivers, hydrological measurements such as the flow discharge and depth are
essential for the design and implementation of river training works and for water
resources management. Manning equation is the simplest computational tool for
changing flow depth to the discharge. This equation gives an adequate estimate
for flow discharge in rivers, provided that no significant flood occurs in such a
way that river overflows its banks. Field and laboratory experiments conducted
by Martin and Myers (1991) and Lai and Bessaih (2004) indicated that the
maximum errors caused by Manning equation are up to 40 and 60 %, respectively.
To overcome this difficulty, various methods have been developed with different
assumptions for compound channels. Off these many approaches, works of Shiono
and Knight (1991), Ackers (1992), and Bousmar and Zech (1999) have good
accuracy and hence, very wide applications in flow discharge computations of
compound channels (Abril and Knight 2004; Unal et al. 2010). However, the above
mentioned approaches are not straightforward to be applied by hydraulic engineers
and also may suffer from long-time computations. Furthermore, efficient solution
of some of these methods mainly depends to numerical solution of differential
equations. For simplifying the computations of conveyance capacity, in this section,
GEP is used for prediction of flow discharge in compound channels.
For training and testing the proposed GEP equation in this research, 394 data sets of
flow hydraulic parameters from 30 different straight laboratory and river compound
sections were selected. Most of these data are gathered form an experimental
program undertaken by HR Wallingford (FCF) in large scale compound channel
flumes (Knight and Sellin 1987). In addition, some extra laboratory data from other
studies were used (Blalock and Sturm 1981; Knight and Demetriou 1983; Lambert
and Sellin 1996; Myers and Lyness 1997; Lambert and Myers 1998; Bousmar
and Zech 1999; Haidera and Valentine 2002; Guan 2003; Lai and Bessaih 2004;
Bousmar et al. 2004). Field data were collected from natural compound rivers of
River Severn at Montford Bridge (Ackers 1992; Knight et al. 1989), River Main
(Martin and Myers 1991) and Rio Colorado (Tarrab and Weber 2004). A typical
geometry for natural compound section having inclined berms is seen in Fig. 4.6.
The domains of main parameters of flow hydraulics and cross section geometry of
compound channels used in this book chapter are mentioned in Table 4.3.
4 Application of Gene-Expression Programming in Hydraulic Engineering 85
bf H
h
bc
Fig. 4.6 Typical river compound channel cross section with berm inclination
For developing a precise explicit equation to obtain total flow discharge in com-
pound channels, GEP has been used. Through following equation, it is assumed that
total flow discharge in compound channels is proportional to three dimensionless
parameters:
Qt QVDCM
D f Dr; COH; (4.24)
Qb Qb
Where Qt is total flow discharge and Qb is bank-full discharge. Of the total data set,
approximately 70 % (272 sets) were selected randomly and used for training. The
remaining 30 % (112 sets) were considered for testing.
Using optimization procedure, following relationship has been obtained for
training data:
Dr
.1COH/2
Qt Dr QVDCM .1COH/
D 3:954Dr 0:457DrC.1COH/ C
Qb Qb
QVDCM
2:462
Dr.1COH/
Qb
.1 COH/Dr (4.25)
In GEP, to avoid overgrowing programs, the maximum size of the program is gen-
erally restricted (Brameier and Banzhaf 2001). This configuration was tested for the
86 A. Zahiri et al.
Predicted Qt/Qb
VDCM
GEP-Train
GEP-Test
Perfect Line
Observed Qt/Qb
Fig. 4.7 Comparison of traditional vertical divided channel method and GEP model for discharge
ratios
proposed GEP model and was found to be sufficient. The best individual (program)
of a trained GEP can be converted into a functional representation by successive
replacements of variables starting with the last effective instruction (Oltean and
Grosan 2003). For developing Eq. (4.25), beside to the basic arithmetic operators
and mathematical functions (C, , , power), a large number of generations (5000)
were used for testing. First, the maximum size of each program was specified as 256,
starting with 64 instructions for the initial program. Table 4.4 shows the operational
parameters and functional set used in the GEP modelling.
The computed discharge ratios (Qt /Qb ) resulted from the GEP model for both
training and testing data as well as the vertical divided method are presented in
Fig. 4.7. It is clearly seen that GEP model in all variable ranges of selected data
(laboratory and field compound sections), has very promised accuracy. Based on
these prediction results, the mean absolute errors of discharge ratios for VDCM and
GEP model have been calculated as 55.2 and 8.5 %, respectively. It indicates that the
GEP model (Eq. 4.25) is highly satisfactory for total flow discharge in compound
open channels.
4 Application of Gene-Expression Programming in Hydraulic Engineering 87
Main channel
PLAN
ym yb
a
b
SECTON
Side sluice gates are underflow diversion structures placed along channels for
spilling part of the liquid through it (Fig. 4.8). These structures are mainly used
in irrigation, land drainage, urban sewage system, sanitary engineering, storm relief
and as head regulators of distributaries (Ghodsian 2003).
The flow through a side sluice gate is a typical case known as spatially varied
flow with decreasing discharge. By considering flow through side sluice gate as an
orifice flow, the flow discharge through side sluice gate under free flow condition
may be written as (Mostkow 1957):
p
Qs D Cd ab 2gym (4.26)
where Cd is discharge coefficient, a is opening height of the side gate, b is the gate
length and ym is upstream flow depth in the main channel.
Review of the literature shows that in spite of the importance of the side sluice
gates, relatively little attention has been given to studying the behavior of flow
through this structure (Panda 1981; Swamee et al. 1993; Ojah and Damireddy 1997;
Ghodsian 2003; Azamathulla et al. 2012). They related the discharge coefficient
of side sluice gates to the approach Froude number and ratio of flow depth to
gate opening. In the recent work, Azamathulla et al. (2012) used GEP technique
for developing a relationship for the discharge coefficient for the computation of
discharge through sluice gates.
88 A. Zahiri et al.
The experimental data of Ghodsian (2003) were used in Azamathulla et al. (2012)
study. The experiments were restricted to subcritical flow in main and side channel.
The range of various parameters used in this study is given in Table 4.5.
Basic arithmetic operators (C, , , /) as well as main basic trigonometric
and mathematical functions (sin, cos, tan, log, power) were used for GEP equation
development. Furthermore, a large number of generations (5000) were tested. The
functional set and operational parameters used in side slice gate flow hydraulic
modelling with GEP during this study are listed in Table 4.6.
The GEP model presented by Azamathulla et al. (2012) for estimating the
discharge coefficient of side sluice gates for free flow conditions is as follows:
h n ym oi
Cd D 0:12574tan1 tan1 sin log C Fr
a
h n y oi
1 m
0:1293 tan sin cos sin
a
" !#
1 Fr
C tan cos 1=3 (4.27)
8:54
Fr
C yam
The correlation coefficient (R) and mean square error (MSE) for model training
(60 data) are 0.976 and 0.0012, respectively, while for testing phase (14 data) are
0.967 and 0.0043, respectively. The performance of the GEP model is shown in
4 Application of Gene-Expression Programming in Hydraulic Engineering 89
1.0
0.9
e
lin
r
ro
er
%
0.8
e
+5
Perfect agreement line
li n
Cd (Computed)
r
ro
er
%
0.7
-5
0.6
0.5
0.4
0.3
0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0
Cd (Observed)
Fig. 4.9 Comparison of computed Cd using GEP with observed ones for training data
Figs. 4.9 and 4.10 for training and testing data, respectively. As can be seen, for both
phases, the computed Cd is within ˙5 % of the observed ones. The mean absolute
percentage error of the computed discharge coefficient by proposed GEP model is
about 2.15. It should be noted that although GEP models are somewhat complicated
in the mathematical form, but are easy to use practically by engineers at the field by
aid of available tools (e.g. spreadsheets).
Bed sills are a common solution to stabilize degrading bed rivers and channels. They
are aimed at preventing excessive channel-bed degradation in alluvial channels by
dividing them into partitions (Zahiri et al. 2014). For practical purposes, designers
and civil engineers are often interested in a short-term local scouring and its extent in
downstream of grade control structures. By this local scour, the structure itself (and
many times other structures in vicinity of it, like bridge piers or abutments, or bank
revetments) might be undermined (Bormann and Julien 1991; Gaudio and Marion
2003). Therefore, most researchers have focused on local scouring at isolated or
series bed sill structures. Summaries of research for the bed sills can be found
in Lenzi et al. (2002). Most of the studies on scouring at bed sills have been
conducted through experimental works (Bormann and Julien 1991; Gaudio et al.
90 A. Zahiri et al.
0.65
0.60
0.50
0.45
e
in
rl
e
ro
in
er
rl
ro
%
er
0.40
+5
%
-5
0.35
0.35 0.40 0.45 0.50 0.55 0.60 0.65
Cd (Observed)
Fig. 4.10 Comparison of computed Cd using GEP with observed ones for testing data
2000; Lenzi et al. 2002, 2003; Lenzi and Comiti 2003; Marion et al. 2004; Tregnaghi
2008; Chinnarasri and Kositgittiwong 2008). In general, in laboratory works a non-
linear regression equation is proposed based on curve fitting of experimental scour
depth data and hydraulic quantities and sediment properties. Some well-known
empirical equations based on regression analysis of experimental data have been
presented in Table 4.7.
These regression equations have one key limitation which mainly originate from
the wide ranges of hydraulic and sediment characteristics of flow in rivers. Owing
to rapid increase in successful applications of artificial intelligence techniques, it is
interesting to explore the applicability of the GEP in prediction of maximum scour
depth at bed sills.
In this section, using the 226 experimental data set of maximum scour depth at
bed sills from literatures in different canal bed slopes, applicability of GEP has been
examined in prediction of relative maximum scour depth at bed sills. These data
4 Application of Gene-Expression Programming in Hydraulic Engineering 91
are collected from Lenzi et al. (2002), Gaudio and Marion (2003), Marion et al.
(2004), Tregnaghi (2008) and Chinnarasri and Kositgittiwong (2006, 2008). Range
of variations as well as the mean values of important flow hydraulic and sediment
characteristics of experimental data are shown in Table 4.8.
Based on dimensional analysis of scour depth downstream bed sills, one can select
the parameters of a/Hs , a/
d50 , L/Hs , d50 /Hs and S0 as input variables and ys /Hs as
output variable. Table 4.9 reports the ranges of input and output parameters, used in
this section.
According to training data (174 data), an explicit equation has been developed GEP
technique as following:
92 A. Zahiri et al.
Hs a Energy line
Bed Sill a
S0
1 h
Maximum
Scour ys
Depth Seq
1
Bed Sill
ld
ls
L
Fig. 4.11 Schematic of scour depth and length downstream of a bed sill (Tregnaghi 2008)
ys a a d50 d50 a S0
D Ln C C 9:8561 C
Hs
d50 Hs Hs Hs
d50 d50 =Hs
d50 =Hs (4.32)
C
Log .a=Hs/ .a= .
d50 //1=3
The prediction results of bed sill scour depth for training and testing (52 data) GEP
model have been showed in Fig. 4.12. Comparison of the GEP model with the
empirical equations of scour depth at bed sills (Eqs. 4.28–4.30) are presented in
Fig. 4.13.
The detailed information of GEP model as well as the empirical equations is
indicated in Table 4.10. Based on this statistical analysis, it is indicated that among
different models considered in this study, Eq. (4.28) (Lenzi et al. 2004) has the
highest errors and therefore, doesn’t recommended for application. On the other
hand, GEP model can be proposed as an option for prediction of maximum scour
depth at bed sills. In addition, the simple equation of Chinnarasri and Kositgittiwong
(2008), with requiring to only two parameters and also having good accuracy, is may
be considered as a suitable approach.
4 Application of Gene-Expression Programming in Hydraulic Engineering 93
12 GEP-
Train
9
Predicted ys/Hs
0
0 3 6 9 12
Measured ys/Hs
Fig. 4.12 GEP model results of relative maximum scour depth for training and testing data
12
Eq.
10
Predicted ys/Hs
0
0 2 4 6 8 10 12
Measured ys/Hs
Fig. 4.13 Comparison of GEP model and well-known empirical equations for prediction of
relative maximum scour depth at bed sills
Table 4.10 Evaluation of empirical equations and GEP model for bed sill scour
depth prediction
Training Testing All data
Method R MSE ¡ R MSE ¡ R MSE ¡
Empirical eqs.
Eq. (4.28) – – – – – – 0.780 402 3.35
Eq. (4.29) – – – – – – 0.956 0.559 0.118
Eq. (4.30) – – – – – – 0.925 1.667 0.200
GEP model 0.976 0.203 0.067 0.986 0.308 0.013 0.979 0.286 0.081
94 A. Zahiri et al.
4.4 Conclusions
References
Abril, J.B., and Knight, D.W. 2004. Stage-discharge prediction for rivers in flood applying a depth-
averaged model. J. Hydraul. Res., IAHR, 122(6): 616–629.
Ackers, P. 1992. Hydraulic design of two-stage channels. Journal of Water and Maritime
Engineering, 96, 247–257.
Ackers, P. 1993. Flow formulae for straight two-stage channels. J. Hydraul. Res., IAHR, 31(4):
509–531.
Atabay, S. 2001. Stage-discharge, resistance and sediment transport relationships for flow in
straight compound channels. PhD thesis, the University of Birmingham, UK.
Azamathulla, H.Md., and Jarrett, R. 2013. Use of Gene-Expression programming to estimate
manning’s roughness coefficient for high gradient streams. Water Resources Management,
27(3): 715–729.
Azamathulla, H.Md., Ahmad, Z., and AbGhani, A. 2012. Computation of discharge through
side sluice gate using gene-expression programming. Irrig. and Drain. Engrg., ASCE, 62(1):
115–119.
Azamathulla, H.Md., and Zahiri, A. 2012. Flow discharge prediction in compound channels using
linear genetic programming. J. Hydrol. Ser. C, 454(455):203–207.
Azamathulla, H.Md., AbGhani, A., Leow, C.S., Chang, C.K., and Zakaria, N.A. 2011. Gene-
expression programming for the development of a stage-discharge curve of the Pahang River.
Water Resources Management, 25(11): 2901–2916.
Azamathulla, H.Md., AbGhani, A., Zakaria, N.A., and Guven, A. 2010. Genetic programming to
predict bridge pier scour. J. Hydraul. Eng., ASCE, 136(3): 165–169.
Blalock, M.E. and Sturm, T.W. 1981. Minimum specific energy in compound channel.
J. Hydraulics Division, ASCE, 107: 699–717.
Bormann, N., and Julien, P.Y. 1991. Scour downstream of grade-control structures, J. Hydraul.
Eng., ASCE, 117(5): 579–594.
Bousmar, D. 2002. Flow modelling in compound channels: momentum transfer between main
channel and prismatic or non-prismatic floodplains. PhD Dissertation, Université catholique de
Louvain, Belgium. http://dial.academielouvain.be/handle/boreal:4996
Bousmar, D. and Zech, Y. 1999. Momentum transfer for practical flow computation in compound
channels. J. Hydraul. Eng., ASCE, 125(7), 696–70.
Bousmar, D., Wilkin, N., Jacquemart, H. and Zech, Y. 2004. Overbank flow in symmetrically
narrowing floodplains. J. Hydraul. Eng., ASCE, 130(4), 305–312.
4 Application of Gene-Expression Programming in Hydraulic Engineering 95
Brameier, M., and Banzhaf, W. 2001. A comparison of linear genetic programming and neural
networks in medical data mining. IEEE Trans. Evol. Comput., 5: 17–26.
Chinnarasri, C., and Kositgittiwong, D. 2006. Experimental investigation of local scour down-
stream of bed sills. J. Research in Engineering and Technology, Kasetsart University, 3(2):131–
143.
Chinnarasri, C., and Kositgittiwong, D. 2008. Laboratory study of maximum scour depth down-
stream of sills. ICE Water Manage., 161(5): 267–275.
Chow, V. T. 1959. Open channel hydraulics, McGraw-Hill, London.
Fernandes, J.N., Leal, J.B. and Cardoso, A.H. 2012. Analysis of flow characteristics in a compound
channel: comparison between experimental data and 1-D numerical simulations. Proceedings
of the 10th Urban Environment Symposium, 19: 249–262.
Ferreira, C. 2001. Gene expression programming: a new adaptive algorithm for solving problems.
Complex Systems, 13(2): 87–129.
Gandomi A.H., and Alavi, A.H., 2011. Multi-stage genetic programming: A new strategy to
nonlinear system modeling. Information Sciences, 181(23): 5227–5239.
Gandomi A.H., and Alavi, A.H., 2012. A new multi-gene genetic programming approach to
nonlinear system modeling. Part II: Geotechnical and Earthquake Engineering Problems.
Neural Computing and Applications, 21(1): 189–201.
Gandomi, A.H., and Roke, D.A. 2013. Intelligent formulation of structural engineering systems.
In: Seventh M.I.T. Conference on Computational Fluid and Solid Mechanics, Massachusetts
Institute of Technology, Cambridge, MA.
Gaudio, R., and Marion, A. 2003. Time evolution of scouring downstream of bed sills. J. Hydraul.
Res., IAHR, 41(3): 271–284.
Gaudio, R., Marion, A., and Bovolin, V. 2000. Morphological effects of bed sills in degrading
rivers. J. Hydraul. Res., IAHR, 38(2): 89–96.
Ghodsian, M. 2003. Flow through side sluice gate. J. Irrig. and Drain. Engrg., ASCE, 129(6):
458–463.
Guven, A., and Azamathulla, H.Md. 2012. Gene-Expression programming for flip-bucket spillway
scour. Water Science & Technology, 65(11): 1982–1987.
Guven, A., and Gunal, M. 2008. Genetic programming approach for prediction of local scour
downstream of hydraulic structures. J. Irrigation and Drainage Engineering, 134(2): 241–249.
Guan, Y. 2003. Simulation of dispersion in compound channels. M.Sc. Thesis in Civil Engineering,
Ecole Polytechnique Federale de Lausanne (EPFL), Switzerland.
Guven, A., Aytek A., 2009. New approach for stage discharge relationship: gene expression
programming. J. Hydrol. Eng., ASCE, 14: 812–820.
Haidera, M.A. and Valentine, E.M. 2002. A practical method for predicting the total discharge
in mobile and rigid boundary compound channels. International Conference on Fluvial
Hydraulics, Belgium, 153–160.
Huthoff, F. 2007. Modeling hydraulic resistance of floodplain vegetation. PhD Thesis, Twente
University, the Netherland.
Huthoff, F., Roose, P.C., Augustijn, D.C.M., Hulscher, S.J.M.H. 2008. Interacting divided channel
method for compound channel flow. J. Hydraul. Eng., ASCE, 134(8): 1158–1165.
Jarrett, R.D. 1984. Hydraulics of high gradient streams. J. Hydraul. Eng., ASCE, 110(1):1519–
1539.
Knight, D.W. and Sellin, R.H.J. 1987. The SERC flood channel facility. Journal of Institution of
Water and Environment Management, 1(2), 198–204.
Knight, D.W. and Demetriou. J.D. 1983. Flood plain and main channel flow interaction.
J. Hydraulic Division, ASCE, 109(8): 1073–1092.
Knight, D.W., Shiono, K., and Pirt, J. 1989. Predictions of depth mean velocity and discharge in
natural rivers with overbank flow. International Conference on Hydraulics and Environmental
Modeling of Coastal, Estuarine and River Waters, England, 419–428.
Lai, S.H. and Bessaih, N. 2004. Flow in compound channels. 1st International Conference on
Managing Rivers in the 21st Century, Malaysia, 275–280.
96 A. Zahiri et al.
Lambert, M.F. and Myers, R.C. 1998. Estimating the discharge capacity in straight compound
channels. Water, Maritime and Energy, 130, 84–94.
Lambert, M.F. and Sellin, R.H.J. 1996. Discharge prediction in straight compound channels using
the mixing length concept. Journal of Hydraulic Research, IAHR, 34, 381–394.
Lenzi, M. A., Marion, A., Comiti, F., and Gaudio, R. 2002. Local scouring in low and high gradient
streams at bed sills. J. Hydraul. Res., IAHR, 40(6): 731–739.
Lenzi, M.A., and Comiti, F. 2003. Local scouring and morphological adjustments in steep channels
with check-dams sequences. Geomorphology, 55: 97–109.
Lenzi, M.A., Marion, A., and Comiti, F. 2003. Interference processes on scouring at bed sills. Earth
Surface Processes and Landforms, 28(1): 99–110.
Lenzi, M.A., Comiti, F., and Marion, A. 2004. Local scouring at bed sills in a mountain river:
Plima river, Italian alps. J. Hydraul. Eng., ASCE, 130(3):267–269.
Liggett, J.A. 2002. What is hydraulic engineering? J. Hydraul. Eng., ASCE, 128(1): 10–19.
Marion, A., Lenzi, M.A., and Comiti, F. 2004. Effect of sill spacing and sediment size grading on
scouring at grade-control structures. Earth Surface Processes and Landforms, 29(8): 983–993.
Martin, L. A. and Myers, R. C. 1991. Measurement of overbank flow in a compound river channel.
Journal of Institution of Water and Environment Management, 645–657.
Mostkow, M. A. 1957. A theoretical study of bottom type water intake. La Houille Blanche, 4:
570–580.
Moussa, Y.A.M. 2013. Modeling of local scour depth downstream hydraulic structures in
trapezoidal channel using GEP and ANNs. Ain Shams Engineering Journal, 4(4): 717–722.
Mujahid, K., Azamathulla, H.Md., Tufail, M. and AbGhani, A. 2012. Bridge pier scour prediction
by gene expression programming. Proceedings of the Institution of Civil Engineers Water
Management, 165 (WM9):481–493.
Myers, R.C. and Lyness. J.F. 1997. Discharge ratios in smooth and rough compound channels.
J. Hydraul. Eng., ASCE, 123(3): 182–188.
Najafzadeh, M., Barani, Gh-A., Hessami Kermani, M.R., 2013. GMDH Network Based Back
Propagation Algorithm to Predict Abutment Scour in Cohesive Soils. Ocean Engineering, 59:
100–106.
Ojah, C.S.P., Damireddy, S. 1997. Analysis of flow through lateral slot. J. Irrig. and Drain. Engrg.,
ASCE, 123(5): 402–405.
Oltean, M., and Grosan, C. 2003. A comparison of several linear genetic programming techniques.
Complex Syst., 14(1):1–29.
Onen, F. 2014. Prediction of scour at a side-weir with GEP, ANN and regression models. Arab J.
Sci Eng, 39: 6031–6041.
Panda, S. 1981. Characteristics of side sluice flow. ME thesis, University of Roorkee, Roorkee,
India.
Prasuhn, A. 1987. Fundamentals of hydraulic engineering. Holt, Rinehart, and Winston, New York,
USA.
Roberson, J.A., Cassidy, J.J., and Chaudhry, M.H. 1998. Hydraulic Engineering. John Wiley and
Sons, Inc., New York, USA.
Sattar, M.A. 2014. Gene expression models for the prediction of longitudinal dispersion coeffi-
cients in transitional and turbulent pipe flow. J. Pipeline Systems Engineering and Practice,
ASCE, 5(1):.
Sharifi, S. 2009. Application of evolutionary computation to open channel flow modeling. PhD
Thesis in Civil Engineering, University of Birmingham, 330 p.
Shiono, K. and Knight, D.W. 1991. Turbulent open-channel flows with variable depth across the
channel. J. Fluid Mechanics, 222: 617–646.
Swamee, P.K., Pathak, S.K., Sabzeh Ali, M. 1993. Analysis of rectangular side sluice gate. J. Irrig.
and Drain. Engrg., ASCE, 119(6): 1026–1035.
Tarrab, L., and Weber, J.F. 2004. Transverse mixing coefficient prediction in natural channels.
Computational Mechanics, 13: 1343–1355.
Tregnaghi, M. 2008. Local scouring at bed sills under steady and unsteady conditions. PhD Thesis,
University of Padova, Italy.
4 Application of Gene-Expression Programming in Hydraulic Engineering 97
Unal, B., Mamak, M., Seckin, G., and Cobaner, M. 2010. Comparison of an ANN approach with 1-
D and 2-D methods for estimating discharge capacity of straight compound channels. Advances
in Engineering Software, 41: 120–129.
Van Prooijen, B.C., Battjes, J.A. And Uijttewaal, W.S.J. 2005. Momentum exchange in straight
uniform compound channel flow. J. Hydraul. Eng., ASCE, 131:175–183.
Wang, C.Y., Shih, H.P, Hong, J.H. and Rajkumar, V.R. 2013. Prediction of bridge pier scour using
genetic programming. J. Marine Science and Technology, 21(4), 483–492.
Whitley, D. 1991. Fundamental principles of description in genetic search. Foundations of genetic
algorithms. G. Rawlings, ed., Morgan Kaufmann.
Wormleaton, P.R. and Merrett, D.J. 1990. An improved method of calculation for steady uniform
flow in prismatic main channel/floodplain sections. J. Hydraul. Res., IAHR, 28: 157–174.
Yang, K., Liu, X., Cao, S., and Huang, E. 2014. Stage-discharge prediction in compound channels.
J. Hydraul. Eng., ASCE, 140(4).
Zahiri, A., Azamathulla, H.Md., and Ghorbani, Kh. 2014. Prediction of local scour depth
downstream of bed sills using soft computing models. T. Islam et al. (eds.), Computational
Intelligence Techniques in Earth and Environmental Sciences, Springer, Chapter 11, 197–208.
Zakaria, N.A, Azamathulla, H.Md, Chang, C.K and AbGhani, A. 2010. Gene-Expression program-
ming for total bed material load estimation – A Case Study. Science of the Total Environment,
408(21): 5078–5085.
Chapter 5
Genetic Programming Applications
in Chemical Sciences and Engineering
5.1 Introduction
networks (ANNs), fuzzy logic (FL), evolutionary algorithms (EA) and support
vector machines/regression (SVM/SVR).
There exists a novel member of the evolutionary algorithms family, namely
genetic programming (GP) (Koza 1992), that addresses the above-stated challenge
by providing a method for automatically creating a computer program that performs
a prespecified task simply from a high-level statement of the problem. Genetic
programming follows Darwin’s theory of biological evolution comprising “survival
of the fittest” and “genetic propagation of characteristics” principles. It addresses
the goal of automatic generation of computer programs by: (i) genetically breeding
a random population of computer programs, and (ii) iteratively transforming the
population into a new generation of computer programs by applying analogs of
nature-inspired genetic operations, namely, selection, crossover and mutation.
The operating mechanisms of GP are similar to that of the genetic algorithms
(GA) (Goldberg 1989; Holland 1975). Though both these formalisms use the same
evolutionary principles, their application domains are very different; while GA
searches and optimizes the decision variables that would maximize/minimize a
specified objective function, GP automatically generates computer codes perform-
ing prespecified tasks. In addition to generating computer programs automatically,
there exist two important data-mining applications, namely, classification and sym-
bolic regression, for which GP has been found to be a suitable methodology. Unlike
the “divide and conquer” approach employed by machine learning algorithms to
perform classification, an evolutionary algorithm such as GP does not directly
construct a solution to a problem (e.g., a decision tree) but rather searches for a
solution in a space of possible solutions (Eggermont et al. 2004). The GP-based
symbolic regression (GPSR) is an extension of the genetic model of learning into the
space of function identification. Here, members of the population are not computer
programs but they represent mathematical models/expressions coded appropriately
using symbols.
As compared to classification, GP has been used extensively to conduct symbolic
regression in chemical sciences and engineering. GPSR possesses several advan-
tages over the two widely employed strategies namely artificial neural networks
(ANNs) and support vector regression (SVR) in developing exclusively data-driven
models. ANNs and SVR construct models in terms of a non-linear transfer function
and a kernel function, respectively. Depending upon the specific application for
which an ANN (SVR) model is being developed and the nature of the nonlinearities
between the corresponding input and output data, the complexity of the model
differs. However, owing to the use of the transfer (kernel) function, the basic
building blocks of the data-driven models fitted by the ANN (SVR) strategy remain
the same irrespective of their application domains. In contrast, GPSR provides a
system-specific linear or a nonlinear model that fits the given input–output data and
that too without making any assumptions regarding the form of the fitting function
(Kotanchek 2006). This is a remarkable feature of GP, which makes it a novel,
ingenious and an effective data-driven modeling formalism. The GP models are
also more compact and utilize less number of parameters than the existing classical
statistical techniques. Some of the comparative studies have indicated GP to be
superior in terms of accuracy of prediction than ANNs (Can and Heavy 2012).
5 Genetic Programming Applications in Chemical Sciences and Engineering 101
Moreover, GP models may enable the user to gain an insight into the fundamental
mechanisms (first principles) underlying the data. In a noteworthy study by Schmidt
and Lipson (2009) GP has been demonstrated to yield a phenomenological model
(natural law) governing the dynamics of a pendulum.
Despite its novelty and potential, GP—unlike ANNs and SVR—has not wit-
nessed an explosive growth for data-driven modeling applications. One possible
reason behind this scenario is that for a long time feature-rich, user-friendly and
efficient GP software packages were not available. The situation has changed in
recent years and a few software packages, both commercial and open source, have
become available for performing GP-based classification and symbolic regression.
These packages have definitely assisted in the development of a large number of
diverse GP applications in various science, engineering and technology disciplines.
In this chapter, GP-based classification and regression applications in chemical
sciences including biochemical sciences and chemical engineering/technology are
reviewed. Owing to the predominance of GPSR over GP-based classification, the
implementation details of the former are presented in greater depth. For an in-
depth generic treatment of the GP-based classification the reader is referred to, for
example, Koza (1991), Bonet and Geffner (1999), Cantú-Paz and Kamath (2003),
Eggermont et al. (2004), and Espejo et al. (2010).
Hereafter, this chapter is structured as follows. Section 5.2 provides a detailed
discussion of symbolic regression, issues involved in conducting GPSR, the step-
wise procedure of GPSR and a short list of GP software packages. In Sect. 5.3,
a review of classification and regression applications of GP in various sub-areas
of chemistry is provided. The specific GP application areas covered in this section
include drug design, environmental chemistry, green technologies, analytical chem-
istry, polymer chemistry, biological chemistry, and proteomics. Section 5.4, presents
GP applications in chemical engineering and technology. Here, the specific GP
application areas that are considered comprise process modeling, energy and fuels,
membrane technology, petroleum processes and heat transfer. Finally, Sect. 5.5
provides concluding remarks.
nonlinear regression analysis, even after expending a major effort in exploring the
function space there is no guarantee that a well-fitting model can indeed be secured
in a finite number of trials.
Most of the above-stated difficulties are overcome by the symbolic regression
(SR). It essentially involves function identification wherein a mathematical model/-
expression coded in a symbolic form is found in a manner such that the model
and the associated parameters provide a good, best, or a perfect fit between a
given finite sampling of values of the independent variables (model inputs) and
the corresponding values of the dependent variable (model output). Notably, SR
does this without making any assumptions about the structure of that model. That
is, it finds an appropriate linear or nonlinear form of the function that fits the data.
A similarity between the original “automatic development of a computer program
doing a specified job” and “symbolic regression” applications of GP is that both
methods take the values of the independent variables as input and produce the
values of the dependent variables as output (Koza 1990). Symbolic regression was
one of the earliest applications of GP (Koza 1992), and continues to be widely
studied (see, for example, Koza and Poli 2005; Poli et al. 2008; Iba et al. 2010;
Keedwell and Narayanan 2005; Sumathi and Surekha 2010; Cartwright 2008; Cai
et al. 2006; Gustafson et al. 2005; Lew et al. 2006). The major drawback of GPSR,
however, is that it is computationally intensive since it searches wide function and
associated parameter spaces. This however does not pose a major difficulty since
GPSR procedure is amenable to parallel processing.
Consider a multiple input—single output (MISO) example data set, D D
f(x1 , y1 ), (x2 , y2 ), : : : , (xN , yN )g, consisting of N patterns, where xn (n D 1, 2, : : : ,N)
denotes an M-dimensional vector of inputs (xn D [xn1 , xn2, : : : , xnM ]T ), and yn
denotes the corresponding scalar output. Using the data set D, the task of GPSR
is to search and optimize the exact form and the associated parameters of that
unknown MISO linear/nonlinear function (f ), which for the given set of inputs
produces the corresponding outputs as closely as possible. The general form of the
function/model to be fitted by GPSR is given as:
y D f .x; ’/ (5.1)
Fig. 5.1 Schematic of genetic programming: (a) basic tree structure, (b) random selection of
branches for reproduction, (c) crossover operation, and (d) mutation operation. Symbols in the
figure denote following operators (function nodes): (“C”) addition, (“”) subtraction, (“*”)
multiplication, (“”) division; xn; n D 1, 2, : : : , N, and numeric values define operands (terminal
nodes)
(x) and parameters (’). The trees in a population are of different sizes and their
maximum size is predefined. An illustrative tree structure representing an expression
“cos2 x1 C xx23 logx4 ” is depicted in Fig. 5.1. Upon forming the initial population
of candidate solutions, following steps are performed: evaluation of fitness scores
of candidate solutions, formation of a mating pool of parents, and actions of the
genetic operations, namely crossover and mutation. An iteration of these steps
produces a new generation of offspring candidate solutions. Several such iterations
are needed before convergence is achieved. The candidate solution possessing the
highest fitness score encountered during the iterative process is chosen as the best-
fitting model.
104 R. Vyas et al.
number of generations have been evolved, and (2) the fitness value of the best
candidate solution in a population no longer increases significantly or remains
constant over a large number of successive generations.
1. Fitness score evaluation: Using the example set inputs fxn g, n D 1,2, : : : N,
compute the output of each candidate solution in the current population. The
computed outputs are utilized to calculate the magnitude of the pre-selected error
metric (e.g., RMSE), which is then used to determine the fitness score of that
solution. This procedure is repeated to calculate fitness scores of all candidate
solutions in the current population.
2. Selection: Create a mating pool of candidate solutions termed “parents” to
undergo crossover operation described in step 3. The members of the mating
pool are selected in a manner such that only those candidate solutions possessing
relatively high fitness scores can enter the pool. There exist a number of methods
for selecting the candidate solutions in a mating pool, such as Roulette-wheel
selection, greedy over-selection, ranking selection, tournament selection and
elite strategy (Iba et al. 2010). Each one of these strategies possesses certain
advantages and limitations.
3. Crossover: This step can be performed multiple ways, for example, single- and
two-point crossover. In the former (see Fig. 5.1b), a pair of parent candidate
solutions is selected randomly from the mating pool and two new candidate
solutions (offspring) are created by slicing each parent tree at a random point
along the tree length and mutually exchanging and recombining the sliced parts
between the parents (see Fig. 5.1c). This crossover operation is conducted with a
pre-specified probability value (termed crossover probability) and repeated with
other randomly chosen parent pairs until Np offspring candidate solutions are
formed. The trees in GP do not have a fixed length and these can grow or shrink
due to the crossover operation.
4. Mutation: In this step, small changes are applied to the operator and operand
nodes of the offspring solutions to produce a new generation of candidate
solutions (see Fig. 5.1d). This step is performed with a small magnitude of the
probability termed “mutation probability.”
Avoiding over-fitting of models An important issue that needs to be addressed
during GPSR is over-fitting of the constructed models. Over-fitting can occur in two
ways, that is, when a model is trained over a large number of iterations (termed
overtraining), and/or the model contains more terms and parameters than necessary
(over-parameterization). Over-parameterization tends to increase the complexity of
the fitted model. Over-fitting results in a model that has learnt even the noise in the
data at the cost of capturing a smooth trend therein. Such a model performs poorly
at generalization, which refers to the model’s ability to accurately predict outputs
corresponding to a new set of inputs. A model incapable of generalization is of no
practical use. To overcome the problem of over-fitting, the available input–output
example set is partitioned into two sets namely training and test sets. While the
former is used to train the model, the test set is used for evaluating the generalization
capability of the model. After each training iteration or convergence the candidate
106 R. Vyas et al.
solutions are assessed for their generalization capability using the test set data
and those predicting the training and test set outputs with high and comparable
accuracies are accepted. Sometimes, a third set known as validation set is formed
from the available example data and a model performing well on all the three i.e.
training, test and validation sets is selected.
GP implementation is computationally very intensive and often the evolved solu-
tion is anything but ideal thus requiring even greater numerical processing to secure
an acceptable solution. The component of the GP algorithm that is computationally
most expensive is fitness evaluation. It is therefore at most important to use fitness
functions that are computationally economical yet efficient.
The most attractive feature of the GP-based symbolic regression is that it searches
and optimizes the “structure” (form) of a suitable linear or nonlinear data-fitting
function. It also obtains values of all the parameters of that function although
relatively this is a less important GP characteristic since several deterministic and
stochastic linear/nonlinear parameter estimation strategies are already available.
Moreover, optimality of the parameter values searched by the GP cannot be guar-
anteed. It is therefore advisable that the parameters of the GP-searched function are
optimized using an appropriate parameter optimization strategy such as Marquardt’s
method (Marquardt 1963).
• RGP (Flasch 2014) is a GP system based on, as well as fully integrated into, the
R environment. It implements classical tree-based GP as well as other variants
including, for example, strongly typed GP and Pareto GP. This package is flexible
enough to be applied in nearly all possible GP application areas, such as symbolic
regression, feature selection, automatic programming and general expression
search.
• ECJ toolkit (ECJ 22 2014) is one of the popular computational tool with full
support for GP. It is a evolutionary computation research system written in
Java and developed at George Mason University’s evolutionary computation
Laboratory. ECJ 22 is the latest version of the toolkit provided with a GUI and
various features such as flexible breeding architecture, differential evolution and
multiple tree forest representation. This toolkit is reviewed by White (2012).
Genetic programming has been used in chemistry with a great success for providing
potential solutions to a variety of classification and data-driven modeling problems
as well as to create new knowledge. It has been also established that GP models
can approximate to a first principles models (Anderson et al. 2000). While GA
has been used extensively for optimization in chemical sciences and chemical
engineering/technology, GP-based applications in these disciplines are relatively
fewer (Aguiar-Pulido et al. 2013). The applications of GP in chemical sciences
have focused mainly on data mining, which can be further broadly categorized into
rule-based classification and symbolic regression based model development (see
Fig. 5.2).
GP is an apt tool for data-mining in chemical sciences. Formally, data mining is
defined as “identification of patterns in large chunks of information” (Cabena et al.
1997). The data could be physicochemical data from small molecule based assays or
spectral data emanating from the analytical instruments used for characterization of
chemical or biological moieties. Several representations of the tree-based GP have
been used exhaustively for the rule-based data classification (Li and Wong 2004).
GP-based intelligent methodologies have been used for developing the rule-based
systems in chemistry and biochemistry domains such as, chemical networks and
reactivity, wherein the computer programs are all functional models of chemical or
biochemical properties (Tsakonas et al. 2004). GP-based regression has been mainly
employed for building quantitative structure—property relationship (QSPR) models
(Barmpalexis et al. 2011). The flexibility of GPSR is at the core of the development
of free form mathematical models from the observed data (Kotanchek 2006). Both
the approaches persist and consequently the subsequent sections are devoted to the
applications of the GPSR and GP-based classification methods, illustrated by using
copious examples from the literature.
108 R. Vyas et al.
and Barrett 2005). The results showed that classification of drugs into ‘high’ and
‘low’ OB classes could be performed on the basis of their molecular structure and
the output given by the developed models would be useful in the pharmaceutical
research. Moreover, the results indicated that the quantitative prediction of the oral
bioavailability is also possible. In a study involving structure-property relationships,
GP was utilized for the prediction of Caco-2 cell permeability (Vyas et al.
2014), which is an important ADMET parameter; the said GP model yielded high
coefficient of correlation (0.85) between the desired and model predicted values
of the permeability and a low RMSE value of 0.4.
Mathematical models have been developed to predict drug release profiles
(Ghosal et al. 2012; Güres et al. 2012). In a related study, GP-based models were
built for the prediction of drug release from the solid-lipid matrices (Costa and
Lobo 2001). Here, GP was used specifically for determining the parameters of
the model—a modified Weibull equation—that is commonly used in the reliability
engineering defined as:
ˇ1
T ˇ
ˇ T
f .T/ D e (5.2)
In the graph theory representation of molecules, atoms are denoted as nodes and
bonds connecting them are represented as edges. In a classic study by Globus et al.
(1999), the molecular design problem was viewed as a search of the space of all
molecules to find a target molecule with the desired properties. The GP fitness
function was used to automatically evolve chains and rings in new molecules by
using the crossover operator to divide trees into fragments (Globus et al. 1999). Ring
evolution was enabled by the mutation operator and the fitness function defined a
distance measure, i.e., Tanimoto coefficient. Likewise, GP has found a potential
use in systems that can be represented using the graph theory such as metabolic
pathways wherein the networks of organic transformations occurring in biological
systems can be represented as program trees (Ivanova and Lykidis 2009).
Whole cell biosensors have become an integral part of the environment monitoring
(Gu et al. 2004). Here, the main task is to detect the substance specific patterns from
the huge biosensor data being monitored continuously. GP has been found to be
a suitable classification technique to handle the stated task. For example, GP has
been employed in the classification of herbicide chemical classes and herbicides
with high sensitivity albeit with a low selectivity (Podola and Melkonian 2012).
Electronic noses are being employed as vapor sensors since they provide rich
information regarding the analyte binding (Persaud and Dodd 1982). GP-based
approaches were able to detect the airborne analytes in real time with a good
sensitivity as also selectivity (Wedge et al. 1999).
Gene expression programming (GEP) is an extension of the genetic programming
(GP) and genetic algorithms (GAs). It is a population-based evolutionary algo-
rithm (Ferreira 2001) wherein a mathematical function defined as a chromosome
consisting of multi-genes is developed using the data presented to it. In GEP, a
mathematical expressions are encoded as simple linear strings of a fixed-length,
which are subsequently expressed as nonlinear entities of different sizes and shapes
(i.e. simple diagram representations or expression trees) (Cevik 2007). Singh
and Gupta (2012) employed GEP for forecasting the formation trihalomethanes
(THMs)—which are toxic to human health—in waters subjected to chlorination.
In this study, five parameters namely dissolved organic carbon normalized chlorine
dose, water pH, temperature, bromide concentration, and contact time, were used as
model inputs. Similar to the GP-based model, ANN and SVM based models were
developed for comparison purposes. The results of this comparison revealed that
the ANN, SVM, and GEP models are capable of capturing the complex nonlinear
relationship between the water disinfection conditions and the corresponding THM
formation in the chlorinated water.
5 Genetic Programming Applications in Chemical Sciences and Engineering 111
1:397
HHV D 0:365 FC C 0:131 VM C
FC
328:568 VM
C (5.3)
10283:138 C 0:531 FC3 ASH 6:893 FC2 ASH
where FC, VM, and ASH are the weight percentages (dry basis) of fixed carbon,
volatile matter, and ash respectively.
• Ultimate analysis based optimal model:
53:883 O C H 115:971
HHV D 0:367 C C 2
C
2:131 C 93:299 10:472 H C 0:129 C O
91:531 232:698
C (5.4)
.35:299 C N/ 77:545 C S
HHVs were high (>0.95) while the corresponding magnitudes of mean absolute
percentage error (MAPE) were low (<4.5 %). The HHV prediction accuracy
and generalization performance of these models were rigorously compared with
the corresponding multilayer perceptron (MLP) neural network based as also
previously available high-performing linear and nonlinear HHV models. This
comparison showed that the HHV prediction and generalization performance
of the GP as also MLP-based models to be consistently better than that of
their linear and/or nonlinear counterparts proposed earlier. The biofuel HHV
prediction models proposed by Ghugare et al. (2014a), due to their excellent
performance, possess a potential of replacing the models proposed earlier. Also,
their GP-based strategy can be extended for developing HHV prediction models
for other types of fuels.
Among the two commonly employed analyses for characterizing biomass fuels,
proximate analysis is relatively easy to perform than the ultimate analysis. Accord-
ingly, GPSR was employed for building non-linear models for the accurate predic-
tion of C, H and O fractions of the solid biomass fuels from the constituents of
the corresponding proximate analysis (Ghugare et al. 2014b). These models were
constructed using a large data set of 830 fuels. For comparison purposes, C, H and
O prediction models were developed using ANN and SVR approaches also. The
results of the comparison of the prediction accuracy and generalization performance
of GP, ANN and SVR based nonlinear models with that of the currently available
linear models indicated that the nonlinear models have consistently and significantly
outperformed their linear counterparts.
GP has been applied for conducting multivariate analysis of the nonlinear dielectric
spectroscopy (NLDS) data of a yeast fermentation process (Woodward et al. 1999).
In this study, GP was found to outperform the conventional methods like partial least
squares (PLS) and ANNs. Genetic programming was also used for recognizing the
bonds taking part in increasing or decreasing the dominant excitation wavelength by
identifying the conjugated … systems for lowest UV transition for a system of 18
anthocyanidins (Alsberg et al. 2000). The model stressed upon the important role
of bond critical point (BCP) characterized by the electron density, the Laplacian
operator and the ellipticity.
The reactivity ratios in free radical copolymerization are routinely estimated using
Alfrey-Price (AP) model (Alfrey and Price 1947). However, the accuracy of
5 Genetic Programming Applications in Chemical Sciences and Engineering 113
Cancer is a major disease for which GP has found numerous applications ranging
from the classification models of cancer tumors to mechanistic understanding of the
underlying pathogenesis. Easy interpretability of these GP models greatly enhances
our understanding of the underlying cellular and disease pathway dynamics at
the systems biology level (Finley et al. 2014). Interested readers are referred to a
comprehensive review of GP applications in cancer research (Worzel et al. 2009).
As of today, for most of the neurodegenerative diseases there is no cure; however
an early detection can provide a better life for the patients suffering from them.
An example of the use of GP in the clinical time series data involves inducing
classifiers capable of recognizing the movements characteristic of patients afflicted
with a disease like Parkinson’s wherein a diagnostic accuracy of 97 % was achieved
(Castelli et al. 2014). Here, GP was used to identify patterns in the slow motor
movements (Bradykinesia) related clinical data. Similar applications are found in
the context of visuo-spatial diseases also where a graph based GP system termed
Implicit Context representation Cartesian Genetic Programming (IRCGP), which
functions similarly to the well known crossover operator was devised (Smith and
Lones 2009).
NMT plcC
Phosphomethyl
ethylamine 1,2-Diacylglycerol
NMT plcC
Phospho- Phosphatidyl-
Choline CDP-choline
choline choline
Fig. 5.4 A schematic representation of reactions involved in the human phospholipid pathway
last product of the predicted network model matched with a high accuracy with
the experimental data; GP could successfully construct two metabolic pathways
viz. phospholipids cycle and degradation of ketone bodies. The eukaryotic phos-
pholipids biosynthetic pathway and its role in the cellular biology has been well
studied (Vamce and Vance 2008) (see Fig. 5.4). Here, the researchers chose four
enzymatic reactions from the said pathway with glycerol and fatty acid as inputs
and diacyl-glycerol as the end product. A tree was constructed programmatically
to represent the chemical reaction functions and selector functions as nodes, and
reaction rates, substrates, products and enzymes, as leaves. The results of the GP run
could be corroborated with the observed experimental data. Thus, GP could create
metabolic pathways that included topological features such as an internal feedback
loop, bifurcation point, accumulation point and rates for all reactions using the time
domain concentration values.
Complex chemical processes are modeled using the input–output data from the
experimental tests. In one of the early significant contributions of GP, McKay et al.
(1997) developed data driven steady-state models for two processes, namely, a
binary vacuum distillation column and a chemical reactor system. For the vacuum
distillation unit a model was developed to infer the bottom product composition.
The vacuum distillation column was equipped with 48 trays, a steam reboiler
and a total condenser. The feed was split into two product streams i.e., the
distillate and bottoms. McKay et al. (1997) obtained the input–output data from a
phenomenological model of the column wherein a set consisting of 150 data points
of the steady-state composition estimates from three trays (numbered 12, 27, and
42 from the top) was considered along with the corresponding values of the bottom
composition. A set of 50 data points was used in the model validation. The models
were accepted only if the validation set RMSE was less than 0.02. An F-test was then
performed to find the best model. The overall best model had an RMSE of 0.011
on the training set data and 0.015 on the validation set data. McKay et al. (1997)
applied a similar method to obtain a functional relationship, the data for which was
generated from an assumed relationship with three inputs (u1 , u2 , and u3 ) and a
single output, y:
y D 1000 u1 exp .5=u2/ C u3 (5.5)
They also modeled a continuous stirred tank reactor (CSTR) system for the
prediction of the product composition.
In chemical processes, operating conditions need to be optimized for a variety
of reasons such as maximization of conversion, profit and selectivity of desirable
products, and minimization of cost and selectivity of undesirable products. For
conducting such an optimization, it is necessary that a representative and accurate
process model is available. Often, process behavior is nonlinear and complex and
therefore developing phenomenological (also termed “first principles” or “mech-
anistic”) process models becomes tedious, costly and difficult. In such instances,
data-driven process models can be constructed. Cheema et al. (2001) presented
GP-assisted stochastic optimization strategies for the optimization of glucose to
gluconic acid bioprocess wherein Aspergillus niger strain was used for producing
gluconic acid. Their study utilized two hybrid process modeling-optimization
approaches wherein a GP-based model was first developed from the process data,
following which the input space of the GP model was separately optimized using
two stochastic optimization (SO) formalisms, namely, genetic algorithms (GA)
and simultaneous perturbation stochastic approximation (SPSA) (Spall 1998). A
schematic of the GP-based process modeling and GA-based optimization strategy
is shown in Fig. 5.5. Cheema et al. (2002) used process data from 46 batch
fermentation experiments conducted by them in building the GP-based model. The
gluconic acid concentration (y)(g/L) which formed the output of the GP model
(see Eq. 5.6) was predicted as a function of three process parameters, namely,
5 Genetic Programming Applications in Chemical Sciences and Engineering 117
Compute
Fitness of
Candidate
New Generation of Solutions
Candidate Solutions
Return Best
Perform Convergence Yes
Crossover &
Candidate
Mutation
Achieved? Solution
No
Select Candidate Solutions
Optimize Input Space of
in a Mating Pool based
Gp-based Overall Best
on their Fitness Values
Process Model Using a
Suitable Optimization
Optimized Algorithm ( e.g. GA)
Process
Operating
Conditions
glucose concentration (x1 )(g/L), biomass concentration (x2 )(g/L), and dissolved
oxygen (x3 ) (mg/L).
ˇ1 x1 1 1
yD (5.6)
.x1 ˇ2 /4 C ˇ3 2
x2 ˇ4 x2 C ˇ5 2
ˇ6 x3 ˇ7 x3 C ˇ8
The values of the eight model parameters fitted by the GPSR are: “1 D 3.1911
1010 , “2 D 158.219, “3 D 2.974 106 , “4 D 5.421, “5 D 107.15, “6 D 0.116,
“7 D 12.752, and “8 D 448.112. The magnitude of the variance (R2 ) pertaining to
the training (test) set output predictions made using Eq. (5.6) was 0.987 (0.986).
Since the form of the model was known (determined by GP), Cheema et al.
(2002) subjected the GP-based model to Marquardt’s non-linear regression analysis
(Marquardt 1963) to explore whether the model’s eight parameters could be fine-
tuned further to improve its prediction accuracy. This parameter fine-tuning indeed
led to a better R2 value of 0.9984 (0. 9979) for the training (test) set. The GP based
118 R. Vyas et al.
where,
s
dp d
log2 Dpc
Dc dp p g 0:3
AD C (5.8)
H0 Dc g
Dc
s 0:8 !2
dp dp Di H0
BD C C (5.9)
Dc Dc Dc Dc
The results obtained from the above GP-based correlation are in good agreement
with the experimental values.
5 Genetic Programming Applications in Chemical Sciences and Engineering 119
Grosman and Lewin (2004) developed a GP model to predict the reaction rate of
toluene. The developed model showed that the GP-based model possesses a good
prediction accuracy (RMSE D 0.0038).
In semiconductor manufacturing, rapid thermal processing (RTP) has gained
importance in recent years. In the processes using RTP, the principal issue is
temperature regulation. Dassau et al. (2006) presented GP for the development of
steady-state and dynamic temperature control models. To improve RTP an NMPC
system was developed using GP, as done by Grosman and Lewin (2002). The models
were based on the mathematical representation of the Steag RTP system and these
were subsequently used in controlling the temperature at three different locations on
the wafer, namely, centre of the wafer, 5 cm from the centre, and edge of the wafer
(9.5 cm from the centre). The advantage of this NMPC system is that the process
approaches the set point easily.
Hinchcliffe and Willis (2003) developed GP-based dynamic process models in
two case studies—a system with time delay (Narendra and Parthasarathy 1990)
and a cooking extruder (Elsey et al. 1997)—which were used to compare the
120 R. Vyas et al.
performance of the GP algorithm with the filter based neural networks (FBNNs).
It was observed that the two approaches exhibit comparable performances although
GP has a potential advantage over FBNNs in that it can measure the performance
of the model during its development. In another study involving development
of GP based steady-state and dynamic input–output models, Willis et al. (1997)
constructed models for a vacuum distillation column and a twin screw cooking
extruder.
Coal gasification is more environment-friendly and efficient process for energy
generation than coal combustion. The performance of a coal gasification process
is significantly dependent on the quality of coal. For instance, process efficiency
is adversely affected by the high ash content in a coal. Such coals are found
in many countries, and form an important raw material for the coal gasification
and combustion. Despite the wide-spread availability of high ash coals, modeling
studies on gasification using these coals are much less in number when compared
with the studies performed using low ash coals. Coal gasification is a nonlinear
and complex process and its phenomenological modeling is a difficult, tedious
and expensive task. In such circumstances, data-driven modeling provides a low-
cost and relatively easier alternative to conduct process modeling if representative,
statistically well-distributed and sufficient process data are available. Accordingly,
Patil-Shinde et al. (2014) developed data-driven steady-state models for a pilot plant
scale fluidized bed coal gasifier (FBCG) utilizing Indian coals with a high ash
content. These models were constructed using process data from 36 experiments
conducted in the FBCG. Specifically, four models predicting gasification related
performance variables, namely, COCH2 generation rate, syngas production rate,
carbon conversion and heating value of syngas, were developed using GP and multi-
layer perceptron (MLP) neural network formalisms, separately. The input space of
these models consisted of eight coal and gasifier process related parameters, namely
fuel ratio, ash content of coal, specific surface area of coal, activation energy of
gasification, coal feed rate, gasifier bed temperature, ash discharge rate and air/coal
ratio. A comparison of the GP and MLP-based models revealed that their output
prediction accuracies and generalization performance vary from good to excellent
as indicated by the high training and test set correlation coefficient magnitudes lying
between 0.920 and 0.996.
Gandomi and Alavi (2011) developed a new strategy for non-linear system
modeling. They proposed a multistage genetic programming (MSGP) formulation
to provide accurate predictions by incorporating the individual effects of predictor
variables and interactions among them. The initial stage of MSGP formulates the
output variable in terms of an influencing variable. Thereafter, a new variable
is defined by considering the difference between the actual and predicted value.
Finally, an interaction term is derived by considering the difference between the
desired output values and those predicted by the individually developed terms.
Gandomi and Alavi (2011) applied this strategy to various engineering problems
such as simulation of pH neutralization process. The results yielded by the MSGP
strategy were observed to be more accurate than that by the standard GP.
5 Genetic Programming Applications in Chemical Sciences and Engineering 121
Oils, fuels, solvents, paints, detergents, organic matter, and rusts, are a few typical
contaminants present in the wastewater. With ever growing need for a high quality
water, the need for treating the waste-water has also increased. In recent years,
membrane technology has assumed an important role in the wastewater treatment.
A significant difficulty with this technology is fouling of the membrane, which
leads to a decline in the permeation flux. Accordingly, Lee et al. (2009) utilized
genetic programming for the prediction of membrane fouling in a microfiltration
(MF) system. The model was developed to predict the membrane fouling rate in a
pilot scale drinking water production system consisting of a hollow fiber membrane
of polyvinylidene fluoride (PVDF). The model was developed using the following
input variables (predictors): operating conditions (flow rate and filtration time) and
feed water quality (turbidity, temperature and algae pH). Lee et al. (2009) collected
data from a membrane filtration system for 470 days, during which chemical
washing was done three times. The operating conditions and water quality data used
in the model development were analyzed separately during the three system runs.
5 Genetic Programming Applications in Chemical Sciences and Engineering 123
The resultant GP model predicted accurately the membrane resistance and yielded
very low RMSE values (ranging from 0.06 to 0.12). Shokrkar et al. (2012) also
developed a GP model for an accurate quantification of the permeation flux decline
during cross-flow membrane filtration of oily wastewater. A total of 327 data points
covering the effects of individual variations in the five predictor variables, namely,
temperature, cross flow velocity, trans-membrane pressure, oil concentration, and
filtration time, were used in the model development. Simulations conducted using
the GP model suggested that by increasing the filtration time, oil layer on the
membrane surface thickens and permeation flux decreases. The permeation flux falls
rapidly in the beginning of the filtration. The GP-model predicted results agreed
within 95 % of the experimental data. This study clearly gives an idea of how each
parameter affects the flux.
Reverse osmosis (RO) has been proved to be a promising technology for desali-
nation. It has been observed that the performance of RO is negatively affected by the
formation of scales of soluble salts. Cho et al. (2010) developed a GP model to study
the effect of CaSO4 scale formation on the RO membrane. The extent of RO fouling
(permeation flux decline) is dependent on the applied pressure, time, volumetric
concentration factor (VCF), stirring speed, and humic acid concentration. After
training and validation, the correlation coefficient for the model predictions was
0.832. It was observed that the applied pressure and VCF have higher impact on
the RO fouling. Park et al. (2012) have also developed a GP model for the analysis
of the performance of an RO process. The input parameters considered by them
are pH, oxidation reduction potential (ORP), conductivity, temperature, flux, TMP,
and recovery time. The GP models were developed for trans-membrane pressure
and membrane permeability separately for an early stage data (0–5 days) and the
late stage data (20–24 days). The models matched the trend of the pilot plant data
well. The sensitivity analysis of the GP models showed that the model-fitted early
stage data exhibit higher sensitivity towards conductivity, flux and ORP, whereas
the late stage data show higher sensitivity to temperature, recovery time, and ORP.
It was also observed that the GP model predictions of membrane permeability
were more accurate than the predictions of transmembrane pressure. In a recent
study, Meighani et al. (2013) have conducted a thorough comparative analysis
of the three modeling techniques (pore-blocking model, ANN and GP) for the
prediction of permeate flux decline. The permeate flux is modeled as a function of
the transmembrane pressure, feed temperature, cross flow velocity, pH and filtration
time. Eight sets of experimental data were compiled from the literature to investigate
the accuracy of the models. The correlation coefficients for the pore blocking, ANN,
and GP models were found to be 0.9799, 0.9999 and 0.9723, respectively.
Okhovat and Mousavi (2012) predicted the performance of a nanofiltration
process for the removal of heavy metals such as arsenic, chromium and cadmium.
GP-based models were developed for studying the membrane rejection of arsenic,
chromium and cadmium ions. Specifically, ions rejection (%) was considered as
the model output while feed concentration and transmembrane pressure, formed the
model inputs. The models showed satisfactory prediction accuracies with RMSE
magnitudes ranging between 0.005 and 0.02. Suh et al. (2011) utilized GP for
124 R. Vyas et al.
estimating the membrane damage during the membrane integrity test of a silica
fluorescent nanoparticle microfiltration membrane. This model predicts the area
of membrane damage for the experimental input parameters (concentration of
fluorescent nanoparticles, permeate water flux, and transmembrane pressure). The
GP model yielded good prediction results with mean absolute error (MAE) of 0.83.
In one of its early applications, Greeff and Aldrich (1998) employed GP in various
leaching experiments as described below.
1. Acid pressure leaching of nickeliferous chromites: Based on the data by Das et al.
(1995) GP models were developed for the dissolution of nickel, cobalt and iron
from the beneficiated lateritic chromite samples as a function of temperature,
ammonium sulfate concentration, and acid concentration at different time inter-
vals.
126 R. Vyas et al.
2. Leaching of uranium and radium: Here, models were developed for the co-
extraction of uranium and radium from a high grade arseniferous uranium ore
(Kondos and Demopoulos 1993), where the percentage of radium and uranium
was modeled as a function of pulp density, concentration of hydrochloric
leaching acid, concentration of the calcium chloride additive and time. Cross-
validation (using leave-k out method, k D3) was performed to evaluate the
prediction performance of the models. All the GP-based models were found to
be comparable and significantly more accurate than those developed by means
of the standard least-squares methods. In another leaching related study, Biswas
et al. (2011) analyzed the leaching of manganese from low grade sources,
using genetic algorithms, GP and other strategies where they made an extensive
comparison of the data-driven modeling techniques.
Wang et al. (2008a) applied GP to a complex heat-integrated distillation system
to synthesize a flow-sheet for separating a multicomponent mixture into pure
components at a minimum total annual cost. Both sharp and non-sharp distillations
were considered in the modeling. Based on the knowledge of chemical engineering,
unique solution encoding methods and solution strategies were proposed in this
study. The GP-based synthesis algorithm automatically optimizes the problem of
complex distillation systems. In related studies, Wang et al. (2008b) and Wang
and Li (2008, 2010) made use of GP for the synthesis of non-sharp distillation
sequences, synthesis of multicomponent products separation sequences, and syn-
thesis of heat integrated non-sharp distillation.
where k denotes the discrete time, ykC1 refers to the one-time-step-ahead process
output, u is the manipulated variable, f refers to the linear/nonlinear functional
5 Genetic Programming Applications in Chemical Sciences and Engineering 127
relationship (to be identified) between ykC1 and the current (kth) and the lagged
values of the input–output variables, and m and n refer to the number of lags
in the output and input variables, respectively. Traditionally, system identification
employs statistical methods for constructing models of dynamical systems from
their experimental (measured) data consisting of yk and uk values. In one of the early
studies on the applications of GP technique in chemical engineering, Kulkarni et al.
(1999) utilized the methodology for system identification by conducting two case
studies involving nonlinear pH and heat exchanger control systems. The objective
in both case studies was to obtain an appropriate non-linear form of f given the time
series values of the process input (u) and the corresponding output (y). To derive GP-
based models, Kulkarni et al. (1999) used synthetic process data. In actual practice,
these data are collected by conducting open-loop tests wherein manipulated variable
(u) is varied randomly and its effect on y is monitored. The specific systems
considered in two case studies were: (1) continuous stirred tank reactor (CSTR)
wherein hydrochloric acid and sodium hydroxide streams are mixed and effluent
stream’s pH is measured and controlled using a model-based control strategy, and
(2) a nonlinear heat exchanger control system wherein heater voltage and exchanger
outlet temperature are the manipulated and controlled variables, respectively. The
CC values in respect of the training and test set output predictions by both the
models were greater than 0.99 indicating an excellent prediction and generalization
performance by the models identified by GP. In a similar work, Nandi et al. (2000)
performed GP-based system identification of a fluidizedcatalytic cracking (FCC)
unit, wherein an exothermic reaction (A ! B ! C takes place. Here, two
GP models each possessing an excellent prediction accuracy and generalization
capability were developed for the prediction of one-time-step-ahead and three-time-
steps-ahead outlet concentrations of species B.
Sankpal et al. (2001) utilized a GP-based model for the monitoring of a process
involving continuous production of gluconic acid by the fermentation of sucrose
and glucose solution in the presence of aspergillus niger immobilized on cellulose
fabric. During the continuous conversion of glucose and sucrose the rate of gluconic
acid formation drops as fermentation progresses. To compensate for this loss in
efficiency the residence time needs a suitable adjustment. As online determination
of the reaction rates and substrate concentration is cumbersome, a GP-based model
was developed to predict the o conversion (z) as a function of the time. For a given
time series fzt ; zt1 ; : : : ; zT , of length T, a model was developed to compute ztC1 ,
where ztC1 D f .zt ; zt1 ; : : : ; zt˛ / I ˛ t .T 1/ I 0 ˛ L; where t refers
to the discrete time, ˛ refers to the number of lags and L denotes the maximum
permissible lags. The expression for the one-time-step-ahead prediction gave high
prediction accuracies with correlation coefficient magnitudes 1.
Timely and efficient process fault detection and diagnosis (FDD) is of critical
importance since it helps in, for example, energy savings, reduction in operating and
maintenance costs, curbing damage to the equipment, avoiding economic losses due
to process down-time, and most importantly preventing mishaps and injuries to plant
personnel. Process identification and FDD are related since a good process model
128 R. Vyas et al.
is needed for conducting the latter. Witczak et al. (2002) used a GP-based approach
for process identification and fault diagnosis of non-linear dynamical systems. They
proposed a new fault detection observer and also demonstrated the use of GP for
increasing the convergence rate of the observer. The reliability and effectiveness
of the identification network proposed by Witczak et al. (2002) were checked by
constructing models for a few individual parts of the evaporation section at Lublin
Sugar Factory S.A. and for an induction motor.
Madar et al. (2005) applied GP to develop nonlinear input–output models for
dynamic systems. They hybridized GP and the orthogonal least squares (OLS)
method for the selection of a model structure using a tree representation based
symbolic optimization. The strategy was implemented using MATLAB GP-OLS
Toolbox—a rapid prototyping system—for predicting (a) the structure of a known
model, (b) the model order for a continuous polymerization reaction, and (c) both
order and structure of the model for Van der Vusse reaction. The results of this
modeling study indicated that the proposed strategy provides an efficient strategy for
the selection of the model order and identification of the model structure. Recently,
Faris and Sheta (2013) also adopted GP for the system identification of Tennessee
Eastman Chemical process reactor.
Other than the main areas of GP applications covered in the preceding subsections,
there exist a number of studies wherein the formalism has been employed to address
diverse problems in chemical engineering.
Genetic programming based methods are frequently employed in the multi-scale
modeling of process and product data (Seavey et al. 2010). A number of properties
such as critical flux heat prediction, viscosity, dew point pressure, compressibility,
permeation flux, solubility and gas consumption have been modeled using the
GP approach. It is often applied in conjunction with the genetic algorithms for
automatically building the kinetic models in terms of ordinary differential equations
(ODEs) to model complex systems of chemical reactions (Cao et al. 1999).
The proportional-integral-derivative (PID) controllers are the most commonly
used industrial process control strategies. Implementation of a PID controller
requires knowledge of three parameters, namely, the proportional gain (Kp ), the
integral time (Ti ) and the derivative time (Td ). Ziegler-Nichols (ZN) proposed
a method to determine these parameters. However, parameters evaluated in this
manner usually have an overshoot of 25 % thus making their fine tuning essential.
Almeida et al. (2005) used GP for fine tuning PID controller parameters designed
via ZN technique. The GP algorithm was programmed to create an initial population
of 500 individuals (candidate solutions), which evolved over 30 generations. The
GP-based fine tuning of PID parameters is a simple and an efficient method and
it improved the settling time of the system with a minimum overshoot and with a
null steady-state error. This performance was clearly seen in the three case studies
5 Genetic Programming Applications in Chemical Sciences and Engineering 129
that Almeida et al. (2005) performed, namely, a high order process, a process
with a large time delay, and a highly non-minimum phase process. Their GP-
based approach when compared with four other fine tuning techniques, was found
to exhibit superior performance. GP has been also used in the implementation of
nonlinear model predictive control strategies for the rapid acquisition of efficient
models to accurately predict the process trajectories (Tun and Lakshminarayanan
2004).
Often a situation arises, in which an appropriate hardware-based sensor for
measuring a process variable is either unavailable or the alternative analytical
procedure for its determination is time-consuming, expensive and tedious. In
such cases, a suitably developed soft-sensor can be employed for estimating the
magnitude of the “tricky-to-measure” process variable/parameter. Soft-sensor is a
software module consisting of a mathematical model that utilizes the available
quantitative information of other process variables and parameters for estimating
the magnitude of the chosen variable/parameter. Recently, Sharma and Tambe
(2014) demonstrated that GP can be effectively used to develop soft-sensors
models for biochemical systems. Specifically, they developed the GP-based soft-
sensors possessing excellent prediction accuracy and generalization capability for
two biochemical processes, namely, extracellular production of lipase enzyme and
bacterial production of poly(3-hydroxybutyrate-co-3-hydroxyvalerate) copolymer.
The strategy developed by Sharma and Tambe (2014) is generic and can be extended
to develop soft-sensors for various other types of processes.
Despite the fact that most industrial reactions employ heterogeneous catalysis,
GP has received little attention in this area. Baumes et al. (2009), however, used
GP for an advanced performance assessment of industrially relevant heterogeneous
catalysts. Epoxides are cyclic ethers with three ring atoms. Their structure is highly
strained, which makes them more reactive than other ethers. Epoxidation of double
bonds to obtain epoxides, is carried out with micro- and meso-porous titanosilicates
(Ti-MCM-41 and Ti-ITQ-2) as catalysts. The catalytic activity of these materials can
be improved by controlling their surface properties. Baumes et al. (2009) achieved
this control by anchoring alkyl-silylated agents onto the catalyst surface, which
modifies the hydrophilic nature of the catalyst. In the absence of a rigorous kinetic
study of the synthesized catalysts, they used GP for obtaining a model for the
conversion of reactant in presence of a catalyst as a function of the reaction time.
The catalyst activity was assessed via the conversion versus reaction time curve.
Catalyst performance evaluation by this method is based on the total reaction time.
The catalyst activity was monitored during 16 h of reaction in a batch reactor and
the GP model constructed thereby resulted in an adjusted fitness of 0.93. Baumes
et al. (2009) also presented a GP algorithm with the context aware crossover (CAX)
operator, which did not perform better than the ordinary crossover operator.
The chemical industry requires reliable and accurate thermodynamic data for
different fluids, covering a wide range of temperature, pressure and composition
(Hendriks et al. 2010). The knowledge of thermodynamic properties of fluids plays
a critical role in the design and operation of chemical processes. A large number
of phenomenological and empirical models have been developed for the prediction
130 R. Vyas et al.
‰; xi ; yi D f1 .T; P; zi / (5.13)
yi D Ki xi (5.15)
X zi .Ki 1/
D0 (5.16)
i
1 C ‰ .Ki 1/
i Psat
Ki D i
(5.17)
P
5 Genetic Programming Applications in Chemical Sciences and Engineering 131
5.5 Conclusion
methodologies such as ANNs and SVR/SVM. Often it has been found that
for non-linear systems, no single modeling/classification formalism consistently
outperforms the other methods. Accordingly, it is advisable to explore multiple
approaches such as GP, ANN and SVR/SVM for conducting modeling/classification
of nonlinear systems and choose the best performing one.
Although there are a few studies demonstrating its capability to obtain phe-
nomenological models (see for example, Schmidt and Lipson 2014) this fascinating
and un-matched feature of GP has been largely ignored. It is thus necessary to
exploit the stated GP characteristic extensively for developing first principles models
in chemistry and chemical engineering.
A large number of semi-quantitative and purely empirical correlations are
routinely used in chemistry and chemical engineering/technology. Prediction accu-
racies of many of these correlations are far from satisfactory. Also, in a number of
instances, linear correlations—since being easy to develop—are utilized although
the underlying phenomena being modelled are nonlinear. It is possible to construct
these correlations freshly using GP to improve their prediction accuracy and
generalization performance. The notable feature of GP that it is by itself capable
of arriving at an appropriate linear or a nonlinear model, can be gainfully exploited
for the development of the correlations alluded to above.
GP possesses certain limitations such as it is computationally demanding, the
solutions provided by it may be over fitted and an extensive heuristics is involved
in obtaining the best possible solution. Despite these limitations it is envisaged
that owing to its several attractive and unique features together with the advent of
user-friendly software such as Eureqa Formulize, Discipilus and ECJ tool kit, GP
will be extensively and fruitfully employed for providing meaningful relationships
and insights into the vast data available in the domain of chemical sciences and
engineering.
References
Anderson AC (2003) The process of structure-based drug design. Chem Biol 10:787–797
Anderson B, Svensson P, Nordahl M et al (2000) On-line evolution of control for a four-legged
robot using genetic programming. In: Cagnoni S et al (eds) Real World Applications of
Evolutionary Computing, Springer, Berlin, p 319-326
Atkinson AJ Jr., Huang S-M, Lertora J, Markey SP (eds) (2012) Principles of Clinical Pharmacol-
ogy. Elsevier, San Diego, USA
Bagheri M, Bagheri M, Gandomi AH, Golbraikh A (2012) Simple yet accurate prediction
method for sublimation enthalpies of organic contaminants using their molecular structure.
Thermochimica Acta 543: 96-106
Bagheri M, Gandomi AH, Bagheri M, Shahbaznezhad M (2013) Multi-expression programming
based model for prediction of formation enthalpies of nitro-energetic materials. Expert Systems
30: 66–78
Bagheri M, Borhani TNG, Gandomi AH et al (2014) A simple modelling approach for prediction
of standard state real gas entropy of pure materials. SAR QSAR Environ Res 25:695-710
Barmpalexis P, Kachrimanis K, Tsakonas A et al (2011) Symbolic regression via genetic
programming in the optimization of a controlled release pharmaceutical formulation. Chemom
Intell Lab Syst 107:75-82
Baumes LA, Blansché A, Serna P et al (2009) Using genetic programming for an advanced per-
formance assessment of industrially relevant heterogeneous catalysts. Mater Manuf Processes
24:282-292
Biswas A, Maitre O, Mondal DN et al (2011) Data-driven multiobjective analysis of manganese
leaching from low grade sources using genetic algorithms, genetic programming, and other
allied strategies. Mater Manuf Processes 26:415-430
Bonet B, Geffner H (1999) Planning as heuristic search: New results. In: Proceedings of the
5th European Conference on Planning (ECP-99). Springer-Verlag, Heidelberg, Germany,
p 360–372
Cabena P, Hadjnian P, Stadler R et al (1997) Discovering data mining: from concept to implemen-
tation, Prentice Hall, USA
Cai W, Pacheco-Vega A, Sen M et al (2006) Heat transfer correlations by symbolic regression. Int
J Heat Mass Transfer 49(23-24):4352–4359
Can B, Heavy C (2012) A comparison of genetic programming and artificial neural networks in
metamodeling of discrete event simulation models. Comput Oper Res 39(2):424-436
Cantú-Paz E, Kamath C (2003) Inducing oblique decision trees with evolutionary algorithms. IEEE
Trans Evol Comput 7(1):54-68
Cao H, Yu J, Kang L et al (1999) The kinetic evolutionary modeling of complex systems of
chemical reactions. Comput Chem 23(2):143-151
Cartwright H (2008) Using artificial intelligence in chemistry and biology. Taylor & Francis, Boca
Raton, FL
Castelli M, Vanneschi L, Silva S (2014) Prediction of the unified Parkinson’s disease rating scale
assessment using a genetic programming system with geometric semantic genetic operators.
Expert Syst Appl 41(10): 4608–4616
Cevik A (2007) Genetic programming based formulation of rotation capacity of wide flange beams.
J Constr Steel Res 63: 884–893
Cheema JJS, Sankpal NV, Tambe SS et al (2002) Genetic programming assisted stochastic
optimization strategies for optimization of glucose to gluconic acid fermentation. Biotechnol
Progr 18:1356-1365
Cho J-S, Kim H, Choi J-S et al (2010) Prediction of reverse osmosis membrane fouling due to scale
formation in the presence of dissolved organic matters using genetic programming. Desalin
Water Treat 15:121-128
Christin C, Hoefsloot HC, Smilde AK et al (2013) A critical assessment of feature selection
methods for biomarker discovery in clinical proteomics. Mol Cell Proteomics 12(1):263-76
Costa P, Lobo JMS (2001) Modeling and comparison of dissolution profile. Eur J Pharm Sci
13:123–133
5 Genetic Programming Applications in Chemical Sciences and Engineering 135
Das GK, Acharya S, Anand S et al (1995) Acid pressure leaching of nickel-containing chromite
overburden in the presence of additives. Hydrometallurgy 39:117-128
Dassau E, Grosman B, Lewin DR (2006) Modeling and temperature control of rapid thermal
processing. Comput Chem Eng 30:686-697
ECJ 22- A java-based evolutionary computation research system (2014) http://cs.gmu.edu/~eclab/
projects/ecj/. (Accessed 3 Dec 2014)
Eggermont J, Kok JN, Kosters WA (2004) Genetic programming for data classification: partition-
ing the search space. Presented at the 19th Annual ACM Symposium on Applied Computing
(SAC’04), Nicosia, Cyprus, 14-17 March 2004, p 1001-1005
Elsey J, Riepenhausen J, McKay B et al (1997) Modelling and control of a food extrusion process.
Comput Chem Eng 21 (supplementary): 5361-5366. Proceedings of the European Symposium
on Computer Aided Process Engineering , ESCAPE-7 , Trondheim, Norway
Elsharkawy A (2002) Predicting the dewpoint pressure for gas condensate reservoir: Empirical
models and equations of state. Fluid Phase Equilib 193:147–165
Espejo P, Ventura S, Herrera F (2010) A Survey on the Application of Genetic Programming to
Classification. IEEE Trans Syst Man and Cybern 40(2): 121-144
Estrada-Gil JK, Fernández-López JC, Hernández-Lemus E, Silva-Zolezzi I, Hidalgo-Miranda A,
Jiménez-Sánchez G, Vallejo-Clemente EE (2007) GPDTI: A Genetic Programming Decision
Tree Induction method to find epistatic effects in common complex diseases. Bioinformatics
23(13):i167-i174
Faris H, Sheta AF (2013) Identification of the Tennessee Eastman chemical process reactor using
genetic programming. Int J Adv Sci Technol 50:121-139
Fattah KA (2012) K-value program for crude oil components at high pressures based on PVT
laboratory data and genetic programming. J King Saud University—Eng Sci 24:141-149
Fattah KA (2014) Gas-oil ratio correlation (RS) for gas condensate using genetic programming.
J Petrol Explor Prod Technol 4:291-299
Ferreira C (2001) Gene Expression Programming in Problem Solving. In: Roy R, Köppen M,
Ovaska S, Furuhashi T, Hoffmann F (eds) Soft computing and industry-recent applications,
Springer-Verlag, Heidelberg, p 635-654
Finley SD, Chu LH, Popel AS (2014) Computational systems biology approaches to anti-
angiogenic cancer therapeutics. Drug Discov Today. doi: 10.1016/j.drudis.2014.09.026.
Flasch O (2014) A Friendly Introduction to RGP, RGP release 0.4-1. Available online: http://cran.
r-project.org/web/packages/rgp/vignettes/rgp_introduction.pdf. Accessed 6 Nov 2014
Gandomi AH, Alavi AH (2011) Multi-stage genetic programming: A new strategy to nonlinear
system modeling. Information Sciences 181: 5227–5239
Garcia B, Aler R, Ledeezma A, Sanchis (2008) Genetic Programming for Predicting Protein
Networks. In: Geffner H et al (eds) Proceedings of the 11th Ibero-American Conference on AI,
IBERAMIA 2008, Lisbon, Portugal. Lecture notes in computer science, vol 5290. Springer-
Verlag, Germany, p 432-441
Garg A, Vijayaraghavan V, Mahapatra SS, Tai K, Wong CH (2014) Performance evaluation of
microbial fuel cell by artificial intelligence methods. Expert Syst Appl 41(4):1389-1399
Ghosal K, Chandra A, Rajabalaya R et al (2012) Mathematical modeling of drug release profiles
for modified hydrophobic HPMC based gels. Pharmazie 67(2):147-55
Ghugare SB, Tiwary S, Elangovan V et al (2014a) Prediction of higher heating value of solid
biomass fuels using artificial intelligence formalisms. Bioenergy Res 7:681-692
Ghugare SB, Tiwary S, Tambe SS (2014b) Computational intelligence based models for prediction
of elemental composition of solid biomass fuels from proximate Analysis. Int J Syst Assur Eng
Manag DOI 10.1007/s13198-014-0324-4 Published online December 7, 2014
Giri BK, Hakanen J, Miettinen K et al (2012) Genetic programming through bi-objective genetic
algorithms with a study of a simulated moving bed process involving multiple objectives. Appl
Soft Comput 13:2613-2623
Globus AI, Lawton J, Wipke T (1999) Automatic molecular design using evolutionary techniques.
Nanotechnology 10: 290-299
136 R. Vyas et al.
Goldberg DE (1989) Genetic algorithms in search, optimization, and machine learning. Addison-
Wesley, Reading, MA
Greeff DJ, Aldrich C (1998) Empirical modelling of chemical process systems with evolutionary
programming. Comput Chem Eng 22:995-1005
Grosman B, Lewin DR (2002) Automated nonlinear model predictive control using genetic
programming. Comput Chem Eng 26:631-640
Grosman B, Lewin DR (2004) Adaptive genetic programming for steady-state process modeling.
Comput Chem Eng 28:2779-2790
Gu MB, Mitchell RJ, Kin BC (2004) Whole cell based biosensors for environmental monitoring
and application. Adv Biochem Eng/Biotechnol 87: 269-305
Guido RV, Oliva G, Andricopulo AD (2008) Virtual screening and its integration with modern drug
design technologies. Curr Med Chem 15(1):37-46
Güres S, Mendyk A, Jachowicz R et al (2012) Application of artificial neural networks (ANNs)
and genetic programming (GP) for prediction of drug release from solid lipid matrices. Int. J.
Pharm. 436(1–2): 877–879
Gustafson SM, Burke EK, Krasnogor N (2005) On improving genetic programming for symbolic
regression. In: Proceedings of the 2005 IEEE Congress on Evolutionary Computation, Edin-
burgh, UK, 2-5 September 2005. Volume 1, IEEE Press, Los Alamitos, p 912–919
Hendriks E, Kontogeorgis GM, Dohrn R et al (2010) Industrial Requirements for Thermodynamics
and Transport Properties. Ind Eng Chem Res 49:11131–11141
Hiden HG, Willis MJ, Tham MT et al (1999) Non-linear principal components analysis using
genetic programming. Comput Chem Eng 23: 413-425
Hinchcliffe MP, Willis MJ (2003) Dynamic systems modeling using genetic programming. Comput
Chem Eng 27:1841-1854
Holland JH (1975) Adaptation in natural and artificial systems. University of Michigan Press, Ann
Arbor
Hosseini SH, Karami M, Olazar M et al (2014) Prediction of the minimum spouting velocity by
genetic programming approach. Ind Eng Chem Res 53:12639-12643
Iba H (1996) Random tree generation for genetic programming. In: Voigt H-M, Ebeling W,
Rechenberg I, Schwefel, H.-P (eds) 4th International Conference on Parallel Problem Solving
from Nature (PPSN IV), Berlin, Germany, 22 - 26 September 1996. Lecture Notes in Computer
Science, vol 1411. Springer-Verlag, Germany, p 144–153
Iba H, Paul TK, Hasegawa Y (2010) Applied genetic programming and machine learning. Taylor
and Francis, Boca Raton, FL
Ivanova N, Lykidis A (2009) Metabolic Reconstruction. Encyclopedia of Microbiology (Third
Edition) p 607–621
Keedwell E, Narayanan A (2005) Intelligent bioinformatics. John Wiley & Sons, Chichester,
England
Khan MW, Alam M (2012) A survey of application: Genomics and genetic programming, a new
frontier. Genomics 100(2): 65-71
Kondos PD, Demopoulos GP (1993) Statistical modelling of O2 -CaCl2 -HCl leaching of a complex
U/Ra/Ni/As ore. Hydrometallurgy 32: 287-315
Kotanchek M (2006) Symbolic regression via genetic programming for nonlinear data modeling.
In: Abstracts, 38th Central Regional Meeting of the American Chemical Society, Frankenmuth,
MI, United States, 16-20 May 2006. CRM-160
Koza JR (1990) Genetic programming: a paradigm for genetically breeding populations of
computer programs to solve problems. Stanford University, Stanford
Koza JR (1992) Genetic Programming: On the programming of computers by means of natural
selection. MIT Press, Cambridge, MA, USA
Koza JR (1991) Concept formation and decision tree induction using the genetic programming
paradigm. In: Proceedings of the 1st Workshop on Parallel Problem Solving from Nature
(PPSN-1), London, UK. Springer-Verlag, Heidelberg, Germany, p 124–128
Koza JR, Mydlowee W, Lanza G, Yu J, Keane MA (2001) Reverse engineering of metabolic
pathways from observed data using genetic programming. Pac Symp Biocomput 6: 434-445
5 Genetic Programming Applications in Chemical Sciences and Engineering 137
Koza JR, Poli R (2005) Genetic programming. In: Burke EK, Kendall G (eds), Search Method-
ologies: Introductory Tutorials in Optimization and Decision Support Techniques, Springer,
Boston, p 127-164
Koza JR, Rice JP (1991). Genetic generation of both the weights and architecture for a neural
network. Neural Networks 2:397-404
Kulkarni BD, Tambe SS, Dahule RK et al (1999) Consider genetic programming for process
identification. Hydrocarbon Process 78:89-97
Langdon WB, Barrett SJ (2005) Evolutionary computing in data mining In: Ghosh A, Jain LC (eds)
Studies in Fuzziness and Soft Computing, vol 163. Springer-Verlag, Heidelberg, p 211-235
Lanza G, Mydlowec W, Koza JR (2000) Automatic creation of a genetic network for the lac operon
from observed data by means of genetic programming. Paper presented at the First International
Conference on Systems Biology (ICSB), Tokyo, 14-16 November 2000
Lee DG, Kim H-G, Baek W-P et al (1997) Critical heat flux prediction using genetic programming
for water flow in vertical round tubes. Int Commun Heat Mass Transfer 24: 919–929
Lee TM, Oh H, Choung YK, et al (2009) Prediction of membrane fouling in the pilot-scale
microfiltration system using genetic programming. Desalination 247:285–294
Lennartsson D, Nordin P (2004) A genetic programming method for the identification of signal
peptides and prediction of their cleavage sites. EURASIP Journal on Applied Signal Processing
1:138-145
Lew TL, Spencer AB, Scarpa F et al (2006) Identification of response surface models using genetic
programming. Mech Syst Sig Process 20(8):1819–1831
Li J, Wong L (2004) Rule-Based Data Mining Methods for Classification Problems in Biomedical
Domains In: A tutorial note for the 15 the European conference on machine learning (ECML)
and the 8th European conference on principles and practice of knowledge discovery in
databases (PKDD), Pisa, Italy, September 2004
Madar J, Abonyi J, Szeifert F (2004) Genetic programming for system identification. In: 4th
International conference on intelligent systems design and application. Available via CiteSeerx.
http://conf.uni-obuda.hu/isda2004/7_ISDA2004.pdf. Accessed 18 Nov 2014
Madár J, Abonyi J, Szeifert F (2005) Genetic programming for the identification of nonlinear input
output models. Ind Eng Chem Res 44:3178-3186
Madar J, Abonyi J, Szeifert F (2005) Genetic programming for system identification. Available
online: http://conf.uni-obuda.hu/isda2004/7_ISDA2004.pdf. Accessed 3 Dec 2014
Manshad AK, Ashoori S, Manshad MK et al (2012a) The prediction of wax precipitation by neural
network and genetic algorithm and comparison with a multisolid model in crude oil systems.
Pet Sci Technol 30:1369-1378
Manshad AK, Manshad MK, Ashoori S (2012b) The application of an artificial neural network
(ANN) and a genetic programming neural network (GPNN) for the modeling of experimental
data of slim tube permeability reduction by asphaltene precipitation in Iranian crude oil
reservoirs. Pet Sci Technol 30:2450-2459
Marquardt DW (1963) An algorithm for least-squares estimation of nonlinear parameters. J Soc
Ind Appl Math 11:431-441
Marref A, Basalamah S, Al-Ghamdi R (2013) Evolutionary computation techniques for predicting
atmospheric corrosion. Int J Corros. doi: http://dx.doi.org/10.1155/2013/805167
Martin MA (2010) First generation biofuels compete. N Biotechnol 27 (5): 596–608
McKay B, Willis M, Barton G (1997) Steady-state modeling of chemical systems using genetic
programming. Comput Chem Eng 21:981-996
McQuiston FC (1978) Heat, mass and momentum transfer data for five plate-fin-tube heat transfer
surfaces. ASHRAE Trans 84:266–293
Meighani HM, Dehghani A, Rekabdar F et al (2013) Artificial intelligence vs. classical approaches
: a new look at the prediction of flux decline in wastewater treatment. Desalin Water Treat
51:7476-7489
Motsinger AA, Lee SL, Mellick G, Ritchie MD (2006) GPNN: power studies and applications of a
neural network method for detecting gene-gene interactions in studies of human disease. BMC
Bioinformatics 7:39
138 R. Vyas et al.
Muzny CD, Huber ML, Kazakov AF (2013) Correlation for the viscosity of normal hydrogen
obtained from symbolic regression. J Chem Eng Data 58:969-979
Nandi S, Rahman I, Tambe SS, Sonolikar RL, Kulkarni BD (2000) Process identification using
genetic programming: A case study involving Fluidized catalytic cracking (FCC) unit. In: Saha
RK, Maity BR, Bhattacharyya D, Ray S, Ganguly S, Chakraborty SL (eds) Petroleum refining
and petrochemicals based industries in Eastern India, Allied Publishing Ltd., New Delhi,
p 195-201
Narendra K, Parthasarathy K (1990) Identification and control of dynamic systems using neural
networks. IEEE Trans Neural Network 1: 4-27
Nemeth LK, Kennedy HTA (1967) Correlation of dewpoint pressure with fluid composition and
temperature. SPE J 7:99–104
Okhovat A, Mousavi SM (2012) Modeling of arsenic, chromium and cadmium removal by
nanofiltration process using genetic programming. Appl Soft Comput 12:793-799
Pacheco-Vega A, Cai W, Sen M et al (2003) Heat transfer correlations in an air–water fin-
tube compact heat exchanger by symbolic regression. In: Proceedings of the 2003 ASME
International Mechanical Engineering Congress and Exposition, Washington, DC, 15-21
November 2003
Pandey DS, Pan I, Das S, Leahy JJ, Kwapinski W (2015) Multi gene genetic programming based
predictive models for municipal solid waste gasification in a fluidized bed gasifier. Bioresource
Technology 179:524-533
Parikh J, Channiwala SA, Ghosal GK (2005) A correlation for calculating HHV from proximate
analysis of solid fuels. Fuel 84:487–494
Park S-M, Han J, Lee S et al (2012) Analysis of reverse osmosis system performance using a
genetic programming technique. Desalin Water Treat 43:281-290
Patil-Shinde V, Kulkarni T, Kulkarni R, Chavan PD, Sharma T, Sharma BK, Tambe SS, Kulkarni
BD (2014) Artificial intelligence-based modeling of high ash coal gasification in a pilot plant
scale fluidized bed gasifier. Ind Eng Chem Res. 53:18678–18689; doi: 10.1021/ie500593j
Paul TK, Hasegawa Y, Iba H (2006) Classification of gene expression data by majority voting
genetic programming classifier. In: 2006 IEEE congress on evolutionary computation, Vancou-
ver, 16-21 July 2006
Peng DY, Robinson DB (1976) A new two constant equation of state. Ind Eng Chem Fundam
15:59–64
Persaud K, Dodd G (1982) Analysis of discrimination mechanisms in the mammalian olfactory
system using a model nose. Nature 299 (5881): 352–355
Podola B, Melkonian M (2012) Genetic Programming as a tool for identification of analyte-
specificity from complex response patterns using a non-specific whole-cell biosensor. Biosens
Bioelectron 33(1):254-259
Poli R, Langdon WB, McPhee NF (2008) A field guide to genetic programming. Available via
lulu: http://lulu.com, http://www.gp-field-guide.org.uk . Accessed on 3 Dec 2014
Porter M, Willis M, Hiden HG (1996) Computer-aided polymer design using genetic programming.
M Eng. Research Project, Dept. Chemical and Process Engng, University of Newcastle, UK.
Register Machine Learning Technologies Inc (2002) Discipulus 5 genetic programming predictive
modelling. http://www.rmltech.com/. Accessed 7 Nov 2014
Ritchie MD, Motsinger AA, Bush WS (2007) Genetic programming neural networks: A powerful
bioinformatics tool for human genetics. Appl Soft Comput 7:471–479
Ritchie MD, White BC, Parker JS et al (2003) Optimization of neural network architecture
using genetic programming improves detection of gene–gene interactions in studies of human
diseases. BMC Bioinf 4:28
Sajedian A, Ebrahimi M, Jamialahmadi M (2012) Two-phase inflow performance relationship
prediction using two artificial intelligence techniques: Multi-layer perceptron versus genetic
programming. Pet Sci Technol 30:1725-1736
Samuel, AL (1983) AI, Where it has been and where it is going. In: Proceedings of the 8th
International Joint Conference on AI (IJCAI-83), Karlsruhe, Germany. p 1152–1157
5 Genetic Programming Applications in Chemical Sciences and Engineering 139
Sankpal NV, Cheema JJS, Tambe SS et al (2001) An artificial intelligence tool for bioprocess
monitoring: application to continuous production of gluconic acid by immobilized aspergillus
niger. Biotechnol Lett 23: 911-916
Schmidt M, Lipson H (2009) Distilling free-form natural laws from experimental data. Science
324:81-85
Schmidt M, Lipson H (2014) Eureqa (Version 0.98 beta) [Software]. Available online: www.
nutonian.com. Accessed 3 Dec 2014
Schneider MV, Orchard S (2011) Omics technologies, data and bioinformatics principles. Methods
Mol Biol 719: 3-30
Searson DP, Leahy DE, Willis MJ (2010) GPTIPS: an open source genetic programming toolbox
for multi-gene symbolic regression. In: International Multi-Conference of Engineers and
Computer Scientists 2010 (IMECS 2010), Kowloon, Hong Kong, 17-19 March 2010. Volume
1. Available online: http://www.iaeng.org/publication/IMECS2010/IMECS2010_pp77-80.pdf
Seavey KC, Jones AT, Kordon AK (2010) Hybrid genetic programming - first principles approach
to process and product modelling. Ind Eng Chem Res 49(5):2273-2285
Sharma S, Tambe SS (2014) Soft-sensor development for biochemical systems using genetic
programming. Biochem Eng J 85:89-100
Shokir EMEl-M (2008) Dewpoint pressure model for gas condensate reservoirs based on genetic
programming. Energy Fuels 22:3194-3200
Shokir EMEl-M, Dmour HN (2009) Genetic programming (GP)-based model for the viscosity of
pure and hydrocarbon gas mixtures. Energy Fuels 23:3632-3636
Shokir EMEl-M, El-Awad MN, Al-Quraishi AA et al (2012) Compressibility factor model of
sweet, sour, and condensate gases using genetic programming. Chem Eng Res Des 90:785-792
Shokrkar H, Salahi A, Kasiri N et al (2012) Prediction of permeation flux decline during MF of
oily wastewater using genetic programming. Chem Eng Res Des 90:846-853
Shrinivas K, Kulkarni RP, Ghorpade RV et al (2015) Prediction of reactivity ratios in free-radical
copolymerization from monomer resonance-polarity (Q-e) parameters: Genetic programming
based models. International Journal of Chemical Reactor Engineering Published Online on
04/03/2015, DOI 10.1515/ijcre-2014-0039.
Singh KP, Gupta S (2012) Artificial intelligence based modeling for predicting the disinfection
by-products in water. Chemom Intell Lab Syst 114:122–131
Smith SL, Lones MA (2009) Implicit Context Representation Cartesian Genetic Programming for
the Assessment of Visuo-spatial Ability Evolutionary Computation. In: 2009 IEEE Congress
on Evolutionary Computation, Trondheim, Norway, 18-21 May 2009
Spall JC (1998) Implementation of the Simultaneous Perturbation Algorithm for Stochastic
Optimization. IEEE Trans Aerosp Electron Syst 34:817-822
Stein L (2001) Genome annotation: from sequence to biology. Nature Reviews Genetics 2:493-503
Sugimoto M, Kikuchi S, Tomita M (2005) Reverse engineering of biochemical equations from
time-course data by means of genetic programming. Biosystems 80(2):155–164
Suh C, Choi B, Lee S et al (2011) Application of genetic programming to develop the model for
estimating membrane damage in the membrane integrity test using fluorescent nanoparticle.
Desalination 281:80-87
Sumathi S, Surekha P (2010) Computational intelligence paradigms theory and applications using
MATLAB. Taylor and Francis, Boca Raton, FL
Tomita Y, Kato R, Okochi M et al (2014) A motif detection and classification method for peptide
sequences using genetic programming. J Biosci Bioeng 85: 89–100
Tsakonas A, Dounias G, Jabtzen J et al (2004) Evolving rule-based systems in two medical domain
using genetic programming. Artif Intell Med 32(3):195-216
Tun K, Lakshminarayanan S (2004) Identification of algebraic and static space models using
genetic programming. Dyn Control Process Syst 1:311-328
Vamce DE, Vance JE (eds) (2008) Biochemistry of lipids, lipoproteins and membranes. Elsevier,
Amsterdam
140 R. Vyas et al.
Venkatraman V, Dalby AR, Yang ZR (2004) Evaluation of mutual information and genetic
programming for feature selection in QSAR. J Chem Inf Comput Sci 44(5):1686–1692.
doi:10.1021/ci049933v
Vyas R, Goel P, Karthikeyan M, Tambe SS, Kulkarni BD (2014) Pharmacokinetic modeling of
Caco-2 cell permeability using genetic programming (GP) method. Lett Drug Des Discovery
11(9): 1112-1118
Wagner S (2009) Heuristic optimization software systems - modeling of heuristic optimization
algorithms in the HeuristicLab software environment. PhD Dissertation, Johannes Kepler
University, Linz, Austria
Wang X-H, Hu Y-D, Li Y-G (2008b) Synthesis of nonsharp distillation sequences via genetic
programming. Korean J Chem Eng 25:402-408
Wang X-H, Li Y-G (2008) Synthesis of multicomponent products separation sequences via
stochastic GP method. Ind Eng Chem Res 47:8815-8822
Wang X-H, Li Y-G (2010) Stochastic GP synthesis of heat integrated nonsharp distillation
sequences. Chem Eng Res Des 88:45-54
Wang X-H, Li Y-G, Hu Y-D et al (2008a) Synthesis of heat-integrated complex distillation systems
via genetic programming. Comput Chem Eng 32:1908-1917
Wedge DC, Das A, Dost R et al (1999) Real-time vapour sensing using an OFET-based electronic
nose and genetic programming. Bioelectrochem Bioenerg 48(2):389-96
White DR (2012) Software review: the ECJ toolkit. Genet. Prog. Evol. Mach 13: 65-67
Willis M, Hiden H, Hinchcliffe M et al (1997) Systems modeling using genetic programming.
Comput Chem Eng 21:1161-1166
Witczak M, Obuchowicz A, Korbics J (2002) Genetic programming based approaches to identifi-
cation and fault diagnosis of non linear dynamic systems. Int J Control 75:1012-1031
Woodward AM, Gilbert RJ, Kell DB (1999) Genetic programming as an analytical tool for non-
linear dielectric spectroscopy. Bioelectrochem Bioenerg 48(2):389-396
Worzel WP, Yu J, Almal AA et al (2009) Applications of Genetic Programming in Cancer
Research. Int J Biochem Cell Biol 41(2): 405–413
Xu C, Rangaiah GP, Zhao XS (2014) Application of neural network and genetic programming
in modeling and optimization of ultra violet water disinfection reactors, Chem Eng Commun
doi: 10.1080/00986445.2014.952813
Chapter 6
Application of Genetic Programming
for Electrical Engineering Predictive
Modeling: A Review
6.1 Introduction
Over the last decade, GP has received the interest of streams of researchers around
the globe. First, we wanted to provide an outline of the basics of GP, to sum up
valuable tasks that gave impetus and direction to research in GP as well as to
discuss some interesting applications and directions. Things change fast in this area,
as researchers discover new paths of doing things, and new things to do with GP.
It is not possible to cover all phases of this field, even within the generous page
limits of this chapter.
GP produces computer models to solve a problem utilizing the principle of
Darwinian natural selection. GP results are computer programs that are represented
as tree structures and shown in a functional programming language (such as LISP)
(Koza 1992; Alavi et al. 2011). In other words, programs evolved by genetic
programming are parse trees whose length can change throughout the run (Hosseini
et al. 2012; Gandomi et al. 2012). GP provides the architecture of the approximation
model together with the values of its parameters (Zhang et al. 2011; Gandomi et al.
2011). It optimizes a population of programs based on a fitness landscape specified
by a program capability to perform a given task. The fitness of each program is
assessed utilized an objective function. Therefore, a fitness function is the objective
function that GP optimizes (Gandomi et al. 2010; Javadi and Rezania 2009; Torres
et al. 2009). GP and other evolutionary methods have been successfully applied
to different supervised learning work like regression (Oltean and Dioan 2009), and
unsupervised learning work like clustering (Bezdek et al. 1994; Jie et al. 2004; Falco
et al. 2006; Liu et al. 2005; Alhajj and Kaya 2008) and association discovery (Lyman
and Lewandowski 2005).
and mutates it (see Fig. 6.3). GEP is a linear variant of GP. The linear variants
of GP make a clear distinction between the genotype and the phenotype of an
individual. Thus, the individuals are represented as linear strings that are decoded
and expressed like nonlinear entities (trees) (Yaghouby et al. 2010; Baykasoglu et al.
2008; Gandomi et al. 2008).
6.3 GP Applications
The following section shows a review of engineering applications of GP. The results
of the literature survey have been organized into the following broad groups:
• Control
• Optimization and scheduling
• Signal processing
• Classification
• Power System Operation
6 Application of Genetic Programming for Electrical Engineering Predictive. . . 145
6.3.1 Control
input and create two outputs values. The output values has been transmitted to the
robot as motor speeds. the population of each individuals is processed by the GP
system.
Figure 6.6 shows a schematic view of the system. This schematic has been
captured from Nordin et al. (1997).
Grimes (1995) were used genetic algorithm (GA) and genetic programming (GP)
methods for track maintenance work with profit as the optimization criteria. The
results were compared with an existing method. It was shown that the GP algorithm
6 Application of Genetic Programming for Electrical Engineering Predictive. . . 147
provided the best results, with the GA approach providing good results for a short
section and poor results for a long section of track. Genetic programming was
used in Stephenson et al. (2003) to optimize the priority functions associated with
two well-known compiler heuristics: predicted hyperblock formation and register
allocation. Their system achieved remarkable speedups over a standard baseline
for both problems. Vanneschi and Cuccu (2009) presented a new model of genetic
programming with variable size population in this paper and applied to the
reconstruction of target functions in dynamic environments. This models suitability
was tested on a set of benchmarks based on some well-known symbolic regression
problems.
Experimental results confirmed that their variable size population model found
solutions of similar quality to the ones found by genetic programming, but with a
smaller amount of computational effort. Ho et al. (2009) developed an algorithm
to derive a distributed method automatically dynamically to optimize the coverage
of a femtocell group utilizing genetic programming. The resulting evolved method
showed the capability to optimize the coverage well. Also, this algorithm was able to
offer increased overall network capacity compared with a fixed coverage femtocell
deployment. The evolution of the best-known schedule illustrated in Langdon and
Treleaven (1997) for the base South Wales problem utilizing genetic programming
starting from the hand coded heuristics. Montana and Czerwinski (1996) applied
a hybrid of a genetic algorithm and strongly typed genetic programming (STGP)
to the problem of controlling the timings of traffic signals that optimize aggregate
performance. STGP learns the single basic decision tree to be executed by all the
intersections when determining whether to change the phase of the traffic signal.
148 S.S. Sadat Hosseini and A. Nemati
6.3.4 Classification
Decision
Trees
6.4 Conclusions
This paper has provided us with the background context required to understand the
reviewed documents and use as a guideline to categorize and sort relevant literature.
While computer scientists have focused on gaining a significant understanding of
the methods the engineering community is solving practical problems, frequently
by introducing accepted systems engineering methodologies and concepts. The
combination of different methods permits us to make the most of several algorithms,
using their strengths and preventing their drawbacks. The flexibility of GP makes it
possible to combine it with very various algorithms. But the combination of GP with
some other methods is not the only option; GP can be employed as a mechanism to
integrate different techniques. It is stressed that GP is a young area of research,
whose practitioners are still exploring its abilities and drawbacks. Therefore, it is
the authors’ belief that the future holds much promise.
References
Ahmad AM, Khan GM (2012) Bio-signal processing using cartesian genetic programming evolved
artificial neural network (cgpann). In: Proceedings of the 10th International Conference on
Frontiers of Information Technology, 261–268
Holladay K, Robbins K (2007) Evolution of signal processing algorithms using vector based
genetic programming. 15th International Conference in Digital Signal Processing, Cardiff,
Wales, UK, 503–506
Harding S, Leitner J, Schmidhuber J (2013) Cartesian genetic programming for image processing.
In Riolo, R., Vladislavleva, E., Ritchie, M. D., and Moore, J. H., editors, Genetic Programming
Theory and Practice X, Genetic and Evolutionary Computation, 31–44. Springer New York
Sharman KC, Alcazar AIE, Li Y (1995) Evolving signal processing algorithms by genetic
programming. First International Conference on Genetic Algorithms in Engineering Systems:
Innovations and Applications, GALESIA, IEE, 414, 473–480
Esparcia Alcázar AI (1998) Genetic programming for adaptive digital signal processing. PhD
thesis, University of Glasgow, Scotland, UK
Esparcia-Alcázar A, Sharman K (1999) Genetic Programming for channel equalisation. In R. Poli,
H. M. Voigt, S. Cagnoni, D. Corne, G. D. Smith, and T. C. Fogarty, editors, Evolutionary
Image Analysis, Signal Processing and Telecommunications: First European Workshop, 1596,
126–137, Goteborg, Sweden, Springer-Verlag
Alcázar, Anna I. Esparcia, and Ken C. Sharman. “Some applications of genetic programming
in digital signal processing.” In Late Breaking Papers at the Genetic Programming 1996
Conference Stanford University, pp. 24–31. 1996
Smart W, Zhang M (2003) Classification strategies for image classification in genetic program-
ming. In: Proceeding of Image and Vision Computing Conference, 402–407, New Zealand
Li J, Li X, Yao X (2005) Cost-sensitive classification with genetic programming. The IEEE
Congress on.Evolutionary Computation
Escalante HJ, Acosta-Mendoza N, Morales-Reyes A, Gago-Alonso A (2009) Genetic Program-
ming of Heterogeneous Ensembles for Classification. in Progress in Pattern Recognition, Image
Analysis, Computer Vision and Applications, Springer, 9–16
Liu KH, Xu CG (2009) A genetic programming-based approach to the classification of multiclass
microarray datasets. Bioinformatics 25(3):331–337
Zhang L, Nandi AK (2007) Fault classification using genetic programming. Mech Syst Signal Pr
21(3):1273–1284
Chaturvedi DK, Mishra RK, Agarwal A (1995) Load Forecasting Using Genetic Algorithms
Journal of The Institution of Engineers (India), EL 76, 161–165
Dr. Hanan Ahmad Kamal (2002) Solving Curve Fitting problems using Genetic Programming
IEEE MELECON May, 7–9
Farahat MA (2010) A New Approach for Short-Term Load Forecasting Using Curve Fitting
Prediction Optimized by Genetic Algorithms 14th International Middle East Power Systems
Conference (MEPCON10)19–21
Chapter 7
Mate Choice in Evolutionary Computation
7.1 Introduction
Darwin’s theory of Natural Selection (Darwin 1859) has been widely accepted and
endorsed by the scientific community since its early years. Described as the result of
competition within or between species affecting its individuals rate of survival, it has
had a deep impact on multiple research field and is at the source of the ideas behind
EC. The theory of Sexual Selection (Darwin 1906) was later developed by Darwin to
account for a number of traits that were observed in various species, which seemed
to have no place in his Natural Selection theory. Darwin described Sexual Selection
as the result of the competition between individuals of the same species affecting
their relative rate of reproduction, a force capable of shaping traits to high degrees of
complexity and responsible for the emergence of rich ornamentation and complex
courtship behaviour.
Despite having been discredited by the scientific community at the time, it is
now widely regarded as a major influence on evolution theory. Interest arose in the
1970s through the works of Fisher (1915, 1930) and Zahavi (1975) and since then
the community as gradually embraced it, having found its place in various research
fields. While it has come a long way, Sexual Selection is still far from understood in
EC, both regarding possible benefits and behaviour.
Electronic supplementary material: The online version of this chapter (doi: 10.1007/978-3-319-
20883-1_7) contains supplementary material, which is available to authorized users.
A. Leitão () • P. Machado
CISUC, Department of Informatics Engineering, University of Coimbra, 3030-290 Coimbra,
Portugal
e-mail: [email protected]; [email protected]
Mate Choice was one of the processes of Sexual Selection described by Darwin
and that mostly attracted his followers. This chapter describes a nature-inspired
self-adaptive Mate Choice setup and covers the design steps necessary for applying
it. A two-chromosome scheme where the first chromosome represents a candidate
solution and the second chromosome represents mating preferences is proposed.
Two approaches for encoding mating preferences are presented and differences
discussed. Details on how to apply each of them to problems with different
characteristics are given and design choices are discussed. The application of both
approaches on different problems is reviewed and the observed behaviour discussed.
Section 7.2 introduces Mate Choice as well as its background, Sect. 7.3 gives
a general overview of Mate Choice in Evolutionary Computation and covers the
state of the art through a classification based on adaptation of mating preferences,
popular preference choices and the role of genders. The section finally introduces
the proposed setup, giving specific details on the ideas behind it and how to
apply it. Section 7.4 describes the application of both the proposed approaches
to multiple problems and discusses the obtained results and behaviours. Finally,
Sect. 7.5 presents a summary.
Since his journey on the Beagle, Darwin has thoroughly studied the forces responsi-
ble for the evolution of species. The result of competition within or between species
affecting their individuals relative rate of survival was named Natural Selection.
Since the publication of Darwin’s On the Origin of Species by Means of Natural
Selection, or the Preservation of Favoured Races in the Struggle for Life in 1859
(Darwin 1859), the theory has become widely accepted by the scientific community.
This was achieved thanks to the evidence gathered by Darwin, its co-discoverer
Alfred Russel Wallace (1858) as well as multiple following researchers, ultimately
overcoming other competitive ideas (Cronin 1993).
Despite such a success, Darwin battled with gaps in its theory. For instance,
Darwin questioned how was it that Natural Selection could account for animal
ornamentation or courtship behaviour. He observed a large number of species, where
individuals displayed rich and costly ornamentations or complex and risky courtship
behaviours that seemed to serve no purpose in survival, sometimes even risking it.
These characteristics challenged the theory of Natural Selection and the idea that
traits adapted to the environment in a purposeful way. Individuals carrying aimless
and costly features should be unfavored in competition, making such features bound
to face extinction.
However, as Darwin observed, that was not the case. Ornamentation and
courtship behaviour were spread across populations and species, although Natural
Selection could not explain their origin. He figured, however, that for these features
to emerge they had to bring some kind of competitive advantage. As they didn’t fit
in Natural Selection, he envisioned the existence of another trait-shaping selection
force in nature, one capable of shaping species in complex and diverse ways, by
7 Mate Choice in Evolutionary Computation 157
causing traits that help in competing for mates to spread through future generations.
These traits were linked to reproduction, as he observed in nature, and brought
evolutionary advantages, even when risking survivability. Darwin developed the
theory of Sexual Selection (Darwin 1906) to explain this phenomena and described
it as the result of competition within species affecting its individuals relative rate of
reproduction.
Darwin therefore saw evolution as the interplay between two major forces,
Natural Selection as the adaptation of species to their environment, and Sexual
Selection as the adaptation of each sex in relation to the other, in a struggle of
individuals of one sex for the possession of individuals of the other in order to
maximize their reproductive advantages. While the outcome of failing in Natural
Selection would be low survivability, the outcome of failing in Sexual Selection
would be a low number or no offspring. From an evolutionary perspective, they
reach the same outcome with competition in reproductive rates between individuals
leading to evolutionary changes across populations.
Unlike his theory of Natural Selection, which easily found support on the
scientific community, his theory of Sexual Selection was mostly rebuffed. The
scientific community was not keen on Darwin’s ideas regarding Sexual Selection,
specially his ideas on Female Mate Choice and the impact it could have on evolution.
It was clear for them that Natural Selection was the only force capable of adapting
species and so a number of theories emerged in order to explain rich displays or
courtship behaviour. One of the most avid opponents of Darwin’s ideas was Wallace
who came up with various reasons for the emergence of traits such as protection
through dull colors in females, recognition of individuals of the same species, usage
of surplus energy on courtship behaviour or non-selective side effects (Cronin 1993).
The community was better prepared to understand such ideas, which were
embraced by various renowned researchers such as Huxley (Cronin 1993). This lead
to a time where Darwin’s ideas on Sexual Selection and specially Mate Choice were
dismissed as non-important, and its impact to be regarded as a small part of Natural
Selection. These ideas remained for over a century, with the exception of a few
works by a select few researchers who made important contributions, such as Fisher
(1915) who explored the origin of mating preferences and runaway sexual selection,
Williams and Burt (1997) who discussed how ornaments should be considered
as important as other adaptations or Zahavi (1975) who expanded on the role of
displays as fitness indicators.
Overtime, the work of these researchers was able to gain some space in the
community, eventually reaching more open-minded generations who were also
better equipped to test and understand the workings of Sexual Selection and Mate
Choice. The resulting discussion has for the past few decades attracted experimental
biologists, psychologists and anthropologists that since then have put Darwin’s
ideas as well as those promoted by his followers to the test, contributing with
increasing evidence to back the ideas behind Sexual Selection through Mate Choice.
Nowadays, there is active research on various fields and the theory has been widely
accepted by the community. Two extensive reviews on Sexual Selection have been
published by Helena Cronin (1993) and Malte Andersson (1994).
158 A. Leitão and P. Machado
Mate Choice is one of the main Sexual Selection processes described by Darwin
and where he put much of his effort, as did most of his followers. They aimed
to explain the emergence of aesthetic features such as ornamentation or courtship
behaviour. One of the pillars of Darwin’s theories on evolution was that species
went through adaptations because they brought some kind of advantage over time.
In this case, these traits emerged because they brought reproductive advantage as a
result of preferences when selecting mating partners (Darwin 1906). Fisher helped
explain the relation between mating preferences and traits, and the genetic link
between them through his theory of runaway sexual selection. He addressed how
displays may arise as a result of positive reinforcement between mental mating
preferences and physical traits, through a feedback loop that can lead to extravagant
adaptations such as the peacock’s tail, colorful appearance or complex courtship
behaviour (Fisher 1915, 1930).
Fisher’s work suggests the inheritance of mating preferences much like any
other trait, therefore adapting throughout the generations. This process can be better
understood if mate choice is thought of like any other adaptive choice such as food
choice (Miller 1994). Still, criticism of theses ideas remained, since the evolution of
traits through such a runaway process with increasing speed could drastically risk
the survival ability of individuals. Zahavi later expanded on this subject, suggesting
that aesthetic displays can act as indicators of fitness, health, energy, reproductive
potential etc. He argued through his handicap principal (Zahavi 1975) that even
in the case of costly displays and behaviour, which seemed to have no purpose
in Natural Selection, it was in fact their high cost that made them good fitness
indicators. As these traits were handicaps, they couldn’t be maintained by weak,
unfit individuals and that only strong healthy individuals would be able to maintain
them and survive. Therefore reinforcement of mating preferences for these traits
would be beneficial for females, which would in turn reinforce such physical traits
in males as suggested by Fisher (1915).
These ideas were explored and discussed by many other researchers, who finally
brought Sexual Selection through Mate Choice into the spotlight. Their work
corroborated Darwin’s ideas and brought new evidence allowing Sexual Selection to
be seen as an important force in Evolutionary Theory. The interplay between Natural
Selection and Sexual Selection was found to have a deep impact on various traits
on many different species, especially among those equipped with complex sensory
systems (Cronin 1993). During the last few decades, Sexual Selection has found
its place on various research fields such as Evolutionary Biology, Evolutionary
Psychology and Evolutionary Anthropology. On the other hand, it is yet to attract the
full attention of the Evolutionary Computation community, despite the publication
of several papers over the last couple of decades.
The possible advantages that Sexual Selection, particularly through Mate Choice
can bring to the field of Evolutionary Computation have been previously discussed
by several researchers, and an extensive discussion on arguably the most relevant
ones has been published by Miller and Todd (1993). They find that the addition
of Mate Choice to Natural Selection can bring advantages such as (1) increased
accuracy when mapping from phenotype to fitness, therefore reducing the “error”
7 Mate Choice in Evolutionary Computation 159
Selection Selection
Parent 1 Parent 2
160 A. Leitão and P. Machado
Selection Selection
Parent 1 Random
candidates
Selection
Most
attractive
Parent 2
Mating preferences may remain static over the generations or undergo adaptation.
In order to address preferences and adaptation mechanisms we will rely on the
classification of adaptation of parameters and operators by Hinterding et al. (1997).
Afterwards, we will look on different preferences and how they can be used to assess
genotypes or phenotypes. Finally we will discuss the use of genders in Mate Choice
strategies.
7 Mate Choice in Evolutionary Computation 161
Various authors rely on guiding the choice of mating partners using pre-established
preferences, which remain static over the evolution process. The most widely known
strategies are probably those that mate individuals based on similarity measures.
Examples found in the literature include using Hamming distances (De et al. 1998;
Fernandes and Rosa 2001; Galan et al. 2013; Ochoa et al. 2005; Ratford et al. 1996,
1997; Varnamkhasti and Lee 2012), Euclidean distances (Galan et al. 2013; Ratford
et al. 1996), number of common building blocks (Ratford et al. 1996) or even
simply using the fitness value as a distance measure (De et al. 1998; Galan et al.
2013; Goh et al. 2003; Ratford et al. 1996; Varnamkhasti and Lee 2012). In some
implementations, the first selected parent chooses the candidate that maximizes
similarity (Fernandes et al. 2001; Hinterding and Michalewicz 1998) while in other
cases the candidate that minimises the measure is considered the best (Varnamkhasti
and Lee 2012). Other approaches attribute a probability of selection proportional or
inversely proportional to the distance measures (Ratford et al. 1996).
Moreover, some authors don’t want to maximize or minimize distances but rather
consider an ideal distance and favour mating candidates that have distances closer to
that pre-established value. In this case the attractiveness of a candidate is established
using bell curves or other functions with a predefined center and width parameters
(Ratford et al. 1997).
Other metrics have been applied such as in Hinterding and Michaelwicz’s study
on constrained optimization (Hinterding and Michalewicz 1998). They suggest
having the first parent select the mating candidate that in conjunction with itself
maximizes the number of constraints satisfied. A second example is the study by
Fernandes et al. on vector quantization problems where a problem specific metric
is used (Fernandes et al. 2001). Sometimes different metrics are combined in a
seduction function. This can be accomplished through the use of rules such as
choosing the fittest candidate if two candidates both maximize or minimize the
similarity measure or choosing between them randomly if they both also share the
same fitness value (Varnamkhasti and Lee 2012). A different approach is to combine
different metrics using different functions (such as weighted functions) which has
been done for instance by Ratford et al. (1996).
While fixed parameters and mating preferences can often achieve competitive
results and reproduce desired behaviours in Mate Choice algorithms, allowing their
online control may be extremely valuable. Such approaches, which allow the Mate
Choice strategy to change online without external control, therefore turning it into a
dynamic process, can be subdivided into three groups: deterministic, adaptative and
self-adaptative.
Deterministic approaches are the least common among the literature. Still,
Ratford et al. give a good example (Ratford et al. 1997). On an aforementioned
study they use a function to calculate a candidate’s attractiveness which has two
variables, centre and depth. It values individuals whose hamming distance from the
first parent is closer to the centre of the function. However they complement this
study by testing an approach where the centre of the function is adjusted at each
162 A. Leitão and P. Machado
generation, so that dissimilar individuals are favoured at the beginning of the run
but the opposite happens at the end.
Adaptative approaches are more common and rely on information about the
evolution process to adapt parameters or preferences. Fry et al. (2005) have applied
such a strategy to control which operator is used to select the second parent, either a
regular tournament selection or a Mate Choice operator. In their study they choose
mating partners by combining fitness with a penalization for similar candidates
(different similarity measures are applied). However they use this operator with a
given probability which is increased or reduced depending on its relative success
in producing enhanced offspring in previous generations. A second example can
be found on studies by Sánchez-Velazco and Bullinaria (2003, 2013). Mating
candidates in this case are evaluated based on a weighted function that combines
three metrics: fitness, likelihood of producing enhanced offspring, and age. While
age is adapted deterministically at each generation, the second factor represents a
feedback on each individual’s ability to produce fit offspring in the past.
Self-adaptive approaches better resemble the workings of Mate Choice in nature.
By allowing preferences and parameters to be encoded in each individual, self-
adaptation allows them to take part in the evolution process and to impact not
only the individuals that encode them but the population as a whole. Possibly the
simplest example of such an approach relies on encoding an index as an extra gene
in each individual. When evaluating its mating candidates, each individual will order
them from best to worst according to a given metric and select the candidate at the
encoded position. Galan et al. (2013) have experimented with this approach using
Euclidean distance and fitness to order mating candidates. On the aforementioned
study by Fry et al. (2005), a second approach was tested, where each individual
encodes its own probability of selecting a mating partner using a mate choice
operator rather than regular tournament selection. The probability is inherited by
the offspring and adapted by comparing their fitness with that of their parents.
More complex mating preferences can be found on self-adaptive approaches.
Miller and Todd (1993) and Todd and Miller (1997) suggest encoding a reference
position on the phenotype space marking each individual’s ideal position for a
mating candidate. When assessing potential mating partners, the probability of
mating varies according to their distance to the reference position. New offspring
inherit genetic material from both parents through two-point crossover thus allowing
for its self-adaptation throughout the evolutionary process.
Holdener and Tauritz (2008) relied on an extra chromosome to encode a list
of desired features to look for on mating candidates. They tackle a problem
using a binary representation on the first chromosome but rely on a real value
representation on the preferences chromosome. This chromosome has the same
size as the first chromosome with each gene representing how much an individual
wants the corresponding gene to be set to 1. This information is used to evaluate
mating candidates by comparing the preferences chromosome with each candidate’s
potential solution, favouring desired genes. Preference genes are inherited from
parents to offspring so that they match the genes they influence and adapt to match
the offspring’s relative success. On a related study, Guntly and Tauritz (2011)
7 Mate Choice in Evolutionary Computation 163
In nature, Mate Choice is almost absolutely on the side of females. Due to their high
reproductive investment, they are picky when selecting a mating partner, looking for
a fit male that can provide good genes. On the other hand males are more willing
164 A. Leitão and P. Machado
When designing Mate Choice approaches we feel that in order to best resemble the
natural process, models should follow three nature-inspired rules:
1. individuals must choose who they mate with based on their own mating
preferences
2. mating preferences, as mental traits, should be inherited the same way as physical
ones
3. mate selection introduces its own selection pressure but is subject to selection
pressure itself
We see the evaluation of mating candidates as a complex process, where the
relation between observed traits, or their weight on each individual’s mate choice
mechanism is difficult to establish beforehand. While some traits could be seen
as valuable on a mating candidate, others could turn out to be irrelevant or even
harmful. The relation between them is also certainly not straightforward as they
could be connected in unforeseeable ways. Moreover, certain displayed traits may
be very important for survival purposes but hold little value for mate choice, or
the other way around. This value can also vary on each selection step, depending
7 Mate Choice in Evolutionary Computation 165
on the characteristics of parent1, its mating preferences and the mating candidates
involved.
This discussion is particularly relevant when we recall the following aspects of
Sexual Selection through Mate Choice:
1. each individual has its own characteristics and may benefit differently from
reproducing with different mates. Each individual also has its own distinct
mating preferences that may value different characteristics in mating candidates.
The reproductive success of individuals depends on choosing appropriate mating
partners the same way that it depends on how attractive they are to others.
2. the paradigm may result on cases where individuals with poor survival abilities
attain a high reproductive success because they display characteristics that are
favoured by mating preferences. Their offspring may achieve low fitness values
but may contribute to exploration and the emergence of innovation which may
eventually turn into ecological opportunities.
3. the handicap principle shows that certain traits may risk the survival ability of
individuals while in fact being indicative of good gene quality, thus reducing the
accuracy of fitness values. Mate Choice mechanisms may be able to help increase
the accuracy of mapping between phenotypes and fitness values and translate that
into reproductive success.
4. mating preferences and evolved physical traits have an intrinsic and deep
dependence between them which results from the feedback loop described in
the theory of runaway sexual selection. The resulting arms race causes mating
preferences to evolve in relation to displayed traits and physical traits to adapt
in relation to enforced mating preferences. This process can lead traits to a high
degree of elaboration.
With these ideas in mind and with the goal of designing mate choice mechanisms
that best resemble the natural process, we refrain from establishing what are good
or bad mating preferences. Also, we avoid linking each individual’s candidate
solution with its mating preferences using any pre-established method, such as
inheritance rules. We therefore treat genetic material regarding physical traits and
mating preferences equally and leave the responsibility of adapting individuals up
to the evolutionary process. Inheritance and selection pressure should be able to
bring reproductive and survival advantages to individuals carrying genes linked
to appropriate phenotypes and mating preferences through the intrinsic relation
between Natural Selection and Sexual Selection through Mate Choice.
The following subsections detail how we design Mate Choice approaches that
meet the presented rules and aspects.
7.3.2.1 Representation
assess potential mating partners. The first chromosome could use any representation
wanted for the problem at hand, be it a Genetic Algorithm (GA) vector, a GP
tree or others. We propose representing mating preferences, therefore the second
chromosome, as a GP tree.
Mating preferences can be represented using two possible approaches: (1)
representing an ideal mating partner; (2) representing an evaluation function. The
first approach requires that we are able to map GP representations to the phenotype
space of the problem at hand. In th