Data and code from: Speciation by local adaptation and isolation by distance in extended environments
Data files
Sep 22, 2026 version files 281.57 MB
-
auxiliar_functions.h
7.13 KB
-
DATA.zip
281.43 MB
-
graph.h
8.46 KB
-
main_text.ipynb
45.67 KB
-
main.cpp
3.97 KB
-
node.h
3.09 KB
-
README.md
14.75 KB
-
speciation_model.h
18.43 KB
-
supplementary_material.ipynb
40.52 KB
Abstract
Speciation is often associated with geographical barriers that limit gene flow. However, species can also emerge in parapatry, even in homogeneous environments, through spatially restricted mating and limited dispersal. When the environment is not homogeneous, natural selection contributes to differentiation by local adaptation and tends to facilitate speciation. To explore how isolation by distance and adaptation combine to determine species diversity, we propose a model regulated by these two components. The former is implemented via mating restrictions on spatial proximity and genetic similarity, whereas the latter is realized by an ecological phenotype subjected to adaptation by natural selection. We consider a scenario where the environment has two distinct optima, and compare the resulting diversity patterns and phenotypic distributions with those of a homogeneous environment, with a single ecological optimum. We show that the interplay between selection and isolation by distance affects not only species formation but also phenotypic distributions and the timing of diversification. Simulating individuals with either restrictive or permissive mating regimes, combined with strong or weak selection, we show that: (i) environmental selection can accelerate diversification but is not always necessary for reproductive isolation; (ii) bimodal phenotypic patterns can arise through different evolutionary pathways and are therefore not necessarily signatures of ecological speciation; and (iii) when selection is weak and mating is restrictive, parapatric speciation begins before pronounced ecological differentiation, while the phenotypic peaks can oscillate and fail to reach a stationary state.
Dataset DOI: 10.5061/dryad.s1rn8pkqm
Speciation Dynamics Model
This dataset contains a C++ simulator implementing an individual-based model (IBM) of speciation in sexually reproducing populations. Individuals occupy a discrete two-dimensional lattice and carry reproductive and phenotypic genomes. The model describes reproduction, mutation, spatial dispersal, environmental selection, genetic similarity, and species formation.
The simulations are used to investigate how spatial structure, environmental heterogeneity, and genetic differentiation contribute to the emergence and maintenance of species. The model is also used to examine how environmental differences and spatial structure influence phenotypic differentiation, genetic divergence, and species formation.
Two environmental scenarios are considered. In the heterogeneous scenario, the landscape is divided into two regions with different phenotypic optima. In the homogeneous scenario, the same phenotypic optimum is applied across the entire landscape.
Project Structure
The dataset is organized into C++ simulation code, simulation data, and analysis notebooks:
- C++ source files: The C++ source files (
.hand.cpp) are located in the root directory:node.h: Defines individual nodes in the spatial population, including their identifiers, spatial coordinates, and spatial neighbors.graph.h: Defines the spatial structure of the simulation, including the two-dimensional lattice, spatial neighborhoods, and offspring dispersal region.speciation_model.h: Implements the evolutionary dynamics of the model, including reproduction, mutation, environmental selection, genetic similarity, and species identification.auxiliar_functions.h: Provides functions for exporting simulation results and other simulation outputs.main.cpp: Initializes and runs the simulations and writes the simulation outputs.
DATA.zip: Contains the simulation outputs and processed data used to reproduce the analyses and figures presented in the manuscript.main_text.ipynb: Contains the analysis code used to generate the figures presented in the main text. The notebook does not contain embedded figure outputs.supplementary_material.ipynb: Contains the analysis code used to generate the figures presented in the supplementary material. The notebook does not contain embedded figure outputs.
1. Simulation of Speciation
During the simulation, individuals reproduce through mating with genetically compatible individuals. Offspring inherit genomic loci from both parents and may undergo mutation. Offspring are spatially dispersed within a defined radius, and their phenotype determines their compatibility with the local environment.
Genetic dissimilarity is used to define genetic connectivity between individuals and to identify species. The model therefore combines spatial structure, environmental selection, reproduction, mutation, and genetic differentiation to simulate the emergence and dynamics of species.
Main simulation parameters
The main simulation parameters are:
| Parameter | Description | Value / input |
|---|---|---|
InitialPop |
Initial number of individuals in the population | Command-line input |
G |
Maximum normalized genetic dissimilarity used to define genetic connectivity and species | Command-line input; fraction between 0 and 1 |
sigma |
Parameter controlling the strength and breadth of environmental selection | Command-line input |
S |
Spatial dispersal radius of offspring | 6.0 in the main simulations; varied in the dispersal-radius analysis |
SizeX |
Width of the lattice | 100 |
SizeY |
Height of the lattice | 100 |
MutationRate |
Probability of mutation at each genomic locus during reproduction | 0.00025 |
qmating |
Mating parameter affecting reproduction probability | 0.37 |
t |
Number of generations simulated | 5000 |
SizeRepGenome |
Number of loci in the reproductive genome | 1500 |
SizePheGenome |
Number of loci in the phenotypic genome | 100 |
Phenotype0 |
Reference phenotype used in environmental selection | 0.50 |
How to Run
Compile directly from main.cpp using:
g++ -O1 main.cpp -o speciation
Run a simulation using:
./speciation <InitialPop> <G> <sigma> <output_directory> [run_id]
2. Data and Figure Generation
The DATA.zip contains DATA/ directory contains the simulation outputs used directly by main_text.ipynb and supplementary_material.ipynb to reproduce the analyses and figures presented in the manuscript.
Not every parameter folder contains every file type described below. Only the files required for the corresponding analysis were retained.
Directory organization
The main directories and files are:
heterogeneous/: Data from simulations performed in heterogeneous environments, in which the landscape contains two regions with different phenotypic optima.homogeneous/: Data from simulations performed in homogeneous environments, in which the same phenotypic optimum is applied throughout the landscape.local/: Simulations used to examine species richness and local adaptation over a broader range ofGandsigmavalues.radiusS/: Simulations used to examine the effect of offspring dispersal radiusS, for both homogeneous and heterogeneous environments.sigma30/: Simulations withsigma = 30and different values ofG, for homogeneous and heterogeneous environments.binomial.dat: Theoretical binomial reference curve used in the phenotype-distribution analysis.
Folder naming convention
The .dat files are plain-text whitespace-delimited files, and the .csv files are standard comma-separated text files. They can be read with standard text-processing software or with Python packages such as NumPy and pandas.
For the main parameter-series folders, fractional parameter values are abbreviated in the folder names by removing the decimal point. For example:
G5corresponds toG = 0.05G10corresponds toG = 0.10G15corresponds toG = 0.15G20corresponds toG = 0.20G25corresponds toG = 0.25G30corresponds toG = 0.30
Similarly:
sigma5corresponds tosigma = 0.05sigma10corresponds tosigma = 0.10sigma20corresponds tosigma = 0.20sigma30corresponds tosigma = 0.30sigma40corresponds tosigma = 0.40
Combined folder names identify both parameters. For example:
G5sigma5corresponds toG = 0.05andsigma = 0.05G5sigma40corresponds toG = 0.05andsigma = 0.40G30sigma5corresponds toG = 0.30andsigma = 0.05G30sigma40corresponds toG = 0.30andsigma = 0.40
The top-level directory DATA/sigma30/ is an exception to this naming convention. In this specific directory name, sigma30 refers to the parameter value sigma = 30, rather than sigma = 0.30.
For dispersal-radius analyses, folders named S6, S12, S18, and S24 correspond to offspring dispersal radius of 6, 12, 18, and 24 lattice units, respectively.
Data files and variables
evolution_global*.dat
Whitespace-delimited time series from simulation replicates.
When a numeric suffix is present, for example evolution_global1.dat, evolution_global2.dat, etc., the suffix identifies an independent simulation replicate.
Each row corresponds to one simulation generation.
The columns are:
| Variable | Description | Units / interpretation |
|---|---|---|
Generation |
Simulation generation | Integer generation |
SpeciesCount |
Number of species identified in the population | Integer count |
PopulationSize |
Total number of individuals in the population | Integer count |
seed_final_*.dat
Whitespace-delimited files storing the final internal state of the std::mt19937 pseudorandom number generator for each independent simulation replicate. The * in the filename identifies the replicate number (for example, seed_final_34.dat corresponds to replicate 34). The integer values represent the serialized state of the random-number generator at the end of the simulation. This state can be reloaded to continue the pseudorandom number sequence from the point at which it was saved.
demographic.csv
Comma-separated time series used for demographic analyses.
The variables are:
| Variable | Description | Units / interpretation |
|---|---|---|
t |
Simulation generation | Integer generation |
varP |
Global variance of phenotypic values in the population | Dimensionless |
mean_fitness |
Mean environmental fitness of the population | Dimensionless, range 0–1 |
All/data[t].dat
Population snapshots at selected generations, where [t] indicates the simulation generation. For example, data5000.dat contains the population at generation 5000.
Each row represents one individual.
The columns, in order, are:
| Variable | Description | Units / interpretation |
|---|---|---|
id |
Individual identifier | Integer |
x |
x-coordinate on the spatial lattice | Lattice position |
y |
y-coordinate on the spatial lattice | Lattice position |
species |
Species identifier based on genetic connectivity | Integer categorical identifier |
degree |
Number of spatial neighbors of the individual | Integer count |
phenotype |
Fraction of phenotypic loci equal to 1 | Dimensionless, range 0–1 |
fitness |
Environmental fitness of the individual | Dimensionless, range 0–1 |
diffphen/
Contains data used for analyses of phenotypic-genome differentiation among species.
data5000.dat: Population snapshot at generation 5000. It contains the same seven variables described above forAll/data[t].dat.phegenome.dat: Binary phenotypic genomes corresponding to the individuals indata5000.dat. Each row represents one individual and contains 100 binary loci (0or1).
phenotype/
Contains data used for the analysis of geographic versus genetic distance.
data.dat: Population data containing the same seven variables described above:id,x,y,species,degree,phenotype, andfitness.genome.dat: Binary reproductive genomes corresponding to the individuals indata.dat. Each row represents one individual and contains 1500 binary loci (0or1).
binomial.dat
Whitespace-delimited two-column file containing the theoretical binomial reference curve used in the phenotype-distribution analysis.
The columns are:
| Column | Description | Units / interpretation |
|---|---|---|
| 1 | Phenotypic value P = k/100, where k is the number of phenotypic loci equal to 1 |
Dimensionless, range 0–1 |
| 2 | Scaled binomial reference value, calculated as 130 × Binomial(k; n = 100, p = 0.5) |
Dimensionless |
The binomial distribution represents the expected phenotype distribution for 100 independent binary phenotypic loci with equal probabilities of 0 and 1.
The factor 130 scales the theoretical curve for comparison with the simulation histogram.
diff_boundary_[t].txt
One-column text files containing phenotypic differences between mating pairs located near the boundary between the two environmental regions. The t indicates the simulation generation at which the data were saved.
Each row corresponds to one reproductive event for which both parents are located within a distance S of the environmental boundary. The recorded value is the signed phenotypic difference between the two parents, calculated as:
Phenotype_father - Phenotype_mother
Phenotypic differences are dimensionless and range from -1 to 1.
3. Software Requirements
The analysis notebooks require Python 3 and Jupyter Notebook.
The main Python packages used in the analysis notebooks are:
- NumPy
- pandas
- Matplotlib
- SciPy
- seaborn
Some figure-generation code uses LaTeX for text rendering.
4. Code Availability
The source code is also available at:
