How morphological development can guide evolution

Organisms result from adaptive processes interacting across different time scales. One such interaction is that between development and evolution. Models have shown that development sweeps over several traits in a single agent, sometimes exposing promising static traits. Subsequent evolution can then canalize these rare traits. Thus, development can, under the right conditions, increase evolvability. Here, we report on a previously unknown phenomenon when embodied agents are allowed to develop and evolve: Evolution discovers body plans robust to control changes, these body plans become genetically assimilated, yet controllers for these agents are not assimilated. This allows evolution to continue climbing fitness gradients by tinkering with the developmental programs for controllers within these permissive body plans. This exposes a previously unknown detail about the Baldwin effect: instead of all useful traits becoming genetically assimilated, only traits that render the agent robust to changes in other traits become assimilated. We refer to this as differential canalization. This finding also has implications for the evolutionary design of artificial and embodied agents such as robots: robots robust to internal changes in their controllers may also be robust to external changes in their environment, such as transferal from simulation to reality or deployment in novel environments.


Rigid-bodied robots.
Rigid-bodied robots and their environment were simulated using Pyrosim (ccappelle.github.io/pyrosim). The robot is a quadruped with a large, spherical abdomen; each leg is attached by a single degree-of-freedom hinge joint.
Morphological development was approximated in rigid bodies using linear actuators to slowly lengthen or shorten the length of each leg, from an evolved starting value (between 0 and 1) to an evolved final value (between 0 and 1). The controller is a simple neural net: two central pattern generators are fully connected to four motor neurons, each of which innervate a separate hinge joint. Controller development was approximated in neural networks through ballistic change to each synaptic weight: As the simulation proceeds, each weight develops linearly, from an evolved starting value (between -1 and +1) to an evolved final value (between -1 and +1).
The genotype spans two arrays: one for initial and final synaptic weights (controller), and another for initial and final leg lengths (morphology). Mutations affect, on average, a single element in each array. Apart from the genotype and its mutations, the evolutionary algorithm is identical to that of the soft robots. However, the task environment now consists of a sloped floor, declined toward a light source; and performance is measured by the average light intensity recorded by a light sensor embedded in the center of the agent's abdomen, according to the inverse square law of light propagation, at each time step in its life. Occlusion of the light caused by interference of the robot's own body parts was not simulated.
The results are presented below in Supplementary Fig. S5.

Mutations for soft robots.
The following derivation shows that there is a negligible difference in the mutations produced by the Evo and Evo-Devo treatments, in terms of the number of voxels modified (Fig. S4).
Each voxel cell of a soft robot has its own material properties that can be changed by the evolutionary algorithm. Evo voxels have two material properties: (1) resting length and (2) phase offset. Evo-Devo voxels have four material properties: (1) initial resting length, (2) final resting length, (3) initial phase offset, and (4) final phase offset.
Mutations are applied by first choosing which material properties to mutate, and then choosing, separately for each property, which voxels to modify. For each of the n material properties, we select it with independent probability p = 1/n. If none are selected, we randomly choose one. This occurs with probability (1 − p) n . Hence the number of selected material properties for mutation is a random variable S which follows a truncated binomial distribution, The expected number of selected material properties is then: For a selected material property, each voxel is mutated independently with probability λ , a hyperparameter we call the mutation rate. The expected number of genotype elements mutated given K total voxels is thus: Dividing by the length of the genome, nK, we get the expected proportion of genotype elements mutated: Note that imposing bilateral symmetry does not change these expected values.
We have K = 48 total voxels, and n = {2, 4} material properties for our two main experimental treatments {Evo, Evo-Devo}, respectively. The expected difference between a robot and its offspring, in terms of genotype elements, is summarized in the following table. n = 2 n = 4 δ gene 60λ 63.8175λ π gene 0.625λ 0.3291λ However, because multiple material properties can be mutated within a single voxel, the expected number of voxels mutated is lower than the expected number of genotype elements mutated. To calculate the average number of voxels mutated we need to consider a hierarchy of binomial distributions. Given that S material properties were selected for mutation, the number of material properties mutated within a single voxel, M follows a binomial distribution, For brevity, let's denote the probability that at least one mutation occurs within the voxel as θ , Then the number of voxels mutated, V , across a total of K voxels and S selected material properties, also follows a binomial distribution:

3/9
And the expected number of voxels mutated (out of K total) is There is an extremely tight bound on the proportion of voxels mutated, π vox = δ vox /K, for any n > 1 (Fig. S4). Thus mutations in Evo (n = 2) and Evo-Devo (n = 4) have practically the same impact in terms of the number of voxels modified. For completeness, the following table displays δ vox for the specific values of λ considered by our hyperparameter sweep (K = 48) (Fig. S3) Figure S2. Experimental treatments. The phase of an oscillating global temperature (A, F) is offset in the k-th voxel by a linear function from φ k to φ * k (B, G). The resting length of the k-th voxel is a linear function from k to * k (C, H). For Evo, there is no development, so φ k = φ * k and k = * k . The offset actuation is added on top of the resting length to give the current length of the k-th voxel (D, I). These example voxel-level changes occur across ontogenetic time (t), independently in each of the 48 voxels, and together interact with the environment to generate robot-level velocity (E, J). To see this, we averaged displacement across intervals of two actuation cycles (0.5 sec) and plotted the smoothed velocities for two Evo-Devo run champions with minimal canalization (J) alongside their non-developmental counterparts (E). Also note that the frequency of actuation is plotted here at 1.4 Hz; but, in our experiments, we used a frequency of 4 Hz. Figure S3. Mutation rate sweep. Median fitness (with 95% bootstrapped confidence intervals) under various mutation rates, λ , a hyperparameter defined in Supplementary Methods which affects the probability a voxel is mutated. In the main experiment of this paper, the mutation rate is evolved for each voxel independently, and is constantly changing. In this mutation rate sweep, λ is held uniform across all voxels. There were two observed basins of attraction in terms of fitness: a slower design that trots/gallops 5-15 body-lengths during the evaluation period, and a faster design type that rolls at 50-70 body-lengths. Although higher mutation rates facilitate the discovery of the superior phenotype, once found, lower mutation rates tend to produce more refined and faster robots. Without development, the search space has a single spike of high fitness. One can not do better than random search in such a space. When λ = 1, optimizing Evo morphologies reduces to random search, and this is the only mutation rate where Evo does not require significantly more generations than Evo-Devo to find the faster design type. This can be observed for λ ∈ {1/6, 1/3, 1/2, 2/3, 1}, by comparing the generation at which the slopes of the fitness curves increase dramatically. However, the best two treatments (Evo-Devo at λ = 1/2 and λ = 2/3), as measured by the highest median speed at the end of optimization, have development, and the robots they produced are significantly faster than those produced by random search (Evo with the highest mutation rate) (p < 0.01). 7/9 Figure S4. Mutational impact. The expected proportion of voxels modified, π vox , where n is the number of material properties that can be mutated, and λ is the mutation rate. A derivation is provided in Supplementary Methods. Regardless of n, when λ = 1, every voxel must be mutated, and when λ = 0, no voxels can be mutated. Between these two points, there is an extremely tight bound on the proportion of voxels mutated for any n > 1. In this paper, we have treatments Evo, with n = 2, and Evo-Devo, with n = 4. 8/9 Figure S5. Differential canalization in rigid bodies. This environment consists of a sloped floor, declined toward a light source. Rigid-bodied robots, which are described in Supplementary Methods, perform phototaxis. Although longer legs produce faster walking behaviors, the shortest leg setting places the robot's large, spherical abdomen onto the ground, allowing the robot to simply roll down the ramp toward the light. An ancestor discovers rolling over toward the end of its evaluation period through ontogenetic morphological change that compresses leg lengths. This rollable morphology is assimilated to the start of development through differential canalization: The sooner a robot shrinks its legs, the sooner it can begin rolling; eventually descendants start with compressed legs, are able to immediately roll, exhibit little to no morphological change, but continue to sweep over many different synaptic weights as they behave. Rolling down the slope is not entirely passive, however, since protruding legs, even at their shortest setting, can slow down or even stop forward movement. The best controllers not only avoid such interference, but also guide rolling toward the light source. However, development in this particular case did not affect evolvability. These results are consistent with predictions made by other quantitative models used in theoretical biology that suggest that plasticity only expedites evolution under restrictive conditions 1 . We hypothesize that, in our case, this is because the search space is not deceptive enough: Once the rigid-bodied robot evolves to compress its legs, and touch its abdomen to the sloped floor, it will tend to roll for the remainder of its evaluation period. That is, in contrast with the soft robots, there is no intermediate stage between walking and rolling that suffers the fitness penalty of no longer being able to move.