PIKAIA

Pikaia Gracilens, a little worm-like beast that crawled in the mud of a long gone seafloor of the Cambrian era, 530 million years ago. While not particularly impressive in the tooth and claw department, Pikaia is believed to be the founder of the phylum Chordata, whose subsequent evolution had consequences still very much felt today by the rest of the ecosystem. Image digitized from the excellent book The Rise of Fishes, by John A. Long (1995, The Johns Hopkins University Press).

Welcome to the PIKAIA Homepage. PIKAIA (pronounced "pee-kah-yah") is a general purpose function optimization FORTRAN-77 subroutine based on a genetic algorithm. PIKAIA is a public domain software available electronically from the anonymous ftp archive of the High Altitude Observatory. The subroutine is particularly useful (and robust) in treating multimodal optimization problems.

The development of genetic algorithm-based inversion methods is but one aspect of the research in helioseismology carried out in the Long-term Solar Variability Section (LSV) of the High Altitude Observatory (HAO), a scientific division of the National Center for Atmospheric Research (NCAR) in Boulder, Colorado. PIKAIA 1.0 was written in 1995 by Paul Charbonneau and Barry Knapp, both then at HAO/NCAR. Version 1.2 was released in April 2002.

Genetic Algorithms
PIKAIA subroutine
PIKAIA usage example
Obtaining a copy of PIKAIA
PIKAIA User Registration
PIKAIA User's List
Bug Reports
PIKAIA in IDL and f90, and parallel PIKAIA
References and further readings
Oslo Tutorial (a.k.a. NCAR TN-450+IA)

1. Genetic Algorithms

Genetic algorithms (hereafter GAs) are a class of search techniques inspired from the biological process of evolution by means of natural selection. They can be used to construct numerical optimization techniques that perform robustly on problem characterized by ill-behaved search spaces.

Consider the following generic modeling task: a model that depends on a set of adjustable parameters is used to fit a given dataset; the task consists in finding the single parameter set that minimizes the difference between the model's predictions and the data. A top-level view of a canonical genetic algorithm for this task could be as follows: start by generating a set ("population") of trial solutions, usually by choosing random values for all model parameters; then:

Evaluate the goodness of fit ("fitness") of each member of the current population (through a chi square measure with the data, for example).
Select pairs of solutions ("parents") from the current population, with the probability of a given solution being selected made proportional to that solution's fitness.
Breed the two solutions selected in (2) and produce two new solutions ("offspring").
Repeat steps (2)–(3) until the number of offspring produced equals the number of individuals in the current population.
Use the new population of offspring to replace the old population.
Repeat steps (1) through (5) until some termination criterion is satisfied (e.g., the best solution of the current population reaches a goodness of fit exceeding some preset value).

Superficially, this may look like some peculiar variation on the Monte Carlo theme. There are two crucial differences: First, the probability of a given solution being selected to participate in a breeding event is made proportional to that solution's fitness (step 2); better trial solutions breed more often, the computational equivalent of natural selection. Second, the production of new trial solutions from existing ones occurs through breeding. This involves encoding the parameters defining each solution as a string-like structure ("chromosome"), and performing genetically inspired operations of crossover and mutation to the pair of chromosomes encoding the two parents, the end result of these operations being two new chromosomes defining the two offspring. Applying the reverse process of decoding those strings into solution parameters completes the breeding process and yields two new offspring solutions that incorporate information from both parents.

Technical details:

A genetic-algorithm based approach to a given optimization task, as defined above, amounts to a form of forward modeling. Generally speaking, adopting a forward modeling approach has both advantages and drawbacks;

Advantages:

No derivatives of the goodness of fit function with respect to model parameters need be computed; it matters little whether the relationship between the model and its parameters is linear or nonlinear.
Nothing in the procedure outlined above depends critically on using a least-squares statistical estimator; any other robust estimator can be substituted, with little or no changes to the overall procedure.

Drawbacks:

In most real applications, the model will need to be evaluated (i.e., given a parameter set, compute a synthetic dataset and its associated goodness of fit) a great many times; if this evaluation is computationally expensive, the forward modeling approach can become impractical.

GA-based optimization retains the advantageous features of forward modeling, while reducing the number of required function evaluations to a level that is often much more computationally manageable.

2. The PIKAIA subroutine

PIKAIA is a fully self-contained, general purpose optimization subroutine. The routine maximizes a user-supplied FORTRAN function, the name of which is passed in as an argument. Generally, we have tried to emulate the general feel and ease of use of the Numerical Recipes optimization subroutine (for example amobea). A call to PIKAIA looks like:

call pikaia(ff,n,ctrl,xb,fb,status)

where ff is the name of the FORTRAN function to be maximized, n the dimension of the search space, and ctrl is a control vector whose elements determine the behavior of the algorithm (robust default values are built in). The output vector xb contains the optimal solution parameters, namely the best solution of the final population at the end of the evolution. The scalar fb is the fitness of the solution (i.e., the value returned by ff with xb as input). The output variable status codes successful termination, or error conditions associated with illegal input parameter values.

The (real) function ff must be declared as:

function  ff(n,x)

dimension x(n)

where the floating-point array x defines one point in the search space, and the function returns a (positive-definite) measure of fitness (or goodness of fit, or whatever) associated with that point. And that is all there is to it.

Technical details:

PIKAIA incorporates only the two basic genetic operators: uniform one-point crossover, and uniform one-point mutation. Unlike many GA packages available commercially or in the public domain, the encoding within PIKAIA is based on a decimal alphabet made of the 10 simple integers (0 through 9); this is because binary operations are usually carried out through platform-dependent functions in FORTRAN. Three reproduction plans are available: Full generational replacement, Steady-State-Delete-Random, and Steady-State-Delete-Worst. Elitism is available and is a default option. The mutation rate can be dynamically controlled by monitoring the difference in fitness between the current best and median in the population (also a default option). Selection is rank-based and stochastic, making use of the Roulette Wheel Algorithm.

PIKAIA is supplied with a ranking subroutine based on the Quicksort algorithm, and a random number generator based on the minimal standard Lehmer multiplicative linear congruential generator of Park and Miller (1988, Communications of the ACM, 31, 1192).

3. An example of PIKAIA usage

Consider an optimization problem that consists of maximizing a function of two variables f(x,y), with x and y bounded between 0 and 1. The function defines a 2-D landscape, in which one is seeking the highest elevation point. If the landscape is smooth and simple this problem is readily treated with conventional hill-climbing methods. However a landscape such as the following would be a much harder task:

This is a surface plot of the function f(x,y), and the inset in the upper right is a color-coded version of the same function. The global maximum (indicated by the arrow and located at (x,y)=(0.5,0.5), where f=1) is surrounded by concentric rings of secondary maxima, where a simple hill-climbing method would most likely get stuck.

This problem is easily solved with PIKAIA. An individual is an (x,y) pair, and fitness can be directly defined as the altitude in the 2-D landscape, i.e., the value returned by f(x,y). Examination of the corresponding fitness function and driver code shows how simple the use of PIKAIA is for such a problem. The following animation illustrates the evolution of the population's distribution in parameter space. Each individual is shown as a solid green dot, and the best of each generation as a larger, yellow dot. Observe how the "best solution remains stuck, for a little while, on the innermost ring of secondary extrema, but eventually "finds" the thin, central global maximum. The small plot in the upper right shows the variation with generation count of the error (defined as log(1-f) for the best (yellow) and median (purple) individuals in the population. The sudden drops in the error reflect the appearance of a favorable mutation, and its subsequent spread in the population. This solution was obtained by running PIKAIA for 50 generations under its default settings.

Click here to view animation [Size: 574 KB]

4. Obtaining a copy of PIKAIA

PIKAIA is public domain software and is available electronically on the High Altitude Observatory anonymous ftp archive. To get there click on the preceding anchor; save all files as Plain text, except the User's Guide, which should be saved as Postscript. Alternately, type on your own system:

ftp ftp.hao.ucar.edu 122

Use anonymous as a username and your e-mail address as a password. Once logged in type:

cd archive/pikaia

Depending on how your computer system is firewalled, you might experience problems with ftp; if so point your browser to:


http://download.hao.ucar.edu/archive/pikaia/

One way or the other, you are now in the main PIKAIA directory. The files in this directory are:

README.txt

This file contains a description of the directory's content and last minute information. Please do read it.

pikaia.f

This file contains the PIKAIA subroutine itself, together with an appropriate random number generator, ranking routine, and driver and fitness function for an installation check problem. This file is completely self-contained, and can be compiled/linked immediately.

userguide.ps

This is a Postscript file containing the User's Guide to PIKAIA 1.0, otherwise known as NCAR Technical Note 418+IA [Size: about 2MB].

userguide_errata.txt

Some minor errors and typos discovered in the User's Guide. These have been corrected in the current userguide.ps, but may be present in earlier versions of the Guide, depending on the date at which the typos were reported. Note that the typos and errors are also present in the hard copy of the Guide published by NCAR in December 1995.

examples

Subdirectory containing drivers and fitness functions for the example problems discussed in chapter 5 of the User's Guide.

examples.tar

UNIX tar file of the examples subdirectory. A convenient way to grab all the examples files in one shot if you are working on, or have access to, a UNIX system.

In case of difficulties installing or running the code, send e-mail to pikaia@hao.ucar.edu.

5. Registering as a PIKAIA user

We invite prospective users of PIKAIA to send a brief e-mail message to pikaia@hao.ucar.edu (if possible, please include a few words on what problem is to be tackled with the subroutine; we're curious). This will allow us to maintain an e-mail list, to be used to send bug reports and/or news/upgrade announcements.

Note that by registering as a PIKAIA user, you will not be added to some open e-mailing list, with the usual consequence of being submerged under mostly unwanted e-mail messages. Again, the list is used exclusively to send bug reports and other news/upgrade announcements.

6. List of registered PIKAIA Users

Current and Past Users:

Frank Crary        CU/LASP
Eugene Lavely      Blackhawk Geosc.    Geoacoustic waveform inversion
Guillermo Torres   Harvard/CfA           
Alex Razoumov      UBC/Vancouver       Nonlinear function minimization
Mark Fardal        CU/CASA             Game theory
Mihaly Horanyi     CU/LASP             Ozone inverse problem
Isodoros Doxas     CU/APAS             Magnetotail modeling
Kelsey Johnson     CU/APAS             Solar magnetogram analysis
Dieter Hartmann    Clemson U.             
Hal Levison        SWRI/Boulder         
Deborah Haber      CU/Boulder          Sunspot helioseismology
Adrian Webster     Roy.Obs. Edinburgh  Data Modeling
Martin Sperl       U. Vienna           Period fitting
Alberto Cappi      U. Bologna              
Ian Howarth        Univ. Coll. London
Richard Boivin     CERCA/Montreal      Engineering design
Jacques Richer     CERCA/Montreal      Nonlinear semi-empirical fits
Gilles Fontaine    U. de Montreal      Statistical Thermodynamics
Claude Carignan    U. de Montreal      Nonlinear least squares
Robert Lamontagne  U. de Montreal      Nonlinear function minimization
Paul Seagraves     HAO/NCAR            Stokes profiles fitting
Tim Brown          HAO/NCAR            Orbital Parameter fits
Karel Schrijver    Lockheed/Palo Alto  Emission measure analysis
Ignaz Wanders      Ohio State U.       Reverberation Mapping
Burkhart Fuchs     U. Heidelberg       Dynamic Modelling
Russeil Delphine   CNRS/Paris          Data Modeling
Luciano Mantegazza INFN/Pavia              
Tuck Stebbins      CU/Boulder
Andreas Bobinger   U. Muenchen         Inverse Modeling
Paul Cally         Monash U.                
Th. Appourchaux    ESA                 Inverse Modeling
Sarah Gibson       NASA/Goddard        Coronal Structure Inversion
Stuart Jeffries    NOAO
Ronnie Hoogerwerf  Leiden Observatory     
Reyco Henning      U. of Denver        Helioseismic inversion
Sami Solanki       ETH/Zuerich         Stokes profiles fitting
Barry Lapthorn     QMW College, U.K.   Helioseismic inversion
Daniel Carpintero  U. La Plata, Arg.
Misha Haywood      Obs. Paris
Steve Arendt       CORA/Boulder
Davis Valls-Gabaud Obs. Strasbourg     fitting Color-Mag diagrams
Andreas Lagg       Johns Hopkins       Spectral Data analysis
Bill Chaplin       ESTEC/Noordwwijk    Helioseismic inversions
Travis Metcalfe    U. Texas/Austin     White dwarf seismology
Christian Theis    U. Kiel, Ger.       Modeling galaxy encounters
Hardi Peter        Kiepenheuer Inst.   Spectral line analysis
Ladislav Hejna     Czech Acad. Sc.     Orbital fitting
Jose Adolfo        Brazil U.?          Thermo-resistive sensorsa
R. Lantosoa        Madagascar Obs.     Gravity/magnetism inversions
Olof Friberg       Chalmers U., Sweden vehicule design (acoustics)
Thierry Nieus      Belgium             biophysical modeling
Antonio Emolo      U. Naples           seismic analysisa
Luca Teriaca       Armagh Obs.         Spectral line analysis
Katherine Manson   Australian Nat. U.  magnetic spectra analysis
S.L. Hidalgo Rguez IAC/Spain           modeling interactive galaxies
Steven Hale        U. Birmingham, UK   helioseismic data processing
Jackie Schoendorf                      Data modeling
Lidong Xia         MPI Aeronomie       SUMER spectral analysis
Simon Casassus     U. Chile
Gordon Chin        NASA Goddard
Andrea Fontana     CERN                spectroscopy of antihydrogen
Bastian Liebrecht  RWTH-Aachen         Aircraft/spacecraft design
Waleed Hamdy       NRIAG/Egypt 
Stephane Courteau  UBC/Canada          galactic light curve decomp.
Charlie Meyer      STMicroelectronics  engineering design
Kevin Flood        SLAC/Stanford       particle data analysis
Ian Robinson       U. Birmingham/UK    solar wind energetic events
Eric Depagne       DASGAL/Obs. Meudon
Scott McIntosh     ESA                 Emission measure modeling
Alan Miller        CSIRO Math.& Inf.
Frank Timmes       U. Chicago          Stellar Opacities
Cesar Aaron Moya   Kyoto U.            Geoseismology
Mohamed Abdelwahed NRIAG/Egypt         Earthquakes
Mark Alston        U. Strathclyde      Spectral Modelling
I.R. De Souza      Unicamp/Brazil      Chemical engineering design
Bev Smith          East Tennessee S.U. 
Lee Rottler        UC/Santa Cruz       Period fitting & spectral analysis
Roman Petryk       UBC/Physics & Astr.
Daniel Kubas       Univ. Potsdam       Gravitational microlensing
Lapo Boschi        U. Federico II      Seismic Imaging
Antonio Piersanti  INGV/Rome           Mantle Rheology
Volker Rath        U. Aachen/RWTH      Geothermal Inverse modeling
Ivan Zevallos      U. Sao Paulo        Geophysical Inversion
Luciano Lamberti   Politecnico di Bari Structural optimization
Kiran Solanki      Mississippi State U Constrained optimization
Paola Rogata       GMV S.A./Madrid     Interplanetary trajectories
Omur Cakirli       Ege U. Observastory Asteroseismology
Jos de Bruijne     ESA/ESTEC           Engineering design
K. Gozdziewski     Torun CfA           Orbital Element Fitting
Liz Humphreys      Harvard CfA
S.A. Klein         U. Wisconsin        Equation-solving
Michiel Min        Astro. Inst.        IR spectral analysis

PIKAIA Publications in refereed Journals:

Kennelly, E.J., and 12 co-authors, 1996, The oscillation modes of theta^2 Tauri, Astronomy & Astrophysics 313, 571-580.
Kaastra, J.S., et al. 1996, Emission measure analysis methods: the corona of AR Lacertae revisited, Astronomy & Astrophysics, 314, 547-557.
Noyes, R.W., et al. 1997, A planet orbiting the star rho Coronae Borealis, The Astrophysical Journal, 483, L111-L114
Charbonneau, P., Tomczyk, S., Schou, J., and Thompson, M.J., 1998 The Rotation of the Solar Core inferred by genetic forward modeling, The Astrophysical Journal, 496, 1015-1030.
Gibson, S., and Charbonneau, P. 1998, Empirical Modeling of the Solar Corona using Genetic Algorithms, Journal of Geophysical Research, 103(A7), 14511-14521.
McIntosh, S.W., Diver, D.A., Judge, P.G., Charbonneau, P., Ireland, J., and Brown, J.C. 1998, Spectral Decomposition by Genetic Forward Modeling Astronomy & Astrophysics (Supplements), 132, 145-153.
Theis, C. 1998, Modeling Encounters of Galaxies: the case of NGC4449, Review of modern astronomy, 12, in press.
Peter, H. 1999, Analysis of transition-region emission-line profiles from full-disk scans of the sun using the SUMER instrument on SOHO, The Astrophysical Journal, 516, 490-504.
Charbonneau, P., Christensen-Dalsgaard, J., Henning, R., Larsen, R.M., Schou, J., Thompson, M.J., and Tomczyk, S. 1999 Helioseismic constraints on the structure of the solar tachocline, The Astrophysical Journal, 527, 445-460.
Peter, H. 1999, The chromosphere in coronal holes and the quiet-sun network: an HeI (584 A) full-disk scan by SUMER/SOHO, The Astrophysical Journal, 522, L77-L80.
Sevenster, M., Saha, P., Valls-Gabaud, D., and Fux, R. 1999, New constraints on a triaxial model of the Galaxy, Monthly Notices of the Royal Astronomical Society, 307, 584-594.
Peter, H., and Judge, P.G. 1999, On the Doppler shifts of solar ultraviolet emission lines, The Astrophysical Journal, 522, 1148-1166
Teriaca, L., Banerjee, D., and Doyle, J.G. 1999, SUMER observations of Doppler shift in the quite sun and in active region, Astronomy and Astrophysics, 349, 636-648
Doyle, J.G., Teriaca, L., and Banerjee, D. 1999, Coronal hole diagnostics out to 8R, Astronomy and Astrophysics, 349, 956-960.
McIntosh, S.W., Charbonneau, P., and Brown, J.C. 2000, The Astrophysical Journal, 529, 1115-1130.
Billieres, M., Fontaine, G., Brassard, P., Charpinet, S., Liebert, J., and Saffer, R.A. 2000, Detection of p-mode pulsations and possible ellipsoidal luminosity variations in the hot subdwarf B star KPD 1930+2752, The Astrophysical Journal, 530, 441-453.
Doyle, J.G., Teriaca, L., and Banerjee, D. 2000, Solar transition region line broadening: limb to limb measurements, Astronomy and Astrophysics, 356, 335-338.
Bobinger, A. 2000, Genetic eclipse mapping and the advantage of black sheep, Astronomy and Astrophysics, 357, 1170-1180.
McIntosh, S.W. 2000, On the inference of differential emission measures using diagnostic line ratios, The Astrophysical Journal, 533, 1043-1052.
Skartlien, R., and Rast, M.P. 2000, p-mode intensity-velocity phase differences and convective sources, The Astrophysical Journal, 535, 464-472.
Peter, H. 2000, Multi-component structure of solar and stellar transition regions, Astronomy and Astrophysics, 360, 761-776.
Harries, T.J., and Howarth, I.D. 2000, Spectropolarimetric orbits of symbiotic stars: five S-type systems, Astronomy and Astrophysics, 361, 139-152.
Banerjee, D., Teriaca, L., Doyle, J.G., and Lemaire, P. 2000, Polar plumes and inter-plume regions as observed by SUMER on SOHO, Solar Physics, 194, 43-58.
Metcalfe, T.S., Nather, R.E., and Winget, D.E. 2000, Genetic-algorithm-based asteroseismological analysis of the DBV white dwarf GD358, The Astrophysical Journal, 545, 974-981
Metcalfe, T.S., Winget, D.E., and Charbonneau, P. 2001, Preliminary constraints on 12C(alpha,gamma)16O from white dwarf seismology, The Astrophysical Journal, 557, 1021-1027
Theis, Ch., and Kohle, S. 2001, Multi-method-modeling for interacting galaxies. I. A unique scenario for NGC 4449? Astronomy and Astrophysics, 370, 365-383.
Torres, G. 2001, The change in the inclination angle of the noneclipsing binary SS Lacertae: future eclipses, The Astronomical Journal, 121, 2227-2238.
Peter, H. 2001, On the nature of the transition region from the chromosphere to the corona of the Sun, Astronomy and Astrophysics, 374, 1108-1120.
Bartzakos, P., Moffat, A.F.J., and Miemela, V.S. 2001, Magellanic Cloud WC.WO Wolf-Rayet stars. II. Colliding winds in binaries, Monthly Notices of the Royal Astronomical Society, 324, 33-50.
Severino, G., Magri, M., Oliviero, M., Straus, Th., and Jefferies, S.M. 2001, The solar intensity-velocity cross spectrum: a powerful diagnostic for helioseismology, The Astrophysical Journal, 561, 444-449.
Gelino, D.M., Harrison, T.E., and Orosz, J.A. 2001, A multiwavelength, multiepoch study of the soft X-ray transient prototype, V616 Monocerotis, The Astronomical Journal, 122, 2668-2678.
Metcalfe, T.S. 2003, White dwarf asteroseismology and the 12C(alpha,gamma)16O rate, The Astrophysical Journal (Letters), in press.

Please send preprints/reprints of any work you have carried out where PIKAIA was used in a non-trivial way. Like we said, we're curious.

7. Bug Reports

Since the official release of PIKAIA 1.0 in January 1996 several minor bugs and a few not-so-minor ones have been reported by users. All users who were registered at the time were notified of the corrective actions to be taken. The following give access to the original texts of the bug reports issued up to the present time. All bugs have been corrected as per the issue date of the corresponding bug report, as listed below:

REPORT ISSUED ON CONCERNING
No. 1     01/20/96    Incorrect declaration ordering in 
                      subroutine report

No. 2     02/14/96    Minor bug in random number initialization

No. 3     02/14/96    Assorted minor and one not-so-minor bugs
		      in some of the example drivers and
                      fitness functions

No. 4     05/20/96    Additional safeguard needed in
                      random number generator routine

8. PIKAIA in other computing languages

Thanks to the kind efforts and generosity of various PIKAIA users, versions in IDL and f90 are now available publicly for general use. An IDL transliteration of the original FORTRAN-77 has been produced by Sarah Gibson, and is available directly by clicking here. A more streamlined and efficient IDL version has also been produced by Scott McIntosh, and can be accessed by clicking here. Scott is currently (April 2003) fixing some last remaining odds and ends, but the code seems to be running OK.

For a FORTRAN-90 version, go to Alan Miller's software page and scroll down to Miscellaneous Software, looking for pikaia.f90

A parallel version of PIKAIA is now available, thanks to the efforts of Travis Metcalfe. Versions exist in both MPI and PVM, and can be found at Travis' Parallel PIKAIA Web Page. [Technical detail: only the Full-Generational-Replacement reproduction plan can be used in these parallel versions, as the steady-state plans are much harder to parallelize in a robust and efficient manner].

9. References and further readings

On the subroutine PIKAIA itself see:

Charbonneau, P., 1995, Genetic Algorithms in Astronomy and Astrophysics, The Astrophysical Journal (Supplements), 101, 309 [abstract]
Charbonneau, P., and Knapp, B. 1995, A User's guide to PIKAIA 1.0, NCAR Technical Note 418+IA (Boulder: National Center for Atmospheric Research) [Table of Content]
Charbonneau, P., 2002, Release Notes for PIKAIA 1.2, NCAR Technical Note 451+STR (Boulder: National Center for Atmospheric Research) [full text, postscript]
Charbonneau, P., 2002, An Introduction to Genetic Algorithms for Numerical Optimization, NCAR Technical Note 450+IA (Boulder: National Center for Atmospheric Research) [full text, postscript]; previously known as the Oslo Tutorial

On genetic algorithms in general, see:

Davis, L. 1991, Handbook of Genetic Algorithms (New York: Van Nostrand Reinhold).
Goldberg, D.E. 1989, Genetic Algorithms in Search, Optimization and Machine Learning (Reading: Addison-Wesley)
Holland, J.H. 1975, Adaptation in Natural and Artificial Systems (Ann Arbor: The University of Michigan Press; Second Edition 1992, MIT Press)
Mitchell, M. 1996, An Introduction to Genetic Algorithms (Cambridge: The MIT Press)
Baeck, T. 1996, Genetic Algorithms in Theory and Practice (Oxford: The Oxford University Press)
Michalewicz, Z. 1994, Genetic Algorithms + Data Structures = Evolution Programs (Berlin: Springer)

And finally, our not particularly well-informed choice for the top five books dealing with relevant aspects of evolutionary biology:

Dawkins, R. 1986, The Blind Watchmaker (New York: W.W. Norton)
Eldredge, N. 1985, Time Frames (New York: Simon and Schuster)
Gould, S. J. 1989, Wonderful Life (New York: W.W. Norton)
Maynard Smith, J. 1989, Evolutionary Genetics (Oxford: The Oxford University Press)
Stebbins, G. L. 1966, Processes of organic Evolution (Englewood Cliffs: Prentice-Hall)

Please address questions or comments to:

paulchar@astro.umontreal.ca.