https://orcid.org/0000-0001-7387-5455
Figures
Abstract
Arvicolinae is one of the most impressive placental radiations with over 150 extant and numerous extinct species that emerged since the Miocene in the Northern Hemisphere. The phylogeny of Arvicolinae has been studied intensively for several decades using morphological and genetic methods. Here, we sequenced 30 new mitochondrial genomes to better understand the evolutionary relationships among the major tribes and genera within the subfamily. The phylogenetic and molecular dating analyses based on 11,391 bp concatenated alignment of protein-coding mitochondrial genes confirmed the monophyly of the subfamily. While Bayesian analysis provided a high resolution across the entire tree, Maximum Likelihood tree reconstruction showed weak support for the ordering of divergence and interrelationships of tribal level taxa within the most ancient radiation. Both the interrelationships among tribes Lagurini, Ellobiusini and Arvicolini, comprising the largest radiation and the position of the genus Dinaromys within it also remained unresolved. For the first time complex relationships between genus level taxa within the species-rich tribe Arvicolini received full resolution. Particularly Lemmiscus was robustly placed as sister to the snow voles Chionomys in the tribe Arvicolini in contrast with a long-held belief of its affinity with Lagurini. Molecular dating of the origin of Arvicolinae and early divergences obtained from the mitogenome data were consistent with fossil records. The mtDNA estimates for putative ancestors of the most genera within Arvicolini appeared to be much older than it was previously proposed in paleontological studies.
Citation: Abramson NI, Bodrov SY, Bondareva OV, Genelt-Yanovskiy EA, Petrova TV (2021) A mitochondrial genome phylogeny of voles and lemmings (Rodentia: Arvicolinae): Evolutionary and taxonomic implications. PLoS ONE 16(11): e0248198. https://doi.org/10.1371/journal.pone.0248198
Editor: Laurent Viriot, Ecole normale superieure de Lyon, FRANCE
Received: February 19, 2021; Accepted: November 3, 2021; Published: November 19, 2021
Copyright: © 2021 Abramson et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.
Data Availability: The data are held in a public repository NCBI (https://www.ncbi.nlm.nih.gov/nucleotide/). All sequences numbers obtained in the current study are given in the supporting materials S1 Table.
Funding: This study was conducted in the Zoological Institute RAS with financial support from the Russian Science Foundation grant N19-74-20110 for all authors.
Competing interests: The authors have declared that no competing interests exist.
Introduction
The subfamily Arvicolinae Gray, 1821 (Rodentia: Cricetidae), voles, lemmings and muskrats, is a highly diverse, young and fast-evolving group within the order Rodentia. Currently, representatives of the subfamily occupy most of the temperate and cold-climate terrestrial habitats across the Northern Hemisphere. Modern global fauna of Arvicolinae consists of more than 150 recent species grouped into 28 genera, with new species being constantly discovered and described [1,2].
Phylogeny of Arvicolinae has been explored using both morphological and genetic methods, allowing comparisons of reconstruction from different datasets for further cross-validation. Application of molecular phylogenetic methods resulted in a series of revisions of phylogenetic relationships and taxonomic structure of several genera and species [3–12], generally supporting the hypothesis of three successive waves of adaptive radiations in the evolutionary history of the group [7]. According to the paleontological records, the ancestors of Arvicolinae were specialised “microtoid cricetids’’, sharing the dental adaptations of increased hypsodonty. Rich fossil records of these lineages in the late Miocene sites of Eurasia and their absence in North America indicates that the origin and primary evolutionary centre for arvicolids was situated in the northern parts of Asia [13,14].
Arvicolinae sensu stricto presumably emerged during the rapid species diversification in late Miocene—early Pliocene, when global climate became cooler and drier and woodland habitats were largely replaced by open grassy landscapes. The members of this first radiation were muskrats (Ondatrini), lemmings (Lemmini and Dicrostonychini) and long-clawed mole voles (Prometheomyini). Second radiation wave is characterized by a divergence of the ancestors of modern Clethrionomyini, the compact tribe of red-backed and mountain voles [15], currently abundant in boreal forests and highlands. The third radiation wave resulted in a formation of the steppe lemmings (Lagurini), mole voles (Ellobiusini) and the richest species group—Arvicolini [7]. The branching order both within the first and the last radiation waves also remained unresolved since all attempts to untangle complex phylogenetic relationships within subfamily were made with the use of only a few mitochondrial or nuclear markers, or the study included an insufficient number of taxa in the analysis [3,5,7,9,10,12,16–21].
The recent review of Muroid rodent phylogeny based on mitochondrial cytochrome b and five nuclear genes confirmed [21] findings of previous studies [7] that the proper taxonomic sampling is crucial to build the large-scale tree reconstructions. The broad introduction of next-generation sequencing (NGS) resulted in constant accumulation of published mitochondrial genomes of Arvicolinae [22–37]. As a result, mitochondrial genome sequences of nearly one-fifth of the total species diversity of the subfamily (ca. 30 species) is already available.
Mitochondrial DNA has many advantages for the phylogenetic studies, as it possesses a strict maternal inheritance [38] with high mutation rate due to a limited repair system (5–10 times that of nuclear DNA) [39] and simple conserved structure. Due to the high number of mitochondria in the cell almost the complete mitochondrial genome may be sequenced from old historical samples from museum collections thus the currently hardly available taxa may be included in the data set. Complete mitochondrial genome datasets provide a robust phylogeny with highly resolved trees having strong branch support [40]. However, the mitochondrial genomes can conflict with a phylogeny estimated from the nuclear genome due to the replacement of organelle genomes from one species or populations with those of another mediated by hybridization and introgression [41–44]. The use of mitochondrial genome datasets may also provide evidences of adaptive radiations, as signatures of selection were detected during reconstructing mitochondrial genome phylogenies for the number of vertebrate taxa including deep-sea fishes [45,46], marine turtles [47,48], and mammals [49–54]. Yet, mitochondrial genomes remain an important tool both for phylogenies resolution and for the species identification by barcoding using individual mitochondrial genes or complete genomes. Consecutive sequencing of new mitochondrial genomes recently provided new insights into many animal groups [55–59], several mammal phylogenies including e.g. fruit bats [58], carnivores [41], wood mice [60] and zokors [61].
Though molecular studies during the last decades have considerably extended and refined our knowledge of the pattern and timescale of Arvicolinae phylogeny, there are important issues that remain to be elucidated. In this study, we were aimed at estimating the phylogeny of voles and lemmings using complete mitogenomes generated using high-throughput sequencing. By significantly increasing the number of newly sequenced mitogenomes representing major tribes of voles and lemmings, we implement the phylogenetic and molecular dating analysis on a dataset consisting of almost all living genera within the subfamily. The following questions were specifically addressed during the study: (1) the order of divergence and interrelationships of taxa within the first, most ancient radiation, (2) the interrelationships of three tribe level taxa Lagurini, Ellobiusini and Arvicolini (3) relative phylogenetic placement of genera Dinaromys Kretzoi, 1955 and Lemmiscus Thomas, 1912, (4) tangled interrelationships of genera and subgenera in the most speciose tribe Arvicolini, (5) the position of Agricola agrestis Linnaeus, 1761 and Iberomys cabrerae Thomas, 1906, and (6) the timing of arvicoline divergences.
Material and methods
Ethics statement
Our study was performed using tissue and vaucher collection of the Zoological Institute Russian Academy of Sciences and the research did not require field work or live animal experimentation. DNA samples for this study were obtained following the guidelines of the Animal Care and Research Committee of the Zoological Institute RAS, Saint Petersburg and in accordance with respective national and international animal care and use policies Permission N 3-7/15-06-2021 (See S1 Table). Tissues of specimens used in the study are publicly deposited and accessible by others in a permanent repository of Zoological Institute Russian Academy of Sciences, laboratory of molecular systematics.
Taxonomic sampling
Fifty-eight species of Arvicolinae, belonging to 27 genera and all the tribe level taxa, as well as six outgroup taxa were used in this study. Complete mitochondrial genomes for 30 Arvicolinae species were sequenced in the current study (including 15 species belonging to the Arvicolini tribe, one for Dicrostonychini, two for Lagurini, three for Lemmini, six species belonging to Clethrionomyini and three crucial species without stable taxonomic position: Prometheomys schaposchnikowi Satunin, 1901 (Prometheomyini), Dinaromys bogdanovi Martino, 1922 and Lemmiscus curtatus Cope, 1868). For 28 species belonging to Arvicolini, Ellobiusini, Clethrionomyini, Dicrostonychini and Ondatrini tribes, sequences were available in the NCBI nucleotide database (https://www.ncbi.nlm.nih.gov/nuccore). The detailed information, GenBank accession numbers, and the voucher IDs for new sequences are given in S1 Table. Hereinafter, we use the taxonomic classification following Gromov & Polyakov [62], Musser & Carleton [1], Abramson & Lissovsky [63] with amendments made in result of the current study.
DNA isolation, NGS library preparation and sequencing
Muscle tissue samples of fresh specimens were collected between 1996–2019 years and stored in 96% ethanol at -20 degrees Celcius in a tissue and DNA collection of the Group of molecular systematics of mammals (Zoological Institute RAS). Historic specimen of Lemmiscus curtatus (sampled in 1927) was obtained from the collection of the laboratory of theriology (Zoological Institute RAS), see S1 Table for details.
Homogenization of tissues was performed using the Qiagen TissueLyser LT (Qiagen). For most samples, genomic DNA was extracted using the Diatom DNA Prep 200 (Isogen, Russia) except for the L. curtatus museum specimen.
To reduce the potential contamination, all manipulations with the L. curtatus were carried out using a PCR workstation (LAMSYSTEMS CC, Miass, Russia) in a separate laboratory room isolated from post-PCR facilities, being used exclusively for studies of historic samples from the collection of Zoological Institute. All the working surfaces, instruments and plastics were sterilized with UV light and chloramine-T. DNA from the museum skin sample (2 × 2 mm piece from the inner side of the lip, dissected by a sterilized surgical blade) was isolated using the phenol-chloroform extraction method according to a standard protocol [64].
The following ultrasound fragmentation of the total genomic DNA was implemented using Covaris S220 focused ultrasonicator instrument (Covaris). The resulting fragmented DNA was purified and concentrated using paramagnetic bead-based chemistry AMPure XP (Beckman-Coulter) using standard workflow. DNA concentration was evaluated using a Qubit fluorometer (Thermo Fisher Scientific).
NGS libraries were prepared using the NEBNext Ultra II DNA Library Prep Kit for Illumina (New England Biolabs). The resulting PCR products were purified and concentrated using AMPure XP beads (Beckman-Coulter). The concentration of samples was measured using a Qubit fluorometer, and quality control of the libraries was implemented using Bioanalyzer 2100 instrument and the DNA High Sensitivity kit (Agilent). Sequencing was performed on an Illumina HiSeq 4000 system, resulting in pair-end reads of 75bp. DNA quality was checked with Qubit, the final distribution of lengths of the libraries adapter content checking was conducted using Bioanalyzer2100 (Agilent). DNA extraction (except the museum specimen of L. curtatus), library preparation and sequencing were performed using resources of the Skoltech Genomics Core Facility (https://www.skoltech.ru/research/en/shared-resources/gcf-2/). Standard precautions were taken to avoid contamination according to the Illumina’s recommendations at all stages of library preparation and sequencing [65].
Read processing, mitogenome assembly and annotation
The quality of raw reads was evaluated using FastQC [66], and parts with the quality score below 20 were trimmed using Trimmomatic-0.32 [67]. Bowtie 2.3.5.1 [68] was used to filter reads with contamination. Complete mitochondrial genomes of other Arvicolinae were used as reference sequences. Also, this was made for the museum specimen to enrich reads with mitochondrial DNA.
Nucleotide misincorporation patterns that can often be observed during the studies of ancient or old museum sample DNA as a result of post-mortem DNA damage in reads from L. curtatus were achieved using mapDamage 2.0 [69].
Complete mitochondrial genome was assembled using plasmidSPAdes [70] with default settings. The resulting contigs were filtered by length, the most similar in size to mitochondrial DNA were selected (size about 16 kb for mammals). Raw reads of L. curtatus were mapped on the mitochondrial genome of Mynomes ochrogaster Wagner, 1842 with manual settings using Geneious Prime 2019.1 (Biomatters Ltd.) due to the fragmentation of the DNA extracted from the museum sample (S1 Fig). The contigs were annotated using the MITOS web server [71] (http://mitos2.bioinf.uni-leipzig.de/index.py), with default settings and the vertebrate mitochondrial genetic code.
Gene boundaries were checked and refined by alignment against 28 published mitogenome sequences of Arvicolinae (see details in S1 Table). All positions of low quality, low coverage, as well as fragments that greatly differed from the reference Arvicolinae mitochondrial genomes, were replaced by N manually. Assembled sequences of protein coding genes (PCGs) were checked for internal stops manually. Raw read data are available from the a SRA database (PRJNA590630), and all assembled and annotated mitogenomes were deposited in NCBI GenBank (S1 Table).
Sequence alignment
The 30 newly obtained mitochondrial genomes were compared with 28 earlier published Arvicolinae mitogenomes mined from NCBI (see accession numbers in S1 Table), including and five mitogenomes obtained by us earlier [72–74]. All mitochondrial genomes were aligned with Mauve (http://darlinglab.org/mauve/mauve.html) in Geneious Prime 2019.1 (Biomatters Ltd.).
In several studies, it has been convincingly shown that protein-coding sequences may have a strong resolving power for inferring phylogenetic interrelationships and divergence time estimates derived from PCG may be quite accurate [58,75,76]. We used this approach, however complete mt genomes will serve as the starting point for further analyses. For the subsequent analyses, the concatenated alignment of 13 PCGs using MAFFT version 7.222 [77] was produced.
Third codon position has previously been shown to bias phylogenetic reconstructions [78]. The phylogeny on a smaller dataset of Arvicolinae turned out to be very poorly resolved with the exception of the third codon position. So we masked transitions in 3rd codon position by RY-coding (R for purines and Y for pyrimidines) as described in Abramson et al. [74].
Thus, two datasets were subsequently analyzed—total alignment of 13 PCGs, where all three codon positions were considered (with a length of 11,391 bp) and RY-coded alignment with transitions in third codon position masked.
Analysis of base composition
The base composition was calculated in Geneious Prime 2019.1 (Biomatters Ltd.). The strand bias in nucleotide composition was studied by calculating the relative frequencies of C and G nucleotides (CG3 skew = [C − G]/[C + G]) [58,79,80]. Both analyses were calculated using full-length mitogenomes.
The PCG-alignment of 64 mitochondrial genomes was used to calculate relative frequencies of four bases (A, C, G and T) at each of three codon positions in MEGA X [81]. The 12 variables, each representing base frequency in first, second or third position, were then summarized by a Principal Component Analysis (PCA) using the PAST v.4.04 [82].
Saturation tests
The presence of phylogenetic signal was assessed with a substitution saturation analysis using the Xia test [83] in the DAMBE 7.2.1 software [84] for the whole alignment of the PCG dataset and 13 separate genes following the procedure described by Xia & Lemey [85], particularly when (a) 1st and 2nd codon position considered and 3rd position is masked from the alignment, and (b) when only 3rd codon position is included in the analysis. The analysis is based on Index of substitution saturation—Iss, and Iss.c is the critical value at which the sequences begin to fail to recover the true tree). Once Iss.c is known for a set of sequences, then we can calculate the Iss value from the sequences and compare it against the Iss.c. If Iss value exceeds the Iss.c, we can conclude that the sequence dataset consists of substitution saturation and cannot be used for further phylogenetic reconstruction.
The proportion of invariant sites was specified for tests considered 1st and 2nd codon positions. The analysis was performed on a complete alignment with all sites considered. Additional analysis of saturation for each of the PCG was estimated using R packages seqinr [86] and ape [87]. P-distances were plotted against K81 distances for transitions and transversions of each codon position.
Phylogenetic analyses
We used PartitionFinder 2.1.1 [88] applying AICc and “greedy” algorithm, when an analysis is based on the a priori features of the alignment, to select the optimal partitioning scheme for each dataset. Our analysis started with the partitioning by codon positions within PCG fragments, each treated as a unique partition. For the complete 13 PCG alignment, GTR+I+G model was suggested for almost all the partitions except ND6 3rd codon position, for which the TRN+I+G model was selected. For the alignment with RY-coded 3rd codon position, two partitions were suggested—1+2nd and 3rd codon positions with GTR+I+G and GTR+G models respectively (S2 Table).
Phylogenetic reconstructions using Maximum Likelihood (ML) and Bayesian Inference (BI) analyses were performed on both complete and RY-masked datasets. Trees were rooted by six species: Akodon montensis Thomas, 1913, Peromyscus megalops Merriam, 1898, Urocricetus kamensis Satunin, 1903 and three species of hamsters from genus Cricetulus Milne-Edwards, 1867 (S1 Table).
Maximum Likelihood (ML) analysis was performed using IQ-TREE web server [89] with 10,000 ultrafast bootstrap replicates [90]. Bayesian Inference (BI) analysis was performed in MrBayes 3.2.6 [91] with next settings: nst = mixed, rates = invgamma and partitions as suggested with PartitionFinder results (S2 Table). Each analysis started with random trees and performed two independent runs with four independent Markov Chain Monte Carlo (MCMC) for 10 million generations with sampling every 1,000th generation, the standard deviations of split frequencies were below 0.01; potential scale reduction factors were equal to 1.0; stationarity was examined in Tracer v1.7 [92]. A consensus tree was constructed based on the trees sampled after the 25% burn-in.
We also conducted ML analysis for each PCG separately with partitions by codon positions and models supposed automatically in IQ-TREE. Hyperacrius fertilis True, 1894 sequence was excluded from the ND4 alignment since this gene was highly fragmented [74]. The mitogenome of Craseomys rufocanus Sundevall, 1846 (accessed from GenBank) completely lacked the ND6 sequence (S2 Table), so this species was excluded from the analysis for this gene.
Divergence dating
Divergence times were estimated on the CIPRES Science Gateway [93] with Bayesian approach implemented in BEAST v.2.6.2 [94] using both complete PCG dataset and the one in which the transitions in the third codon position were masked with RY-coding. Datasets were partitioned according to the recommendations of PartitionFinder. All site model parameters were chosen for separate partitions with corrected Akaike’s information criterion (AICc) in JMODELTEST 2.1.1 [95]. Eight fossil calibrations were used (S3 Table). Lognormal prior distributions were applied to all the calibrations with offset values and 95% HPD intervals based on first appearance data (FAD) and stratigraphic sampling downloaded from the Paleobiology Database on 01.12.2020 using the parameters “Taxon = fossil species, Timescale = FAD” (S3 Table). To check the stability of the result, we performed an analysis with the alternate exclusion of each of the eight calibrations used (S4 Table).
BEAST analyses under the birth-death process used a relaxed lognormal clock model and the program’s default prior distributions of model parameters. Each analysis was run for 100 million generations and sampled every 10 000 generations. The convergence of two independent runs was examined using Tracer v1.7 [92], and combined using LogCombiner, discarding the first 25% as burn-in. Trees were then summarized with TreeAnnotator using the maximum clade credibility tree option and fixing node heights as mean heights. Divergence time bars were obtained automatically in FigTree v1.4.3 (http://tree.bio.ed.ac.uk/software/figtree/) from the output using the 95% highest posterior density (HPD) of the ages for each node. The comparison with an empty ‘prior run’ showed that the data were informative for estimating the divergence dates.
Results
Mitochondrial genome assembly and annotation
We sequenced, assembled and annotated mitochondrial genomes for the 30 new taxa of Arvicolinae. The mapDamage analysis implemented on the raw reads of Lemmiscus curtatus (S2 Fig) showed a low variation of deamination misincorporations values. C to T misincorporations varied from 10 to 15%, G to A from 10 to 12% and were equal to the results of similar studies [96]. Since the relative level of observed misincorporations was not significantly different from the other substitution variants, the mitogenome of L. curtatus was assembled using the same pipeline as for the rest of taxa.
The assembled mitogenomes, circular double-stranded DNA of the same organization as in other mammals, contained 13 PCGs, 22 transfer RNAs (tRNA), two ribosomal RNAs (rRNAs), and a non-coding region corresponding to the control region (D-loop). Nine genes (ND6 and eight tRNAs) were oriented in the reverse direction, whereas the others were transcribed in the forward direction. All the assembled mitogenomes contained all the genes listed above, but in some species demonstrated incomplete gene sequences (see S5 Table for details).
Mitochondrial genome sequences were deposited in GenBank under accession numbers indicated in S1 Table. In the subsequent analyses, the PCG dataset, containing 11,391 bp was used.
Variation in base composition
Comparison of base composition calculated using the alignment of full-length mitogenome sequences of Arvicolinae showed that mitochondrial genomes of taxa from tribes Clethrionomyini (28.36% С) and Ellobiusini (28.7% C) have the highest GC-content. Arvicolini, Dicrostonychini and Lemmini had slightly smaller values: 27.69, 28.03 and 27.77% C respectively. Lagurini were found to have the most AT-skewed base composition of mitogenomes: 31.20 and 31.30% of adenine, respectively (S1 Table). Lagurini and Arvicolini also demonstrated the highest GC-skew values (-0.32 in both cases). Ellobiusini and Clethrionomyini occupied an intermediate position in terms of it with -0.33 and -0.34, respectively. Dicrostonychini and Lemmini with equal value -0.35 have the smallest GC-skew values.
The base composition (frequency of the nucleotides A, C, G, and T) was further analyzed at the three codon positions in the concatenated alignment of PCGs for each species separately (S1 Table). The 12 variables measured for 64 taxa were summarized by a PCA, based on the first two components, which contributed 73.7% and 18.8% of the total variance, respectively (Fig 1). Most of the observed variation was related to the percentage of base composition in a third codon position (S6 Table). The first component (PC1) demonstrated a high positive correlation (0.98) with the percentage of C3 (percentage of cytosine in the third position) and a high negative correlation (-0.93) with T3 (thymine in the third position). The second component (PC2) positively correlated with G3 (0.67) and negatively correlated with A3 (-0.88). Most of the Arvicolinae formed a compact group on the PCA graph. Three species of genus Cricetulus, C.longicaudatus, C. migratorius and C. griseus clustered separately from the main group of Arvicolinae along the PC1. Both species demonstrated approximately 10% lower percentage of C3 and 10% higher percentage in T3 (S6 Table). The other outgroup taxa, Urocricetus kamensis and Akodon montensis grouped with Arvicolinae, and A. montensis showed similar base composition to Hyperacrius fertilis. Among Arvicolinae, the most dissimilar base composition was observed in Mynomes longicaudus Merriam, 1888, Chionomys gud Satunin, 1909 and Arvicola amphibius Linnaeus, 1758, showing higher than group average percentage of T3 and lower than group average percentage of C3. The mitochondrial genome of Ondatra zibethicus Linnaeus, 1766 was characterized by the highest percentage of adenine in the third position (45.7%) compared to other Arvicolinae. Ellobius lutescens Thomas, 1897 demonstrated the highest percentage of С3 (37.1%) among the complete PCG dataset (Fig 1). Since most of the variation has been observed in the third codon position (S6 Table), the complete dataset was subsequently used in the following phylogenetic reconstructions.
The frequency of the four bases (A, C, G, and T) at each codon position (first, second and third) in concatenated alignment was used as 12 variables for PCA. Tribes are indicated by colours. Species names are given only for taxa that differ in their nucleotide composition. Details with percentages for each position are given in S6 Table.
Substitution saturation analysis
Substitution saturation decreases phylogenetic information contained in the sequences and plagues the phylogenetic analysis involving deep branches.
According to the analysis implemented in DAMBE software (S7 Table), the observed Iss saturation index was significantly (P<0.0001) lower than critical Iss.c value for both symmetrical and asymmetrical topology tests indicating the lack of saturation in the studied Arvicolinae dataset.
The results of saturation plots for separate genes (S1 File) show the same pattern of negligible saturation. As a result for all 13 PCGs, no significant saturation for the 1st and 2nd codon position, and they are all suitable for the phylogenetic inference. CYTB, ND1 and ND6 show the same for 3rd codon position in contrast with ND2, ND3, ND4, ND4L. For other genes there was no significant saturation even for the 3rd codon position considering symmetrical topology for more than 32 OTUs (number of operational taxonomic units, S7 Table).
Time-calibrated mitochondrial genome phylogeny of Arvicolinae
The maximum-likelihood (ML) and Bayesian inference (BI) trees reconstructed using complete and RY-coded alignments of PCGs had similar topology (Fig 2, S1 File). Overall, ML analysis demonstrated lower node supports compared to BI analysis. In total 70% of the nodes were highly supported by ML and BI, with Bayesian probabilities BP>0.95 and ML bootstrap support BS>95 (Fig 2, nodes with a black dot).
Node labels display the following supports: BI complete / BI RY-coded 3rd codon position / ML complete / ML RY-coded 3rd codon position. Black circles show nodes with 0.95–1.0 BI and 95–100 ML support. All letters at nodes correspond to fossil constraints in S3 Table. Traditional tribal designations are also given above the branches and corresponding branches distinguished by different colors.
The monophyly of subfamily Arvicolinae was strongly supported by BI and ML analyses. However, several nodes, predominantly the internal nodes representing deeper phylogenies, which were highly supported by Bayesian analysis, did not receive high BS values. The divergence time between Arvicolinae and Cricetinae was estimated as Late Miocene, ca. 11.31 (9.48–13.3) / 10.7 (8.42–13.31) Ma, based on the complete and RY-coded alignments respectively (Table 1). When one of the eight calibration was discarded by turn, the results generally remained stable (S4 Table), with the exception of the analysis with excluded calibration of MRCA for the subfamily Arvicolinae (node“A”, Fig 2). This run yielded in younger dates for most of the nodes (S4 Table)).
The earliest radiation of the proper arvicolines (tribes Lemmini, Prometheomyini, Ondatrini and Dicrostonychini) dates back to the Late Miocene with mean at 7.36 (7.04–7.78) / 7.33 (7.05–7.73) Ma. Despite the high node support for the nodes marking Lemmini and Dicrostonychini, the basal part of the phylogenetic tree remains unresolved and represents a polytomy with several nodes not receiving significant BI and ML support. Analysis based on the PCG dataset where the third codon position was not masked with RY-coding, indicated significantly high Bayesian support for node C, combining Ondatrini and Dicrostonychini (Fig 2). The time to MRCA of node C is about 5.85 (5.15–6.56) / 5.54 (4.41–6.59) Ma and the MRCA of proper Dicrostonychini at 4.89 (4.08–5.7) / 4.49 (3.23–5.86) Ma.
The tribe Clethrionomyini representing second radiation of Arvicolinae received high BI and ML support, and nodes within the clade were also highly supported. The MRCA for Clethrionomyini dates back to 4.02 (3.33–4.72) / 4.46 (3.35–5.64) Ma.
The cluster containing tribes Ellobiusini, Lagurini, Arvicolini and genera Dinaromys, Arvicola Lacepede, 1799 and Hyperacrius Miller, 1896, i.e. the third radiation of Arvicolinae, was robustly supported by BI using both alignments and received reliable support by ML only with RY-masked alignment (Fig 2). Within this cluster, nodes marking tribes were highly supported by BI and ML. At the level of terminal branches within this cluster Dinaromys bogdanovi, Hyperacrius fertilis and Arvicola amphibius were the only to lack a certain phylogenetic position. D. bogdanovi grouped with Ellobiusini showing high BI support and no ML support. The water vole, Arvicola amphibius, clustered with Lagurini (high BP and no BS support) thus being paraphyletic to Arvicolini. The sagebrush vole Lemmiscus curtatus was sister to snow voles, Chionomys gud and C. nivalis Martins, 1842 with a robust support obtained in all analyses. The cluster of Chionomys Miller, 1908 and Lemmiscus was the earliest derivative in the highly supported group uniting all known vole genera of Arvicolini tribe. Arvicolini sensu stricto (excluding Arvicola) was fully resolved: both node H marking the whole tribe and all nodes within the tribe received robust support in ML and BI analyses.
The estimated time of the largest radiation event within the subfamily and TMRCA for the trichotomy Arvicolini—Ellobiusini—Lagurini dates back to 6.2 (5.65–6.76) / 6.11 (5.17–6.92) Ma. All following divergence events within this radiation according to the obtained estimates took place very close to each other, the 95% HPD of the diverging branches leading to MRCA of existing tribes are highly overlapping. Thus the date estimate to the MRCA of Ellobiusini was 4.97 (4.21–5.69) / 4.58 (3.42–5.68) Ma. The MRCA of the Lagurini tribe is around 3 Ma: 3.1 (2.59–3.75) / 3.05 (2.56–3.75) (Fig 2, Table 1). The estimate for the earliest split within the Arvicolini tribe radiation (not including Arvicola and Hyperacrius) with all recent genera is about 4.9 (4.33–5.47) / 5.02 (4.12–5.89) Ma, that coincides with the onset of Pliocene period, whereas the major part of recent genera, excluding early derivating Chionomys, Lemmiscus and Proedromys Thomas, 1911, according to obtained estimates appear either in the middle of the Pliocene or close to the boundary of Late Pliocene-Early Pleistocene (Fig 2, Table 1).
Gene trees
The topology of the Arvicolinae phylogeny varied between the 13 PCG trees (S1 File). While tribe level nodes received good support at most of the trees, the phylogenetic relationships between the taxa remained unresolved. The ATP8- and COX2-based trees lacked resolution at both deep and shallow nodes, and therefore, these trees resulted in a complete polytomy. The only node at the ATP8-based tree that retained its integrity with high support was the tribe Clethrionomyini. Noteworthy that this Clethrionomyini node had high support and was consistent at the majority of the gene trees, except for the ND3. The node containing the taxa of the Arvicolini tribe (excluding Arvicola amphibius) received high support on the COX1, ATP6, ND3, ND5, ND6, ND1 and CYTB gene trees. The ND4 gene tree yielded in a highly supported node grouping the semiaquatic species—Ondatra zibethicus and Arvicola amphibius, the result was not supported by any other gene tree and mitogenome BI and ML phylogenetic reconstructions (Fig 2). Positions of these two species, as well as Dinaromys bogdanovi were very unstable across the individual gene trees.
These taxa often occupied different positions and clustered with other species randomly. Remarkable that even in the case when tribal support and content was consistent across various trees and with the mitogenome tree, the interrelationships between tribes at the individual gene trees were unresolved. The lack of resolution especially at the deep nodes may be related to high saturation that is demonstrated with some genes, particularly ATP6, ATP8, ND1, ND2, ND3, ND5 and ND4 (maximum saturation, S7 Table), since phylogenetic signal disappears when divergence is over 10%.
Discussion
Our phylogenetic reconstruction of the subfamily Arvicolinae is based on a PCG dataset of mitochondrial genome sequences of 58 species of voles and lemmings with the outgroup of six hamsters. The dataset included 30 original sequences, and for 10 genera the mitogenomic sequencing was implemented for the first time. To date, this is the most comprehensive mitogenome dataset aimed at the revision of the Arvicolinae phylogeny considering almost all recent genera represented by nominal species. While the monophyletic origin of Arvicolinae has always been considered indisputable, previous attempts to resolve phylogenetic relationships within the subfamily using either morphological analysis or combinations of mitochondrial and nuclear markers yielded in several hard polytomies [16] or conflicting topologies [3,5–7,9,12,21,97,98]. The more taxa and more markers were considered in the analysis, the better resolution for the nodes marking major tribes within Arvicolinae has been obtained [7,20,21]. However, the diversification events within major radiation waves remained unresolved. The phylogenetic position of the genera Prometheomys Satunin, 1901, Arvicola, Ondatra Link, 1795 and Dinaromys in reconstructions performed with mitochondrial and nuclear markers was controversial [3,5,7,20,21] and genera Hyperacrius and Lemmiscus received little attention, their phylogenetic position was arguable.
Mismatches between the Bayesian and Maximum Likelihood support for the tribes and three waves of radiation within Arvicolini
The topology of the mitochondrial genome tree of Arvicolinae obtained in this study, in general, was in good agreement with previous large-scale phylogenetic reconstructions of the group based on mitochondrial and nuclear genes [7,20,21,97,98]. Using the concatenated alignment of 13 PCG, the present reconstruction resulted in high support for the nodes marking tribes in both Bayesian and Maximum Likelihood analyses. While Bayesian analysis also provided high posterior probability support for the basal nodes, ML approach failed to recover relationships and order of divergence between the basal branches. Previously, these deep divergences were identified as three waves of rapid radiations [7].
The first radiation within the subfamily is represented by four tribes—Lemmini, Prometheomyini, Dicrostonychini (including Phenacomys Merriam, 1889) and Ondatrini. The order of divergence and phylogenetic relationships between these ancient tribes remains unresolved using mitochondrial genome data. In particular, the absence of resolution is highlighted by the unstable position of the Prometheomys or Lemmini as the basal lineages. The basal position of Lemmini and Prometheomys with Ondatrini/Dicrostonichini both postdating them was also obtained earlier [20] on a very large set of Arvicolinae taxa analyzed using several mitochondrial and nuclear markers. It is important to underline that in each case this clustering has no support. In a number of other reconstructions [3,6,21,23] Prometheomys is the earliest split within the subfamily. To check if this inconsistency may be related to nucleotide composition bias we withdrew three species from the outgroup of the genus Cricetulus, (C.longicaudatus, C. migratorius and C. griseus) that demonstrated approximately 10% lower percentage of C3 and 10% higher percentage in T3 than most Arvicolinae and left in the outgroup only Urocricetus kamensis and Akodon montensis that showed similar base composition with species of Arvicolinae and carried out 4 variants Bayesian (on full and RY-coded alignment) and ML, respectively. Only Bayesian inference on full alignment resulted in the basal position of Prometeomys, but with minimal support (S3 Fig). With both BI and ML and RY masked alignment the same topology as on the full set of taxa (Fig 2), with basal Lemmini (but again without support) was obtained.
The second radiation is represented exclusively by the large monophyletic tribe Сlethrionomyini (Fig 2). These are predominantly forest-dwelling taxa that originated in Eurasia with only a few species penetrating North America during the Pleistocene. According to our data, the monophyly of Clethrionomyini was supported in analyses of either concatenated alignment or individual mitochondrial genes except for the short ND4L (S1 File). With all the nodes receiving high BP and ML support, the internal topology of branches within Clethrionomyini obtained in this study was similar to previous reconstructions of this tribe based on one mitochondrial and three nuclear loci [11].
The third radiation comprises three tribes Arvicolini, Ellobiusini and Lagurini. While these tribes, as well as most genera within the tribes, received strong support, and our reconstruction demonstrates that all the taxa of the third radiation share the same putative common ancestor, their interrelationships within this large clade also were not recovered, actually representing polytomy. According to our data, the genera Dinaromys, Hyperacrius and Lemmiscus whose assignment to certain tribes has previously been doubtful (Fig 2) also belong to the third radiation. Their taxonomic position, as well as the position of the genus Arvicola that suddenly appeared to be paraphyletic to other Arvicolini, are discussed below.
Phylogenetic relationships of the genus level taxa. Monotypic and low-diverse genera of uncertain position
The subfamily Arvicolinae includes several seriously understudied genera of unclear taxonomic position. For these genera, molecular data include either only mitochondrial CYTB sequences [5,16–18] or several additional mitochondrial and nuclear markers [3,7,9–11,21]. These genera are often the orphan genera, i.e. being represented by a single extant species. Considering such taxa is of remarkable importance for the reconstruction of high-level phylogenies, but their position on a tree can often be contradictory due to long-branch attraction [99,100]. While the resolving power of the phylogenetic reconstruction increases with the number of genes in analysis, several studies of rapid radiations based on organellar genomes pointed out the effect of long branch attraction [56,101,102] and references therein. Our study, among other, considers five genera of the unclear position either within the first (Prometheomys and Ondatra) or third (Dinaromys, Hyperacrius and Lemmiscus) radiation waves.
The Balkan vole, Dinaromys bogdanovi, endemic species from Balkan Peninsula was attributed to either Ondatrini [103] or Prometheomyini [104], but conventionally to Clethrionomyini [62,105,106]. Morphologically Dinaromys is mostly close to the extinct Pliocene genus Pliomys Méhely, 1914 [62,107–109], which distinguishes it from the rest of extant vole taxa. The genus Pliomys, in turn, has generally been considered the ancestral form for the whole Clethrionomyini tribe. That was the main reason [62] to distinguish a separate subtribe Pliomyi within the latter consisting of the two genera—extant Dinaromys and extinct Pliomys. Until recently, CYTB was the only studied locus for Dinaromys, and it was placed as sister to Prometheomys, another monotypic genus, and both were close to the Ellobiusini—Arvicolini—Lagurini group [5]. This grouping was strongly rejected by the following attempts to build molecular phylogeny of Arvicolinae showing the position of Prometheomys as the earliest derivative within the subfamily [3,7,9,21], and Dinaromys within the clade uniting Ellobiusini, Lagurini and Arvicola, i.e. the third radiation [9,20,21]. According to mitochondrial genome phylogeny (Fig 2), Dinaromys does not have putative MRCA with monophyletic Clethrionomyini tribe and most likely belongs to the third radiation, yet the certain position of this genus within this large group remains unclear.
The analysis of partial mitochondrial CYTB sequence [11] demonstrated that genus Hyperacrius does not seem to belong to the Clethrionomyini tribe. By analysing the set of mitogenomes of Clethrionomyini and Arvicolini it was recently suggested that Hyperacrius has the basal position within the tribe Arvicolini [74]. Here, using the broader taxonomic sampling, we confirm these previous findings showing that Hyperacrius predates the diversification of all main genera of Arvicolini.
Reconstructions performed using the individual mitochondrial genes often placed genus Ondatra as sister to Arvicola [5,16] or Clethrionomyini tribe [97] with low support. In all studies involving varying sets of nuclear genes Ondatra was among early diverging lineages [7] and sister to Neofiber True, 1884 if it was included in the analysis [21,97]. Such a position better corresponds to conventional taxonomy and paleontological data. Our results placed Ondatra sister to the Dicrostonychini tribe, hence with low support (except BI with transitions in the 3rd position included). Similar topology was observed by Lv et al. [20].
Lemmiscus curtatus—the sagebrush vole—is the only extant representative of the genus, it inhabits semi-arid prairies on the western coast of North America. For a long time, Lemmiscus was considered as closely related to the steppe voles Lagurini of the Old World and even as a subgenus within Lagurus Gloger, 1841 [62,110–112]. The close affinity between Lemmiscus and the Palearctic Lagurini was then seriously criticized from the paleontological perspective. Morphological similarities among the two groups were interpreted as a result of the parallel evolution at open, steppe-like landscapes, and Lemmiscus was proposed to be close to the tribe Arvicolini, particularly the genus Microtus Schrank, 1798 [113]. These data corroborated the previous grouping of Microtus and Lemmiscus in phylogenetic reconstruction based on restriction fragment LINE-1 [114], yet their taxonomic sampling did not include Lagurini and most genera of the Arvicolini tribe. In a recent reconstruction using mitochondrial CYTB and the only nuclear gene Lemmiscus clustered with Arvicola amphibius, yet with no support [21].
Using mitochondrial genomes to reconstruct Arvicolinae phylogeny, we sensationally show that Lemmiscus appears to be sister to the snow voles genus Chionomys. This clustering was obtained in all variants of the analysis, and node support values were significant. The snow voles unite three species occurring only in the Old World, particularly mountain systems of Southwestern, Central and Southeastern Europe and Southwestern Asia. Snow voles inhabit rocky patches of a subalpine and alpine belt from 500 up to 3500 m above the sea level [62,115]. Reliable pre-Pleistocene fossil remains of Chionomys are unknown, and the origin of the genus was previously attributed to the mid-Pleistocene [116]. Our data strongly contradicts this conventional view, and both Lemmiscus and Chionomys probably are more ancient taxa. Also, Lemmiscus and Chionomys occur at different continents and occupy contrasting ecological niches; they are also very dissimilar morphologically. These findings, broadly discussed below, are important for the understanding of the migration events of Arvicolinae from Eurasia to North America.
Phylogenetic relationships within the tribe Arvicolini sensu stricto
By using the mitochondrial genome data, we obtained good support for the nodes within the tribe Arvicolini except for the Arvicola amphibius that clustered with Lagurini (Fig 2). This evident artefact may be related to nucleotide composition bias (see Fig 1, S6 Table). Note that with the third position masked as RY this clustering has weak or no support. On the other hand the evidently wrongposition of A. amphibius can be a consequence of the long-branch attraction effect, and further studies should consider including sequences of the e.g. southern water wole, A. sapidus Miller, 1908 and nuclear genome data for better phylogenetic position resolution of the genus. The other nodes within Arvicolini, hereafter called Arvicolini sensu stricto marking the genus and subgenus level taxa, were recovered as monophyletic and clearly resolved.
The phylogenetic pattern indicates two major migration waves of voles to the Nearctic. The earliest derivative from the MRCA is a branch leading to Chionomys—Lemmiscus node and this gives clear indication on the first dispersal of common ancestors of the group from Palearctic to Nearctic. The only recent descendant of this lineage in North America is Lemmiscus.
The next split of ancestral lineage evidently took place in Asia and is represented by poorly diversified genus Proedromys and highly diversified cluster, uniting all the rest recent vole genera. This latter cluster further splits into highly supported clade of Asian voles showing sister relationships of genus Neodon Horsfield, 1841 and genera Alexandromys Ognev, 1914 and Lasiopodomys Lataste, 1887 and a cluster uniting two sister clades: Nearctic voles with following fast radiation resulting in nearly 20 recent species (here named Mynomes Rafinesque, 1817 after the earliest name of the generic group level), and a clade that further splits into Western Palearctic (Microtus s.str., Terricola Fatio, 1867 and Sumeriomys Argyropulo, 1933) and one containing taxa distributed in Central Asia (Blanfordimys Argyropulo, 1933), Westernmost Europe (Iberomys Chaline, 1972) and wide-ranged Agricola Blasius, 1857 (from Western Europe to Siberia). It is important to note that trees uniting Nearctic “Microtus’’ species in one cluster were obtained in various studies [12,18,19,21], but for the first time this cluster receives robust support, justifying the genus level status under the name Mynomes.
Another significant finding is the more clear assignment of Iberomys cabrerae and Agricola agrestis, both species conventionally assigned to Microtus [62], but whоse position at the molecular trees within the Arvicolini tribe was always uncertain. The tendency for clusterization of A. agrestis and Blanfordimys, though without support, was shown earlier [9,12,18,20,21]. In the paper where both I. cabrerae and A. agrestis were analysed in a comprehensive dataset with CYTB [18] these species appeared in different clusters: A. agrestis with Blanfordimys, while I. cabrerae within the cluster of Nearctic voles, however later on [10] in a detailed study of Asian voles came up with analogous to reported here clustering of I. cabrerae and A. agrestis with Blanfordimys. An important contribution was recently made by Barbosa et al. [34] who used a genomic approach for resolving phylogeny of speciose Microtus voles. According to their results both species appear to be monophyletic, however this study was based only on eight species and lacked most of the genera of the group. According to our results the cluster showing close relationships of these species with Blanfordimys is highly supported.
Divergence time estimates for the major Arvicolinae lineages in the context of fossil record
Dating the origin of Arvicolinae.
Our data estimated the time of radiation from the MRCA of all Arvicolinae as ca. 7.3 Ma (Table 1), i.e. the Late Miocene and divergence time of Cricetinae and Arvicolinae from common ancestors around 11 Ma. These dating estimates correspond with fossil records [117] and molecular dating obtained by previous studies [21]. Between the 11.1 and 7.75 Ma (from Early Vallesian to Late Turolian) in Eurasia appeared many taxa, conventionally referred to as microtoid cricetids. These forms were characterised by the arvicoline−like prismatic dental pattern with variously pronounced hypsodonty [117–120] and are generally considered as the ancestors for arvicolids [62,117–119,121]. There is no single opinion concerning the earliest known form assigned to the true Arvicolinae. Several authors attribute Pannonicola (= Ishymomys) sp. As the first fossil representative of Arvicolinae, dated as ca. 7.3 Ma from Middle Turolian, Hungary [122], and Asia [123]. Fejfar and coauthors [117] considered Pannonicola as the starting point of arvicolid evolution, though it is far from the generally accepted point of view. According to other, more widely distributed opinions, the first fossil Arvicolinae is Promimomys Kretzoi 1955 known from the Miocene-Pliocene boundary, dated around 5.3 Ma [14,124,125]. While calibrating the offset of Arvicolinae by Pannonicola at 7 Ma, we received the mean divergence date for the clade Arvicolinae as 7.36 Ma, which is unrealistically too close to FAD of Pannonicola. By excluding this calibration point from analysis, we obtained the 1 Ma younger estimate for the Arvicolinae node (S4 Table). However, during the consistent removal of one calibration point in turn, the estimate of the Arvicolinae origin between 7.36–7.43 remained stable.
Ancient radiation of Arvicolinae and the first migration event from Palearctic to Nearctic.
It is generally agreed that Arvicolinae originated from arvicolid-like cricetids that first appeared in Eurasia but not in North America [13,14,124]. Chronology of immigration events supported by fossil evidence was detaily considered earlier [14,126,127 and references therein]. Below we examine how divergence dating obtained in this study is consistent with fossil evidence.
The divergence between Dicrostonychini and Ondatrini took place ca. 6 Ma according to our data, indicating that the ancestors of this group were very closely related to the first arvicolines. Feifar et al. [117] pointed that the molar pattern of Pannonicola Kretzoi, 1965, according to the opinion of the authors, the oldest known fossil Arvicolinae, show similarity with Dolomys Nehring, 1898 the putative ancestor of Ondatrini, and possibly Dicrostonychini, indicating their closer relationships. Our data corroborate this grouping. The conventional idea of Promimomys as the first true arvicoline taxon indicates the first immigration timing around 5.5 Ma (Martin, 2010, 2015) as there is no record of this taxon in the older sediments.
From the mitochondrial genome data, the date for the MRCA of Lemmini was estimated as 4.81–4.37 Ma. This molecular dating consider ancestors of Lemmini almost a million years older than the fossil remains reliably attributed to Lemmini in Europe [128,129] and Asia [130] dated as the Early Villaniyan (Mammal Neogene zone MN16, 3.2 Ma), while North American fossils of Lemmini were dated at ca. 3.9 Ma or Late Blancan according to Ruez & Gensler [131].
However, these lemming fossils are characterized by very advanced unrooted teeth and masticatory patterns, close to the recent forms of lemmings. Among the potential ancestors here can be mentioned Tobenia kretzoi Fejfar, Repenning, 1998, a species with rooted molars known from the early Pliocene of Wölfersheim, Germany [132]. This record refers to MN15 that is ca. 4 Ma, close to our molecular dating.
The molecular estimate of the two major radiation waves of Arvicolinae, leading to the tribes Clethrionomyini, Arvicolini, Ellobiusini and Lagurini, dates back to 6 Ma, Late Miocene, MN 13 (Late Turolian). These findings correspond to paleontological data and confirm the estimates received previously in the study based on nuclear genes [7]. Slightly later, firstly in North America [14] and then in Europe [117,125] and Western Asia [133], the most plesiomorphic arvicoline genus Promimomys Kretzoi, 1955 appeared. This form is considered ancestral to numerous species emerged in the Early Pliocene and conventionally assigned to highly mixed and species-rich genus Mimomys Forsyth-Major, 1902. According to the generally accepted view, different forms of this complex “Mimomys’’ group represented the starting point for all subsequent lineages of Arvicolines. The concept of common ancestry for these forms within this geological period does not contradict the data obtained in the present study and hypothesis proposed by paleontologists [117]. The radiation of common ancestors for all Clethrionomyini species starts later, since Late Ruscinian (MN 15), around 4 Ma.
The origin of the tribe Arvicolini sensu stricto: Second trans-Beringian dispersal.
The molecular estimate for the MRCA of node H, Arvicolini s.str. (Fig 2) is ca. 5 Ma. (4.33–5.47). Considering that most primitive forms of the genus Mimomys are among the MRCA candidates for all main genera within the tribe Arvicolini s.str, the obtained time estimate i.e. the very beginning of the Pliocene, may be considered as consistent with the fossil record. However, there is a certain probability of the overestimation of the age of this node during the divergence dating analysis.
One of the earliest records of Mimomys in North America was dated as 4.75 Ma [134], while fossil remains from Asia are slightly older [135]. This is the time of the second dispersal of arvicolids from Asia to the Nearctic. The only recent descendant of these immigrants in North America is Lemmiscus curtatus. According to our data, the starting point of evolutionary history for this lineage is around 4 Ma. The earliest remains assigned to the genus are known from the end of Early Pleistocene from the SAM Cave in New Mexico [136] in the sediments according to paleomagnetic and faunistic data that may be dated as 1.8 Ma. Repenning [136] deduced Lemmiscus from plesiomorphic Allophaiomys Kormos, 1932. Remains assigned to the latter taxon are widely distributed among the Early Pleistocene sites dated between 2.2 and 1.6 Ma in both the Palearctic and Nearctic. Yet, Allophaiomys is a rather lumped taxon [137] presumably accepted as ancestral to most Microtus species and associated genera. Tesakov and Kolfschoten [113] suggested the hypothesis of a Mimomys–Lemmiscus phyletic lineage. However, their hypothesis also presumes that ancestral Mimomys (Cromeromys Zazhigin, 1980), a form having rooted molars and inhabiting vast areas from Western Europe to Beringia dispersed southwards across North America in the late Early Pleistocene and evolved there into rootless Lemmiscus. Thus, our dating conflicts with both views and supports the idea of dispersal and further evolution from “Mimomys” stage [113,136] in the middle of the Pliocene, ca. 4 Ma. The fossil remains of Chionomys are known only from the Pleistocene sediments [116]. According to our dating based on mtDNA, the diversification of ancestral lineage may have started in Western Palearctic as early as in the Middle Pliocene, i.e. either late Zanclean or early Piacenzian periods.
Diversification within Arvicolini sensu stricto: Late Pliocene exchange between Palearctic and Nearctic faunas.
The other genera within Arvicolini were monophyletic according to Bayesian and ML analyses with high node support. According to conventional view, this group originated from Allophaiomys [62,136], a highly complex taxon common in the Early Pleistocene (ca. 2 Ma) faunas of the Nearctic and Palearctic. Our results on divergence dating raise another hypothesis on the starting point for the group taking place at the level of “Mimomys’’ stage, i.e. in the Pliocene. The genus Proedromys is the first derivative from this common stem, most likely in the middle of the Pliocene (approx. 4 Ma). The standalone position of this genus among other genera of Arvicolini that plausibly derived from Allophaiomys has been earlier underlined by Gromov, Polyakov [62] and Repenning [136].
A further split within Arvicolini took place in the late Pliocene and resulted in the entirely Asian lineage which is currently represented by genera Neodon, Alexandromys and Lasiopodomys. The other, sister lineage emerged in the Early Pliocene, around 3.8–4 Ma, also from the pre-Allophaiomys stage and diverged into two branches. Ancestors of the first branch (Mynomes) penetrated the Nearctic during the third Nearctic immigration event, where they diversified into 20 species. Descendants of the other, Palearctic branch, further produced two lineages. Among them, the first apparently originated in Central Asia and dispersed westwards during the Late Pliocene. Iberomys cabrerae, inhabiting the Iberian Peninsula and the foothills of Pyrenees is a relict descendant of this line [138] and references therein. Our mitogenomic data supports the hypothesis of long independent evolution of Iberomys cabrerae lineage previously confirmed by several unique morphological, biological and ecological traits [138]. The first fossil remains of Iberomys were found in the Early Pleistocene [139] sediments in Spain predating the Jaramillo reversal event (approx. 1.2 Ma). According to the scenario set by Cuenca-Bescós et al. [138], the genus Iberomys is a basal sister group of Arvicolini evolved since its origins in the Early Pleistocene in the western Mediterranean region. The most probable origin and vicariant speciation according to these authors was linked to the stock of Allophaiomys taxa with plesiomorphic character states. The results reported here, in a whole are in a good agreement with this scenario although indicate an earlier time of origin and speciation from the stock predating Allophaiomys-stage and going far back to the Mimomys stage of the Late Pliocene.
The other descendants of the same stock in the recent fauna are represented by Agricola and Blanfordimys. The range of Agricola covers whole Europe and stretches to the east up to the Lake Baikal and watershed between the Yenisey and Lena Rivers [62,63]. Recent studies showed that Agricola is represented by three highly divergent lineages, possibly a species level taxa [140]. Three species of Blanfordimys occur in high mountain forests and steppes of Central Asia and are characterized by a very primitive molar pattern, similar to Allophaiomys. The idea that Agricola and Iberomys represent relicts of a very early colonization of Arvicolini to Western Europe was earlier suggested by Martinkova and Moravec [19] and well agrees with the given data.
The second lineage of Palearctic branch has evolved in Western Palearctic and in the modern fauna is represented by species-rich genus Microtus (with subgenera Microtus s.str and Sumeriomys) and Terricola (around 14 species, mainly found in South and Southwestern Europe). The divergence between these lineages corresponds to the Late Pliocene, however the speciation events coincided with Early Pleistocene for genus Terricola and early Middle Pleistocene for subgenera Sumeriomys and Microtus. The latter dating matches with known fossil records [109,141].
Summing up the comparison between the molecular estimates of divergence times reported here and known paleontological data, it is curious to note that while the dates for MRCA for most genera within Arvicolini s.str. are significantly older than was previously supposed [109,117,142], dating of speciation events within the genera (Lasiopodomys, Alexandromys, Terricola) are consistent with fossil records [109,134,136,141–143] and others.
Although our results are consistent both with fossil record and previous molecular dating results it is worth to mention, that the lack of calibration points for terminal nodes in our analysis may be the source of underestimation of basal node ages, while ages of several youngest nodes in Arvicolini could be overestimated.
Systematic remarks
While systematic relationships of higher taxa within Arvicolinae undoubtedly require further studies involving genomic approaches, some amendments to the current taxonomic system could be made already at this step of the research. Our study provided significant input for the potential review of taxonomic structure and composition of the tribe Arvicolini (S2 File). Our data shows that the position of the genus Arvicola is still unresolved. On the contrary, genera Lemmiscus and Hyperacrius certainly should be considered as members of the tribe Arvicolini. The further grouping of species into genera and subgenera within this highly diverse tribe always was very subjective and debatable. Most arguable was the composition of the genus Microtus. The current system [62] where Blanfordimys, Neodon and Lasiopodomys have generic status, while Alexandromys, Stenocranius and Terricola are referred as subgenera within the genus Microtus is strongly outdated and contradicts the data of recent phylogenetic studies. The last checklist [2] partly modified this scheme and following Abramson and Lissovsky [63] elevated Alexandromys to full genus and Stenocranius considered as subgenus of Lasiopodomys, and Neodon as a genus. However, despite the accumulated evidence from several previous papers [9,19,144], in this reference book without substantiation the status of Blanfordimys was downgraded [145] while three well differentiated lineages (Blanfordimys, all Nearctic microtines, Terricola, Microtus and Sumeriomys) were illogically united in one genus Microtus. These well-differentiated lineages together form the sister branch to one with similar branching pattern and recognized three genera: Alexandromys, Lasiopodomys and Neodon. It is widely known, that the better is the phylogenetic resolution of any species-rich group the more complicated it matches the conventional hierarchical categories of Linnean system. Trying to keep as much stability of nomenclature retaining the already commonly used names that correspond to certain lineages from one hand and to reflect robust phylogenetic nodes in a formal classification from the other, here we suggest the following system of generic group taxa within the tribe Arvicolini sensu stricto which is given in S2 File.
Conclusions
Our phylogenetic analysis based on a complete mitochondrial genomes confirmed the monophyly of the subfamily, monophyly of the most tribes originated during the three subsequent radiation events. While order of divergence between ancient genera belonging to the first radiation were not uniformly supported by Bayesian and Maximum Likelihood analyses, our study reports the high node statistical support for the groups of genera within the highly diverse tribe Arvicolini. Mitogenome phylogeny resolved several previously reported polytomies and also revealed unexpected relationships between taxa. The robust placement of Lemmiscus as sister to the snow voles, Chionomys in the tribe Arvicolini, in contrast with a long-held belief of its affinity with Lagurini, is an essential novelty of our phylogenetic analyses. Our results resolve some of the ambiguous issues in phylogeny of Arvicolinae, but some phylogenetic relationships require further genomic studies, e.g. the evaluation of the precise positions of Arvicola, Dinaromys and Hyperacrius.
Here, we provide the evidence of high informativeness of the mitogenomic data for phylogenetic reconstruction and divergence time estimation within Arvicolinae, and suggest that mitogenomes can be highly informative, when the number of extant and extinct forms are comparable (the case of Arvicolini) and insufficient when extant forms represent single lineages of once rich taxon (most cases of early radiation in the subfamily).
The accuracy and precision of previous divergence time estimates derived from multigene nDNA and nDNA–mtDNA datasets are here refined and improved. The estimates for subfamily origin and early divergence are consistent with fossil record, however mtDNA estimates for putative ancestors of most genera within Arvicolini appeared to be much older than it was supposed from paleontological studies.
Supporting information
S1 Fig. Mapping of raw Illumina reads of Lemmiscus curtatus against the reference mitochondrial genome of Mynomes ochrogaster using Geneious Prime.
https://doi.org/10.1371/journal.pone.0248198.s001
(PDF)
S2 Fig.
Nucleotide misincorporations at 5’-termini (A) and 3’-termini (B) of the Lemmiscus curtatus calculated using mapDamage. All possible misincorporations are plotted in gray, except for guanine to adenine (G>A, blue lines) and cytosine to thymine (C>T, red lines).
https://doi.org/10.1371/journal.pone.0248198.s002
(PDF)
S3 Fig. A molecular phylogeny for the subfamily Arvicolinae reconstructed using the complete PCG dataset and Bayesian algorithm without three Cricetidae species.
Major Arvicolinae tribes are indicated by color coding (Arvicolini—light blue, Lagurini—blue, Ellobiusini—purple, Clethrionomyini—magenta, Dicrostonychini—dark green, Ondatrini—light green, Prometheomyini—yellow, Lemmini—red, nomen nudum species—black). Black circles show nodes with 0.95–1.0 BI and 95–100 ML support.
https://doi.org/10.1371/journal.pone.0248198.s003
(JPG)
S1 Table. Material, GenBank accession numbers and mitogenome characteristics.
https://doi.org/10.1371/journal.pone.0248198.s004
(XLSX)
S2 Table. Optimal partitioning schemes of the concatenated full and RY-coding alignments for mitochondrial protein-coding genes obtained using PartitionFinder for the subsequent Bayesian phylogenetic analysis.
https://doi.org/10.1371/journal.pone.0248198.s005
(XLSX)
S3 Table. Fossil calibrations used in dating analysis.
Node labels correspond to Fig 2, ages are given in million years ago (Ma). FAD—age of the first appearance of the taxon in the fossil record (first appearance date).
https://doi.org/10.1371/journal.pone.0248198.s006
(XLSX)
S4 Table. Additional divergence time estimates for the major lineages within the subfamily Arvicolinae.
Complete—the main variant of divergence dating analysis (Table 1); A-H—results of the analyses with alternate exclusion of calibrations (labels show what calibration was excluded).
https://doi.org/10.1371/journal.pone.0248198.s007
(XLSX)
S5 Table. Completeness of analyzed mitogenomes.
Mitogenomes obtained in the current study are marked in bold, partial genes colored by yellow. Absent genes indicated with red color. The length of protein-coding genes is given in nucleotides.
https://doi.org/10.1371/journal.pone.0248198.s008
(XLSX)
S6 Table. Source table for the PCA analysis.
Each codon position is indicated by numbers. For nucleotides, percentages are given.
https://doi.org/10.1371/journal.pone.0248198.s009
(XLSX)
S7 Table. Test of substitution saturation.
Analysis performed on all sites for 1&2nd and 3rd codon position separately. Iss—index of substitution saturationIssSym is Iss.c assuming a symmetrical topology, IssAsym is Iss.c assuming an asymmetrical topology, NumOTU—number of operation taxonomic units. Red color indicates P-value < 0.05.
https://doi.org/10.1371/journal.pone.0248198.s010
(XLSX)
S1 File. A molecular phylogeny for the subfamily Arvicolinae reconstructed using the complete PCG dataset and each of the 13 individual genes with the corresponding saturation plots indicated.
Major Arvicolinae tribes are indicated by color coding (Arvicolini—light blue, Lagurini—blue, Ellobiusini—purple, Clethrionomyini—magenta, Dicrostonychini—dark green, Ondatrini—light green, Prometheomyini—yellow, Lemmini—red, nomen nudum species—black). Bayesian topology was used to plot the tree for the complete 13 PCGs dataset. Node labels display the following supports: BI complete / BI RY-coded 3rd codon position / ML complete / ML RY-coded 3rd codon position. Black circles show nodes with 0.95–1.0 BI and 95–100 ML support. For each of the PCGs, maximum likelihood topology is given, node labels display ultrafast ML bootstrap above 50%. Saturation plots are indicated on the side insets, where colors mark the following partitions: 1st codon position transitions (ts)—brown, 1st transversions (tv)—red, 2nd ts—blue, 2nd tv—green, 3rd ts—pink, 3rd tv—black.
https://doi.org/10.1371/journal.pone.0248198.s011
(PDF)
S2 File. The proposed system of generic group taxa within the tribe Arvicolini sensu stricto.
https://doi.org/10.1371/journal.pone.0248198.s012
(DOCX)
Acknowledgments
We are grateful to colleagues who shared the material for the study—Abramov A.V., Golenishchev F.N., Stekolnikov A.A., Bannikova A.A., Chabovskiy, A.V., Kowalskaya Yu.M., Smorkatcheva A.V., Dokuchaev N.E., Bogdanov A.S., Grafodatsky A.S., Buzan E. We would like to thank Margarita Ezhova and Maria Logacheva from the Genomics Core Facility of Skolkovo Institute of Science and Technology for NGS library preparation. Special thanks to Dr. Rudolf Haslauer, Zoological Society for the Conservation of Species and Populations (ZGAP), Poernbach, Bavaria, for fruitful discussion on taxonomy issues. Two anonymous reviewers provided valuable comments and raised important issues that helped to improve the manuscript. We are very grateful to R.A. Martin, Department of Biological Sciences Murray State University, for providing updated literature on fossil Arvicolinae and critical comments.
References
- 1.
Musser GG, Carleton MD. Superfamily Muroidea. Mammal species of the world: a taxonomic and geographic reference. JHU Press; 2005. pp. 894–1531.
- 2. Pardiñas UFJ, Myers P, León-Paniagua L, Ordóñez Garza N, Cook JA, Kryštufek B, et al. Family Cricetidae (true hamsters, voles, lemmings and new world rats and mice). Handb Mamm World. 2017;7: 204–279.
- 3. Galewski T, Tilak M, Sanchez S, Chevret P, Paradis E, Douzery EJ. The evolutionary radiation of Arvicolinae rodents (voles and lemmings): relative contribution of nuclear and mitochondrial DNA phylogenies. BMC Evol Biol. 2006;6: 80. pmid:17029633
- 4. Lebedev VS, Bannikova AA, Tesakov AS, Abramson NI. Molecular phylogeny of the genus Alticola (Cricetidae, Rodentia) as inferred from the sequence of the cytochrome b gene. Zool Scr. 2007;36: 547–563.
- 5. Buzan EV, Krystufek B, Hänfling B, Hutchinson WF. Mitochondrial phylogeny of Arvicolinae using comprehensive taxonomic sampling yields new insights. Biol J Linn Soc. 2008;94: 825–835.
- 6. Robovský J, ŘIčánková V, Zrzavý J. Phylogeny of Arvicolinae (Mammalia, Cricetidae): utility of morphological and molecular data sets in a recently radiating clade. Zool Scr. 2008;37: 571–590.
- 7. Abramson NI, Lebedev VS, Tesakov AS, Bannikova AA. Supraspecies relationships in the subfamily Arvicolinae (Rodentia, Cricetidae): An unexpected result of nuclear gene analysis. Mol Biol. 2009;43: 834.
- 8. Bannikova AA, Lebedev VS, Lissovsky AA, Matrosova V, Abramson NI, Obolenskaya EV, et al. Molecular phylogeny and evolution of the Asian lineage of vole genus Microtus (Rodentia: Arvicolinae) inferred from mitochondrial cytochrome b sequence. Biol J Linn Soc. 2010;99: 595–613.
- 9. Fabre P-H, Hautier L, Dimitrov D, P Douzery EJ. A glimpse on the pattern of rodent diversification: a phylogenetic approach. BMC Evol Biol. 2012;12: 88. pmid:22697210
- 10. Liu SY, Sun ZY, Liu Y, Wang H, Guo P, Murphy RW. A new vole from Xizang, China and the molecular phylogeny of the genus Neodon (Cricetidae: Arvicolinae). Zootaxa. 2012;3235: 1–22.
- 11. Kohli BA, Speer KA, Kilpatrick CW, Batsaikhan N, Damdinbaza D, Cook JA. Multilocus systematics and non-punctuated evolution of Holarctic Myodini (Rodentia: Arvicolinae). Mol Phylogenet Evol. 2014;76: 18–29. pmid:24594062
- 12. Pradhan N, Sharma AN, Sherchan AM, Chhetri S, Shrestha P, Kilpatrick CW. Further assessment of the Genus Neodon and the description of a new species from Nepal. PLOS ONE. 2019;14: e0219157. pmid:31314770
- 13. Fejfar O, Heinrich W-D, Kordos L, Maul L. Microtoid cricetids and the Early history of arvicolids (Mammalia, Rodentia). Palaeontol Electron. 2011;14.
- 14. Martin RA. THE NORTH AMERICAN PROMIMOMYSIMMIGRATION EVENT. Paludicola. 2010; 8(1):14–21.
- 15. Kryštufek B, Tesakov AS, Lebedev VS, Bannikova AA, Abramson NI, Shenbrot G. Back to the future: the proper name for red-backed voles is Clethrionomys Tilesius and not Myodes Pallas. Mammalia. 2020; 84(2):214–217.
- 16. Conroy CJ, Cook JA. MtDNA Evidence for Repeated Pulses of Speciation Within Arvicoline and Murid Rodents. J Mamm Evol. 1999;6: 221–245. pmid:1020561623890.
- 17. Conroy CJ, Cook JA. Molecular Systematics of a Holarctic Rodent (Microtus: Muridae). J Mammal. 2000;81: 344–359.
- 18. Jaarola M, Martínková N, Gündüz İ, Brunhoff C, Zima J, Nadachowski A, et al. Molecular phylogeny of the speciose vole genus Microtus (Arvicolinae, Rodentia) inferred from mitochondrial DNA sequences. Mol Phylogenet Evol. 2004;33: 647–663. pmid:15522793
- 19. Martínková N, Moravec J. Multilocus phylogeny of arvicoline voles (Arvicolini, Rodentia) shows small tree terrace size. J Vertebr Biol. 2012;61: 254–267.
- 20. Lv X, Xia L, Ge D, Wu Y, Yang Q. Climatic niche conservatism and ecological opportunity in the explosive radiation of arvicoline rodents (Arvicolinae, Cricetidae). Evolution. 2016;70: 1094–1104. pmid:27061935
- 21. Steppan SJ, Schenk JJ. Muroid rodent phylogenetics: 900-species tree reveals increasing diversification rates. PLOS ONE. 2017;12: e0183070. pmid:28813483
- 22. Bendová K, Marková S, Searle JB, Kotlík P. The complete mitochondrial genome of the bank vole Clethrionomys glareolus (Rodentia: Arvicolinae). Mitochondrial DNA Part A. 2016;27: 111–112. pmid:24438307
- 23. İbiş O, Selçuk AY, Sacks BN, Yıldız B, Özcan S, Kefelioğlu H, Tez C. Whole mitochondrial genome of long-clawed mole vole (Prometheomys schaposchnikowi) from Turkey, with its phylogenetic relationships. Genomics. 2020; 112(5): 3247–3255. pmid:32512144
- 24. Zhao H, Qi X, Li C. Complete mitochondrial genome of the muskrat (Ondatra zibethicus) and its unique phylogenetic position estimated in Cricetidae. Mitochondrial DNA Part B. 2018;3: 296–298. pmid:33474150
- 25. Folkertsma R, Westbury MV, Eccard JA, Hofreiter M. The complete mitochondrial genome of the common vole, Microtus arvalis (Rodentia: Arvicolinae). Mitochondrial DNA Part B. 2018;3: 446–447. pmid:33474199
- 26. Kim HR, Park YC. The complete mitochondrial genome of the Korean red-backed vole, Myodes regulus (Rodentia, Murinae) from Korea. Mitochondrial DNA. 2012;23: 148–150. pmid:22409761
- 27. Yang C, Hao H, Liu S, Liu Y, Yue B, Zhang X. Complete mitochondrial genome of the Chinese oriental vole Eothenomys chinensis (Rodentia: Arvicolinae). Mitochondrial DNA. 2012;23: 131–133. pmid:22397384
- 28. Cao W, Xia Y, Dang X, Xu Q. The first complete mitochondrial genome of the Microtus ochrogaster. Mitochondrial DNA Part A. 2016;27: 3682–3683. pmid:26305486
- 29. Chen S, Chen G, Wei H, Wang Q. Complete mitochondrial genome of the Père David’s Vole, Eothenomys melanogaster (Rodentia: Arvicolinae). Mitochondrial DNA Part DNA Mapp Seq Anal. 2016;27: 2496–2497. pmid:26024146
- 30. Cong H, Kong L, Liu Z, Wu Y, Motokawa M, Harada M, et al. Complete mitochondrial genome of the mandarin vole Lasiopodomys mandarinus (Rodentia: Cricetidae). Mitochondrial DNA Part A. 2016;27: 760–761. pmid:24845440
- 31. Fedorov VB, Goropashnaya AV. Complete mitochondrial genome of the Eurasian collared lemming Dicrostonyx torquatus Pallas, 1779 (Rodentia: Arvicolinae). Mitochondrial DNA Part B. 2016;1: 824–825. pmid:28670624
- 32. Fedorov VB, Goropashnaya AV. Complete mitochondrial genomes of the North American collared lemmings Dicrostonyx groenlandicus Traill, 1823 and Dicrostonyx hudsonius Pallas, 1778 (Rodentia: arvicolinae). Mitochondrial DNA Part B. 2016;1: 878–879. pmid:28642935
- 33. Yu P, Kong L, Li Y, Cong H, Li Y. Analysis of complete mitochondrial genome and its application to phylogeny of Caryomys inez (Rodentia: Cricetidae: Arvicolinae). Mitochondrial DNA Part B. 2016;1: 343–344. pmid:33644377
- 34. Barbosa S, Paupério J, Pavlova SV, Alves PC, Searle JB. The Microtus voles: Resolving the phylogeny of one of the most speciose mammalian genera using genomics. Mol Phylogenet Evol. 2018;125: 85–92. pmid:29574272
- 35. Jiang J-Q, Wu S-X, Chen J-J, Liu C-Z. Characterization of the complete mitochondrial genome of short-tailed field vole, Microtus agrestis. Mitochondrial DNA Part B. 2018;3: 845–846. pmid:33474341
- 36. Li J-Q, Li L, Fu B-Q, Yan H-B, Jia W-Z. Complete mitochondrial genomes confirm the generic placement of the plateau vole, Neodon fuscus. Biosci Rep. 2019;39. pmid:31262975
- 37. Mu Y, Duan Y, Di Z, Wang Z, Wanlong Z. The complete mitochondrial genome of the Yunnan red-backed vole Eothenomys miletus (Rodentia: Cricetidae) and its phylogeny. Mitochondrial DNA Part B. 2019;4: 1424–1425.
- 38. San Mauro D, Gower DJ, Zardoya R, Wilkinson M. A hotspot of gene order rearrangement by tandem duplication and random loss in the vertebrate mitochondrial genome. Molecular Biology and Evolution. 2006; 23(1): 227–234. pmid:16177229
- 39. Brown WM, George M, Wilson AC. Rapid evolution of animal mitochondrial DNA. Proceedings of the National Academy of Sciences. 1979; 76(4): 1967–1971. pmid:109836
- 40. De Mandal S, Chhakchhuak L, Gurusubramanian G, Kumar NS. Mitochondrial markers for identification and phylogenetic studies in insects–A Review. DNA Barcodes. 2014; 2(1): 1–9.
- 41. Hassanin A, Veron G, Ropiquet A, Jansen van Vuuren B, Lécu A, Goodman SM, et al. Evolutionary history of Carnivora (Mammalia, Laurasiatheria) inferred from mitochondrial genomes. PloS one. 2021; 16(2): e0240770. pmid:33591975
- 42. Govindarajulu R, Parks M, Tennessen JA, Liston A, Ashman TL. Comparison of nuclear, plastid, and mitochondrial phylogenies and the origin of wild octoploid strawberry species. American Journal of Botany. 2015; 102(4): 544–554. pmid:25878088
- 43. Abramson NI, Rodchenkova EN, Kostygov AY. Genetic variation and phylogeography of the bank vole (Clethrionomys glareolus, Arvicolinae, Rodentia) in Russia with special reference to the introgression of the mtDNA of a closely related species, red-backed vole (Cl. rutilus). Russian Journal of Genetics. 2009; 45(5): 533–545.
- 44. Bodrov SY, Vasiljeva VK, Okhlopkov IM, Mamayev NV, Zakharov ES, Oleinikov AY, et al. Evolutionary history of mountain voles of the subgenus Aschizomys (Cricetidae, Rodentia), inferred from mitochondrial and nuclear markers. Integrative zoology. 2020; 15(3): 187–201. pmid:31631516
- 45. Shen X, Pu Z, Chen X, Murphy RW, Shen Y. Convergent evolution of mitochondrial genes in deep-sea fishes. Frontiers in genetics. 2019; 10: 925. pmid:31632444
- 46. Strohm JH, Gwiazdowski RA, Hanner R. Fast fish face fewer mitochondrial mutations: Patterns of dN/dS across fish mitogenomes. Gene. 2015; 572(1): 27–34. pmid:26149654
- 47. Duchene S, Frey A, Alfaro-Núñez A, Dutton PH, Gilbert MTP, Morin PA. Marine turtle mitogenome phylogenetics and evolution. Molecular phylogenetics and evolution. 2012; 65(1): 241–250. pmid:22750111
- 48. Ramos EKS, Freitas L, Nery MF. The role of selection in the evolution of marine turtles mitogenomes. Sci Rep. 2020; 10: 16953. pmid:33046778
- 49. da Fonseca RR, Johnson WE, O’Brien SJ, Maria JR, Agostinho A. The adaptive evolution of the mammalian mitochondrial genome. BMC Genomics. 2008; 9: 119. pmid:18318906
- 50. Tomasco IH, Lessa EP. The evolution of mitochondrial genomes in subterranean caviomorph rodents: adaptation against a background of purifying selection. Molecular Phylogenetics and Evolution. 2011; 61(1): 64–70. pmid:21723951
- 51. Luo Y, Gao W, Gao Y, Tang S, Huang Q, Tan X, et al. Mitochondrial genome analysis of Ochotona curzoniae and implication of cytochrome c oxidase in hypoxic adaptation. Mitochondrion. 2008; 8(5–6): 352–357. pmid:18722554
- 52. Shen YY, Liang L, Zhu ZH, Zhou WP, Irwin DM, Zhang YP. Adaptive evolution of energy metabolism genes and the origin of flight in bats. Proceedings of the National Academy of Sciences. 2010; 107(19): 8666–8671. pmid:20421465
- 53. Blier PU, Dufresne F, Burton RS. Natural selection and the evolution of mtDNA-encoded peptides: evidence for intergenomic co-adaptation. TRENDS in Genetics. 2001; 17(7): 400–406. pmid:11418221
- 54. Mishmar D, Ruiz-Pesini E, Golik P, Macaulay V, Clark AG, Hosseini S, et al. Natural selection shaped regional mtDNA variation in humans. Proceedings of the National Academy of Sciences. 2003; 100(1): 171–176. pmid:12509511
- 55. Vilstrup JT, Ho SY, Foote AD, Morin PA, Kreb D, Krützen M, et al. Mitogenomic phylogenetic analyses of the Delphinidae with an emphasis on the Globicephalinae. BMC Evol Biol. 2011;11: 65. pmid:21392378
- 56. Li T, Hua J, Wright AM, Cui Y, Xie Q, Bu W, et al. Long-branch attraction and the phylogeny of true water bugs (Hemiptera: Nepomorpha) as estimated from mitochondrial genomes. BMC Evol Biol. 2014;14: 99. pmid:24884699
- 57. Lee Y, Kwak H, Shin J, Kim S-C, Kim T, Park J-K. A mitochondrial genome phylogeny of Mytilidae (Bivalvia: Mytilida). Mol Phylogenet Evol. 2019;139: 106533. pmid:31185299
- 58. Hassanin A, Bonillo C, Tshikung D, Shongo CP, Pourrut X, Kadjo B, et al. Phylogeny of African fruit bats (Chiroptera, Pteropodidae) based on complete mitochondrial genomes. J Zool Syst Evol Res. 2020;58: 1395–1410.
- 59. da Silva FS, Cruz ACR, de Almeida Medeiros DB, da Silva SP, Nunes MRT, Martins LC, et al. Mitochondrial genome sequencing and phylogeny of Haemagogus albomaculatus, Haemagogus leucocelaenus, Haemagogus spegazzinii, and Haemagogus tropicalis (Diptera: Culicidae). Sci Rep. 2020;10: 16948. pmid:33046768
- 60. Nicolas V, Fabre PH, Bryja J, Denys C, Verheyen E, Missoup AD, et al. The phylogeny of the African wood mice (Muridae, Hylomyscus) based on complete mitochondrial genomes and five nuclear genes reveals their evolutionary history and undescribed diversity. Molecular phylogenetics and evolution. 2020; 144: 106703. pmid:31816395
- 61. Zou Y, Xu M, Ren S, Liang N, Han C, Nan X, et al. Taxonomy and phylogenetic relationship of zokors. Journal of Genetics. 2020; 99(1): 1–10. pmid:32482927
- 62.
Gromov IM, Polyakov IY. Voles (Microtinae). Brill; 1992.
- 63.
Abramson NI, Lissovsky AA. Subfamily Arvicolinae. In: Pavlinov IY, Lissovsky AA, editors. The mammals of Russia: A taxonomic and geographic reference. KMK Scientific Press; 2012. pp. 127–141.
- 64.
Köchl S, Niederstätter H, Parson W. DNA Extraction and Quantitation of Forensic Samples Using the Phenol-Chloroform Method and Real-Time PCR. In: Carracedo A, editor. Forensic DNA Typing Protocols. Totowa, NJ: Humana Press; 2005. pp. 13–29. https://doi.org/10.1385/1-59259-867-6:013
- 65.
Effects of Index Misassignment on Multiplexing and Downstream Analysis. Illumina. Available: https://www.illumina.com/content/dam/illumina-marketing/documents/products/whitepapers/index-hopping-white-paper-770-2017-004.pdf.
- 66.
Andrews S. FastQC: A Quality Control Tool for High Throughput Sequence Data. 2016. Available: http://www.bioinformatics.babraham.ac.uk/projects/fastqc/.
- 67. Bolger AM, Lohse M, Usadel B. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics. 2014;30: 2114–2120. pmid:24695404
- 68. Langmead B, Trapnell C, Pop M, Salzberg SL. Ultrafast and memory-efficient alignment of short DNA sequences to the human genome. Genome Biol. 2009;10: R25. pmid:19261174
- 69. Jónsson H, Ginolhac A, Schubert M, Johnson PLF, Orlando L. mapDamage2.0: fast approximate Bayesian estimates of ancient DNA damage parameters. Bioinformatics. 2013;29: 1682–1684. pmid:23613487
- 70. Bankevich A, Nurk S, Antipov D, Gurevich AA, Dvorkin M, Kulikov AS, et al. SPAdes: A New Genome Assembly Algorithm and Its Applications to Single-Cell Sequencing. J Comput Biol. 2012;19: 455–477. pmid:22506599
- 71. Bernt M, Donath A, Jühling F, Externbrink F, Florentz C, Fritzsch G, et al. MITOS: Improved de novo metazoan mitochondrial genome annotation. Mol Phylogenet Evol. 2013;69: 313–319. pmid:22982435
- 72. Bondareva OV, Abramson NI. The complete mitochondrial genome of the common pine vole Terricola subterraneus (Arvicolinae, Rodentia). Mitochondrial DNA Part B. 2019;4: 3925–3926. pmid:33366254
- 73. Bondareva OV, Mahmoudi A, Bodrov SY, Genelt-Yanovskiy EA, Petrova TV, Abramson NI. The complete mitochondrial genomes of three Ellobius mole vole species (Rodentia: Arvicolinae). Mitochondrial DNA Part B. 2020;5: 2485–2487. pmid:33457837
- 74. Abramson NI, Golenishchev FN, Bodrov SY, Bondareva OV, Genelt-Yanovskiy EA, Petrova TV. Phylogenetic relationships and taxonomic position of genus Hyperacrius (Rodentia: Arvicolinae) from Kashmir based on evidences from analysis of mitochondrial genome and study of skull morphology. PeerJ. 2020;8: e10364. pmid:33240667
- 75. Li B, Wolsan M, Wu D, Zhang W, Xu Y, Zeng Z. Mitochondrial genomes reveal the pattern and timing of marten (Martes), wolverine (Gulo), and fisher (Pekania) diversification. Mol Phylogenet Evol. 2014;80: 156–164. pmid:25132128
- 76. Song N, Li X, Na R. Mitochondrial genomes of stick insects (Phasmatodea) and phylogenetic considerations. PLOS ONE. 2020;15: e0240186. pmid:33021991
- 77. Katoh K, Misawa K, Kuma K, Miyata T. MAFFT: a novel method for rapid multiple sequence alignment based on fast Fourier transform. Nucleic Acids Res. 2002;30: 3059–3066. pmid:12136088
- 78. Breinholt JW, Kawahara AY. Phylotranscriptomics: Saturated Third Codon Positions Radically Influence the Estimation of Trees Based on Next-Gen Data. Genome Biol Evol. 2013;5: 2082–2092. pmid:24148944
- 79. Arabi J, Cruaud C, Couloux A, Hassanin A. Studying sources of incongruence in arthropod molecular phylogenies: Sea spiders (Pycnogonida) as a case study. C R Biol. 2010;333: 438–453. pmid:20451886
- 80. Lobry JR. Properties of a general model of DNA evolution under no-strand-bias conditions. J Mol Evol. 1995;40: 326–330. pmid:7723059
- 81. Kumar S, Stecher G, Li M, Knyaz C, Tamura K. MEGA X: Molecular Evolutionary Genetics Analysis across Computing Platforms. Mol Biol Evol. 2018;35: 1547–1549. pmid:29722887
- 82. Hammer O, Harper D, Ryan P. PAST: Paleontological Statistics Software Package for Education and Data Analysis. Palaeontol Electron. 2001;4: 1–9.
- 83. Xia X, Xie Z, Salemi M, Chen L, Wang Y. An index of substitution saturation and its application. Mol Phylogenet Evol. 2003;26: 1–7. pmid:12470932
- 84. Xia X. DAMBE7: New and Improved Tools for Data Analysis in Molecular Biology and Evolution. Mol Biol Evol. 2018;35: 1550–1552. pmid:29669107
- 85.
Xia X, Lemey P. Assessing substitution saturation with DAMBE. 2nd ed. In: Vandamme A-M, Salemi M, Lemey P, editors. The Phylogenetic Handbook: A Practical Approach to Phylogenetic Analysis and Hypothesis Testing. 2nd ed. Cambridge: Cambridge University Press; 2009. pp. 615–630. https://doi.org/10.1017/CBO9780511819049.022
- 86.
Charif D, Lobry JR. SeqinR 1.0–2: A Contributed Package to the R Project for Statistical Computing Devoted to Biological Sequences Retrieval and Analysis. In: Bastolla U, Porto M, Roman HE, Vendruscolo M, editors. Structural Approaches to Sequence Evolution: Molecules, Networks, Populations. Berlin, Heidelberg: Springer; 2007. pp. 207–232. https://doi.org/10.1007/978-3-540-35306-5_10
- 87. Paradis E, Schliep K. ape 5.0: an environment for modern phylogenetics and evolutionary analyses in R. Bioinformatics. 2019;35: 526–528. pmid:30016406
- 88. Lanfear R, Frandsen PB, Wright AM, Senfeld T, Calcott B. PartitionFinder 2: New Methods for Selecting Partitioned Models of Evolution for Molecular and Morphological Phylogenetic Analyses. Mol Biol Evol. 2017;34: 772–773. pmid:28013191
- 89. Trifinopoulos J, Nguyen L-T, von Haeseler A, Minh BQ. W-IQ-TREE: a fast online phylogenetic tool for maximum likelihood analysis. Nucleic Acids Res. 2016;44: W232–W235. pmid:27084950
- 90. Hoang DT, Chernomor O, von Haeseler A, Minh BQ, Vinh LS. UFBoot2: Improving the Ultrafast Bootstrap Approximation. Mol Biol Evol. 2018;35: 518–522. pmid:29077904
- 91. Ronquist F, Teslenko M, van der Mark P, Ayres DL, Darling A, Höhna S, et al. MrBayes 3.2: Efficient Bayesian Phylogenetic Inference and Model Choice Across a Large Model Space. Syst Biol. 2012;61: 539–542. pmid:22357727
- 92. Rambaut A, Drummond AJ, Xie D, Baele G, Suchard MA. Posterior Summarization in Bayesian Phylogenetics Using Tracer 1.7. Syst Biol. 2018;67: 901–904. pmid:29718447
- 93. Miller MA, Pfeiffer W, Schwartz T. Creating the CIPRES Science Gateway for inference of large phylogenetic trees. 2010 Gateway Computing Environments Workshop (GCE). 2010. pp. 1–8.
- 94. Bouckaert R, Vaughan TG, Barido-Sottani J, Duchêne S, Fourment M, Gavryushkina A, et al. BEAST 2.5: An advanced software platform for Bayesian evolutionary analysis. PLOS Comput Biol. 2019;15: e1006650. pmid:30958812
- 95. Darriba D, Taboada GL, Doallo R, Posada D. jModelTest 2: more models, new heuristics and parallel computing. Nat Methods. 2012;9: 772–772. pmid:22847109
- 96. Molto JE, Loreille O, Mallott EK, Malhi RS, Fast S, Daniels-Higginbotham J, et al. Complete Mitochondrial Genome Sequencing of a Burial from a Romano–Christian Cemetery in the Dakhleh Oasis, Egypt: Preliminary Indications. Genes. 2017;8: 262. pmid:28984839
- 97. Chen W-C, Hao H-B, Sun Z-Y, Liu Y, Liu S-Y, Yue B-S. Phylogenetic position of the genus Proedromys (Arvicolinae, Rodentia): Evidence from nuclear and mitochondrial DNA. Biochem Syst Ecol. 2012;42: 59–68.
- 98. Schenk JJ, Rowe KC, Steppan SJ. Ecological Opportunity and Incumbency in the Diversification of Repeated Continental Colonizations by Muroid Rodents. Syst Biol. 2013;62: 837–864. pmid:23925508
- 99. Wiens JJ. Can Incomplete Taxa Rescue Phylogenetic Analyses from Long-Branch Attraction? Syst Biol. 2005;54: 731–742. pmid:16243761
- 100. Brinkmann H, van der Giezen M, Zhou Y, de Raucourt GP, Philippe H. An Empirical Assessment of Long-Branch Attraction Artefacts in Deep Eukaryotic Phylogenomics. Syst Biol. 2005;54: 743–757. pmid:16243762
- 101. Straub SCK, Moore MJ, Soltis PS, Soltis DE, Liston A, Livshultz T. Phylogenetic signal detection from an ancient rapid radiation: Effects of noise reduction, long-branch attraction, and model selection in crown clade Apocynaceae. Mol Phylogenet Evol. 2014;80: 169–185. pmid:25109653
- 102. Zou H, Jakovlić I, Zhang D, Hua C-J, Chen R, Li W-X, et al. Architectural instability, inverted skews and mitochondrial phylogenomics of Isopoda: outgroup choice affects the long-branch attraction artefacts. R Soc Open Sci. 7: 191887. pmid:32257344
- 103.
Corbet G. B. The mammals of the palaearctic region. a taxonomic review. London and Ithaca (NY): British Museum (Natural History) and Cornell University Press; 1978.
- 104.
Pavlinov IY. Systematics of Recent mammals. Moscow: Moscow University Press; 2003.
- 105.
Hooper ET, Hart BS. A synopsis of Recent North American microtine rodents. Miscellaneous Publications of the Museum of Zoology, University of Michigan; 1962.
- 106. Mckenna MC, Bell SK. Review of Classification of Mammals above the Species Level. J Vertebr Paleontol. 1997;19: 191–195.
- 107. Kretzoi M. 1969. Skizze einer Arvicoliden-Phylogenie–Stand 1969. Vertebr Hung Musei Hist Hung. 1969;11: 155–193.
- 108.
von Koenigswald W, KOENIGSWALD V. Schmelzstruktur und Morphologie in den Molaren der Arvicolidae (Rodentia). 1980.
- 109. Chaline J, Brunet-Lecomte P, Montuire S, Viriot L, Courant F. Anatomy of the arvicoline radiation (Rodentia): palaeogeographical, palaeoecological history and evolutionary data. Annales Zoologici Fennici. JSTOR; 1999. pp. 239–267.
- 110. Carroll LE, Genoways HH. Lagurus curtatus. Mamm Species. 1980; 1–6.
- 111.
Hall ER. The mammals of North America. New York; 1981.
- 112. Honacki JH, Kinman KE, Koeppl JW. Mammals species of the world; a taxonomic and geographic reference. 1982.
- 113. Tesakov AS, van Kolfschoten T. The Early Pleistocene Mimomys hordijki (Arvicolinae, Rodentia) from Europe and the origin of modern nearctic sagebrush voles (Lemmiscus). Palaeontol Electron. 2011;14: 1–11.
- 114. Modi WS. Phylogenetic history of LINE-1 among arvicolid rodents. Mol Biol Evol. 1996;13: 633–641. pmid:8676737
- 115. Barros P, Vale-Gonçalves HM, Paupério J, Cabral JA, Rosa G. Confirmation of European snow vole Chionomys nivalis (Mammalia: Rodentia: Cricetidae) occurrence in Portugal. Ital J Zool. 2016;83: 139–145.
- 116. Yannic G, Burri R, Malikov VG, Vogel P. Systematics of snow voles (Chionomys, Arvicolinae) revisited. Mol Phylogenet Evol. 2012;62: 806–815. pmid:22182990
- 117. Fejfar O, Heinrich W-D, Kordos L, Maul L. Microtoid cricetids and the Early history of arvicolids (Mammalia, Rodentia). Palaeontol Electron. 2011;14.
- 118. Schaub S. Über einige Simplicidentaten ausChina und Mongolei. Abh Schweiz-Schen Paläontol Ges. 1934;54: 1–40.
- 119.
Fejfar O. Microtoid Cricetids. The Miocene Land Mammals of Europe. Rössner GE, Heissig K, editors. München, Germany.: Verlag Dr. Friedrich Pfeil; 1999.
- 120. Maridet O. An Early Miocene microtoid cricetid (Rodentia, Mammalia) from the Junggar Basin of Xinjiang, China. Acta Palaeontol Pol. 2014 [cited 17 Feb 2021]. pmid:25908897
- 121. Kretzoi M. Acta Geologica Hungarica. Dolomys Ondatra. 1955;3: 347–355.
- 122. Kretzoi M. Pannonicola brevidens n.g n.sp., ein echter Arvicolide aus dem ungarischen Unterpliozän. Vertebr Hung Musei Hist-Nat Hung. 1965;7: 131–139.
- 123.
Zazhigin VS. Order Rodentia-rodents. In: Shanzer EV, editor. Stratigraphy of the USSR: The Quaternary system. Moscow: Nedra; 1982. pp. 294–305.
- 124.
Martin RA. In: Evolution of Tertiary Mammals of North America, Vol. 2. ed. Janis C. M., Gunnell G. F., and Uhen M. D. Published by Cambridge University Press; 2007. pp. 480–497.
- 125. Hordijk K, de Bruijn H. The succession of rodent faunas from the Mio/Pliocene lacustrine deposits of the Florina-Ptolemais-Servia Basin (Greece). Hellenic Journal of Geosciences. 2009; 44: 21–103.
- 126.
Repenning CA. Biochronology of the microtine rodents of the United States. Cenozoic Mammals of North America: Geochronology and Biostratigraphy. University of California Press, Berkeley, California; 1987. pp 236–268.
- 127. Martin RA. A preliminary diversity-divergence model for North American arvicolid rodents. Palaeobiodiversity and Palaeoenvironments. 2015; 95(3): 253–256.
- 128. Kowalski K. Fossil lemmings (Mammalia, Rodentia) from the Pliocene and early Pleistocene of Poland. 1977;22: 297–317.
- 129. Sukhov VP. Remains of lemmings in the Pliocene deposits of Bashkiria. Proc Zool Inst. 1976;66: 117–121.
- 130. Devyatkin EV, Zazhigin VS. Eopleistocene deposits and new mammalian find-localities of Northern Mongolia. 1974;1: 357–363.
- 131. Ruez DR, Gensler PA. An unexpected early record of Mictomys Vetus (Arvicolinae, Rodentia) from the Blancan (Pliocene) Glenns Ferry Formation, Hagerman Fossil Beds National Monument, Idaho. J Paleontol. 2008;82: 638–642.
- 132. Fejfar O, Repenning CA. The ancestors of the lemmings (Lemmini, Arvicolinae, Cricetidae, Rodentia) in the early Pliocene of Wölfersheim near Frankfurt am Main; Germany. Senckenberg Lethaea. 1998;77: 161–193.
- 133.
Zazhigin VS, Zykin VS. Novye dannye po stratigrafii pliotsena yuga Zapadno-Sibirskoi ravniny. Stratigrafiya pogranichnykh otlozhenii neogena i antropogena Sibiri. Novosibirsk: Inst. geologii i geofiziki SO AN SSSR; 1984. pp. 29–53.
- 134. Repenning CA. Beringian climate during intercontinental dispersal: a mouse eye view. Quat Sci Rev. 2001;20: 25–40.
- 135. Zazhigin V. Late Pliocene and Anthropogene rodents of the south of Western Siberia. Acad Sci USSR Trans. 1980;339: 1–156.
- 136.
Repenning CA. Allophaiomys and the age of the Olyor Suite, Krestovka sections, Yakutia. US Government Printing Office; 1992. pmid:1391607
- 137. Martin RA, Teskov AS. Introductory remarks: does Allophajomys exist? Paludicola.1998; 2(1): 1–7.
- 138. Cuenca‐Bescós G, López‐García JM, Galindo-Pellicena MA, García‐Perea R, Gisbert J, Rofes J, et al. Pleistocene history of Iberomys, an endangered endemic rodent from southwestern Europe. Integr Zool. 2014;9: 481–497. pmid:25236417
- 139. Cuenca Bescós G, Laplana C. Evolución de Iberomys (Arvicolidae, Rodentia, Mammalia) durante el Cuaternario español. XI Jorn Paleontol. 1995; 69–72.
- 140. Paupério J, Herman JS, Melo‐Ferreira J, Jaarola M, Alves PC, Searle JB. Cryptic speciation in the field vole: a multilocus approach confirms three highly divergent lineages in Eurasia. Mol Ecol. 2012;21: 6015–6032. pmid:23163319
- 141.
Tesakov AS. Biostratigraphy of Middle Pliocene–Early Pleistocene of southern part of Eastern Europe (based on small mammals). Moscow: Nauka; 2004.
- 142.
Agadjanian AK. Pliocene-Pleistocene Small Mammals of the Russian Plain. Nauka; 2009.
- 143. Ivanova VV, Erbajeva MA, Shchetnikov AA, Kazansky AY, Matasova GG, Alexeeva NV, et al. Tologoi key section: A unique archive for pliocene-pleistocene paleoenvironment dynamics of Transbaikalia, Bikal rift zone. Quat Int. 2019;519: 58–73.
- 144. Liu S, Jin W, Liu Y, Murphy RW, Lv B, Hao H, et al. Taxonomic position of Chinese voles of the tribe Arvicolini and the description of 2 new species from Xizang, China. J Mammal. 2017;98: 166–182. pmid:29674783
- 145.
Haslauer R. Microtus juldaschi. In: Wilson DE, Lacher TE Jr, Mittermeier RA, editors. Handbook of the Mammals of the World. Barcelona: Lynx Edicions; 2019. p. 332.