Skip to main content
Download PDF

Population structure of Apodemus flavicollis and comparison to Apodemus sylvaticus in northern Poland based on RAD-seq

BMC Genomics volume 21, Article number: 241 (2020) Cite this article

Abstract

Background

Mice of the genus Apodemus are one the most common mammals in the Palaearctic region. Despite their broad range and long history of ecological observations, there are no whole-genome data available for Apodemus, hindering our ability to further exploit the genus in evolutionary and ecological genomics context.

Results

Here we present results from the double-digest restriction site-associated DNA sequencing (ddRAD-seq) on 72 individuals of A. flavicollis and 10 A. sylvaticus from four populations, sampled across 500 km distance in northern Poland. Our data present clear genetic divergence of the two species, with average p-distance, based on 21377 common loci, of 1.51% and a mutation rate of 0.0011 - 0.0019 substitutions per site per million years. We provide a catalogue of 117 highly divergent loci that enable genetic differentiation of the two species in Poland and to a large degree of 20 unrelated samples from several European countries and Tunisia. We also show evidence of admixture between the three A. flavicollis populations but demonstrate that they have negligible average population structure, with largest pairwise FST<0.086.

Conclusion

Our study demonstrates the feasibility of ddRAD-seq in Apodemus and provides the first insights into the population genomics of the species.

Background

Mice of the genus Apodemus (Kaup, 1829) (Rodentia: Muridae) are one the most common mammals in the Palaearctic region [39]. The genus comprises of three subgenera (Sylvaemus, Apodemus and Karstomys) [39], however the systematic classification of the 20 species belonging to the genus [17] is not fully settled [33]. In the Western Palearctic, the yellow-necked mice A. flavicollis (Melchior, 1934) and the woodmice A. sylvaticus (Linnaeus, 1758) are widespread, sympatric and occasionally syntopic species. They are often difficult to distinguish morphologically in their southern range [28], but in the Central and Northern Europe both are easily recognisable by the full yellow collar around the neck of A. flavicollis, which only forms a narrow elongated spot on the breast in A. sylvaticus [52].

Their prevalence in Western Palearctic and common status in Western and Central Europe made them one of the model organisms to study post-glacial movement of mammals [22, 41]. Both species have traditionally been studied in a parasitological context, as one of the vectors of Borellia-carrying ticks Ixodes ricinus, who often feed on Apodemus [43, 58], tick-borne encephalitis virus [14] and hantaviruses [31, 46] and have been used as markers for environmental quality [36, 63]. Lastly, they have extra-autosomal chromosomes, called B chromosomes, with varied distribution among the populations [56] and suggested involvement in a variety of physiological phenomena, from cell division and development to immune response [64].

Previous studies on Apodemus typically employed a small number of microsatellite [59] and mtDNA markers [22, 38, 40, 41], which are insufficient to learn about the species’ population structure and admixture patterns in detail, or to identify loci under selection. In the absence of high-quality reference genome, which remains cost-prohibitive for complex genomes, whole-genome marker discovery enabled by restriction site-associated DNA sequencing presents a cost-effective method to study species on a population scale even with no previous genetic and genomic resources available [5].

Here we employ the double-digest restriction site-associated DNA sequencing (ddRAD-seq) to elucidate the genetic structure and connectivity of three populations of A. flavicollis and compare it to a population of A. sylvaticus in Poland. We demonstrate clear divergence between the two species and very low differentiation between populations of A. flavicollis. Our results provide the first estimates of population parameters in A. flavicollis based on thousands of loci, calculation of p-distance between the two Apodemus species, as well as a selection of loci enabling their accurate identification.

Results

Sequencing and variant calling

The sequencing produced a total of 92741120 reads. The number of reads per individual varied from 346810 to 4157586, with an average of 1078385 reads per individual and median of 905786,5 (Supplementary Table S2). The best parameters for calling the stacks and variants for the entire dataset were: minimum number of identical, raw reads required to create a stack m = 2, number of mismatches allowed between loci for each individual M = 4 and number of mismatches allowed between loci when building the catalogue n = 5 (Supplementary Figure S1). The best parameters calculated for A. flavicollis samples only were: m = 2, M = 4 and n = 3 (Supplementary Figure S3). The coverage per sample ranged from 4.95x to 26.20x with an average of 10.13x and median of 9.32x for the entire dataset (Supplementary Figures S2 and S4).

SNPs and loci co-identification rates

Analysis of the duplicated samples showed that loci and allele misassignment rates were of similar magnitude, on average, between all pairs of duplicates. The duplicate pair F06-B02 showed the highest discrepancy between loci, of 10%, and also between alleles, of 8%. When only shared loci were included in the comparisons, all four sets of duplicates showed on average 0.5% ±0.2% SNPs called differently (Table 1).

Table 1 Error rates calculated by comparing four sets of duplicated samples. D1/D2: ratio of reads from Duplicate 1 to Duplicate 2. Locus misassignment rate: the percentage of unidentified loci, calculated by dividing the number of loci found only in one of the duplicates by the total number of loci in each sample. Allele misassignment rate: the percentage of mismmatches between the IUPAC consensus sequences between homologous loci from each pair of duplicates. SNP error rate 1: the percentage of different SNPs called in each of the duplicated samples using either 10178 SNPs. Shared SNP error rate: the percentage of different SNPs called in each of the duplicated samples after excluding missing data between duplicate samples
Full size table

Comparison of A. flavicollis and A. sylvaticus

The number of assembled loci per individual ranged from 46286 to 117366 (mean: 73711, median: 71395, standard deviation: 29917). 52494 loci passed the population filters established for species differentiation (see Methods, section "Variant calling and filtering"), representing 8,3% of the total 632063 loci included in the catalogue. Out of 158144 SNPs called, 60366 (38.1%) were removed after filtering for minor allele frequency (MAF) and 52298 (33%) were removed after failing the HWE test at p<0.05; further 35302 (22.3%) were removed due to a minimum mean depth lower than 20, leaving 10178 SNPs (6.6%) to be used in the downstream analyses (Fig. 1). PCA plot of the first two components (Fig. 2), accounting for 13.13% of the total variance, shows differentiation of the two species but also distinguish different populations of A. flavicollis.

Fig. 1
figure1

Summary of cataloque construction and SNP filtering steps for the complete dataset (left) and Apodemus flavicollis dataset. The graphic includes: Stacks parameters values (m, M, n), number of loci in the catalogue, number of SNPs filtered by minor allele frequency (MAF), which failed the Hardy-Weinberg equilibrium test at p<0.05 (HWE), SNPs removed due to an average depth, across individuals, lower than 20 (min-meanDP) and the total number of SNPs retained for further analysis

Full size image
Fig. 2
figure2

Principal Component Analysis of all samples analysed in the study. Each point represents one sample; the shape of the point represents the species (circles: Apodemus flavicollis (n = 72), triangles: Apodemus sylvaticus (n = 10), whereas the colour represents the location where the samples were collected: Bial - Białowieża, Kadz - Kadzidło, Hack - Haćki, Bory - Bory Tucholskie

Full size image

Similarly, the phylogenetic tree shows A. sylvaticus as a separate clade to the three populations of A. flavicollis, with A. flavicollis from geographically closer regions (Białowieża and Haćki, 50 km) grouped closer than a population from Bory Tucholskie, 450 km away from Białowieża (Fig. 3). The A. sylvaticus and A. flavicollis clusters have high bootstrap value support (100% and 99% respectively).

Fig. 3
figure3

Maximum likelihood phylogenetic tree of all the samples analysed in the study. Colour represents the species: A. sylvaticus (n=10) in orange and A. flavicollis (n=72) in black. Duplicates samples are included: F06-B02 from Bory Tucholskie, F12-A12 and H11-G06 from Białowieża and G02-D01 from Haćki. Bootstrap support values from 100 replicates are indicated at the nodes of the tree. Bial - Białowieża, Kadz - Kadzidło, Hack - Haćki, Bory - Bory Tucholskie

Full size image

We then investigated the suitability of the loci we identified on Polish populations to distinguish A. sylvaticus and A. flavicollis from other European populations. The genotyping of the extra 10 samples from each species (see Methods) produced 179763 SNPs. 62158 (34.58%) were removed after filtering for MAF and 69125 (38.45%) were removed after failing the HWE test at p<0.05; further 42054 (23.39%) were removed due to a minimum mean depth lower than 20 and 5203 (2.89%) were removed due to more than 5% missing data, leaving 1223 SNPs (0.68%) to be used in the downstream analyses.

The first axis of the PCA plot (Fig. 4) constructed from this data accounts for the 65.73% of the total variance and shows clear differentiation between the two species. All the A. flavicollis samples cluster with the Polish A. flavicollis samples, while all but Tunisian samples of A. sylvaticus cluster with the Polish samples of the same species. Tunisian A. sylvaticus appear as a separate cluster but still closer to the A. sylvaticus group. The catalogue of loci used for species identification is included in the Supplementary Materials, Section 6.

Fig. 4
figure4

Species identification through Principal Component Analysis using a catalogue of 632060 loci and 1223 final SNPs. Light colours represent samples from Poland while dark colours represent samples from other European regions and Tunisia (collectively named “Europe”; Tunisian samples are marked with a circle). Green: A. sylvaticus, blue: A. flavicollis

Full size image

Genetic diversity and population structure of A. flavicollis

The number of assembled loci per individual in the Polish populations ranged from 46286 to 117366 (mean: 72738, median: 70592, stdev: 12575). 30722 loci passed the population filters established for population differentiation, representing and 4,43% of the total 691960 loci included in the catalog. Out of 63742 SNPs called, 31401 (49.26%) were removed after filtering for MAF and 10034 (15.74%) were removed after failing the HWE test at p<0.05. Further 9653 (15.14%) were removed due to a minimum mean depth lower than 20, leaving 12654 (19.85%) SNPs to be used in the downstream analyses (Fig. 1).

PCA plot (Fig. 5) shows differentiation between the three Polish A. flavicollis populations, with PC1 and PC2 cumulatively explaining 10.47% of the total variance. Haćki population shows larger diversity than the other populations, with some Haćki individuals closer to Białowieża individuals than to others from this location. Phylogenetic tree (Fig. 6) supports this pattern of differentiation. Bory Tucholskie and Haćki populations each form a cluster with a 100% of bootstrap support value, whereas Białowieża forms a third cluster with an 95% of bootstrap support. Białowieża and Bory Tucholskie population together form a large cluster with a 100% bootstrap support.

Fig. 5
figure5

PCA plot showing Polish samples of A. flavicollis from Białowieża (red) (n=35), Haćki (blue) (n=14) and Bory Tucholskie (green) (n=23). Bial - Białowieża, Kadz - Kadzidło, Hack - Haćki

Full size image
Fig. 6
figure6

Maximum ilkelihood phylogenetic tree of n = 72 A. flavicollis samples from Bialowieżdot;a (red, n = 35), Haćki (blue, n = 14) and Bory Tucholskie (green, n = 23). Bootstrap support values from 100 replicates are indicated at the nodes of the tree

Full size image

In the ADMIXTURE analysis, the lowest cross-validation errors [2] were always found for K = 3, indicating contribution of three ancestral populations (Fig. 7). Majority of samples from each of the populations show a single dominant component of ancestry with little contribution from other populations, with the exception of four individuals from Haćki, which show clear admixture of the Białowieża population.

Fig. 7
figure7

Maximum likelihood Admixture analysis of all A. flavicollis samples for the optimal K = 3. Each bar represents an individual and each colour represents its ancestry component (red: Białowieża, blue: Haćki, green: Bory Tucholskie)

Full size image

Recognising that STRUCTURE-type analyses (on which ADMIXTURE is based) may be sensitive to the effects of uneven number of samples in compared groups [54], we repeated the ADMIXTURE analysis 10 times, each time randomly drawing the same number of individuals (n = 15) from each population. In all cases, the lowest cross-validation errors were found for K = 2, followed by K = 3 (Supplementary Figure S5). At even sampling, ADMIXTURE pattern found for K = 3 was the closest to the observed ecological and geographical distribution of the samples and closely matched our results when all samples were included (Supplementary Figure S6).

The patterns of heterozygosity highlight Haćki as the only population where the values of Ho is higher than He, where the FIS is negative (Table 2). As parameters such as number of private alleles, nucleotide diversity and heterozygosity can vary with sample size, we performed 100 calculations of the above parameters using random sampling of the same number of individuals (n = 15) from each population. The parameters showed similar relationships except for the number of private alleles (data not shown).

Table 2 Population genetic parameters calculated based on 12654 SNPs from all 72 individuals of A. flavicollis. N, number of individuals; Npa, number of private alleles; Ind per loci, Mean number of individuals per locus in this population; Ho observed and He expected heterozygosity;, average nucleotide diversity; FIS inbreeding coefficient
Full size table

Fst values are consistently very low between all the populations, even though populations from Haćki and Bory Tucholskie show three-fold higher Fst values that for the other two pairs of populations (Table 3).

Table 3 Pairwise FST values for the three populations of A. flavicollis
Full size table

Species divergence

Finally, we calculated that the average p-distance between A. flavicollis and A. sylvaticus, based on 21377 shared loci, is 1.51% (standard deviation = 1.11%).

We then identified the top 117 most divergent loci between the species, which all had the divergence larger than 4.9% (The loci ID are provided in the Supplementary Table S3), and checked whether these loci alone allow for accurate assignment of samples to the two species. We constructed PCA plots from the Polish samples only and from the Polish, other European and Tunisian samples together. They demonstrate that while the 117 loci are sufficient to clearly assign Polish samples to the two species (Supplementary Figure S8), some uncertainty remains when we use these loci for the broader set of samples. Whereas all A. flavicollis samples do cluster together, A. sylvaticus samples do not form a clearly differentiated group (Supplementary Figure S9).

We also identified fixed loci, where all individuals within each species have identical sequences. There were 3526 such fixed loci for A. flavicollis and 5843 for A. sylvaticus. We then used 1273 of those loci that were shared among the two species and calculated that the average p-distance based on fixed differences is 0.97% (standard deviation = 0.94%).

Discussion

RAD-sequencing approaches, including double-digest RAD-seq and its variants [6, 19, 42, 49, 50], have allowed a cost-effective discovery of thousands of genetic markers in both model and non-model organisms [21, 60], proving to be a transformative research tool in population genetics [8, 13, 24], phylogeography and phylogenetics [4, 23, 27, 57], marker development [48], linkage mapping studies [7], species differentiation [45] and detecting selection [62]. However, despite the widespread use of this approach to marker discovery, only few studies have used RAD-seq in mammals [18, 30, 32, 44, 61]. Here, we have identified over 10000 markers in two closely related and common species of Apodemus in Western Palearctic, characterised the population structure of A. flavicollis and compared it to A. sylvaticus, for the first time providing estimates of the species divergence and population genetic parameters based on thousands of SNPs.

Technical considerations

We have used four pairs of technical duplicates to check the accuracy of the RAD-seq genotyping based on the Poland protocol [51]. The largest source of discrepancy in SNP calls between the duplicates is caused by unequal identification of loci: the difference in our case averaged approximately 10% (Table 1) and was similar to allele misindentification rates. However, when considering only shared loci between the duplicates, the discrepancy in SNP calls fell by over an order of magnitude to an average of 0.5%, indicating high accuracy and reliability of calls in once-defined shared loci. Our finding of loci calls being the major source of genotyping variability agrees with Mastretta et al. (2015), although our discrepancies are almost an order of magnitude smaller. Moreover, despite the differences in number of loci included in the analysis, each duplicated pair of samples clustered together with a 100% bootstrap values support and branch length equal to 0 on the phylogenetic tree (Fig. 6), indicating that the samples were identical. Overall, our finding reiterates the importance of the influence of stochastic events and imprecise size selection in the library preparations on genotyping calls [37]. We note that some of these variables could be better controlled with automated size-selection approaches [49]. Our findings also illustrate the usefulness of including technical replicates during library preparation.

Effect of group size

Permutations performed for the calculations of genetic diversity parameters (Table 4) have shown that with the exception of the number of private alleles, the results are comparable, regardless of the number of samples included per each population.

Table 4 Average genetic diversity parameters for Apodemus flavicollis calculated from 100 permutations of 45 individuals (15 samples per population, 12654 SNPs). N, number of individuals; Npa, number of private alleles; Ind per loci, Mean number of individuals per locus in this population; Ho observed and He expected heterozygosity;, average nucleotide diversity; FIS inbreeding coefficient
Full size table

In the ADMIXTURE, we observed different optimal K depending on whether all samples were included in the analysis (K = 3) or a set of 15 randomly-chosen set of samples from each population (K = 2, although closely followed by K = 3). While the previously reported tendency of STRUCTURE-like analyses to produce ΔK = 2 does not apply in our case due to different method to select optimal number of clusters [26], we chose to use K = 3 for our analyses due to close match to the spatial and ecological locations from which our populations were sampled. The results obtained for K = 3 in the evenly-sampled dataset were similar to the clusters obtained for K = 3 with the complete dataset. ADMIXTURE plots for K = 2 to K = 5 are shown as Supplementary Figure S7.

Population structure

The FST values calculated in this study between all three pairs of populations of A. flavicollis, based on 12654 SNPs, are consistently low and are not affected when we randomly draw the same number of individuals from each population to compute pairwise FST (Table 5). Previous studies of A. flavicollis populations in north-eastern Poland based on a small number of microsatellites showed similarly and consistently low values [15, 20], even though [20] suggested some population structure based on statistically significant differences between very low pairwise FST values. [15] also suggest large, broadly geographically defined clusters of A. flavicollis in north-eastern Poland that are separated by highly admixed individuals, but, again, FST between those clusters are as low as those reported by [20] and this study.

Table 5 Average pairwise FST values ± standard deviation for the three populations of A. flavicollis calculated from 100 permutations of 45 individuals (15 samples per population, 12654 SNPs)
Full size table

We would argue, based on a much larger set of markers reported here, that A. flavicollis has a negligible population structure across the entire area studied. Large number of markers nevertheless allows us to discover evidence for admixture of Białowieża population and Haćki (Fig. 7), further indicated by relatively high heterozygosity and negative FIS in this population. It is therefore intriguing that such a low differentiation occurs across hundreds of kilometres of varying landscape in a species that typically has a limited range of about 4 km2 and that suffers up to 86% winter mortality rate [53], which would lead to multiple bottlenecks and drift-driven population differentiation. With this in mind, our data suggests a much larger dispersal ability of the species, a much better connectivity between populations, or both.

Both low overall FST and moderate heterozygosity suggest it would be worthwhile to conduct a genome-wide scan for selection using FST as a metrics of local genomic differentiation to identify geographically local regions under selection. This, however, is not yet possible given the lack of high-quality reference genome for Apodemus and unknown synteny to the available genome of Mus musculus.

Divergence and differentiation of A. flavicollis and A. sylvaticus

Given that accurate identification of the two species using morphological characters is problematic, especially in their southern range [10], a large collection of markers identified in this study allowed us to create a catalogue of 632060 loci that allow clear differentiation between species. This identification is somewhat biased, as the catalogue was built using many more samples of A. flavicollis than A. sylvaticus (72 vs 10) and both from a relatively limited geographical range. Nevertheless, it allowed for accurate assignment of A. flavicollis samples and to a large degree of A. sylvaticus, as we demonstrated on a set of 20 independent samples from other European countries and Tunisia (Figure 4). The imperfect clustering of some European A. sylvaticus samples on Figure 4 may reflect the higher genetic differentiation of population from the edge of the range of the species: samples near the centre of the PC1 axis on Figure 4 (green dots in a circle) all come from Tunisia. Given the wide distribution of both species in Western Palearctic, a more representative sample from both species from a broader geographic range would likely provide more accurate set of markers for their identification.

Finally, we calculated the nucleotide divergence between the two species, based on 21377 shared loci, which is 1.51%. Considering a divergence time between A. flavicollis and A. sylvaticus estimated from archeological data of 4 Mya [38], the evolutionary rate is 0.0019 substitutions per site per million of years. While this estimate of sequence divergence level is in broad agreement with calculations based on mitochondrial 12S rRNA, IRBP and Cytochrome b genes [39], we note that Aghova et al. (2018) [1] provided an updated timing of a split between A. sylvaemus and A. mystacinus, which, at 9.6Mya, is 2Mya older than previously suggested by Michaux et al. (2002) [39]. If we use that date as a reference and move the presumed split between A. flavicollis and A. sylvaticus by 2.6Mya, the estimated evolutionary rate from our data would be 0.0011. In both cases, our calculation is likely an underestimate, as we only used shared loci to calculate divergence and did not include the potential impact of insertion/deletion events, which can significantly affect the total genomic divergence between species [9, 35]. Highly divergent sequences would have been identified as different loci, and would not be compared to their true homologous sequences.

Conclusions

We have successfully applied the ddRAD-seq approach to discover tens of thousands of SNPs in wide-spread and common mammalian species of A. flavicollis and A. sylvaticus. The high resolution data obtained here allowed us to delineate geographically close populations, including identifying admixture between them, but suggest that A. flavicollis effectively forms a single population in an entire sampling area that spans 500 km in the W-E direction. Comparing A. flavicollis and A. sylvaticus, we have calculated their divergence and identified a set of genomic loci that enable effective molecular identification of the species. We anticipate that with the development of further whole-genome resources, Apodemus, thanks to its common status, broad geographic range and long history of ecological observations, will become an excellent model species for evolutionary and ecological research in the genomic era.

Methods

Sample collection and dNA extraction

Eighty two individuals (10 Apodemus sylvaticus and 72 Apodemus flavicollis) from four locations in northern Poland spanning 500 km were trapped in 2015 (Fig. 8). A. flavicollis were collected in Białowieża (E23.8345814, N52.7231935), an oak-lime-hornbeam forest (n = 35), Bory Tucholskie (E17.5160265, N53.7797608), in an oak-lime-hornbeam and pine forest (n = 23) and Haćki (E23.1793284, N52.834369), in a xerothermic meadow (n = 14). A. sylvaticus were trapped in Białowieża (E21.3778496, N53.2089113) in a dry pine forest (n = 5) and in Bory Tucholskie, mainly in a pine forest (n = 5) (Supplementary Table S1). While A. flavicollis are present in all sampled locations, there have been no trappings of A. sylvaticus in Białowieża for the last 20 years, despite Białowieża being within the European range of this species (Dr Karol Zub, personal communication). The sampling procedures were approved by the Local Ethical Commission on Experimentation on Animals in Białystok, Poland, under permission number 2015/99. All animals were released after sample collection.

Fig. 8
figure8

Locations of the Polish samples used in this study. Red circles represent samples from Apodemus flavicollis while blue dots represent samples from Apodemus sylvaticus. The number inside the circles are the number of samples from each locality. Bial - Białowieża, Kadz - Kadzidło, Hack - Haćki, Bory - Bory Tucholskie

Full size image

Tail clippings were collected, preserved in ≥95% ethanol and stored at −20∘C until DNA extraction. The tissues were digested by incubating at 55∘C overnight with lysis buffer (10mM Tris, 100mM NaCl, 10mM EDTA, 0.5% SDS) and proteinase K (20mg/ml). Subsequently, potassium acetate and RNAse A were used to remove protein and RNA contamination. Three ethanol washes were performed using Sera-Mag SpeedBeads solution (GElifesciences, Marlborough, MA, USA). The quality and integrity of the DNA was tested in a 2% agarose gel. Twenty-fold dilutions of the samples were used to measure the DNA concentration using Quant-iT PicoGreen dsDNA assay kit (Invitrogen, Carlsbad, CA, USA) and concentration of each sample was then normalised to 10 ng/ μl in 20 μl volume. Four samples were used as technical duplicates (F06-B02, G02-D01, H11-G06, F12-A12). Technical duplicates had the same DNA but were digested and ligated to barcodes independently.

ddRAD-seq library preparation

ddRAD-seq library was prepared following the protocol from [50], adapted to a different combination of enzymes. Briefly, genomic DNA was digested in a 20 μl reaction with CutSmart® buffer, 8 units of SbfI and 8 units of HF-MseI (New England Biolabs, Frankfurt am Main, Germany). Digestion was performed at 37∘C for 2 hours. Enzymes were inactivated at 65∘C for 20 minutes and the reactions were kept at 8∘C. Adapter ligation was performed at 22∘C for 2 hours and the ligase was inactivated by incubating the samples at 65∘C for 20 minutes. Samples were cooled down to 8∘C and multiplexed by combining 5 μl of each sample. P1 adapters contained barcodes with a length between 5 and 10 bp.

PCR amplification was conducted in 25 μl with 1 μl of each primer (IlluminaF_PE: AATGATACGGCGACCACCGAGATCTACACTCTTTCCCTACACGACGCTCTTCCGATCT and IlluminaR_PE: CAAGCAGAAGACGGCATACGAGATCGGTCTCGGCATTCCTGCTGAA) at 10mM, 0.5 μl of 10 mM dNTPs, 13.25 μl of PCR-grade water, 5 μl of 5x Phusion HF Buffer, 0.25 μl of Phusion DNA Polymerase (New England Biolabs, Frankfurt am Main, Germany) and 4 μl of the multiplexed DNA. After an initial denaturation step of 30s at 98∘C, PCR reaction was carried out for 12 cycles (10s at 98∘C, 20s at 58∘C and 15s at 72∘C). Final elongation step was performed at 72∘C for 5 minutes.

PCR products were loaded into a single lane on a 1% agarose gel with 100 bp DNA ladder (New England Biolabs, Frankfurt am Main, Germany). Fragments between 200 and 500 bp were cut from the gel with a scalpel and purified using the QIAquick gel extraction kit (QIAgen, Hilden, Germany), followed by the second cleanup step with Sera-Mag SpeedBeads (GElifesciences, Marlborough, MA. USA). Sizing, quantification and quality control of the DNA was performed using Bioanalyzer 2100 (Agilent, Santa Clara, CA, USA) before paired-end sequencing on an Illumina HiSeq 3500 with 150 bp read length.

Processing of RAD-tags

Sequences were analysed with Stacks version 1.48 [11]. Samples were demultiplexed using process_radtags allowing no mismatches in barcodes and cutting sites. Sequences with uncalled bases and low quality scores were removed and all reads were trimmed to 141 bp. The four files generated per sample by process_radtags were concatenated using a custom bash script. The best parameters for building and calling SNPs de novo, using denovo_map, were calculated following [47] approach, using either samples from both species or only from A. flavicollis. Secondary reads were not used to call haplotypes in denovo_map (option -H).

SNPs and loci co-identification rates

We estimated the loci and SNP co-identification rates by analysing a set of four samples that were prepared and sequenced in duplicates. Sequences for 52494 loci from both species, were extracted using –fasta_samples option from the population package in Stacks. We extracted sequences for each of the duplicated samples with a custom script and calculated co-identification rates as described by [37]. Briefly, the locus misassignment rate is the percentage of unidentified loci, calculated by dividing the number of loci found only in one of the duplicates by the total number of loci in each sample. The allele misassignment rate is the percentage of mismmatches between the IUPAC consensus sequences between homologous loci from each pair of duplicates. Finally, the two SNP error rates: the percentage of different SNPs called in each of the duplicated samples using either all 10178 SNPs or using the SNPs called without missing data between duplicate samples excluded (see Table 1).

Variant calling and filtering

We combined the data from A. sylvaticus and A. flavicollis to establish species differentiation and then filtered the SNPs using the population package from Stacks [11] and VCFtools [16]. We kept SNPs from the loci present in the 80% of the individuals in each species (p=1, r=0.8) and excluded SNPs with minor allele frequencies MAF<0.05 and which deviated from the Hardy-Weinberg equilibrium (HWE) at P<0.05. We also removed sites with mean depth values lower than 20. We manually modified the chromosome numbers in the vcf file to input it into SNPhylo [34], which we used to build the tree. We set a missing rate (-M) of 1, minor allele frequencies (-m) of 0, linkage disequilibrium threshold (-l) of 1 and the -r option to skip the step of removing low quality data. Confidence values were estimated using 1000 bootstrap replicates. The root was manually fixed to separate both species. Principal Component Analysis (PCA) was performed using the R package Adegenet [29] (Fig. 3).

The set of divergent loci identified between the two species in Polish samples was tested for its ability to differentiate an extra set of samples from other locations in Europe and Tunisia. Ten A. flavicollis (2 samples from Austria, 5 from Lithuania and 3 from Romania) and 10 samples of A. sylvaticus (4 samples from Wales, 3 from Tunisia and 3 from Scotland) were kindly provided by Dr Jeremy Herman, National Museums Scotland, Dr Johan Michaux, University of Liege and Dr Karol Zub, Mammal Research Institute of the Polish Academy of Sciences (MRI) (Supplementary Table S4). We considered all 20 test samples as a different group from Polish A. sylvaticus and A. flavicollis for SNP calling. We kept SNPs from the loci present in the 80% of the individuals in each group (p=1, r=0.8) and excluded SNPs with minor allele frequencies MAF<0.05, SNPs which deviated from the Hardy-Weinberg equilibrium (HWE) at P<0.05, sites with mean depth values lower than 20 and with more than 5% of missing data.

Population divergence

To analyse genetic diversity and population connectivity within A. flavicollis, we analysed the three populations (Bory Tucholskie, Białowieża and Haćki) separately (p=3, r=0.8), while keeping the other parameters as described above. We note that absence of filtering out SNPs in Hardy-Weinberg disequlibrium did not markedly affect the reported pairwise Fst values (data not shown). Due to the lack of outgroup, a mid-point root was chosen in the phylogenetic tree. Individual ancestries were estimated following a maximum likelihood approach with ADMIXTURE [3], after conversion of the VCF file to ped with plink version 1.9 [12, 55]. ADMIXTURE analysis was run for each of K=1 to K=5, each using 10 different seeds. Weighted (Weir-Cockerham) Fst was calculated with VCFtools v0.1.13. Heterozygosity, Pi and Fis were calculated with the population package from Stacks [11].

Species divergence

To calculate the p-distance between the two species, first a set of common loci was extracted with a custom script and the strict consensus sequences for each species were calculated with Consensus.pl script [25] with no threshold parameter set, such that the most common nucleotide was set as a consensus. p-distance was then calculated using a custom R script that counts the number of differences between pairs of sequences (one from each species, for each of the 21377 loci) (Supplementary Materials, Section 9). The same set of scripts was also used to select loci with fixed differences within each species and then to calculate the p-distance between them.

Availability of data and materials

The sequence data, in the format of the output from the Stacks 1.48 proces_radtags (four files per sample), have been deposited in EBI SRA repository under accession number PRJNA554851. The scripts used in the analysis are deposited at GitHub.com, with links available in the Supplementary Materials.

Abbreviations

ddRAD-seq:

double-digest restriction site-associated DNA sequencing

Fst:

fixation index

HWE:

Hardy-Weinberg equilibrium

MAF:

minor allele frequency

mtDNA:

mitochondrial DNA

PCA:

principal component analysis

SNP:

single nucleotide polymorphism

VCF:

variant call format

References

  1. 1

    Aghová T, Kimura Y, Bryja J, Dobigny G, Granjon L, Kergoat GJ. Fossils know it best: Using a new set of fossil calibrations to improve the temporal phylogenetic framework of murid rodents (Rodentia: Muridae). Mol Phylogenet Evol. 2018; 128:98–111.

    PubMed  Google Scholar 

  2. 2

    Alexander DH, Lange K. Enhancements to the admixture algorithm for individual ancestry estimation. BMC bioinformatics. 2011; 12(1):246.

    PubMed  PubMed Central  Google Scholar 

  3. 3

    Alexander DH, Novembre J, Lange K. Fast model-based estimation of ancestry in unrelated individuals. Genome Res. 2009; 19(9):1655–1664.

    CAS  PubMed  PubMed Central  Google Scholar 

  4. 4

    Alter SE, Munshi-South J, Stiassny MLJ. Genome wide snp data reveal cryptic phylogeographic structure and microallopatric divergence in a rapids-adapted clade of cichlids from the congo river. Mol Ecol. 2017; 26(5):1401–1419.

    PubMed  Google Scholar 

  5. 5

    Andrews KR, Good JM, Miller MR, Luikart G, Hohenlohe PA. Harnessing the power of radseq for ecological and evolutionary genomics. Nat Rev Genet. 2016; 17(2):81.

    CAS  PubMed  PubMed Central  Google Scholar 

  6. 6

    Baird NA, Etter PD, Atwood TS, Currey MC, Shiver AL, Lewis ZA, Selker EU, Cresko WA, Johnson EA. Rapid snp discovery and genetic mapping using sequenced rad markers. PloS one. 2008; 3(10):e3376.

    PubMed  PubMed Central  Google Scholar 

  7. 7

    Baxter SW, Davey JW, Johnston JS, Shelton AM, Heckel DG, Jiggins CD, Blaxter. ML. Linkage mapping and comparative genomics using next-generation rad sequencing of a non-model organism. PloS one. 2011; 6(4):e19315.

    CAS  PubMed  PubMed Central  Google Scholar 

  8. 8

    Blanco-Bercial L, Bucklin A. New view of population genetics of zooplankton: Rad-seq analysis reveals population structure of the north atlantic planktonic copepod centropages typicus. Mol Ecol. 2016; 25(7):1566–1580.

    CAS  PubMed  Google Scholar 

  9. 9

    Britten RJ. Divergence between samples of chimpanzee and human dna sequences is 5%, counting indels. Proc Natl Acad Sci. 2002; 99(21):13633–13635.

    CAS  PubMed  Google Scholar 

  10. 10

    Bugarski-Stanojević V, Blagojević J, Adnadjević T, Jovanović V, Vujošević M. Identification of the sibling species apodemus sylvaticus and apodemus flavicollis (rodentia, muridae)—comparison of molecular methods. Zool Anz J Comp Zool. 2013; 252(4):579–587.

    Google Scholar 

  11. 11

    Catchen JM, Amores A, Hohenlohe P, Cresko W, Postlethwait JH. Stacks: building and genotyping loci de novo from short-read sequences. G3: Genes, genomes, genetics. 2011; 1(3):171–182.

    CAS  Google Scholar 

  12. 12

    Chang CC, Chow CC, Tellier LC, Vattikuti s, Purcell sm. Second-generation plink: rising to the challenge of larger and richer datasets. Gigascience. 2015; 4(1):7.

    PubMed  PubMed Central  Google Scholar 

  13. 13

    Cromie GA, Hyma KE, Ludlow CL, Garmendia-Torres C, Gilbert TL, May P, Huang AA, Dudley AM, Fay JC. Genomic sequence diversity and population structure of saccharomyces cerevisiae assessed by rad-seq. G3: Genes, Genomes, Genet. 2013; 3(12):2163–2171.

    Google Scholar 

  14. 14

    Cull B, Vaux AGC, Ottowell LJ, Gillingham EL, Medlock JM. Tick infestation of small mammals in an english woodland. J Vector Ecol. 2017; 42(1):74–83.

    PubMed  Google Scholar 

  15. 15

    Czarnomska SD, Niedziałkowska M, Borowik T, Jędrzejewska B. Regional and local patterns of genetic variation and structure in yellow-necked mice-the roles of geographic distance, population abundance, and winter severity. Ecol Evol. 2018; 8:8171–8186.

    PubMed  PubMed Central  Google Scholar 

  16. 16

    Danecek P, Auton A, Abecasis G, Albers CA, Banks E, DePristo MA, Handsaker RE, Lunter G, Marth GT, Sherry ST, et al.The variant call format and vcftools. Bioinformatics. 2011; 27(15):2156–2158.

    CAS  PubMed  PubMed Central  Google Scholar 

  17. 17

    Darvish J, Mohammadi Z, Ghorbani F, Mahmoudi A, Dubey S. Phylogenetic relationships of apodemus kaup, 1829 (rodentia: Muridae) species in the eastern mediterranean inferred from mitochondrial dna, with emphasis on iranian species. J Mamm Evol. 2015; 22(4):583–595.

    Google Scholar 

  18. 18

    Fernández R, Schubert M, Vargas-Velázquez AM, Brownlow A, Víkingsson GA, Siebert U, Jensen LF, Øien N, Wall D, Rogan E, Mikkelsen B. A genomewide catalogue of single nucleotide polymorphisms in white-beaked and Atlantic white-sided dolphins. Mol Ecol Resour. 2016; 16(1):266–276.

    PubMed  Google Scholar 

  19. 19

    Franchini P, Parera DM, Kautt AF, Meyer A. quaddrad: a new high-multiplexing and pcr duplicate removal ddrad protocol produces novel evolutionary insights in a nonradiating cichlid lineage. Mol Ecol. 2017; 26(10):2783–2795.

    CAS  PubMed  Google Scholar 

  20. 20

    Gortat Tt, Gryczyńska-Siemiątkowska A, Rutkowski R, Kozakiewicz A, Mikoszewski A, Kozakiewicz M. Landscape pattern and genetic structure of a yellow-necked mouse apodemus flavicollis population in north-eastern poland. Acta Theriol. 2010; 55(2):109–121.

    Google Scholar 

  21. 21

    Hammerman NM, Rivera-Vicens RE, Galaska MP, Weil E, Appledoorn RS, Alfaro M, Schizas NV. Population connectivity of the plating coral agaricia lamarcki from southwest puerto rico. Coral Reefs. 2018; 37(1):183–191.

    Google Scholar 

  22. 22

    Herman JS, Jóhannesdóttir FJ, Jones EP, McDevitt AD, Michaux JR, White TA, Wójcik JN, Searle JB. Post-glacial colonization of europe by the wood mouse, apodemus sylvaticus: evidence of a northern refugium and dispersal with humans. Biol J Linn Soc. 2017; 120(2):313–332.

    Google Scholar 

  23. 23

    Hipp AL, Eaton DAR, Cavender-Bares J, Fitzek E, Nipper R, Manos PS. A framework phylogeny of the american oak clade based on sequenced rad data. PLoS One. 2014; 9(4):e93975.

    PubMed  PubMed Central  Google Scholar 

  24. 24

    Hohenlohe PA, Day MD, Amish SJ, Miller MR, Kamps-Hughes N, Boyer MC, Muhlfeld CC, Allendorf FW, Johnson EA, Luikart G. Genomic patterns of introgression in rainbow and westslope cutthroat trout illuminated by overlapping paired-end rad sequencing. Mol Ecol. 2013; 22(11):3002–3013.

    CAS  PubMed  PubMed Central  Google Scholar 

  25. 25

    Hughes J. Sequence-manipulation. GitHub Repository. 2011. https://github.com/josephhughes/Sequencemanipulation.

  26. 26

    Janes JK, Miller JM, Dupuis JR, Malenfant RM, Gorrell JC, Cullingham CI, Andrew RL. The k= 2 conundrum. Mol Ecol. 2017; 26(14):3594–3602.

    PubMed  Google Scholar 

  27. 27

    Jeffries DL, Copp GH, Handley LL, Olsén KH, Sayer CD, Hänfling B. Comparing radseq and microsatellites to infer complex phylogeographic patterns, an empirical perspective in the crucian carp, carassius carassius, l. Mol Ecol. 2016; 25(13):2997–3018.

    PubMed  Google Scholar 

  28. 28

    Jojić V, Bugarski-Stanojević V, Blagojević J, Vujošević M. Discrimination of the sibling species apodemus flavicollis and a. sylvaticus (rodentia, muridae). Zool Anz J Comp Zool. 2014; 253(4):261–269.

    Google Scholar 

  29. 29

    Jombart T. adegenet: a r package for the multivariate analysis of genetic markers. Bioinformatics. 2008; 24(11):1403–1405.

    CAS  PubMed  Google Scholar 

  30. 30

    Knowles LL, Massatti R, He Q, Olson LE, Lanier HC. Quantifying the similarity between genes and geography across alaska’s alpine small mammals. J Biogeogr. 2016; 43(7):1464–1476.

    Google Scholar 

  31. 31

    Kolodziej M, Melgies A, Joniec-Wiechetek J, Michalski A, Nowakowska A, Pitucha G, Niemcewicz M. First molecular characterization of dobrava-belgrade virus found in apodemus flavicollis in poland. Ann Agric Environ Med. 2018; 25(2):368–373.

    CAS  PubMed  Google Scholar 

  32. 32

    Lanier HC, Massatti R, He Q, Olson LE, Knowles. LL. Colonization from divergent ancestors: glaciation signatures on contemporary patterns of genomic variation in collared pikas (ochotona collaris). Mol Ecol. 2015; 24(14):3688–3705.

    PubMed  Google Scholar 

  33. 33

    Lapinski AG, Pavlenko MV, Solovenchuk LL, Gorbachev VV. Some limitations in the use of the mitochondrial dna cytb gene as a molecular marker for phylogenetic and population-genetic studies by the example of the apodemus genus. Russ J Genet Appl Res. 2016; 6(1):84–90.

    Google Scholar 

  34. 34

    Lee T-H, Guo H, Wang X, Kim C, Paterson AH. Snphylo: a pipeline to construct a phylogenetic tree from huge snp data. BMC Genomics. 2014; 15(1):162.

    PubMed  PubMed Central  Google Scholar 

  35. 35

    Li W-H, Tanimura M, Sharp PM. An evaluation of the molecular clock hypothesis using mammalian dna sequences. J Mol Evol. 1987; 25(4):330–342.

    CAS  PubMed  Google Scholar 

  36. 36

    Martiniaková M, Omelka R, Jancová A, Stawarz R, Formicki G. Heavy metal content in the femora of yellow-necked mouse (apodemus flavicollis) and wood mouse (apodemus sylvaticus) from different types of polluted environment in slovakia. Environ Monit Assess. 2010; 171(1-4):651–660.

    PubMed  Google Scholar 

  37. 37

    Mastretta-Yanes A, Arrigo N, Alvarez N, Jorgensen TH, Piñero D, Emerson BC. Restriction site-associated dna sequencing, genotyping error estimation and de novo assembly optimization for population genetic inference. Mol Ecol Resour. 2015; 15(1):28–41.

    CAS  PubMed  Google Scholar 

  38. 38

    René Michaux J, Magnanou E, Paradis E, Nieberding C, Libois R. Mitochondrial phylogeography of the woodmouse (apodemus sylvaticus) in the western palearctic region. Mol Ecol. 2003; 12(3):685–697.

    Google Scholar 

  39. 39

    Michaux JR, Chevret P, Filippucci M-G, Macholan M. Phylogeny of the genus apodemus with a special emphasis on the subgenus sylvaemus using the nuclear irbp gene and two mitochondrial markers: cytochrome b and 12s rrna. Mol Phylogenet Evol. 2002; 23(2):123–136.

    CAS  PubMed  Google Scholar 

  40. 40

    Michaux JR, Libois R, Paradis E, Filippucci M-G. Phylogeographic history of the yellow-necked fieldmouse (apodemus flavicollis) in europe and in the near and middle east. Mol Phylogenet Evol. 2004; 32(3):788–798.

    CAS  PubMed  Google Scholar 

  41. 41

    Michaux JR, Libois R, Filippucci MG. So close and so different: comparative phylogeography of two small mammal species, the yellow-necked fieldmouse (apodemus flavicollis) and the woodmouse (apodemus sylvaticus) in the western palearctic region. Heredity. 2005; 94(1):52–63.

    CAS  PubMed  Google Scholar 

  42. 42

    Miller MR, Dunham JP, Amores A, Cresko WA, Johnson EA. Rapid and cost-effective polymorphism identification and genotyping using restriction site associated dna (rad) markers. Genome Res. 2007; 17(2):240–248.

    CAS  PubMed  PubMed Central  Google Scholar 

  43. 43

    Mlera L, Bloom ME. The role of mammalian reservoir hosts in tick-borne flavivirus biology. Front Cell Infect Microbiol. 2018; 8:298.

    PubMed  PubMed Central  Google Scholar 

  44. 44

    Moura AE, Kenny JG, Chaudhuri R, Hughes MA, Welch AJ, Reisinger RR, de Bruyn PJN, Dahlheim ME, Hall N, Hoelzel AR. Population genomics of the killer whale indicates ecotype evolution in sympatry involving both selection and drift. Mol Ecol. 2014; 23(21):5179–5192.

    PubMed  PubMed Central  Google Scholar 

  45. 45

    Pante E, Abdelkrim J, Viricel A, Gey D, France SC, Boisselier M-C, Samadi S. Use of rad sequencing for delimiting species. Heredity. 2015; 114(5):450.

    CAS  PubMed  Google Scholar 

  46. 46

    Papa A, Rogozi E, Velo E, Papadimitriou E, Bino S. Genetic detection of hantaviruses in rodents, albania. J Med Virol. 2016; 88(8):1309–1313.

    PubMed  Google Scholar 

  47. 47

    Paris JR, Stevens JR, Catchen JM. Lost in parameter space: a road map for stacks. Methods Ecol Evol. 2017; 8(10):1360–1373. https://doi.org/10.1111/2041-210x.12775.

    Google Scholar 

  48. 48

    Pegadaraju V, Nipper R, Hulke B, Qi L, Schultz Q. De novo sequencing of sunflower genome for snp discovery using rad (restriction site associated dna) approach. BMC Genomics. 2013; 14(1):556.

    CAS  PubMed  PubMed Central  Google Scholar 

  49. 49

    Peterson BK, Weber JN, Kay EH, Fisher HS, Hoekstra HE. Double digest radseq: an inexpensive method for de novo snp discovery and genotyping in model and non-model species. PloS One. 2012; 7(5):e37135.

    CAS  PubMed  PubMed Central  Google Scholar 

  50. 50

    Poland JA, Rife TW. Genotyping-by-sequencing for plant breeding and genetics. Plant Genome. 2012; 5(3):92–102.

    CAS  Google Scholar 

  51. 51

    Poland JA, Brown PJ, Sorrells ME, Jannink J-L. Development of high-density genetic maps for barley and wheat using a novel two-enzyme genotyping-by-sequencing approach. PloS One. 2012; 7(2):e32253.

    CAS  PubMed  PubMed Central  Google Scholar 

  52. 52

    Pucek Z, Stachurska J, Bibrich M, et al.Keys to vertebrates of poland: mammals. 1981.

  53. 53

    Pucek Z, Jędrzejewski W, Jędrzejewska B, Pucek M. Rodent population dynamics in a primeval deciduous forest (białowieża national park) in relation to weather, seed crop, and predation. Acta Theriol. 1993; 38(2):199–232.

    Google Scholar 

  54. 54

    Puechmaille SJ. The program structure does not reliably recover the correct population structure when sampling is uneven: subsampling and new estimators alleviate the problem. Mol Ecol Resour. 2016; 16(3):608–627.

    PubMed  Google Scholar 

  55. 55

    Purcell S, Chang CC. Plink 1.9 package. 2018. https://www.cog-genomics.org/plink/1.9/. Accessed 27 Feb 2019.

  56. 56

    Rajičić M, Romanenko SA, Karamysheva TV, Blagojević J, Adnađević T, Budinski I, Bogdanov AS, Trifonov VA, Rubtsov NB, Vujošević M. The origin of b chromosomes in yellow-necked mice (apodemus flavicollis)—break rules but keep playing the game. PloS One. 2017; 12(3):e0172704.

    PubMed  PubMed Central  Google Scholar 

  57. 57

    Reitzel AM, Herrera S, Layden MJ, Martindale MQ, Shank TM. Going where traditional markers have not gone before: utility of and promise for rad sequencing in marine invertebrate phylogeography and population genomics. Mol Ecol. 2013; 22(11):2953–2970.

    CAS  PubMed  PubMed Central  Google Scholar 

  58. 58

    Richter D, Schlee DB, Matuschka F-R. Reservoir competence of various rodents for the lyme disease spirochete borrelia spielmanii. Appl Environ Microbiol. 2011; 77(11):3565–3570. https://doi.org/10.1128/aem.00022-11.

    CAS  PubMed  PubMed Central  Google Scholar 

  59. 59

    Rico A, Kindlmann P, Sedláček F. Can the barrier effect of highways cause genetic subdivision in small mammals?. Acta Theriol. 2009; 54(4):297–310.

    Google Scholar 

  60. 60

    Rodriguez-Ezpeleta N, Álvarez P, Irigoien X. Genetic diversity and connectivity in maurolicus muelleri in the bay of biscay inferred from thousands of snp markers. Front Genet. 2017; 8:195.

    PubMed  PubMed Central  Google Scholar 

  61. 61

    Shafer A, Peart CR, Tusso S, Maayan I, Brelsford A, Wheat CW, Wolf JBW. Bioinformatic processing of rad-seq data dramatically impacts downstream population genetic inference. Methods Ecol Evol. 2017; 8(8):907–917.

    Google Scholar 

  62. 62

    Shultz AJ, Baker AJ, Hill GE, Nolan PM, Edwards SV. Snp s across time and space: population genomic signatures of founder events and epizootics in the house finch (haemorhous mexicanus). Ecol Evol. 2016; 6(20):7475–7489.

    PubMed  PubMed Central  Google Scholar 

  63. 63

    Velickovic M. Measures of the developmental stability, body size and body condition in the black-striped mouse (apodemus agrarius) as indicators of a disturbed environment in northern serbia. Belg J Zool. 2007; 137(2):147.

    Google Scholar 

  64. 64

    Vujošević M, Rajičić M, Blagojević J. B chromosomes in populations of mammals revisited. Genes. 2018; 9(10):487.

    PubMed Central  Google Scholar 

Download references

Acknowledgements

The authors wish to thank Jeremy Herman (National Museums Scotland) and Johan Michaux (Université de Liège) for providing the European/Tunesian samples for testing our catalog of SNPs. MLMC and JB also wish to acknowledge the use of the Orion High Performance Computing cluster at the School of Applied Sciences, University of Huddersfield.

Funding

This project was supported by the University of Huddersfield, the Friedrich Miescher Laboratory of the Max Planck Society and Microsoft Azure Research Award CRM:0518338.

Author information

Affiliations

  1. School of Applied Sciences, University of Huddersfield, Quennsgate, Huddersfield, UK

    Maria Luisa Martin Cerezo & Jarosław Bryk

  2. AVIAN Behavioural Genomics and Physiology Group, IFM Biology, Department of Zoology, Linköping University, Linköping, Sweden

    Maria Luisa Martin Cerezo

  3. Friedrich Miescher Laboratory of the Max Planck Society, Tübingen, Germany

    Marek Kucka & Yingguang Frank Chan

  4. The Mammal Research Institute, Polish Academy of Sciences, Białowieża, Poland

    Karol Zub

Authors
  1. Maria Luisa Martin Cerezo
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Marek Kucka
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. Karol Zub
    View author publications

    You can also search for this author in PubMed Google Scholar

  4. Yingguang Frank Chan
    View author publications

    You can also search for this author in PubMed Google Scholar

  5. Jarosław Bryk
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

MLMC designed the study, prepared the sequencing library, analysed the data and wrote the manuscript. MK contributed to the experimental design and preparation of the sequencing library. KZ contributed to the experimental design, performed the sampling and contributed to the manuscript. YFC contributed to the analysis of the data and the manuscript. JB designed the study, analysed the data and wrote the manuscript. The author(s) read and approved the final manuscript.

Corresponding author

Correspondence to Jarosław Bryk.

Ethics declarations

Ethics approval and consent to participate

The sampling procedures were approved by the Local Ethical Commission on Experimentation on Animals in Białystok, Poland, under permission number 2015/99.

Consent for publication

Not Applicable.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary information

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

Martin Cerezo, M.L., Kucka, M., Zub, K. et al. Population structure of Apodemus flavicollis and comparison to Apodemus sylvaticus in northern Poland based on RAD-seq. BMC Genomics 21, 241 (2020). https://doi.org/10.1186/s12864-020-6603-3

Download citation

Keywords

  • RAD-seq
  • genotyping
  • population structure
  • rodents
  • Apodemus flavicollis
  • Apodemus sylvaticus

BMC Genomics

ISSN: 1471-2164

Contact us