Primary hypotheses:
H1: Epistasis + pleiotropy accelerate coordinated trait emergence relative to additive architecture.
H2: Prey coevolution increases oscillatory dynamics and increases time to persistent coordination, but may raise amplitude of punctuations.
H3: Bottlenecks increase the probability of rapid coordinated fixation (conditional on survival).
Tests & methods:
Compare distributions (e.g., time to coordination) with nonparametric tests (Mann--Whitney U, Kolmogorov--Smirnov) across variants.
Regression / generalized additive models (GAMs) to associate parameter values with outputs (time to emergence, peak trait values).
Survival analysis (Kaplan--Meier / Cox proportional hazards) for time-to-event (emergence) analyses.
Multivariate analysis: PCA on trait space; clustering validation indices; Mantel tests for genotype--phenotype associations.
Model selection: AIC/BIC for simplified ODE approximations fit to ABM averages.
Bootstrap for confidence intervals on ensemble summaries.
Empirical calibration / comparison:
If genomic datasets exist (candidate loci), compute LD and co-selection signatures (haplotype blocks) and compare with simulated LD patterns.
Use Approximate Bayesian Computation (ABC) to estimate parameter posteriors given empirical summaries.
8. Implementation details & reproducibility
Primary platforms:
SLiM 3.x --- recommended for forward-time, genotype-explicit population genetic simulations with selection and epistasis (scriptable Eidos language).
Python (3.8+) --- driver scripts, data processing, plotting. Key packages: NumPy, SciPy, pandas, h5py, scikit-learn, statsmodels, seaborn/matplotlib.
Mesa or NetLogo --- optional ABM frameworks for spatial models if SLiM spatial scripting is insufficient.
R --- statistical analysis (survival, GAMs) and plotting (ggplot2).
Parallelization & HPC: run replicates in parallel (SLURM jobs / GNU parallel) --- store outputs per replicate to avoid IO contention.
Version control & containers: keep SLiM scripts, Python drivers, and analysis code in git. Provide a Docker/Singularity container with exact software versions to ensure reproducibility.
Storage & metadata: archive raw outputs (HDF5) and processed summaries; retain parameter + seed manifest (CSV/JSON) per replicate.
9. Performance & practical notes
Start with reduced settings (L=20, N_P=200) to debug and calibrate.
For final analyses, scale up (L=100, N_P1000) on HPC.
Use checkpointing for long runs.
Log wall-time & memory per run to plan resources.
10. Deliverables from simulation suite
Ensemble of simulation runs with full metadata.
Figures: representative allele-trajectory heatmaps; trait mean/variance trajectories; phase plots predator mean trait vs prey mean trait; distributions of time-to-emergence across variants; LD heatmaps showing co-selection.
Statistical tables summarizing hypothesis tests, sensitivity analyses, and parameter estimates.
Reproducible code package (SLiM scripts + Python analysis notebooks) and container image.
C.2. Â Simulation Design and Parameterization
To evaluate whether coordinated adaptive packages can emerge under a CAS framework, we implemented both replicator--mutator equations and forward-time, agent-based simulations. The simulations were designed to be fully reproducible, with explicit specification of genetic architecture, ecological coupling, and parameterization.
Model Variants
We explored six variants of the model to test distinct hypotheses. Variant A implemented a baseline additive architecture without prey coevolution. Variant B introduced epistasis and pleiotropy while maintaining fixed prey traits. Variant C incorporated full eco-evolutionary coupling, with prey evolving under replicator--mutator dynamics. Variant D introduced demographic bottlenecks to test their effect on the speed of adaptation. Variant E applied an NK-model genotype--fitness mapping with varying ruggedness (K{0,2,4,6}K \in \{0,2,4,6\}K{0,2,4,6}). Variant F added spatial structure by simulating populations on a two-dimensional grid with local dispersal. Each variant was run in ensembles of 30--200 replicates to ensure statistical robustness.
Parameters
The predator genome was modeled with L=20L=20L=20 loci in exploratory runs and L=100L=100L=100 in full analyses, with each locus influencing one or more of four focal traits: vision (T1T_1T1), respiration (T2T_2T2), neuromuscular control (T3T_3T3), and wing morphology (T4T_4T4). Alleles were binary or quantitative, with a per-locus mutation rate =105--103\mu = 10^{-5} \text{--} 10^{-3}=105--103 and recombination rate r=0.01--0.1r=0.01\text{--}0.1r=0.01--0.1. Predator population sizes ranged from 200--5000 individuals (NPN_PNP), with prey populations initialized at Nprey=5000N_{prey}=5000Nprey=5000. Selection strength was scaled by parameter s=0.01--1.0s=0.01\text{--}1.0s=0.01--1.0. Cost coefficients for trait expression (i\lambda_ii) were chosen to ensure trade-offs produced intermediate optima. Epistasis and pleiotropy were encoded by sparse Gaussian matrices, with 10% of locus pairs nonzero for epistasis and 30% of locus--trait links nonzero for pleiotropy. Prey population growth followed logistic dynamics with rprey=1.0r_{prey}=1.0rprey=1.0 and carrying capacity equal to the initial population size.
Simulation Protocol
