Spatial processes and evolutionary models: a critical review

Evolution is a fundamentally population level process in which variation, drift, and selection produce both temporal and spatial patterns of change. Statistical model fitting is now commonly used to estimate which kind of evolutionary process best explains patterns of change through time, using models like Brownian motion, stabilizing selection (Ornstein-Uhlenbeck), and directional selection on traits measured from stratigraphic sequences or on phylogenetic trees. But these models assume that the traits possessed by a species are homogeneous. Spatial processes such as dispersal, gene flow, and geographic range changes can produce patterns of trait evolution that do not fit the expectations of standard models, even when evolution at the local-population level is governed by drift or a typical OU model of selection. The basic properties of population level processes (variation, drift, selection, and population size) are reviewed and the relationship between their spatial and temporal dynamics is discussed. Typical evolutionary models used in palaeontology incorporate the temporal component of these dynamics, but not the spatial. Range expansions and contractions introduce rate variability into drift processes, range expansion under a drift model can drive directional change in trait evolution, and spatial selection gradients can create spatial variation in traits that can produce long-term directional trends and punctuation events depending on the balance between selection strength, gene flow, extirpation probability, and model of speciation. Using computational modelling that spatial processes can create evolutionary outcomes that depart from basic population-level notions from these standard macroevolutionary models.

Methods
These data were produced by a computational agent-based model set on an island platform that was entirely exposed during lowstands with three peaks 2 m above platform height that remained exposed as isolated islands during highstands. The entire island platform was gridded into 5,000 cells (50 x 100), each of which could potentially be occupied by a local population.

Eustasy  (sea level change) was modelled as a sine wave through 2.5 cycles using the equation -5·Sin(0.005π·x – π) -4, where x is time measured in model steps (Fig. 3B). This equation causes sea-level to start 3 m below platform height (thus exposing the entire platform at the beginning of each model run), cresting at 2 m above platform height at highstands (thus inundating everything except the very peaks of the three islands) and dropping to 8 m below platform height during lowstands.

An evolving species was modelled through space and time using a metapopulation concept. At the beginning of each run, a single founder population was randomly placed on the platform with a local population size (N) of 10, a trait value of 0, a phenotypic variance (P) of 0.002, and a heritability (h2) of 0.5. Recalling that the rate of genetic drift is h2·P/N, drift in local populations (demes) was modelled as a random change drawn from a normal distribution with a mean of zero and a variance of 0.001. The heritability value was chosen realistic for morphological traits, but the other two values were arbitrarily chosen (see discussion below about consequences of higher or lower rates of drift). At each model step each local population: (1) reproduced a new generation in its own cell, undergoing a drift event at the same time; (2) had a chance of expanding its range by giving rise to a new local population in any or all of the adjacent cells (P = 0.8 for each cell); and becoming locally extinct (P = 0.2). Each new population, like its parent, had N=10, P=0.2, h2=0.5, and a trait mean equal to the parent’s after a new drift event. If a population expanded into an already occupied cell, the two populations are merged through reproduction by averaging their phenotypes and resetting N=10. The averaging of traits in cohabiting populations creates gene flow since a phenotype in one area can spread via dispersal and reproduction to another. If the cell was covered by water the probability of extirpation rose to 1.0. 

The geographic range of the species expands when populations move into empty cells and contracts when populations are extirpated. Sea-level thus changes the geographic range and the number of extant local populations. At lowstands the species tends to grow to approximately 5,000 local populations that cover the entire platform, and at highstands it contracts to as few as 12 populations subdivided between the three isolated islands. 

The computational model was run in Mathematica (Wolfram Research, 2018) with the aid of the Phylogenetics for Mathematica 5.1 package (Polly, 2018a) and the Quantitative Paleontology for Mathematica 5.0 package (Polly, 2016). Code for the model and full outputs of all nine model runs can be found in associated data S2 of Polly (2019).  The raw output is available in this archive. Models were run on Indiana University’s Karst high-throughput computer cluster. 

A full description of this data set can be found in Polly (2019).

Polly, P.D. 2016. Quantitative Paleontology for Mathematica. Version 5.0. Department of Earth and Atmospheric Sciences, Indiana University: Bloomington, Indiana. https://github.com/pdpolly/Quantitative-Paleontology-for-Mathematica 

Polly, P.D. 2018a. Phylogenetics for Mathematica. Version 5.0. Department of Earth and Atmospheric Sciences, Indiana University: Bloomington, Indiana. https://github.com/pdpolly/Phylogenetics-for-Mathematica 

Polly, P. D. 2019. Spatial processes and evolutionary models: a critical review. Palaeontology, 62: 175-195. (10.1111/pala.12410)