AUTOMATED DISCOVERY OF NUMERICAL APPROXIMATION FORMULAE VIA GENETIC ...

102 Pages·2001·0.68 MB·English

AUTOMATED DISCOVERY OF NUMERICAL APPROXIMATION FORMULAE
VIA GENETIC PROGRAMMING

by

Matthew J. Streeter

A Thesis

Submitted to the Faculty

of the

WORCESTER POLYTECHNIC INSTITUTE

in partial fulfillment of the requirements for the

Degree of Master of Science

in

Computer Science

by

_______________________________
May 2001

APPROVED:

_____________________________________________
Dr. Lee A. Becker, Major Advisor

_____________________________________________
Dr. Micha Hofri, Head of Department
i
Abstract

This thesis describes the use of genetic programming to automate the discovery of numerical approximation formulae.
Results are presented involving rediscovery of known approximations for Harmonic numbers and discovery of rational
polynomial approximations for functions of one or more variables, the latter of which are compared to Padé
approximations obtained through a symbolic mathematics package. For functions of a single variable, it is shown that
evolved solutions can be considered superior to Padé approximations, which represent a powerful technique from
numerical analysis, given certain tradeoffs between approximation cost and accuracy, while for functions of more than
one variable, we are able to evolve rational polynomial approximations where no Padé approximation can be computed.
Furthermore, it is shown that evolved approximations can be iteratively improved through the evolution of
approximations to their error function. Based on these results, we consider genetic programming to be a powerful and
effective technique for the automated discovery of numerical approximation formulae.
ii

ACKNOWLEDGEMENTS

The author wishes to thank Prof. Lee Becker and Prof. Micha Hofri of Worcester Polytechnic Institute for valuable
advice and feedback received during the course of this project.
iii

CONTENTS

1. INTRODUCTION
1.1. Introduction to Genetic Algorithms ........................................................................................................ 1
1.2. Introduction to Genetic Programming ................................................................................................... 3
1.3. Using Genetic Programming to Discover Numerical Approximation Formulae ................................ 5
1.4. Evaluating Approximations ..................................................................................................................... 6
1.5. Related Work ............................................................................................................................................ 6
1.6. Summary of Report .................................................................................................................................. 6

2. OUR GENETIC PROGRAMMING SYSTEM
2.1. GA Framework ......................................................................................................................................... 8
2.2. GP Representation .................................................................................................................................... 8
2.3. Primitive Function Costs........................................................................................................................... 9
2.4. Program Output ........................................................................................................................................ 12
2.5. Consistency with Other Genetic Programming Systems ...................................................................... 17

3. OPTIMIZING GP PARAMETERS
3.1. Experiments with Initial Test Suite ......................................................................................................... 19
3.2. Experiments with Revised Test Suite ...................................................................................................... 21
3.3. Impracticality of Optimizing GP Parameters in this Manner .............................................................. 23

4. REDISCOVERY OF HARMONIC NUMBER APPROXIMATIONS ................................................................ 24

5. DISCOVERY OF RATIONAL POLYNOMIAL APPROXIMATIONS FOR KNOWN FUNCTIONS
5.1. Introduction ............................................................................................................................................... 27
5.2. Comparison with Padé Approximations ................................................................................................. 27
5.3. Avoiding Division by Zero ........................................................................................................................28
5.4. Results ........................................................................................................................................................ 28

6. APPROXIMATING FUNCTIONS OF MORE THAN ONE VARIABLE .......................................................... 34

7. REFINING APPROXIMATIONS
7.1. Approximating Error Function of Evolved Approximations ............................................................... 36
7.2. Other Possible Approaches to Refinement of Approximations ............................................................ 40

8. ATTEMPTED REDISCOVERY OF NEURAL NETWORK ACTIVATION FUNCTIONS ........................... 41
iv
9. ATTEMPTED PIECEWISE APPROXIMATION OF FUNCTIONS
9.1. Introduction and Preliminary Work ....................................................................................................... 46
9.2. Piecewise Rational Polynomial Approximations of Functions of a Single Variable ........................... 47
9.3. Piecewise Rational Polynomial Surface Approximations ..................................................................... 56
9.4. 3-D Surface Generation ............................................................................................................................ 57

10. FUTURE WORK ..................................................................................................................................................... 58

11. SUMMARY AND CONCLUSIONS
11.1 Summary .................................................................................................................................................. 59
11.2 Conclusions ............................................................................................................................................... 59

APPENDIX A: EXTENDED RESULTS FOR RATIONAL POLYNOMIAL APPROXIMATION .................... 60
OF FUNCTIONS

APPENDIX B: CODE DOCUMENTATION .............................................................................................................. 88

REFERENCES ............................................................................................................................................................... 92
v
LIST OF TABLES AND FIGURES

Table 1.1: GP Parameters ................................................................................................................................................. 5
Table 2.1: Primitive Functions ......................................................................................................................................... 8
Table 2.2: Pentium II-233 Timing Data ........................................................................................................................... 10
Table 2.3: Pentium-60 Timing Data ................................................................................................................................ 10
Table 2.4: Final Assigned Costs of Primitive Functions ................................................................................................. 11
Table 2.5: Parameter Settings for Reproduction of Symbolic Regression Experiment ................................................... 18
Table 2.6: Results for Reproduction of Symbolic Regression Experiment ..................................................................... 18
Table 3.1: Experiments with Initial Test Suite ................................................................................................................ 20
Table 3.2: Experiments with Revised Test Suite Using Population Size = 250 .............................................................. 22
Table 3.3: Experiments with Revised Test Suite Using Population Size = 500 .............................................................. 22
Table 3.4: Number of Runs Required to Produce Definitive Results for Experiments with Revised Test Suite ............ 23
Table 4.1: Evolved Harmonic Number Approximations ................................................................................................ 24
Table 4.2: Accuracy of Asymptotic Expansion ............................................................................................................... 26
Table 5.1: Evolved Approximations for ln(x) .................................................................................................................. 29
Table 5.2: Maple Evaluation of Approximations for ln(x) .............................................................................................. 30
Table 5.3: Final Evolved Approximations for ln(x) ........................................................................................................ 30
Table 5.4: Final Evolved Approximations for sqrt(x) ...................................................................................................... 31
Table 5.5: Final Evolved Approximations for arcsinh(x) ................................................................................................ 31
Table 5.6: Final Evolved Approximations for exp(-x) .................................................................................................... 32
Table 5.7: Final Evolved Approximations for tanh(x) ..................................................................................................... 33
Table 6.1: Final Evolved Approximations for xy ............................................................................................................. 35
Table 7.1: Maple Evaluation of Approximations for sin(x) ............................................................................................. 36
Table 7.2: Final Evolved Approximations for sin(x) ....................................................................................................... 37
Table 7.3: Final Evolved Approximations for Refinement of Candidate Approximation 3 for sin(x) ............................ 38
Table 7.4: Final Evolved Approximations for Refinement of Candidate Approximation 7 for sin(x) ............................ 39
Table 7.5: Final Evolved Approximations for Refinement of Candidate Approximation 8 for sin(x) ............................ 39
Table 7.6: Final Refined Approximations for sin(x) ........................................................................................................ 40
Table 8.1. Rational Polynomial Approximations for Perceptron Switching Function ................................................... 41
Table 8.2. Approximations for Perceptron Switching Function Using Function Set {*,+,/,-,EXP} ............................... 43
Table 9.1: Error of Best Evolved Piecewise Approximation to ln(x) Using Various Function Sets ............................... 47
Table 9.2: Evolved Piecewise Rational Polynomial Approximations for Three-Peaks Function .................................. 48
Table 9.3: Evolved Rational Polynomial Approximations for Three-Peaks Function ................................................... 50
Table 9.4: Evolved Piecewise Rational Polynomial Approximations for Two-Peaks Function .................................... 52
Table 9.5: Evolved Rational Polynomial Approximations for Two-Peaks Function ..................................................... 55
Table 9.6: Evolved Rational Polynomial Approximations for Hemispherical Surface .................................................. 56
Table A.1: Evolved Approximations for sqrt(x) .............................................................................................................. 60
vi
Table A.2: Maple Evaluation of Approximations for sqrt(x) .......................................................................................... 62
Table A.3: Final Evolved Approximations for sqrt(x) ..................................................................................................... 63
Table A.4: Evolved Approximations for arcsinh(x) ........................................................................................................ 64
Table A.5: Maple Evaluation of Approximations for arcsinh(x) ..................................................................................... 66
Table A.6: Final Evolved Approximations for arcsinh(x) ............................................................................................... 67
Table A.7: Evolved Approximations for exp(-x) ............................................................................................................. 68
Table A.8: Maple Evaluation of Approximations for exp(-x) ......................................................................................... 69
Table A.9: Final Evolved Approximations for exp(-x) .................................................................................................... 70
Table A.10: Evolved Approximations for tanh(x) ........................................................................................................... 70
Table A.11: Maple Evaluation of Approximations for tanh(x) ........................................................................................71
Table A.12: Final Evolved Approximations for tanh(x) .................................................................................................. 72
Table A.13: Padé Approximations for ln(x) .................................................................................................................... 72
Table A.14: Padé Approximations for sqrt(x) ................................................................................................................. 74
Table A.15: Padé Approximations for arcsinh(x) ............................................................................................................ 76
Table A.16: Padé Approximations for exp(-x) ................................................................................................................ 77
Table A.17: Padé Approximations for tanh(x) ................................................................................................................. 79
Table A.18: Evolved Approximations for xy ................................................................................................................... 80
Table A.19: Maple Evaluation of Approximations for xy ................................................................................................ 81
Table A.20: Final Evolved Approximations for xy .......................................................................................................... 81
Table A.21: Evolved Approximations for sin(x) ............................................................................................................. 82
Table A.22: Evolved Approximations for Refinement of Candidate Approximation 3 for sin(x) ................................. 83
Table A.23: Maple Evaluation of Approximations for Refinement of Candidate Approximation 3 for sin(x) .............. 84
Table A.24: Final Evolved Approximations for Refinement of Candidate Approximation 3 for sin(x) ......................... 84
Table A.25: Evolved Approximations for Refinement of Candidate Approximation 7 for sin(x) ................................. 85
Table A.26: Maple Evaluation of Approximations for Refinement of Candidate Approximation 7 for sin(x) .............. 85
Table A.27: Final Evolved Approximations for Refinement of Candidate Approximation 7 for sin(x) ......................... 86
Table A.28: Evolved Approximations for Refinement of Candidate Approximation 8 for sin(x) ................................. 86
Table A.29: Maple Evaluation of Approximations for Refinement of Candidate Approximation 8 for sin(x) .............. 87
Table A.30: Final Evolved Approximations for Refinement of Candidate Approximation 8 for sin(x) ......................... 87

Figure 1.1: Parental LISP Expressions ............................................................................................................................ 3
Figure 1.2: Child LISP Expression .................................................................................................................................. 4
Figure 2.1: Example HTML Summary ............................................................................................................................ 12
Figure 2.2: Fitness Curve ................................................................................................................................................. 14
Figure 2.3: Error Curve .................................................................................................................................................... 14
Figure 2.4: Cost Curve ..................................................................................................................................................... 14
Figure 2.5: Adjusted Error Curve .................................................................................................................................... 14
Figure 2.6: Adjusted Cost Curve ...................................................................................................................................... 15
vii
Figure 2.7: Candidate Solutions Imported into Maple ..................................................................................................... 15
Figure 2.8: Convergence Probability Curve .................................................................................................................... 16
Figure 2.9: Individual Effort Curve ................................................................................................................................. 17
Figure 2.10: Expected Number of Individuals to be Processed........................................................................................ 17
Figure 5.1: Pareto Fronts for Approximations of ln(x) .................................................................................................... 31
Figure 5.2: Pareto Fronts for Approximations of sqrt(x) ................................................................................................. 31
Figure 5.3: Pareto Fronts for Approximations of arcsinh(x) ........................................................................................... 32
Figure 5.4: Pareto Fronts for Approximations of exp(-x) ................................................................................................ 32
Figure 5.5: Pareto Fronts for Approximations of tanh(x) ................................................................................................ 33
Figure 6.1: f(x)=xy ............................................................................................................................................................ 35
Figure 6.2: x/(y2+x-xy3) ................................................................................................................................................... 35
Figure 7.1: Error Function for Candidate Approximation (3) ......................................................................................... 37
Figure 7.2: Error Function for Candidate Approximation (7) ......................................................................................... 38
Figure 7.3: Error Function for Candidate Approximation (8) ......................................................................................... 38
Figure 8.1: Plot of Rational Polynomial Approximations for Perceptron Switching Function ....................................... 42
Figure 8.2: Plot of Approximations for Perceptron Switching Function Evolved Using Function Set {*,+,/,-,EXP} .... 44
Figure 9.1: Graph of Three-Peaks Function ..................................................................................................................... 48
Figure 9.2: Graph of Two-Peaks Function ....................................................................................................................... 52
1

1 INTRODUCTION
Numerical approximation formulae are useful in two primary areas: firstly, approximation formulae are used in industrial
applications in a wide variety of domains to reduce the amount of time required to compute a function to a certain degree
of accuracy (Burden and Faires 1997), and secondly, approximations are used to facilitate the simplification and
transformation of expressions in formal mathematics. The discovery of approximations used for the latter purpose
typically requires human intuition and insight, while approximations used for the former purpose tend to be polynomials
or rational polynomials obtained by a technique from numerical analysis such as Padé approximants (Baker 1975;
Bender and Orszag 1978) or Taylor series. Genetic programming (Koza 1989, Koza 1990a, Koza 1992) provides a
unified approach to the discovery of approximation formulae which, in addition to having the obvious benefit of
automation, provides a power and flexibility that potentially allows for the evolution of approximations superior to those
obtained using existing techniques from numerical analysis. In this thesis, we discuss a number of experiments and
techniques which demonstrate the ability of genetic programming to successfully discover numerical approximation
formulae, and provide a thorough comparison of these techniques with traditional methods.

1.1 INTRODUCTION TO GENETIC ALGORITHMS

Genetic algorithms represent a search technique in which simulated evolution is performed on a population of entities or
objects, with the goal of ultimately producing an individual or instance which satisfies some specified criterion.
Specifically, genetic algorithms operate on a population of individuals, usually represented as bit strings, and apply
selection, recombination, and mutation operators to evolve an individual with maximal fitness, where "fitness" is
measured in some domain-dependent fashion as the extent to which a given individual represents a solution to the
problem at hand. Genetic algorithms are a powerful and practical method of search, supported both by successful real-
world applications and a solid theoretical foundation. Genetic algorithms have been applied to a wide variety of
problems in many domains, including problems involving robotic motion planning (Eldershaw and Cameron 1999),
modelling of spatial interactions (Diplock 1996), optimization of database queries (Gregory 1998), and optimization of
control systems (Desjarlais 1999).
As an example application for genetic algorithms, suppose we wished to find a real value x which satisfies the quadratic
equation x2 + 2x + 1 = 0, and did not have access to the quadratic formula or any applicable technique from numerical
analysis. We might choose to represent possible solutions to this equation as a set or population (initially generated at
random) of single-precision IEEE floating point numbers, each of which requires 32 bits of storage. The fitness measure
for an individual with x-value x could be defined as the squared difference between x2 + 2x + 1 and the target value of
I I I
0. As in nature, the more fit individuals of the population will be more likely to reproduce and will tend to have a larger
number of children. A "child" can be produced from two IEEE floating point numbers by generating a random bit
position n (0<=n<=31), copying the leftmost n bits from one parent, and taking the remaining 32-n bits from the other, a
process referred to as single-point crossover. As the process continues over many generations, it is likely that an
individual will eventaully be evolved which satisfies the equation exactly (in this case, an individual with x-value x=-1).
I
2
This search method turns out to be highly flexible and efficient. In general, any search problem for which an appropriate
representation and fitness measure can be defined may be attempted by a genetic algorithm.
In Adaptation in Natural and Artificial Systems, John Holland laid the foundation for genetic algorithms. The "Holland
GA" employs fitness-proportionate reproduction, crossover, and (possibly) mutation. Pseudo code for this algorithm is
given below:
1. Initialize a population of randomly created individuals
2. Until an individual is evolved whose fitness meets some pre-established criterion:
2.1. Assign each individual in the population a fitness, based on some domain-specific fitness function.
2.2. Set the "child population" to the empty set.
2.3. Until the size of the child population equals that of the parent population:
2.3.1. Select two members of the parent population, with the probability of a member being
selected being proportionate to its fitness (the same member may be selected twice).
2.3.2. Breed these two members using a crossover operation to produce a child.
2.3.3. (Possibly) mutate the child, according to some pre-specified probability.
2.3.4. Add the new child to the child population.
2.4. Replace the parent population with the child population.

The theoretical underpinnings of this algorithm are given in the Schema theorem (Holland 1975), which establishes the
near mathematical optimality of this algorithm under certain circumstances.
A schema is a set of bit-strings defined by a string of characters, with one character corresponding to each bit. The
characters may be '1', indicating that a 1 must appear in the corresponding bit-position, '0', indicating that a 0 must
appear in the corresponding bit position, or '*', indicating that either value may appear.
For example, the an individual encoded by the bit string:
10101101000010011110101100011111
would be an instance of the schemata 1*******************************,
*01*****************************, and ***0***********************11111, but not of the
schema 0*******************************. Since we can create a schema of which an individual is an
instance using two possible characters for each bit position ('*' and the character which actually occurs in that bit
position), each individual will be a member of 232 different schemata. Following the definition of the fitness of an
individual, the fitness of a schema can be defined as the average fitness of all individuals which are instances of that
schema. For example, the schema 1*. . .* (1 followed by 31 *'s), which in our representation denotes the set of all
negative IEEE floating-point numbers, might be expected to have a different fitness than the schema 0*. . .*, which
denotes the set of all positive numbers.