Although the pervasiveness of intraspecific wing-size polymorphism and transitions to flightlessness have long captivated biologists, the demographic outcomes of shifts in dispersal ability are not yet well understood and have been seldom studied at early stages of diversification. Here, we use genomic data to infer the consequences of dispersal-related trait variation in the taxonomically controversial short-winged (Chorthippus corsicus corsicus) and long-winged (Chorthippus corsicus pascuorum) Corsican grasshoppers. Our analyses revealed lack of contemporary hybridization between sympatric long- and short-winged forms and phylogenomic reconstructions supported their taxonomic distinctiveness, rejecting the hypothesis of intraspecific wing polymorphism. Statistical evaluation of alternative models of speciation strongly supported a scenario of Pleistocene divergence (<1.5 Ma) with ancestral gene flow. According to neutral expectations from differences in dispersal capacity, historical effective migration rates from the long- to the short-winged taxon were three-fold higher than in the opposite direction. Although populations of the two taxa present a marked genetic structure and have experienced parallel demographic histories, our coalescent-based analyses suggest that reduced dispersal has fueled diversification in the short-winged C. c. corsicus. Collectively, our study illustrates how dispersal reduction can speed up geographical diversification and increase the opportunity for allopatric speciation in topographically complex landscapes.

Genomic library preparation

We used NucleoSpin Tissue (Macherey-Nagel, Düren, Germany) kits to extract and purify DNA from a hind leg of each individual. We processed genomic DNA into one genomic library using the double-digestion restriction-site associated DNA sequencing procedure (ddRAD-seq) described in Peterson et al. (2012). In brief, we digested DNA with the restriction enzymes MseI and EcoRI (New England Biolabs, Ipswich, MA, USA) and ligated Illumina adaptors including unique 7-bp barcodes to the digested fragments of each individual. We pooled ligation products and size-selected them between 475-580 bp with a Pippin Prep instrument (Sage Science, Beverly, MA, USA). We amplified the fragments by PCR with 12 cycles using the iProofTM High-Fidelity DNA Polymerase (BIO-RAD, Veenendaal, Netherlands) and sequenced the library in a single-read 150-bp lane on an Illumina HiSeq2500 platform at The Centre for Applied Genomics (Toronto, ON, Canada).

Genomic data assembling and filtering

Raw sequences were demultiplexed and preprocessed using stacks v. 1.35 (Catchen et al., 2013) and assembled using pyrad v. 3.0.66 (Eaton, 2014). Libraries were demultiplexed and filtered for overall quality using process_radtags (Catchen et al., 2011, 2013), retaining reads with a Phred score > 10 (using a sliding window of 15%), no adaptor contamination, and that had an unambiguous barcode and restriction cut site. Raw sequence data quality was checked in fastqc v. 0.11.5 (http://www.bioinformatics.babraham.ac.uk/projects/fastqc/) and sequences were trimmed to 129 bp using seqtk (Heng Li, https://github.com/lh3/seqtk) in order to remove barcodes and low-quality reads near the 3´ ends. We assembled our sequences into de novo loci using pyrad v. 3.0.66 (Eaton, 2014). Briefly, reads retained after process_radtags were further quality-filtered with pyrad to convert base calls with a Phred score <20 into Ns and discard reads with >2 Ns. Retained reads were clustered within- and across samples considering a threshold of sequence similarity (W_clust) of 85% and clusters with a coverage depth <5 were discarded. Loci containing one or more heterozygous sites across >15% of individuals were excluded, as we expect that this represents a fixed difference among clustered paralogs rather than a true polymorphism (Eaton, 2014). Unless otherwise indicated, all downstream analyses were performed using datasets of unlinked SNPs (i.e., a single SNP per RAD locus) obtained with pyrad considering a clustering threshold of sequence similarity of 0.85 (W_clust = 0.85) and discarding loci that were not present in at least 50 % individuals (minCov = 50 %).

References

Catchen, J. M., A. Amores, P. Hohenlohe, W. Cresko, and J. H. Postlethwait. 2011. stacks: Building and genotyping loci de novo from short-read sequences. G3-Genes Genom. Genet. 1:171-182.

Catchen, J., P. A. Hohenlohe, S. Bassham, A. Amores, and W. A. Cresko. 2013. stacks: an analysis tool set for population genomics. Mol. Ecol. 22:3124-3140.

Eaton, D. A. R. 2014. pyrad: assembly of de novo RADseq loci for phylogenetic analyses. Bioinformatics 30:1844-1849.

Peterson, B. K., J. N. Weber, E. H. Kay, H. S. Fisher, and H. E. Hoekstra. 2012. Double digest RADseq: An inexpensive method for de novo SNP discovery and genotyping in model and non-model species. PLoS One 7:e37135.

Samples.xlsx: Description of individual and population codes used in the different analyses and genomic datasets

c85d5m9sh15.nex: Input file (NEXUS format) used for phylogenomic analyses in SNAPP

c85d5m24sh15.unlinked_snps: Input file (NEXUS format) used for phylogenomic analyses in SVDQuartets

c85d5m24sh15.treemix: Input file (TREEMIX format) used for phylogenomic analyses in TREEMIX

STRUCTURE_PCA.zip: This ZIP folder contains the genetic datasets (STR format) used to perform clustering (STRUCTURE) and principal component analyses (PCA)

BPP.zip: This ZIP folder contains the genetic dataset, control file and individual map file used to perform BPP analyses

c85d5m24sh15_Outgroup.loci: Input file used to perform D-statistic (ABBA-BABA) tests in PYRAD

Fastsimcoal2.zip: Input files for the alternative models tested using FASTSIMCOAL2

StairwayPlot.zip: Input files for demographic reconstructions in STAIRWAYPLOT

PYRAD.zip: Output files from PYRAD obtained for populations of the two taxa (c85d5m24sh15) and only considering populations of Chorthippus corsicus corsicus (C_corsicus_corsicus_c85d5m12sh15) and Chorthippus corsicus pascuorum (C_corsicus_pascuorum_c85d5m12sh15)