An elitist algorithm is used. It means that only the best ant is
allowed to add some pheromone on the path that it has built. This
value for the extra amount of pheromone is equal to 1. Once again, to
avoid local minima, the pheromone trail is reinforced only if a new,
better solution is found by the ants. In other words, if the same best
solution is found at iteration k and iteration k +1, then no extra
pheromone is added at iteration k +1.
Our goal is to provide an algorithm that could be used in the
automatic control field. In this domain, people are not very used to
metaheuristic algorithms. Thus, it is very important to have an
algorithm that does not require a fine tuning of parameters. Therefore,
the value of parameters has not been fully optimized. Further, several
versions of the ant colony algorithm (see [DOR 05], for instance)
could be tested. However, it is not in the scope of this study to find the
best algorithm for a given problem, but to prove that the use of the
algorithm can lead to satisfying results.
2.4.2. Experimental results
Table 2.1 presents the results obtained when the depth of the tree is
set to 3. As the algorithm is a stochastic algorithm, it can be only
validated by means of statistical results. Therefore, the following
results are given in Table 2.1:
– The function to be identified. It is used to create the data
( , ), 1, , , i i x y i = N with N =100.
– The threshold. A test is successful when the least square criterion
[2.3] is less than this threshold.
– The number of successes. For each experiment, 100 tests are
performed.
– Mean number of iterations to compute the solution, when the test
is successful.
The maximum allowed number of iterations is set to 500. For some
functions, noise is added. For this purpose, a random value, in a range
of ±10% of the measured data, is added.
As can be seen from this table, results are very satisfactory for
these low-depth functions as the true initial function is always found
by the algorithm. The average number of iterations required for the
result to be found is approximately 125 iterations. The use of a
threshold as a stopping criterion is motivated by the following two
points:
– The approach is a “black box” approach. The interest is to find a
global input/output relation that fits with the data.
– In the real world, data are measured and so measurement noise
has to be taken into account: in this case, the cost can never be zero.
– In the real world, the system can be expressed by symbols that
are not considered in the algorithm.
The value of the threshold is not so crucial for “pure data”,
in other words, without any noise. The last case is a more realistic case,
as measurement noise is added to the data. In such a situation, the value
of the threshold should be adapted to the level of the noise.
Table 2.2 presents some results obtained for higher depths in the
search space. The same kind of resul
ts are given, once again
exhibiting satisfactory results.
As can easily be seen and assumed, the number of required
iterations increases with the depth of trees (the size of the search space
increases exponentially with the depth) and with the decrease in the
stopping threshold.
For the last case, the algorithm almost never finds the true initial
function y = exp(x +1) + sin(cos((x +1) / 3)) , but “only” the main
component y = exp(x +1) . The initial function can be found when the
threshold decreases. Of course, this leads to an increase in the number
of iterations.
The results have been obtained with Matlab 2007b on an Intel Core
Duo CPU 2.5 GHz. Mean computation times are given in Table 2.3,
depending on the chosen depth for the solution tree.
2.5. Discussion
2.5.1. Considering real variables
In the experiments presented in section 2.4, only a finite set of
constants are used (the integer numbers from 1 to 9). This can be seen
as a strong restriction of the proposed algorithm. In this section, an
extension to the identification of functions, which may depend on real
variables, is proposed. Two possibilities can be given for this purpose.
The first possibility uses a post-processing step. In this stage, a
classical nonlinear least square method is used, keeping the structure
found by the ant colony.
In other words, suppose that the ant colony has found the function
y = 2sin(x +1) whose parameters are only integer numbers. Then, the
following optimization problem can be solved in a second stage:
2
, 1
min ( sin( ))
N
a b i i i
y a x b
=
− + [2.10]
It is well known that nonlinear optimization results strongly
depend on the initial point. However, in this case, the constant values
found by the ant colony can be used as the initial point for the
optimization algorithm, which should be a good initialization. This
method leads to a two-step identification procedure, where the ant
colony finds the structure of the model, and the nonlinear least squares
method gives the fine values of parameters.
The second possibility is to use an extension of ant colony
algorithms for continuous problems, see, for instance, [BIL 95] and
[SOC 08].
2.5.2. Local minima
Ant colony optimization belongs to the class of metaheuristic
methods. As mentioned several times earlier, there is no guarantee on
the global optimality of the obtained solution. However, the goal ofthe identification procedure is usually to obtain a representative and
usable model of the plant and this goal is achieved even if the solution
is a local minimum.
2.5.3. Identification of nonlinear dynamical systems
The natural extension of this work deals with the identification of
nonlinear dynamical systems. The use of a state-space representation
of the plant seems natural to extend the results:
1 ( , )
( , )
k k k
k k k
x f x u
y h x u
+ =
=
[2.11]
If some experimental data ( , ) k k u y are available, there is no
theoretical difficulty to extend the results presented in this chapter, as the
problem refers to the identification of nonlinear functions f and h,
which can be represented by tree and define from a set of symbols and
parameters. The method can be used to define reduced-order models as
the number of state variables has to be chosen. However, if the system is
multi-input multi-output (MIMO) and/or has several state variables,
the number of functions to identify increases, and so the complexity of
the problem. Forthcoming works deal with this problem.
2.6. A note on genetic algorithms for symbolic regression
As mentioned in section 2.1, the use of a genetic algorithm is an
interesting alternative to the solution of symbolic regression. In this
case, a population of candidate functions has to evolve. New functions
are created with the help of the genetic operator usually denoted as the
crossing-over and mutation operator.
The crossing-over operator refers to the computation of children
functions from the random choice of two parent functions. Using the
tree representation of functions, the crossing-over can be interpreted
as an exchange of some sub-trees of the parent functions. Figure 2.4
shows an example of such transformations.
Figure 2.4. Crossing-over illustration
In this figure, two parent functions, f1(x) = sin(x + 3)×exp(x) and
2f (x) = log(2x) − cos(x) +8 , are crossed to define two children
functions, 3f (x) = sin(x + 3)×(cos(x) +8) and f4(x) = log(2x)
–exp(x).
The mutation operator consists of computing a child function from
a random choice of a parent. In this case, one node can be randomly
chosen and the corresponding symbol is replaced by another one,
which is also chosen at random.
An example of mutation is shown in Figure 2.5 when the function
1f (x) = sin(x + 3)×exp(x) is changed into 1f (x) = sin(x +3)×log(x) .