Populations at risk: conservation genetics of kangaroo mice (Microdipodops) of the Great Basin Desert

Abstract The Great Basin Desert of western North America has experienced frequent habitat alterations due to a complex biogeographic history and recent anthropogenic impacts, with the more recent alterations likely resulting in the decline of native fauna and flora. Dark (Microdipodops megacephalus) and pallid (M. pallidus) kangaroo mice are ecological specialists found within the Great Basin Desert and are potentially ideal organisms for assessing ecosystem health and inferring the biogeographic history of this vulnerable region. Herein, newly acquired nuclear-encoded microsatellite loci were utilized to assess patterns of variation within and among spatially discrete groups of kangaroo mice and to evaluate gene flow, demographic trends, and genetic integrity. Results confirm that there are at least three genetically distinct units within M. megacephalus and two such units within M. pallidus. The three units of M. megacephalus appear to have different demographic histories, with effectively no gene flow among them since their divergence. Similarly, the two units of M. pallidus also appear to have experienced different demographic histories, with effectively no gene exchange. Contemporary effective population sizes of all groups within Microdipodops appear to be low (<500), suggesting that each genetic lineage may have difficulty coping with changing environmental pressures and hence may be at risk of extirpation. Results of this study indicate that each Microdipodops group should be recognized, and therefore managed, as a separate unit in an effort to conserve these highly specialized taxa that contribute to the diversity of the Great Basin Desert ecosystem. The Great Basin Desert of western North America has experienced frequent habitat alterations due to a complex biogeographic history and recent anthropogenic impacts, with the more recent alterations likely resulting in the decline of native fauna and flora. Herein, newly acquired nuclear-encoded microsatellite loci were utilized to assess patterns of variation within and among spatially discrete groups of the dark (Microdipodops megacephalus) and pallid (M. pallidus) kangaroo mouse, and to evaluate gene flow, demographic trends, and genetic integrity. Results of this study indicate that each Microdipodops group should be recognized, and therefore managed, as a separate unit in an effort to conserve these highly specialized taxa that contribute to the diversity of the Great Basin Desert ecosystem (photo credit J. C. Hafner).


Introduction
The Great Basin Desert of western North America is characterized by a series of alternating islands of mountain ranges and desert basins (Fiero 1986) that formed a backdrop to a dynamic biogeographic history (Davis 2005). The glacial-interglacial cycles of the Pleistocene (Riddle 1995) and the associated rise and fall of pluvial lakes (Benson 1981), shifting climatic patterns (Atvens 1952), and floristic transitions (Reveal 1979) have caused numerous habitat alterations throughout the Great Basin Desert. More recently, anthropogenic habitat alterations (e.g., introduction of nonnative plant species, increased wildfires, and cultivation and irrigation) have also plagued the area (Hafner and Hafner 1998). These alterations have caused a significant loss of available habitat and subsequent reduction in the abundance of native fauna and flora. For example, representatives of the rodent genus Microdipodops (kangaroo mice; family Heteromyidae) have become increasingly rare members of the Great Basin Desert community (Hafner and Upham 2011).
Two species of Microdipodops are currently recognized: the dark kangaroo mouse (M. megacephalus) and the pallid kangaroo mouse (M. pallidus). Both species are sand-obligate endemics to the Great Basin Desert and, as such, are highly specialized to survive in an extreme environment (Hafner 1981). In fact, morphology within the genus is extremely conserved with only slight differences between sibling taxa (Hafner et al. 2008). Given their ecological specializations, these small nocturnal rodents likely serve as indicator species of healthy, sandy desert habitats of the Great Basin. Field observations, however, have concluded that the numbers of both M. megacephalus and M. pallidus are dwindling (Hafner 1981;Hafner and Hafner 1998;Hafner et al. 2008;Hafner and Upham 2011), as is the case for other flora and fauna distributed across the Great Basin Desert (Brussard et al. 1998). However, both Microdipodops species are listed as "Least Concern" by the International Union for Conservation of Nature (IUCN) and are not protected . Given their decreasing numbers, this listing is outmoded and management of kangaroo mice, along with other Great Basin Desert organisms, will be necessary to help preserve this threatened ecosystem.
Microdipodops megacephalus and M. pallidus have unique habitat associations within the Great Basin Desert. Although their distributions overlap (Fig. 1), these species show differential niche specializations. Microdipodops megacephalus is primarily restricted to sandy soils with gravel overlay and found in association with sagebrush and/or rabbit brush (Hafner and Upham 2011;and references therein); whereas M. pallidus prefers greasewood and fine soils with no gravel overlay (Hafner 1981; and references therein). Ancient and current habitat alterations have led to fragmented distributions for both species such that current intraspecific ranges are disjunct (Figs. 1,2), separated either by geological barriers (e.g., mountain ranges) or unsuitable habitat (Hafner et al. 2008;Hafner and Upham 2011).
These unique, fragmented distributions and ecological specializations have made kangaroo mice the recent subjects of several studies that used mitochondrial DNA (mtDNA) gene regions to elucidate the biogeographic history of the Great Basin Desert (Hafner et al. 2006, 2008; Hafner and Upham 2011;Light et al. 2013). These studies identified and supported four distinct mtDNA clades in M. megacephalus (the eastern, central, western, and Idaho clades; Fig. 2A) and two distinct mtDNA clades in M. pallidus (the eastern and western clades; Fig. 2B). While the identification of genetically discrete units within each species is important, additional analyses using fast-evolving nuclear markers, such as microsatellites, are necessary to verify the results of the mtDNA data. These markers also can help to estimate parameters for conservation and management of these specialized taxa; for example, estimates such as rates of gene flow and effective population sizes for Microdipodops are currently unknown. Lastly, examination of multiple markers can facilitate a better understanding of genetic lineages within a species (Avise 1994), especially as these markers may have different evolutionary histories (e.g., Yang and Kenagy 2009).
Herein, we use microsatellite markers to provide an assessment of nuclear variation within each Microdipodops species and to test the findings from previous studies which used mtDNA sequence data (Hafner et al. 2006(Hafner et al. , 2008Hafner and Upham 2011;Light et al. 2013). We hypothesize that microsatellite markers will support discrete genetic units within each Microdipodops species and uncover the same geographic groups found in previous studies. Due to the wealth of information available regarding Microdipodops biogeography, population-level analyses are performed on microsatellite data with sam-ples disaggregated into geographic regions identified in previous studies (Fig. 2) and results are interpreted in reference to Great Basin biogeography. These findings will help to identify evolutionarily significant units and address issues of management, conservation, and desert biogeography that can be applied to other flora and fauna of the threatened Great Basin Desert.

Specimens examined
A total of 184 specimens of M. megacephalus from 46 localities, and a total of 105 specimens of M. pallidus from 27 localities, were used in this study (Table A1; terminology follows Hafner et al. 2008;Hafner and Upham 2011;and Light et al. 2013). The majority of these specimens were collected for use in prior studies: the M. megacephalus specimens were collected between 1975 and 1976, 1999 and 2007, and in 2011. The M. pallidus specimens were collected between 1999 and 2005, with one individual sampled in 1975. Any newly collected specimens were collected according to procedures approved by the Occidental College's Animal Care and Use Committee and the American Society of Mammalogists (Sikes et al. 2011). All tissues were stored in a À80°C freezer.
For many of the analyses, populations were defined by grouping specimens together within each species based on (A) (B) Figure 2. Detailed geographic distributions of mtDNA clades within dark and pallid kangaroo mice. (A) Geographic distribution of Microdipodops megacephalus, with labels corresponding to mtDNA clades (eastern, central, western, and Idaho) discussed in prior studies (Hafner et al. 2006;Hafner and Upham 2011). The genetically distinct Valley Falls subunit (which is nested within the western clade; Hafner and Upham 2011;Light et al. 2013) also is labeled. (B) Geographic distribution of M. pallidus, with labels corresponding to mtDNA clades (eastern and western) from prior studies (Hafner et al. 2008); the isolated Deep Springs locality (which is nested within the western clade; Hafner et al. 2008) also is labeled. Dots indicate exact collecting localities for specimens used in this study and identified in Hafner et al. (2008) and Hafner and Upham (2011); outline of the State of Nevada provides proper orientation. geography. These geographic units correspond to previously identified mtDNA clades (Fig. 2) and subclades. Previous studies recognized four geographic units within M. megacephalus (Hafner and Upham 2011;Light et al. 2013): the eastern geographic unit (n = 49) with two subunits (eastern subunit, n = 25; western subunit, n = 24), the central geographic unit (n = 69) with two subunits (central subunit, n = 19; western subunit, n = 50), the western geographic unit (n = 62) with one subunit (Valley Falls, n = 9), and the Idaho geographic unit (n = 4). Two geographic units were recognized within M. pallidus (Hafner et al. 2008;Light et al. 2013): the eastern geographic unit (n = 42) with two subunits (eastern subunit, n = 18; south-central subunit, n = 24), and the western geographic unit (n = 63) with one subunit (Deep Springs, n = 10).

Laboratory methods
DNA extracts were available from previous studies (Hafner et al. 2006(Hafner et al. , 2008Hafner and Upham 2011). When original extractions were depleted, DNA was extracted from liver or kidney tissues as described by Hafner et al. (2006). Seventeen polymorphic microsatellite loci, developed for Microdipodops by Lance et al. (2010), were tested on preliminary samples and loci that did not reliably amplify were subsequently removed. Polymerase chain reactions (PCR) followed Boutin-Ganache et al. (2001) and contained three primers: a forward primer with an attached 16-bp tail sequence (5'-CAGTCGGGCGTCAT CA-3'), a 6-FAM or HEX (Dye Set D, Applied Biosystems, Foster City, CA) labeled tail sequence (defined above), and an unlabeled reverse primer. Amplified DNA from each PCR reaction was combined with a 400 HD Rox size-standard DNA ladder (Applied Biosystems) and electrophoresed on 6% polyacrylamide gels using an ABI PRISM 377 DNA Sequencer (Applied Biosystems). Sizes of microsatellite fragments were visualized in GENESCAN v. 3.1.2 (Applied Biosystems) and assessed using GENO-TYPER v. 2.5 (Applied Biosystems).

Data analysis
Each microsatellite locus was tested for conformance to the expectations of Hardy-Weinberg equilibrium (HWE) using GENEPOP v. 4.0 (Raymond and Rousset 1995;Rousset 2008). Significance was assessed at the 0.05 level, using exact tests with 20 batches and 5000 iterations per batch, and sequential Bonferroni adjustment was used to correct for multiple testing (Rice 1989). Loci that differed significantly from the expectations of HWE were assessed either by rescoring gels and/or rerunning PCR to determine if genotyping error caused spurious results. GENEPOP also was used to calculate the expected and observed numbers of heterozygotes, test for genotypic disequilibrium, and calculate gene frequencies when null alleles were present. Number of alleles and allelic richness (i.e., number of alleles per locus, averaged over the smallest population) for each locus were calculated with FSTAT v. 2.9.3.2 (Goudet 1995). Population structure within each species was first assessed to test for genetic homogeneity with an analysis of molecular variance (AMOVA; Excoffier et al. 2005) implemented in ARLEQUIN v. 3.5.1.2 (Excoffier and Lischer 2010). AMOVA was performed in a hierarchical fashion with populations grouped a priori by geographic unit, and significance was assessed by 10,000 randomization replicates. Using the same assortment of geographic units, F ST and R ST (a F ST analog assuming a stepwise mutation model; Slatkin 1995) statistics were estimated with ARLE-QUIN, and significance at the 0.05 level was assessed by permuting individuals between samples 10,000 times. Allele size permutation tests were performed to compare F ST and R ST statistics using SPAGeDi 1.4 (Hardy and Vekemans 2002).
The Bayesian multilocus clustering algorithm found in STRUCTURE v. 2.3.3 (Pritchard et al. 2000) was used to examine fine-scale population structure without defining populations a priori. Analyses were run in a hierarchical manner, first within M. megacephalus and M. pallidus, and then within each geographic unit; the eastern, central, western, and Idaho geographic units of M. megacephalus, and the eastern and western geographic units of M. pallidus). The population admixture model was used with 10 replicate runs from K = 1 to K = 10, where K is a userdefined number of clusters. Each run consisted of a burnin of 10,000 steps followed by 100,000 additional steps. To evaluate the most likely K value, STRUCTURE HARVESTER (Earl and vonHoldt 2012) was used to graph both the mean estimated ln Prob (Data) and DK (change in ln Prob (Data) between successive K values) as suggested by Evanno et al. (2005).
MIGRATE-N v. 3.2.1.6 (Beerli and Felsenstein 1999;Beerli 2006) was used to estimate theta (h; h = 4N eLT l, where N eLT is the long-term effective population size and l is the per-generation mutation rate) and M (mutationscaled migration rate) among geographic units within each species using Bayesian inference. Due to the small sample size of the M. megacephalus Idaho geographic unit (n = 4), it was excluded from the analysis. Theta was estimated to detect if there were significant differences in N eLT among geographic units, while M was used to quantify average, long-term gene flow between geographic units. Preliminary runs were performed to estimate priors for M and h for final runs. Final runs were executed in replicate at different starting points and parameter esti-mates were examined to ensure chain mixing and convergence. For M. megacephalus, runs consisted of three long chains, geometric heating, and a burn-in of 100,000 steps followed by 1,000,000 steps with a tree recorded every 100 steps, resulting in 10,000 trees sampled. For M. pallidus, runs consisted of one long chain and a burn-in of 10,000 steps followed by 100,000 steps with a tree recorded every 100 steps, resulting in 1000 trees. In all analyses, effective sample sizes were >50.
IMA (Hey and Nielsen 2007) also was used to determine h and M among geographic units of M. megacephalus and M. pallidus in a pairwise manner. IMA differs from MIGRATE-N in that it takes coancestry into account when looking at migration and it assumes there is one ancestral panmictic population for each extant population. This assumption allows the estimation of the ancestral effective population size and time since divergence (t), where a positive t value would indicate divergence and a value that peaked at or very close to zero would reveal no divergence (Portnoy et al. 2010). Additionally, while both MIGRATE-N and IMA use the Metropolis-Hastings criterion, IMA incorporates a Metropolis-Coupled version of the algorithm which enables multiple heated chains to search the parameter space simultaneously and can provide a more thorough mixing of chains (Hey and Rasmus 2004). Preliminary runs were performed to assess whether the heating conditions and M, h, and t priors were appropriate for the data set. Final runs consisted of 50 chains with geometric heating (to ensure acceptable chain mixing and low autocorrelations) and a burn-in of at least 1,000,000 generations followed by at least 90,000 generations (resulting in effective sample sizes >50). Final runs were executed in replicate with different starting seeds to ensure convergence, and the R package BOA (Smith 2005) was used to visually assess convergence of posterior distributions and examine autocorrelation at different lags to determine appropriate run time.
Bayesian inference of immigration rates (BIMr; Faubet and Gaggiotti 2008) was used to estimate current rates of gene exchange among geographic units, thus facilitating comparison with long-term estimates of M from MIGRATE-N and IMA. Preliminary pilot runs (each at a length of 2000 steps) were executed to provide a rough estimation of starting points for final runs. Replicate runs consisted of a burn-in of 20,000 iterations, followed by an additional 100,000 and 60,000 iterations for M. megacephalus and M. pallidus, respectively. The R package BOA (Smith 2005) was used to examine autocorrelation at different lags to determine appropriate run time and visually assess convergence of posterior distributions. Density functions were analyzed and the mode (point estimate) and 95% highest posterior density interval (HPDI) were noted.
The program LDNE (Waples and Do 2008) was used to estimate contemporary effective population size (N e ) via the modified linkage disequilibrium approach (Hill 1981;Waples 2006) for each geographic unit within M. megacephalus and M. pallidus. Effective population size is a crucial parameter in conservation and wildlife management because of its influence on population viability and ability to predict extinction risk (Luikart et al. 2010). LDNE assumes that the correlation of unlinked alleles at unlinked loci arises from genetic drift in an isolated population (Hill 1981;Wang 2005) with estimates reflecting the number of parents that contributed to the sample (Waples 2005). Microdipodops are semelparous with a generation time of 1 year (Hall 1941), meaning that estimates are of contemporary N e rather than effective number of breeders (N b ; Jorde and Ryman 1995;Waples 2005). Because allele frequencies close to 0 or 1 can skew N e results (Waples 2006;Portnoy et al. 2009), alleles that had a frequency of <2% were omitted from analyses. For all analyses, a random mating model was assumed and 95% jackknife confidence intervals were assessed (Waples 2006).
Extended Bayesian Skyline Plots (EBSP; Heled and Drummond 2008) were used to estimate N e through time within each geographic unit of M. megacephalus and M. pallidus. EBSP differs from other demographic analyses in that it estimates the population function directly from the data. Furthermore, unlike estimates of contemporary N e , EBSP estimates N e through a coalescent approach, and N e estimates can therefore be used to find varying historic demographic changes across lineages. Each geographic unit was analyzed individually (although individuals from the M. megacephalus Idaho geographic unit were excluded due to small sample size). MSVAR v. 1.3 (Beaumont 1999) was used to estimate the average mutation rate for all loci within each geographic unit. Uniform rate analyses were run using a strict molecular clock following a stepwise mutation model. A minimum of two runs of 1 billion generations were performed, with a tree recorded every 40,000 steps after a 10% burn-in. Effective sample sizes and number of population size changes were assessed in Tracer v 1.5 (Rambaut and Drummond 2007). Population size data were plotted using R (R Development Core Team 2011).

Results
Of the 17 initial polymorphic loci screened, 11 and 10 loci amplified successfully and were used in the population genetic analyses of M. megacephalus and M. pallidus, respectively (summary data available in Tables A2, A3) for multiple tests, genotypes at two loci (Mime11 and Mime32) in M. pallidus from the western geographic unit deviated significantly from the expectations of HWE. This was due to the isolated population from Deep Springs (Fig. 2B) where homozygote excess occurred at both loci. When Deep Springs was excluded from analysis, all loci conformed to the expectations of HWE. Less computationally intensive analyses (e.g., AMOVA, STRUCTURE, pairwise R ST ) were run including and excluding the two deviated loci, and there was no difference in the results. Therefore, results reported in this study included all loci that amplified successfully.
Allele size permutation tests indicated that R ST values were consistently significantly greater than F ST values, indicating that F ST may be underestimating actual values of genetic structure (Table A4; Hoffman et al. 2005). Thus, only R ST results are presented here. AMOVA revealed significant population structure among geographic units and among subunits within geographic units in both species (P < 0.001; Table 1 Additional bar plots with increasing number of clusters also were analyzed in case discrete populations could be distinguished. However, population structure became less resolved with K > 3. Results from the DK metric suggested by Evanno et al. (2005) indicated that K = 2 was the most strongly supported (D ln Prob (Data) = 388.02 and 551.32 for M. megacephalus and M. pallidus, respectively), while K = 3 was the next most strongly supported (D ln Prob (Data) = 106.02 and 24.86 for M. megacephalus and M. pallidus, respectively). The 2 clusters were an eastern/ central and western/Idaho group for M. megacephalus and an eastern and western group for M. pallidus. When individual geographic units were analyzed separately, K = 1 was the most likely number of clusters of nuclear variation for all units with the exception of the M. pallidus western geographic unit, where K = 2 was most likely number of clusters (corresponding to Deep Springs and the rest of the western geographic unit).
Gene flow estimates from MIGRATE-N for M. megacephalus were low, with modal values for M ranging from 0 (eastern ↔ western and central ? western) to 0.05 (eastern ? central; Table 2). While estimates of M from central ? eastern were higher (0.18), the 2.5% and 97.5% bounds were 0.01 and 0.92, respectively. With such a large confidence interval it was therefore unclear how much long-term gene exchange was occurring between these two geographic units. Estimates of h were not statistically different among geographic units (Table 3).
Estimates of M from MIGRATE-N within M. pallidus were 0.03 (eastern ? western) and 0.01 (western ? eastern; Table 2). The estimates were not significantly different from zero and the upper bounds were 0.08 and 0.05, respectively. Estimates of h for the M. pallidus eastern and western geographic units were not statistically different (Table 3).
Results from IMA analysis of M. megacephalus revealed that the lower bound of time since divergence (t) did not include zero for the eastern, central, and western geographic units, indicating divergence from an ancestral, panmictic population. Estimates for t were quite low (ranging from ca. 8600-13,900 years before present) with large confidence intervals (often over hundreds of thousands of years; data available upon request). Estimated long-term gene flow (M) from IMA was very small, with a lower confidence interval and mode for all three groups  (Table 3). IMA analysis of M. pallidus revealed that the lower bound of t did not include zero indicating that the eastern and western geographic units had diverged from an ancestral, panmictic population. The estimate for t was low (ca. 9500 years before present), with a rather large confidence interval (spanning nearly 500,000 years). Estimated rates of gene flow for eastern ? western and western ? eastern were 0.09 and 0.04, respectively (Table 2). These estimates of M indicate no or extremely low levels of possible long-term gene exchange. Theta estimates for the eastern and western M. pallidus geographic units were not statistically different (Table 3).
Estimates of current gene flow rates from BIMr for geographic units of M. megacephalus showed modal values from 2.19 9 10 À11 to 3.56 9 10 À16 (Table 4). Such small estimates suggest effectively no gene exchange across geographic units within the last generation. Modal estimates for the geographic units of M. pallidus (Table 4), while larger than those for M. megacephalus, similarly suggest effectively no current gene flow between the eastern and western geographic units within the last generation.
Point estimates of contemporary N e , as well as minimum and maximum estimates (based on 95% confidence intervals obtained by jackknifing), for M. megacephalus and M. pallidus are presented in Table 3. For all populations of both species, point estimates were <500. Minimum estimates of N e (based on 95% confidence intervals), which may serve as a conservative estimate for wildlife management (Waples and Do 2010), ranged from 108.8 individuals (western geographic unit) to 179.8 (central geographic unit) in M. megacephalus, and were 95.1 and 80.5 in the eastern and western geographic units of M. pallidus, respectively. The eastern unit of M. megacephalus and the eastern unit of M. pallidus were the only geographic units with upper limits of infinity (∞).
Mutation rates estimated by MSVAR averaged 2.40 9 10 À4 , 2.75 9 10 À4 , 2.45 9 10 À4 , 3.23 9 10 À4 , and 3.89 9 10 À4 for the M. megacephalus eastern, central,   and western geographic units, and the M. pallidus eastern and western geographic units, respectively. EBSP results showed that the three M. megacephalus geographic units might have undergone a recent population expansion, whereas the M. pallidus eastern unit remained fairly constant and the western unit underwent a recent population contraction (Fig. 4). None of these results, however, were significant (Fig. 4). While these results seem to contradict one another, it is important to note that EBSP is generating long-term estimates of N e while LDN e is generating estimates of contemporary N e . The two estimates may therefore differ because (a) the time periods to which the two effective size estimates apply are not necessarily concordant (Waples 2005) and (b) long-term estimates are more affected by long-term gene flow, even from extinct demes, than contemporary estimates and may reflect global effective size rather than local (Schwartz et al. 1999). One additional parameter that can be estimated from EBSP is the number of population size changes since time of coalescence. When examining the 95% HPD of demographic population size changes, we failed to reject a constant population size in the M. megacephalus eastern and central geographic units and both M. pallidus geographic units (population size changes ranged from 0 to 3 in all units, except the M. pallidus eastern unit which ranged from 0 to 2). We could, however, reject a constant population size in the M. megacephalus western geographic unit (population size changes ranged from 1 to 3).

Discussion
Microsatellite markers reveal that M. megacephalus and M. pallidus are comprised of multiple genetically distinct units within the Great Basin Desert. Primary population genetic analyses (AMOVA and pairwise R ST ) support genetic heterogeneity within each species and STRUCTURE analyses revealed that K = 3 was the most likely number of clusters for both M. megacephalus and M. pallidus (Fig. 3). Although the DK metric indicated that K = 2 was the most likely number of clusters for M. megacephalus and M. pallidus, this method is more conservative and often underestimates K with insufficient sample sizes (Evanno et al. 2005). Thus, the STRUCTURE results are more appropriate for the Microdipodops data examined in this study. The genetic clusters identified here correspond to the mtDNA clades identified in previous studies (Hafner et al. 2006(Hafner et al. , 2008Hafner and Upham 2011;Light et al. 2013). The only exception is the lack of recognition of a cluster corresponding to the Idaho clade within M. megacephalus. Hierarchical STRUCTURE analyses failed to reveal two genetic clusters within the M. megacephalus western/ Idaho cluster; however, increasing K to 3 showed a clearly diverged Idaho cluster with an unresolved western cluster (results available upon request). Increasing the number of individuals from the Idaho clade in future studies would probably tease it apart from the western cluster (Evanno et al. 2005;Hale et al. 2012). In addition to supporting Central ? Eastern unit 4.71 9 10 À9 3.56 9 10 À16 1.15 9 10 À6 Central ? Western unit 5.63 9 10 À9 1.82 9 10 À15 9.18 9 10 À7 Western ? Eastern unit 2.43 9 10 À9 2.19 9 10 À11 5.93 9 10 À10 Western ? Central unit 1.26 9 10 À12 2.19 9 10 À11 5.9 9 10 À10 M. pallidus Eastern ? Western unit 3.4 9 10 À3 2.33 9 10 À3 4.54 9 10 À2 Western ? Eastern unit 2.3 9 10 À4 1.7 9 10 À3 3.24 9 10 À2 distinct genetic lineages within each Microdipodops species, this study provides an in-depth assessment of parameters important for conservation and management, including patterns of current and historical connectivity (gene flow), effective population sizes, and demographic history. Results of this study therefore provide information that can be used for management strategies and conservation efforts specific to each evolutionarily significant unit within the Great Basin Desert. Significant population structure detected within M. megacephalus supports the perspective that kangaroo mice found in the eastern, central, and western units are, at minimum, distinct populations (specimens representing the Idaho clade of M. megacephalus could not be analyzed rigorously due to a small sample size). In agreement with previous phylogenetic analyses, population genetic analyses of microsatellite data reveal a close affinity between the eastern and central populations (Table A4) and a clearly more differentiated western population (Hafner and Upham 2011). Divergence following isolation may partly explain these genetic differences.
Our results indicate that since their divergence there has been effectively no gene flow among the eastern, central, and western populations of M. megacephalus (Table 2). Previous molecular evidence using mtDNA data suggests that lineage divergence within M. megacephalus occurred in the Pliocene,~4 million years ago (Ma; Hafner and Upham 2011;Hafner et al. 2008), and fossil evidence from the late Blancan (1.9-2.9 Ma) supports that kangaroo mice diversified prior to the Pleistocene outside the Great Basin Desert. The significant difference between R ST and F ST values (Table A4) suggests that the populations have been isolated for a sufficiently long period of time such that mutation has played a relatively important role in genetic differentiation. The lack of significant differences among our historical h values, and t parameters significantly larger than zero, suggest that each lineage may have diverged from a single ancestral population (Table 3). These results fail to reject the hypothesis that multiple lineages of M. megacephalus diverged from a common ancestral population and that some or all of these lineages invaded the Great Basin Desert in the early Pleistocene (coincident with the formation of appropriate sandy habitats; see Hafner and Upham 2011 and references therein). Interestingly, this early Pleistocene colonization has been observed in other Great Basin taxa, such as pikas, brown creepers, and mountain chickadees (Grayson 2005;Spellman et al. 2007;Manthey et al. 2011). It is important to note that our IMA estimates of t are significantly more recent than divergence times estimated in previous studies (Hafner et al. 2008;Hafner and Upham 2011). This discrepancy may be due to a complicated biogeographic history of the region making it difficult to track species history, sex-biased dispersal, and associated complications of using different genetic markers. Although our demographic analyses postdate the Pliocene-Pleistocene transition (Fig. 4), we do observe a fairly constant population size over the past 200,000 years followed by possible recent population expansions. The recent expansion within the M. megacephalus central population is strongly supported by previous studies using mtDNA Bayesian Skyline Plots (BSP; Light et al. 2013) and directional analyses of phylogeographic patterns (DAPP; Hafner and Upham 2011). The recent expansions within the M. megacephalus eastern and western populations are not as strongly supported in previous analyses (Hafner and Upham 2011;Light et al. 2013), but this may be attributed to incomplete lineage sorting. Additionally, a constant population size may seem unlikely over a time period filled with climatic oscillations. Therefore, it seems more reasonable that the data do not contain enough demographic signals for these past events. Furthermore, it is important to note that due to excessively long computation times all EBSP analyses were performed using a simple model of evolution. Future studies comparing the results of more complex models of evolution (although previous studies note that skyline plots can be similar regardless of the model used; Allen et al. 2012) and assessing the utility of EBSP analyses on microsatellite data should be performed.
Lack of evidence for current gene flow (Table 4), significant differences in microsatellite allele and genotype distributions, and previously documented reciprocal monophyly among M. megacephalus populations using mtDNA data (Hafner and Upham 2011) support the view that each population, at the very least, should be managed as an evolutionarily significant unit. As noted by Hafner and Upham (2011), the populations are very similar morphologically and are distributed in an allopatric manner (each population is separated by unsuitable habitat or geological barriers). Despite evidence for recent expansions, N e estimates for all populations had lower bounds of confidence intervals and point estimates <500 (Table 3). As an N e > 50 is needed to avoid inbreeding and an N e > 500 to avoid extinction due to an inability to cope with environmental change (Franklin 1980;Jamieson and Allendorf 2012), our results suggest that these populations may be unable to adapt to environmental change and could be at risk for extirpation (Franklin 1980). It is important to note that the exact N e required for both long-and short-term sustainability has been disputed, and the minimum N e may be higher than 50 (Nunney and Campbell 1993;Lande 1995), and the appropriate N e may vary among populations (Flather et al. 2011). Regardless, measures must be taken to conserve each genetically distinct lineage with appropriate management techniques for each population.
Similarly, the eastern and western M. pallidus populations are genetically distinct units that likely diverged 4 Ma (Hafner et al. 2008) with effectively no gene flow (far less than one migrant per generation) between them. Again, it is important to note that our IMA estimates of t are significantly more recent than divergence times estimated in previous studies (Hafner et al. 2008), possibly due to a variety of reasons (see above). Similar to M. megacephalus, it is possible that one panmictic ancestral population (supported by our homogenous h estimates, and positive t estimate) diverged outside of the Great Basin Desert and two independent lineages invaded the region at the beginning of the Pleistocene (supported by the significant difference between R ST and F ST values [ Table A4]). The series of mountain chains that currently serve as a physiographic baffle between the eastern and western populations, likely prevented historic gene flow between these two lineages allowing for further divergence. Demographic analyses also suggest the M. pallidus western population has undergone a recent population contraction while the eastern population has remained constant, or has possibly undergone a population expansion, indicating that these two populations have historically been demographically independent from each other. Demographic results based on EBSP, however, should be interpreted cautiously (see above).
The lower bounds and point estimates of N e of both the eastern and western populations of M. pallidus are well below 500 (Franklin 1980); the western population even has an upper bound below 500 (Table 3). While these estimates may seem low, similar results have been found in other terrestrial vertebrates, some of whom are endangered (Nunney 1993;Nunney and Campbell 1993;Frankham 1995;Phillipsen et al. 2011;Hurtado et al. 2012). Additionally, the small N e of the western population is consistent with results from this study and a previous study using mtDNA, both indicating a recent population contraction (Light et al. 2013). Thus, both the eastern and western populations may be in danger of extirpation and separate management practices for each population should be enforced (Traill et al. 2010). To adequately measure the risk of extirpation, it will be important to further assess census sizes for both M. megacephalus and M. pallidus, which may be 2-10 times larger than these effective population size estimates (Nunney 1993;Nunney and Campbell 1993;Frankham 1995).

Broad implications
The amount of available habitat within the Great Basin Desert is decreasing as a result of a variety of anthropogenic alterations, and future climate change is predicted to reduce available habitat even further. Chaplin et al. (2000) ranked the Great Basin as second in number of imperiled species among ecoregions of the United States. Habitat loss through agricultural practices, wildfires, and invasive plants has devastated the low-elevation areas where kangaroo mice from the eastern and western populations of M. megacephalus are distributed. Recent attempts to trap dark kangaroo mice from northern localities where mice were once abundant have been unsuccessful (J. C. Hafner, unpubl. data). Furthermore, repeated efforts to collect M. pallidus in once fruitful areas have either proven to be increasingly difficult or completely unsuccessful (Hafner et al. 2008;J. C. Hafner, unpubl. data). As rare and highly specialized members of the Great Basin Desert, Microdipodops likely serve as indicator species of a healthy sandy desert ecosystem (Light et al. 2013). Reduction in Microdipodops abundance may signal deterioration of the habitat, and further reduction in their abundance may prove detrimental to the survival of individual populations. The genus Microdipodops is a rare and highly specialized endemic of the Great Basin Desert, and this study provides further support that management and conservation efforts should be applied to each population in an effort to conserve these valuable taxa and the imperiled habitats of the Great Basin Desert.