Evolution repeats itself in replicate long-term studies in the wild

The extent to which evolution is repeatable remains debated. Here we study changes over time in the frequency of cryptic color-pattern morphs in 10 replicate long-term field studies of a stick-insect, each spanning at least a decade (across 30 years of total data). We find predictable ‘up-and-down’ fluctuations in stripe frequency in all populations, representing repeatable evolutionary dynamics based on standing genetic variation. A field experiment demonstrates that these fluctuations involve negative frequency-dependent natural selection (NFDS). These fluctuations rely on demographic and selective variability that pushes populations away from equilibrium, such that they can reliably move back towards it via NFDS. Finally, we show that the origin of new cryptic forms is associated with multiple structural genomic variants such that which mutations arise affects evolution at larger temporal scales. Thus, evolution from existing variation is predictable and repeatable, but mutation adds complexity even for traits evolving deterministically under natural selection.

Methods
We compiled data on morph frequencies in Timema cristinae using samples collected in the spring using sweep nets between 1990 and 2023. All individuals were scored as 'striped', 'unstriped', or 'melanistic', or occasionally when it was difficult to distinguish between the first two categories as 'intermediate-striped'.  Samples from 1990 to 1999 were collected and scored by Cristina Sandoval, who then trained PN in 2000. PN collected and scored most samples from 2000 to 2023. The host-plant collected on was recorded for all individuals. 

Datasets from two field experiments are included as well.

Experiment 1. A total of 513 individuals were collected from populations LA (latitude 34.41, longitude -119.80), PRC (latitude 34.53, longitude -119.86), PRNC (latitude 34.53, longitude -119.85),  VPC (latitude 34.53, longitude -119.85), and WTA (latitude 34.52, longitude -119.78) between April 30th and May 4th, 2021. Numbers were as follows: green-unstriped morphs, LA 17, PRC 64, PRNC 58, VPC 71, WTA 48; striped morphs, LA 136, PRC 0, PRNC 3, VPC 26, WTA 90. Individuals were randomly assigned to one of 21 treatments (n = 20 stick insects per treatment), which varied the initial stripe frequency from 0% (0 striped, 20 unstriped) to 100% (20 striped, 0 unstriped) in 5% increments. Each of these treatments of 20 individuals was then randomly assigned to one of 21 experimental bushes (near site N2, in the general area of latitude 34.51 and longitude -119.80). Each bush was cleared of existing T. cristinae  by sampling it before release. Individuals were released on the morning of May 5th by dumping them from containers; the insects cling to the branches when released in this way because of a strong clinging reflex. Individuals were recaptured using visual surveys and sweep nets on the evening of May 7th and scored as striped or unstriped.

Experiment 2. We conducted a second field transplant experiment to ask whether striped T. cristinae morphs exhibit higher fitness on Adenostoma than melanic morphs. For this experiment, stick insects were collected from a single location, PRNC (latitude 34.53, longitude -119.85). Stick insects for each pair of eight experimental blocks were captured on a single day in 2023: May 9th for blocks 1 and 2, May 11th for blocks 3 and 4, May 13 for blocks 5 and 6, and May 15 for blocks 7 and 8. For each block, 10 striped and 10 melanic stick insects (n = 20 total stick insects per block) were released onto a single Adenostoma bush at the experimental transplant site (latitude 34.51, longitude = -119.80, as in the first experiment) the day after they were collected. Two (blocks 1-6) or three (blocks 7 and 8) days after release, surviving stick insects were recaptured using visual surveys and sweep nets and scored as striped or melanic. 

We used genome-wide association (GWA) mapping to identify regions of the genome associated with color-pattern (striped vs. unstriped) variation in T. cristinae. For this, we used a previously published DNA sequence data set comprising partial genome sequences (genotyping-by-sequencing or GBS data) for 602 T. cristinae from a single population (FHA, latitude 34.52, longitude -119.80).