CORNELL UNIVERSITY

“SYNTHETIC” EVOLUTION

SUBMITTED TO PROFESSOR WILL PROVINE FOR HISTORY 287

DEPARTMENT OF HISTORY

BY WILL REGAN

ITHACA, NEW YORK DECEMBER 2005

INTRODUCTION Despite all the fine evidence Darwin cites in On the Origin of Species, he never gives an example of evolution in action. He could not tell us anywhere to look in nature in order to clearly see an animal evolving. But of course, this is an impossible task. Due to the slow, incremental nature of evolution (in human terms, not geological terms), a human observer would have to live several lifetimes to see a clear change in a species. Thus, even with the large body of biological research available today, we still cannot give empirical “proof” of the action of evolution. Within the framework of the natural world, we can only infer that the agent of change is evolution by common descent through natural selection. But these days, there exist other potential ways of studying evolution. With the great processing power of modern computers, we can make interesting models of previously incalculable physical systems. I don’t feign that we can accurately model even a small organism down to the atomic level. For instance, a typical person composed of about 1027 atoms could only be “completely described” by a computer with more disk space than the combined storage of all existing computers. But we can often make very predictive guesses of natural systems with computer simulation. Thus, it seems natural to reason that we can simulate our paradigm of evolution in a simulated environment, an experiment in which (since we control the rate of time) we can live to see the results. Luckily for me, a recent field of engineering known as “genetic programming” (GP) uses a kind of synthetic evolution to “evolve” optimal solutions to problems. In this project, I will try to tailor the methods used in GP to our model of natural evolution. Then, I will let GP be the agent of change in a very simplistic (but hopefully logical) choice of synthetic “world.” I will perform pseudo-experiments to represent some of Darwin’s inferences about evolution, such as

the success of small and large groups, the interaction between species, and other such inferences. This experiment, a vast oversimplification of nature, might very likely yield complete garbage or results that seem contradictory to what we expect of natural evolution. But if we see a beautiful correlation to Darwin’s predictions, it would lend more feasibility to the theory that natural selection is the agent of evolution. SIMPLISTIC COMPUTER SIMULATION A computer simulation was devised, using a GA strategy proposed in Mitchell.1 This simulation examined the effect of population size on a species’ survival. A population of N individuals, where each individual is defined by a 20 bit binary number, must adapt to an environment where fitness is measured as the number of bits held in common with a fixed 20 bit number (the “perfect individual” given the current environment). Thus, an average individual would have a fitness of 10 and a superb individual would have fitness close to 20. In each generation, individuals produced offspring proportional to their fitness ratings, and all suffered mutations with some finite (small) probability. The bulk of this code is found in the Appendix, though some modifications were made between tests. First, the code was used to observe the average fitness over time of groups subject to mutation and natural pressures to evolve. The following graph shows the rates of evolution of three groups, of sizes N = 10, 50, and 100. Since the creatures were presumed to die at the end of their mating cycle, N remained fixed (roughly) for each group over all time.

1

Mitchell, 31.

Clearly, as Darwin would have reasoned, the initial rates of evolutionary progress went far quicker for larger groups. This can be explained by the fact that all groups suffer mutations with the same probability, and larger groups would tend to accumulate favorable changes more rapidly simply because there are more individuals to act on. But what happens in the long term for these groups? Carrying the simulation out to a large 10,000 generations, I observed the following behavior:

In the long term, we see N = 10 averaging about 13 fitness, N = 50 averaging about 18 fitness, and N = 100 averaging almost perfect fitness. This long term behavior can best be explained by the insulating effects of large populations. Once populations achieve near perfection, they can only be hurt by random mutations. Thus, smaller populations cannot perform nearly as well, since even the most unfit must breed to ensure the species’ survival. Again, Darwin’s deductions prevail. Next, the code was modified to see how fast populations could respond to changes in fitness. The response rate to changes in the environment (and hence changes in the ideal fitness), which I chose as the time required to achieve an average population fitness of 17, was plotted for a large range of N.

This behavior closely resembles an exponential decay. Thus, groups of size N are able to adapt at roughly (where K and c are constants of the system): T = K e – N/c This again clearly favors large groups. But this does not mean all species will want to grow without bound. If we represent natural checks on population size by multiplying this rate of adaptation by N (for instance, if the time required to mate in a large population is proportional to N), we observe the following behavior:

At low values of N, the scaled response time drops dramatically with increasing N. However, as N goes above 80 or so, the scaled response begins to increase again. Thus, even in a very crude model, we see that natural checks will conspire to make the ideal population size below some limit. INTER-SPECIES INTERACTIONS

It became quickly evident that studying inter-species competition and cooperation by computational genetic algorithms would require more time than was available. However, many researchers have done work in this field. For instance, Alexander Bisler has done interesting simulations of the coevolution of insects and flowers. In his symbiotic model, where insects are evolving differing sight abilities and flowers are evolving their shape and color, he finds interesting naturalistic patterns of specialization. He notes, “The specialization process leads to an increased productivity for the insects (more rewards per visits) and to an increased efficiency of pollination which is a benefit for the flowers.” 2 The graph below shows the progression of insect visits for 4 different flower species during coevolution.

3

Other researchers, such as Vladimir Levchenko, have even simulated long-term evolutionary trees. Levchenko produced a program known as Macrophylon, which simulates the evolution and branching of species in environments of varying competition. Below is a diagram

2 3

Bisler, 6. Ibid, 8.

of some evolutionary trees produced by Macrophylon for hard, medium, and soft amounts of competition.4

The work of these and other genetic algorithm researchers often produces realistic results useful for modeling nature. The fact that genetic algorithms are consistently performed with random perturbations, a framework consistent with Darwin’s descent with modification via natural selection, lends more weight to Darwin’s theory being the correct model. NATURAL SELECTION AND HEREDITY As shown above, Darwin’s assumptions held up very well in my simple simulations. But could they be tested on a more complex physical system? It would be interesting, for instance, to see if the inherited flapping-wing flight skills of birds could be generated on a mechanical bird merely by randomly varying the controls and letting natural selection take “control.” And this is exactly what I investigated by performing a genetic algorithm (GA) on an ornithopter (flapping

4

Levchenko.

wing aircraft), which I have been building for a mechanical engineering research project (MAE 490).5 A quick summary of the mechanical bird’s construction will also be given. Genetic Algorithm The physics of controlled flapping flight is quite complicated. The chances of manually designing a successful flapping pattern, not to mention the mechanical design to achieve that motion, are very slim. Yet there are thousands of species of insects, birds, and bats that are all very well adapted to flapping flight. So it is very suitable to use evolution on such a problem. The strategy of the GA for this problem is very simple. Behavior (a cycle of wing movements for flapping) is encoded into a “genome.” A large population of genomes are created and are evaluated (via simulation or testing) for their “fitness,” which would be the maximum average lift or height attained by the craft for the particular flapping genome. After all genomes have been tested for fitness, they undergo a process of mutation (slight random modification) and recombination with other genomes, with the fittest genomes being able to breed the most. In this way, desired behavior is eventually developed. The specifics of our GA schemes (deterministic crowding and elitism) will be omitted for brevity. Simulation To begin the task of evolving wing patterns, we initially performed a series of computer flight simulations, using a C++ physics library known as ODE. We assumed a simple construction of 4 wings flapped symmetrically by three servo motors (connected by cables to three points on each wing, giving the wings well-defined motion in three directions). Aerodynamics was modeled using the following simplified drag equation:

Note: I’ve been working closely with another undergraduate, Floris van Breugel, so any references to “us” or “our” are meant to refer to both him and me. Also, this work is in the process of being published in a conference journal, so please don’t share the results with other people. I trust you, Luana . 5

Drag 

1 2 C d A V 2

To keep the simulation somewhat realistic (so we could use the results for our physical model), full gravity and reasonable weights, sizes, torques and flapping frequencies were used. A snapshot of the simulated model is seen above. A symmetrical 4 wing arrangement was chosen for its stability (ideally, there should be no torques on the system). An arbitrary fixed mechanical setup was used and Bezier curve patterns were evolved to control the servo motor arm positions. Using a standard elitism genetic algorithm (EGA), hovering flight was achieved in simulation. With a population of 100, a total craft weight of 310 grams, the fittest member flew upwards with a maximum total lifting force of 190 N and a maximum upwards velocity of 12 m/s. The wing acceleration was sometimes quite large and the flapping frequency approached a brisk 6Hz, but overall the setup seemed likely to be capable of producing lift with a real model. Physical Model Given the encouraging results from ODE simulation, we set upon designing and constructing a physical ornithopter based on the four-winged cable-controlled simulations. The main difficulty in designing a feasible flapper was finding servo motors with a high enough power to weight (P:W) ratio. Most available servos and DC motors have awfully low P:W ratios, which is why many lightweight electromechanical machines tend to rely on light actuators like shape memory alloys and Piezo-electric elements. Eventually, we did manage to find servos that were both strong and quick enough for our needs and only contributed about 100 g to the total weight. 6

6 As a rough measure of servo fitness, we required that the servos deliver energy to the air (through drag) at twice the rate needed to cover losses in gravitational potential (GPE) for freefall; that is, in one flapping period, the drag’s lifting power would be more than double the GPE lost in one period of freefall.

We found light mylar wings (~ 15 g each) from a biplane kit that provided a large 2200 cm2 of total wing area. We added a 40 g Lithium Polymer (LiPoly) 7.4V battery to provide onboard power and soldered together a 25 g control circuit (Basic Stamp and PWMPAL) to control the servos. The rest of the craft was then designed in SolidWorks CAD and printed with lightweight ABS plastic on an FDM 3000 Rapid Prototyper (~ 80 g). These estimates brought the total weight to about 313 g, making actual flight seem within reach. Design Beyond being strong enough to withstand the forces of flapping and having a small footprint relative to the wing area, the body shape was somewhat irrelevant. The main consideration was the method of actuation. To simplify construction, we only made ½ of the 4wing ornithopter, with the intent of building another half if needed for stabilization. Thus, the physical model was comprised of three servos symmetrically controlling the three degrees of freedom (DoF) of two wings, with the wings arranged about the body's center like a bird. Pushrod Actuation We tried a direct form of actuation faithful to the simulated model. The simplest solution found was rigid pushrod connections from the servo horns to points on the wing frame. To deliver maximal forces for flapping, we placed two strong servos symmetrically about the center of the craft, directly under the wings. The third servo (for yaw control) was set at the

Ornithopter and computer interface

equilibrium height of the wings away from the back and far enough to be out of the way. The wings were connected to a central tower and elevated so that the lifting forces were applied well above the craft’s center of mass for stability. This model was built and can be seen in the above figure. The total weight of the craft was roughly 450 grams, or about 1 pound. Interfacing To evolve flapping on the physical model, we needed an interfacing loop. First, each flapping genome was translated into PBasic servo commands, which were sent to and executed by the onboard Basic Stamp via the serial port. Flapping cycles were repeated up to 8 times to get a good sense of the lifting ability of a particular pattern. During execution of the commands, we measured the vertical forces on the craft with a 1DoF LCEB-10 Omega Load cell (capable of measuring +/- 10 pounds of force). Though more might be learned about the system with a three dimensional multiaxis force and torque cell, or perhaps tilt sensors in all directions, we reasoned that the best hovering individuals must maximize their lift and minimize all other forces and vibrations. These force measurements were buffered and sent back to the computer through a Hytek Automation U120816 12-bit USB data acquisition module. The fitness assigned to each genome was the average lift measured during the flapping. Using these fitness measurements, evolution would then select the best genomes to reproduce. Results Initially, we evolved Bezier curves from which we sampled servo commands, but soon found it was difficult to get smooth flapping. We reasoned that the Bezier curves (which have fixed periods) were simply driving the wings at a frequency not quite at the natural frequency of the wings. We found, however, that genomes composed of random data could be evolved to

allow for smoother motion. Also, random data would be more faithful to Darwin’s natural selection process, which should be randomly determined. We next tested for noise resistance by running a known set of 10 individuals 5 times and comparing the fitness from each run. The variance in these force measurements, seen below, was small enough to convince us that the fitness was indeed meaningful. Fitnes s Vs. Individual for 10 Individuals 0.58 0.56 0.54

Fitness

0.52 0.5 0.48 0.46 0.44 0.42 0.4 0.38 1

2

3

4

5

6

7

8

9

10

Individual

Finally, we began the actual evolution, letting the algorithm run many generations, where each individual in the population is tested and some produce new offspring. The results are shown below, where hovering fitness is measured as the Average Percent of Total Weight Supported (0 = Free Fall, 100 = Hovering). As the data shows, fitness improved dramatically at first, stagnated, and later began to again improve quickly. This is best explained by the fact that

Fitness vs. Generation: Blue dots represent individual fitness, and the red line represents average generational fitness.

hovering flight is a difficult problem, and with a small population size (we used 40 individuals) it can be difficult for evolution to find a good solution. As the fitness approached 20% lifting ability, epoxy on the wing frame broke and forced us to stop the evolution. These results are still quite encouraging. When the evolution terminated, the fitness showed no signs of slowing down. With tweaks to our servos and frame, our GA, and a good bit of epoxy, we expect to achieve at least 50% lifting ability. Furthermore, real birds generate most of their lift from the lift on airfoils resulting from forward motion. So if we were to cut this loose and let it move forward, given more suitable wings, it could possibly already fly like a bird. Most importantly, we observed some very interesting naturalistic flapping patterns. In one such behavior, a small phase difference between the two powerful servos caused the wings to present the full surface area on the down stroke and an angled (and therefore smaller) area on the upstroke, yielding a net lifting force. This pattern can best be described as a figure-8 motion, akin to the flapping patterns exhibited by hummingbirds and dragonflies. In another pattern, the wings performed a sort of swimmer’s frog kick, with a slow up-stroke and quick down-stroke. These types of behaviors are exactly what I was looking for. Given the bird-like construction, a completely random process of evolution (i.e. simulated natural selection in a real environment) resulted in the production of flapping patterns quite similar to real birds and flying insects. In other words, this work shows that random processes could have been used when real birds first began to evolve the ability to fly. Darwin’s assumption that descent with modification via natural selection can evolve the complex behaviors observed in nature has again held its weight. CONCLUSIONS This paper has explored the similarities in the behavior of “synthetic” evolution (genetic algorithms) and the predictions of Darwin’s theory of descent with modification via natural

selection. For the very simplistic genetic algorithms I created, I observed behavior that Darwin would have predicted: the success of large groups adapting quicker than small groups, the greater health of established large groups when subject to mutation, and even the checks on population growth imposed by simple environmental limits. Furthermore, I demonstrated an example of a more complex evolution on my mechanical ornithopter. Through a similar synthetic natural selection process, wing patterns were evolved that directly paralleled those seen in nature, in bird and insect flight. Again in agreement with Darwin’s theory, this parallel suggests that the “natural selection pressure” could have been the agent of change in the original evolution of bird flight.

REFERENCES Augustsson, P., Wolff, K., and Nordin, P. 2002. Creation of a Learning, Flying Robot by Means of Evolution. In Proceedings of the Genetic and Evolutionary Computation Conference (GECCO 02). Bisler, Alexander. 2005. Simulation of the Coevolution of Insects and Flowers. In Proceedings of the VIIIth European Conference on Artificial Life (ECAL 2005). http://www.ecal2005.org/workshopsCD/systemsbiol/pdf/sysbio_bisler_coevolution.pdf http://www.eca

Levchenko, Vladimir F. Computer Simulation of Evolution. Accessed Online 11/30/05. http://www.evol.nw.ru/lew/macrophylon/article.htm http://www.evol.nw.ru/lew/macrophylon/macrophylon.htm

Michelson, R. and Naqvi, M. 2003. Beyond Biologically-Inspired Insect Flight. In Proceedings of von Karman Institute for Fluid Dynamics RTO/AVT Lecture Series on Low Reynolds Number Aerodynamics. Brussels, Belgium. Mitchell, Melanie. 1998. An Introduction to Genetic Algorithms. Cambridge, Mass.: MIT Press. Nolfi, S. and Floreano, D. 2000. Evolutionary Robotics: The Biology, Intelligence, and Technology of SelfOrganizing Machines. Cambridge, Mass.: The MIT Press. Smith R. 2000. Open Dynamics Engine (ODE), http://ode.org/ http://ode.org/. van Breugel, F. 2005. Evolving Buildable Flapping Ornithopters. In Proceedings of the Genetic and Evolutionary Computation Conference (GECCO 05).

APPENDIX //ga.cpp Genetic Algorithm to Simulate Evolution on Population of N Individual #include #include #include #include #include #define ARRAYT_BOUNDS_CHECK //Array class courtesy of Professor Earl Kirkland, A&EP Cornell #include "array.hpp" static int N=100,maxgen=1000; double pc=0.7, pm=0.001;

//force N to be even for ease of reproduction //probs of crossing and mutation

int seed=time(NULL); double ran1 (int &idum) /* Random number generator. Returns random number of double precision between 0 and 1. Source: Numerical Recipes in C++, pg 284 */ { const int IA=16807,IM=2147483647,IQ=127773,IR=2836,NTAB=32; const int NDIV=(1+(IM-1)/NTAB); const double EPS=3.0e-16,AM=1.0/IM,RNMX=(1.0-EPS); static int iy=0; static array iv(NTAB); int j,k; double temp; if (idum <= 0 || !iy) { if (-idum < 1) idum=1; else idum = -idum; for (j=NTAB+7;j>=0;j--){ k=idum/IQ; idum=IA*(idum-k*IQ)-IR*k; if (idum < 0) idum += IM; if (j < NTAB) iv(j) = idum; } iy=iv(0); } k=idum/IQ; idum=IA*(idum-k*IQ)-IR*k; if (idum < 0) idum += IM; j=iy/NDIV; iy=iv(j); iv(j) = idum; if ((temp=AM*iy) > RNMX) return RNMX; else return temp; } //end ran1() void fitness (array& popfit,array& pop,int& fitsum) { int i,j; fitsum=0; for (i=0;i
} } } void init (array& pop) { int i,j; for (i=0;i0.5) pop(i,j)=0; else pop(i,j)=1; } } void weight_probs (array& probdis,array& popfit,int fitsum) { int i,j,curr=0,weight; for (i=0;i& pop,array& newpop,int offsp) { int i,cp; //cp = cross point is in [0,...,19] cp=rand()%20; for (i=0;i& pop,array& newpop,int offsp) { int i; for (i=0;i<20;i++){ newpop(offsp,i)=pop(p1,i); newpop(offsp+1,i)=pop(p2,i); } } void mutate (int offsp,array& newpop) { int i; for (i=0;i<20;i++){ if(ran1(seed)
//invert from 0,1 to 1,0 //invert from 0,1 to 1,0

} int maxft (array& popfit) { int i,cmax=popfit(0); for (i=1;icmax) cmax=popfit(i); return cmax; } int main() { int i, j, fitsum, gen, offsp,p1,p2; do { cout << "\n Enter N, maxgen: "; cin >> N >> maxgen; } while (N<1 || maxgen<1); array pop(N,20), newpop(N,20); array popfit(N); array probdis(N*20); //holds indices which weight occurrence of parents by fitness //N*20=max total fitness array fit_time(maxgen); array max_fit(maxgen); FILE* fp; fp = fopen("avgfit.dat","w+"); if( fp == NULL ) { printf("Can't open file."); exit(0); } FILE* fo; fo = fopen("maxfit.dat","w+"); if( fo == NULL ) { printf("Can't open file."); exit(0); }

/*

srand(seed); do { cout << "\n Enter pc, pm: "; cin >> pc >> pm; } while (pc<0 || pc>1 || pm<0 || pm>1); */ init(pop); gen=0; do { fitness(popfit,pop,fitsum); max_fit(gen)=maxft(popfit); fit_time(gen)=double(fitsum)/double(N); weight_probs(probdis,popfit,fitsum); //weight chance of choosing parent by fitness offsp=0; //start creating offspring while (offsp
//direct copy

mutate(offsp,newpop); offsp+=2; } pop=newpop; gen++; } while (gen<(maxgen-1)); for (i=0;i