Research Article
Print
Research Article
Genetic divergence and demography of pudu deer (Pudu puda) in five provinces of southern Chile, analyzed through latitudinal and longitudinal ranges
expand article infoNelson Colihueque, Javier Cabello§, Andrea Fuentes-Moliz|
‡ Universidad de Los Lagos, Osorno, Chile
§ Centro de Conservación de la Biodiversidad Chiloé-Silvestre, Ancud, Chile
| Universidad de Sevilla, Sevilla, Spain
Open Access

Abstract

Pudu deer (Pudu puda) is endemic to the temperate rainforests of Chile. Genetic studies at different geographic scales for this species are required to better determine the genetic divergence within and among populations and their demography across the distribution range. These data can provide unique insights into the species or population status for conservation plans and decision-makers. We analyzed the mtDNA control region (CR) and cytochrome b (Cyt b) sequences of pudu deer in five provinces of southern Chile located at different latitudinal locations (Cautín, Valdivia, Osorno, Llanquihue and Chiloé Island) and three geographic areas within the studied provinces, representative of different longitudinal sites (Andes range, Central Valley and Coastal Range), to understand their genetic divergence and demography. The haplotype (H) and nucleotide (Π) diversities of CR and Cyt b ranged from 0.64286 to 0.98333 and from 0.00575 to 0.01022, respectively. CR diversity was significantly different among provinces, with Valdivia showing higher values than Llanquihue and Chiloé Island (H = 0.98333 vs. 0.64286–0.92727, P < 0.05). Cyt b variation also showed significant differences among provinces, particularly, among Cautín and Llanquihue (H = 1.000 vs. 0.222, P < 0.05). Genetic structuring among provinces was relatively high, as indicated by the FST index (FST = 0.41905). Clustering analysis indicated the presence of a distinctive cluster for Chiloé Island individuals. Fu’s FS and Tajima’s D based on CR revealed significant, negative deviations from equilibrium for Chiloé Island (D = -1.65898), Valdivia (Fs = -7.75335) and Llanquihue (Fs = -3.93267), suggesting population expansion in these provinces. Analysis at the longitudinal range showed significant differences among areas based on Π (P < 0.05), with the Andes range and Central Valley showing higher diversity than the Coastal Range. Neither population structuring (FST = 0.01360, P > 0.05) nor distinctive clusters in the longitudinal range were observed. Fu’s Fs and Tajima’s D were negative and significant for the Coastal Range based on CR (Fs = -6.64752, P < 0.001) and Cyt b (D = -1.74110, P < 0.05), suggesting the existence of population expansion. Our results suggest that pudu deer in the analyzed provinces is a genetically structured species, which could be associated with reduced panmixia among populations. The genetic divergence pattern and the population expansion recorded are likely to be associated with past processes of recolonization after Pleistocene glaciation events.

Keywords

Control region, cytochrome b, genetic divergence, population structure, pudu deer

Introduction

Pudu deer, Pudu puda (Molina, 1782), is endemic to southern South America and is characterized as being one of the smallest deer in the world due to its short shoulder height (30–40 cm) and small body weight (< 15 kg) (Jiménez 2010). This species is distributed in Chile from 35°10'S to 46°45'S, i.e., from the Maule to Aysen regions (Jiménez 2010), occupying an area of 128,278 km2 according to the International Union for Conservation of Nature (IUCN), and is mainly located in Chile. Pudu deer characteristically inhabit pristine temperate rainforests, particularly in areas of dense understory growth and native bamboo thickets (Eldridge et al. 1987; Meier and Merino 2007), but they can also be found in disturbed and secondary forest habitats (Jiménez 2010). The current conservation status of pudu deer is near threatened according to the IUCN or vulnerable based on the threatened species list of the Chilean Ministry of Environment. This conservation status is related to different threats that appear to have affected the viability of the species. Among them, local threats have been identified to be linked to the expansion of human activities, such as forest loss and fragmentation, predation by domestic dogs, competition with exotic species, and poaching activities (Miller et al. 1973; Wemmer et al. 1998; Silva-Rodríguez et al. 2010; Silva-Rodríguez and Sieving 2012; Jiménez and Ramilo 2013). For example, loss or fragmentation of forest in Osorno and Llanquihue provinces from southern Chile is a matter of particular concern for pudu deer conservation, given that between 2006 and 2013, the area planted with exotic tree species increased significantly (+20.6% and +61.3%, respectively), and a large part of this growth occurred at the expense of native forest (CONAF-UACh (Corporación Nacional Forestal-Universidad Austral de Chile) 2014). This process is more accentuated in the central valleys located between the Andean and Coastal mountain ranges (Miranda et al. 2017). Thus, the conservation of pudu deer in southern Chile is a matter of public interest that requires combined efforts in several research areas, including genetic studies to assess the effect of different threats affecting their survival in their natural environment (Silva-Rodríguez et al. 2011).

Cautín, Valdivia, Osorno, Llanquihue and Chiloé Island are provinces in Southern Chile (39–44°S latitude) that contain a significant remnant of native temperate rainforest that covers a large proportion of each province. For example, Osorno Province contains approximately 42.9% of native forest, and in Llanquihue and Chiloé provinces they cover approximately 54.5% and 68.3%, respectively (CONAF-UACh (Corporación Nacional Forestal-Universidad Austral de Chile) 2014), mostly distributed in the Coastal and Andes mountain ranges in the first two provinces (Miranda et al. 2017). Both mountain ranges present a predominance of vegetational formations comprising temperate laurifoliar rainforests that include the Valdivian, North Patagonian and Subantarctic types (Villagrán and Hinojosa 2005). The temperate rainforest of Chile encompasses the Valdivian rainforest ecoregion, which has been listed among the most endangered ecoregions of the world and has a critical conservation status (Dinerstein et al. 1995; Olson and Dinerstein 1998; Miranda et al. 2017). In addition, the Valdivian rainforest ecosystem is considered a biodiversity hotspot and, therefore, a region of high conservation priority (Ormazábal et al. 1993; Myers et al. 2000; Smith-Ramírez 2004). Previous (Vanoli 1967) and recent (Pavez-Fox and Estay 2016; Colihueque et al. 2020) records support the presence of pudu deer in the aforementioned provinces, with higher probabilities of occurrence in the Coastal mountain range, where extensive coverage of native forest persists (Colihueque et al. 2020).

Inferences based on analyses of the geographic distribution of genetic diversity indicate that the Coastal mountain range in Osorno Province and other neighboring provinces, such as Cautín, Valdivia and Llanquihue, was a lowland refuge for vertebrates during the Pleistocene glaciation in the Last Glacial Maximum (LGM), which occurred approximately 23,000–18,000 years ago (Hulton et al. 2002), since this site is part of a refugial area identified in the Chilean Coastal Range north of 42°S (Sérsic et al. 2011). Paleobotanical evidence regarding the Quaternary history of the Chilean Coastal mountain range also supports this statement (Villagrán et al. 2019). In this area, the vertebrates persisted during this geological period to later recolonize the areas where the ice sheet had retreated (Rodríguez-Serrano et al. 2006; Himes et al. 2008; Victoriano et al. 2008). Genetic analyses of different vertebrate species distributed in this area (e.g., 39°S–42°S) have shown that they present private haplotypes and moderate to high levels of genetic diversity in comparison to those located along the Andean Mountain north–south axis of Patagonia, which comprises previously glaciated areas (Sérsic et al. 2011). Moreover, it should be noted that provinces in southern Chile contain several rivers that originate in the western Andean Mountains, which flow in a relatively straight line until reaching the Pacific Ocean, such as the Cautín, San Pedro and Bueno rivers (Errazuriz et al. 2000). These rivers constitute important geographic barriers due to their constant and large flow; therefore, it is likely that they may have promoted the genetic differentiation among populations of pudu deer. Moreover, it is likely that an important genetic divergence exists between pudu deer from Osorno and those from Chiloé Island because available studies support that continental populations of this cervid present a strong genetic divergence with respect to pudu deer inhabiting this island (ϕst = 0.75, Fuentes-Hurtado et al. 2011). This level of genetic divergence has been attributed to the effect of the Pleistocene glaciations that gave rise to a particular refugia on this island, whose populations were maintained in isolation from the continent by the Chacao channel (Fuentes-Hurtado et al. 2011). In fact, the estimated divergence time among Chiloé Island and continental populations of pudu deer would have occurred at least 217,500 years ago (Fuentes-Hurtado et al. 2011).

Phylogenetic studies performed on pudu deer populations from the northernmost distribution of southern Chile (36–42°S) using mtDNA control region and cytochrome b markers (Fuentes-Hurtado et al. 2011) have indicated a significant genetic structure, with two main divergent clades corresponding to continental and Chiloé Island groups, and with lower genetic divergence among continental populations than between continental and island populations. Thus, these data suggested the existence of two different subspecies across the distribution range of the species; these subspecies have been attributed to ancient geographic isolation of the original population. In addition, a recent study that analyzed several populations that inhabit most areas of distribution of the species, including samples from the southernmost distribution (northern and southern Patagonia, 42–45°S) and Chiloé Island, confirmed this genetic structuring but revealed that both Patagonian groups are more related to the Chiloé Island group than to other continental groups based on the consistent subclusters that they formed (Cabello 2019). Despite the research progress on the phylogenetic relationships of pudu deer populations, until now, no detailed studies have addressed the genetic divergence of this cervid at a narrower geographic scale, for example, within and among provinces in southern Chile. This issue is interesting to address given that it is likely that the genetic structure of pudu deer populations of these sites could be affected by ancient paleoclimatic events or recent habitat modifications. For example, it is unknown whether the consequence of paleoclimatic events, such as the LGM, on the genetic structure of pudu deer populations differs across provinces or geographic areas. It is likely that presumably less affected populations from the Pacific coast in the Coastal Range should show a contrasting pattern of genetic variation. Supporting evidence indicates that these areas were relatively stable during this period (Villagrán and Armesto 2005; Villagrán et al. 2019); therefore, stronger signals of demographic equilibrium with increasing northern latitudes should be expected. In contrast, in areas with a pattern of stronger influences of glaciations (i.e., southernmost), which contain glaciated and recolonized areas, a predominance of derived haplotypes should be expected. Thus, by assessing the extent of gain or loss of genetic variation, the level of genetic divergence within and among populations and the identification of unique haplotypes could provide evidence to clarify the effect of paleoclimatic events on the modification of the genetic structure of pudu deer in southern Chile.

The objective of this study was to assess the level of genetic divergence and the demography of pudu deer in five provinces in southern Chile that are located at different latitudinal sites (from 39°S to 42°S) (Cautín, Valdivia, Osorno, Llanquihue and Chiloé Island) and three characteristic geographic areas within these provinces representative of distinctive longitudinal sites (Andes range, Central Valley (i.e., Longitudinal Valley) and Coastal Range), based on new and publicly available mtDNA control region (CR) and cytochrome b (Cyt b) sequences. We hypothesized that among provinces located at the extreme of the distribution, including those that inhabit an island, pudu deer populations must exhibit marked genetic divergence and reduced gene flow, given the existence of geographic barriers or the effect of Pleistocene glaciations that have modified the genetic structure of the populations. Moreover, pudu deer populations from the Coastal Range, are likely to have experienced an expansion process, since this site was a lowland refuge for vertebrates during the Pleistocene glaciations, to later recolonize the areas where the ice sheet retreated. In addition, populations located in this geographic area should exhibit a higher level of genetic diversity than populations located in other sites, for example, to those located in the Andean mountain range, due to its refugial nature. It is also possible that pudu deer population from some provinces may conform to a distinctive genetic cluster since these populations can present characteristic haplotype frequencies due to local selection pressure or population-specific haplotypes due to its relictual nature. Clarifying this issue is relevant for the sustainability of pudu deer in southern Chile, given that some populations may constitute an evolutionarily significant unit, whose identification is of major importance when conservation plans for threatened species are considered (Höglund 2009; Torres-Florez et al. 2018).

Methods

Sampling and DNA extraction

Pudu deer samples collected from Osorno (n = 6), Llanquihue (n = 8) and Chiloé (n = 6) provinces were analyzed (Table 1). Muscle tissue were obtained from individuals found dead along the highways that were donated by different institutions to our laboratory (Museo de Historia Natural, Purranque; Parque Nacional Puyehue, Osorno; Sitio Paleontológico de Pilauco, Osorno; Centro de Reproducción del Pudu, Osorno). Tissue samples were then fixed in 80% ethanol. Moreover, preserved blood samples taken from live individuals at the Wildlife Rehabilitation Center (Centro de Rehabilitación de Fauna Silvestre) of the Universidad San Sebastián, Puerto Montt, Chile, were also used. DNA was extracted from the fixed muscle or blood tissue using the phenol-chloroform method, as described in Taggart et al. (1992). Extracts were standardized at 100 ng/μL using Tris‒EDTA buffer at pH 8.0. To increase the sample size aimed at robustness of the genetic analysis, CR (n = 27) and Cyt b (n = 11) sequences previously published by Fuentes-Hurtado et al. (2011) were included (Table 1). These sequences represent different haplotypes observed for pudu deer populations in the analyzed areas, particularly, the Y02–Y25 haplotypes for CR and the A–C and E–I haplotypes for Cyt b.

Table 1.

Sampling for genetic analysis of pudu deer from southern Chile. † Coordinates data were recovery from GeoNames (http://www.geonames.org/) according to district location. ‡ Fuentes-Hurtado et al. (2011) § Unknown locality.

Site number, locality and district Province Longitudinal range Coordinates (Lat., Long.) Sampling year Specimen voucher
1. Parque Nacional Puyehue, Puyehue Osorno Andes range 40.7820°S, 72.2114°W 2014 1279ULA
2. Las Cascadas, Puerto Octay Osorno Andes range 41.0869°S, 72.6360°W 2016 1330ULA
3. Colonia Zagal, Purranque Osorno Coastal range 40.9644°S, 73.4288°W 2016 1331ULA
4. Chan Chan, Osorno Osorno Central Valley 40.7244°S, 73.0231°W 2016 1337ULA
5. Parque Nacional Puyehue, Puyehue Osorno Andes range 40.7353°S, 72.3081°W 2016 1505ULA
6. Puaucho, San Juan de la Costa Osorno Coastal range 40.6142°S, 73.3732°W 2017 1515ULA
7. Ensenada, Puerto Varas Llanquihue Andes range 41.2082°S, 72.5383°W 2018 1520ULA
8. Chinquihue, Puerto Montt Llanquihue Coastal range 41.6050°S, 73.2623°W 2019 1522ULA
9. U§, Puerto Montt Llanquihue Coastal range 41.48917°S, 72.79531°W 2018 1524ULA
10.U, Calbuco Llanquihue Coastal range 41.72311°S, 73.19511°W 2018 1525ULA
11. Guayún, Calbuco Llanquihue Coastal range 41.6992°S, 73.2445°W 2018 1526ULA
12. U, Maullín Llanquihue Coastal range 41.63242°S, 73.5083 °W 2018 1527ULA
13. U, Puerto Montt Llanquihue Coastal range 41.48917°S, 72.79531°W 2019 1529ULA
14. Ensenada, Puerto Varas Llanquihue Andes range 41.2082°S, 72.5383°W 2017 1532ULA
15. Chepu, Ancud Chiloé 42.0439°S, 73.9679°W 2019 1516ULA
16. Degañ, Quemchi Chiloé 42.1426°S, 73.4743°W 2019 1517ULA
17. Pauldeo, Ancud Chiloé 41.8991°S, 73.8893°W 2019 1518ULA
18. Chonchi, Chonchi Chiloé 42.7231°S, 73.7855°W 2019 1519ULA
19.Butalcura, Ancud Chiloé 42.2332°S, 73.7499°W 2019 1523ULA
20. Cotao, Quellón Chiloé 43.15178°S, 73.9954°W 2019 1531ULA
21. Pucón, Pucón Cautín Andes range 39.2500°S, 71.9000°W NA Previous study‡
22. Villarrica, Villarrica Cautín Andes range 39.3000°S, 72.2666°W NA Previous study‡
23. Villarrica, Villarrica Cautín Andes range 39.3000°S, 72.2666°W NA Previous study‡
24. Loncoche, Loncoche Cautín Central Valley 39.3666°S, 72.6666°W NA Previous study‡
25. Lican Ray, Villarrica Cautín Andes range 39.4666°S, 72.0833°W NA Previous study‡
26. Coñaripe, Panguipulli Valdivia Andes range 39.6166°S, 72.0500°W NA Previous study‡
27. Caleta Mehuín, San José de la Mariquina Valdivia Coastal range 39.4333°S, 73.2000°W NA Previous study‡
28. San José de la Mariquina, San José de la Mariquina Valdivia Coastal range 39.5333°S, 72.9500°W NA Previous study‡
29. Máfil, Máfil Valdivia Central Valley 39.6666°S, 72.8500°W NA Previous study‡
30. Cayumapu, Valdivia Valdivia Coastal range 39.7166°S, 73.1000°W NA Previous study‡
31. Cayumapu, Valdivia Valdivia Coastal range 39.7166°S, 73.1000°W NA Previous study‡
32. Lago Panguipulli, Panguipulli Valdivia Andes range 39.7666°S, 72.1000°W NA Previous study‡
33. Collico, Valdivia Valdivia Coastal range 39.8000°S, 73.1666°W NA Previous study‡
34. Cutipay, Valdivia Valdivia Coastal range 39.8000°S, 73.2166°W NA Previous study‡
35. Valdivia, Valdivia Valdivia NA 39.8166°S, 73.1833°W NA Previous study‡
36. Valdivia, Valdivia Valdivia NA 39.8166°S, 73.1833°W NA Previous study‡
37. Valdivia, Valdivia Valdivia NA 39.8166°S, 73.1833°W NA Previous study
38. Valdivia, Valdivia Valdivia NA 39.8166°S, 73.1833°W NA Previous study‡
39. Valdivia, Valdivia Valdivia NA 39.8166°S, 73.1833°W NA Previous study‡
40. Valdivia, Valdivia Valdivia NA 39.8166°S, 73.1833°W NA Previous study‡
41. Valdivia, Valdivia Valdivia NA 39.8166°S, 73.1833°W NA Previous study‡
42. Puyehue, Puyehue Osorno Andes range 40.7166°S, 72.3166°W NA Previous study‡
43. Lago Llanquihue, Llanquihue Llanquihue Central Valley 41.2833°S, 73.0166°W NA Previous study‡
44. Lago Llanquihue, Llanquihue Llanquihue Central Valley 41.2833°S, 73.0166°W NA Previous study‡
45. Lago Llanquihue, Llanquihue Llanquihue Central Valley 41.2833°S, 73.0166°W NA Previous study‡
46. Cucao, Chonchi Chiloé 42.6166°S, 74.1000°W NA Previous study‡
47. Cucao, Chonchi Chiloé 42.6166°S, 74.1000°W NA Previous study‡

MtDNA CR and Cyt b amplification

The MtDNA control region was amplified through polymerase chain reaction (PCR) using the primers Pudu-1F (5’ CCCCACTATCAACACCCAAA 3’) and Pudu-2R (5’ ACCATTATGGGGATGCTCAA 3’), according to Fuentes-Hurtado et al. (2011). In the case of Cyt b, PCR amplification was performed using the primers Cytb-PuduF (5’ GACCTCCCAGCTCCATCAAA 3’) and Cytb-PuduR (5’ GCTGCCCTCCGATTCATGTA 3’). These primer sets were designed using PRIMER3 software version 4.0 (Koressaar and Remm 2007) based on the DQ379309 sequence from Pudu puda published in GenBank for the Cyt b gene (Gilbert et al. 2006). For PCRs, we used a 45 µL volume of a reaction mix composed of 9 μL of Taq polymerase buffer A (1x), 0.9 μL of dNTPs (0.2 mM), 1.35 μL of MgCl2 (1.5 mM), 0.9 μL of each primer (0.2 μM), 0.18 μL of Taq DNA polymerase (0.02 U/μL) (Kapa Biosystems), 9 μL of template DNA (20 ng/μL), and 22.77 μL of DNAse/RNAse free distilled water (Gibco). The thermal profile was performed as follows: initial denaturation at 94 °C for 2 min; followed by 25 cycles of 94 °C for 45 s, 62 °C for 45 s with -0.5 °C per cycle and 72 °C for 55 s; 15 cycles of 94 °C for 45 s, 54 °C for 45 s, and 72 °C for 55 s; and a final extension step at 72 °C for 5 min. PCR products were visualized on 2% agarose gels and cleaned prior to sequencing with a QIAquick Gel Extraction Kit (Qiagen). PCR products were bidirectionally sequenced on an Applied Biosystems ABI377 automated sequencer. Sequence records were assembled from forward and reverse reads using GENEIOUS PRIME version 2019.2.3 software (Biomatters Ltd.) to obtain consensus sequences for all individuals. The sequences were deposited in GenBank under accession numbers MK905227MK905230, MZ006215MZ006230 and OM222935OM222953.

Population differentiation and genetic diversity

We used DNASP version 5.1 software to estimate genetic variation within populations based on the number of haplotypes (h), polymorphic sites (s), haplotype diversity (H), nucleotide diversity (Π) and average number of nucleotide differences (k) (Librado and Rozas 2009). Coalescent analysis based on a neutral infinite-site model and assuming a large constant population size (Hudson 1990), with theta per gene being estimated from the data, no recombination and 1,000 replicate settings, was also performed using the same software. This approach allowed us to estimate the average number and variance (95% confidence interval) of expected haplotypes (Neh) to assess whether sample sizes were sufficient to capture the genetic diversity. We estimated pairwise FST values based on the Kimura 2P distance method (Nei and Kumar 2000) and haplotype frequencies that were implemented in ARLEQUIN version 3.1 (Excoffier and Lischer 2010) to quantify the extent of genetic differentiation among populations. For pairwise FST, the significance was tested using 1,023 permutations between populations with a level of significance of α = 0.05. The heterogeneity of haplotype frequencies was assessed using the chi-square homogeneity test (Sokal and Rohlf 2009). In this test the null hypothesis of equality in haplotype frequencies within each geographic area was tested. We performed the Welch t-test for unequal variances to test for statistical significance of comparisons among study areas based on the H and Π indices (Ruxton 2006). Analysis of molecular variance (AMOVA) (Excoffier et al. 1992) implemented in the program ARLEQUIN version 3.1 (Excoffier and Lischer 2010) was used to assess the geographic pattern of population subdivision. In AMOVA the correlation of haplotype frequency was used as an F-statistic analogue at various hierarchical levels. To carry out this analysis, samples were separated into five groups according to their sampling province, which have different latitudinal locations from north to south as follows: Cautín Province (CAU, 39°0'0"S, 72°30'0"W), Valdivia Province (VAL, 39°45'0"S, 72°30'0"W), Osorno Province (OSO, 40°35'0"S, 73°10'0"W), Llanquihue Province (LLA, 41°20'0"S, 72°50'0"W) and Chiloé Island Province (CHI, 42°30'0"S, 74°0'0"W) (Fig. 1a). We selected these provinces because there was a sufficient sample size to carry out the genetic analyses in these geographic areas (from 5 to 16). In addition, to test the possible existence of genetic differentiation between populations located at different longitudinal locations, samples from different provinces, except those of Chiloé Island Province, were grouped as belonging to the Coastal Range (COA), Central Valley or Longitudinal Valley (CEN) and Andes Mountain range (AND) based on their geographic coordinates (Fig. 1b).

Figure 1. 

Sampling sites of pudu deer in southern Chile. A) Sample sites at the latitudinal range (i.e., for provinces) and B) sample sites at the longitudinal range (i.e., for the Andes range, Central Valley and Coastal Range). In A), the blue dots indicate sampling sites within each province, and in B), sampling sites identified by dark dots within specific area of the longitudinal range, with different colors for the Andes range, Central Valley and Coastal Range areas. Sampling site numbers are detailed in Table 1.

Gene flow among provinces

The number of migrants per generation (Nm) among provinces was estimated from FST values using ARLEQUIN version 3.1 (Excoffier and Lischer 2010). This estimation assumes that two populations of size N drawn from a large pool of populations exchange a fraction m of migrants each generation and that the mutation rate u is negligible compared to the migration rate m.

Population clustering

We used a clustering analysis through kMeans version 1.1 (Meirmans and Van Tienderen 2004) to split objects into groups by minimizing the within-group genetic diversity while maximizing it among-groups. For this analysis, the input file comprised a matrix of allele frequencies coded as binary data. An a priori assigned number of groups (k) from 1 to 10 was tested with 50,000 permutations. We used the Caliński–Harabasz pseudo-F statistic (Caliński and Harabasz 1974) to determine the optimal number of clusters. In the pseudo-F statistic, the optimal clustering is the one with the highest value. Simulated annealing as a clustering algorithm was used based on the Monte Carlo Markov chain (MCMC), which prevents local optima becoming stuck.

Population dynamics

We tested for patterns of demographic history with Tajima’s D (Tajima 1989) and Fu’s Fs (Fu 1997) tests of neutrality. Tajima’s D test is based on the infinite-site model without recombination and is appropriate for short DNA sequences. The significance of the D statistic is tested by generating random samples under the hypothesis of selective neutrality and population equilibrium using a coalescent simulation algorithm adapted from Hudson (1990). Significant D values can be due to factors other than selective effects, such as population expansion, bottleneck, or heterogeneity of mutation rates (Tajima 1996). Fu’s Fs test is based on an infinite-site model without recombination. This test evaluates the probability of observing a random neutral sample with a number of alleles that is equal to or smaller than the observed number of alleles. When the Fs statistic is large and negative, one can infer demographic expansion. This test has been shown to be especially sensitive to departure from population equilibrium, as in the case of population expansion. Tajima’s D and Fu’s Fs tests were run in ARLEQUIN version 3.1 (Excoffier and Lischer 2010).

Results

Analysis at the latitudinal range

CR and Cyt b haplotypes and nucleotide polymorphism

Thirty-three CR haplotypes were detected in 47 pudu deer individuals across the five provinces analyzed (Table 2). The coalescent analysis indicated that populations approached the estimated total haplotype number relatively well since, with the exception of Valdivia Province, all provinces fell within a 95% confidence interval (Table 3). Thus, this analysis indicated that most haplotypes were captured. Thirty-six of 626 nucleotide positions were variable (5.75%), including 27 parsimony informative sites and nine singleton variable sites. The number of haplotypes per province ranged from 4 to 14, with a mean of 7.2. Haplotype diversity (H) ranged from 0.64286 to 0.98333 among provinces, and nucleotide diversity (Π) ranged from 0.00575 to 0.01022 (Table 3). Haplotype diversity (H) and nucleotide diversity (Π) were significantly distinct among provinces especially among Valdivia and Llanquihue plus Chiloé provinces in the former (0.98333 versus 0.64286–0.92727, Welch´s t-test, P < 0.05) (Fig. 2a), and between Cautín and the other four provinces in the latter (0.01022 versus 0.00575–0.00797, Welch´s t-test, P < 0.05) (Fig. 2b).

Table 2.

Recorded haplotypes for pudu deer of five provinces from southern Chile based on mtDNA Control region (CR) and Cytochrome b (Cyt b) markers. Bold indicates new recorded haplotypes.

Province Coordinates (Lat., Long.) Sample size (n) No. of haplotypes Haplotypes
CR
Cautín 39°0'0"S, 72°30'0"W 5 4 Y02, Y08, Y16, Y18
Valdivia 39°45'0"S, 72°30'0"W 16 14 Y17, Y23, Y24, Y10, Y14, Y21, Y25, Y04, Y11, Y09, Y19, Y22, Y03, Y20
Osorno 40°35'0"S, 73°10'0"W 7 6 Y17, Y05, Y27, Y28, Y15, Y26
Llanquihue 41°20'0"S, 72°50'0"W 11 8 Y17, Y06, Y36, Y35, Y34, Y32, Y07, Y05
Chiloé 42°30'0"S, 74°0'0"W 8 4 Y12, Y13, Y31, Y29
Overall 47 33
Cyt b
Cautín 39°0'0"S, 72°30'0"W 3 3 A, G, H
Valdivia 39°45'0"S, 72°30'0"W 4 4 A, E, I, F
Osorno 40°35'0"S, 73°10'0"W 7 4 A, I, L, M
Llanquihue 41°20'0"S, 72°50'0"W 9 2 A, C
Chiloé 42°30'0"S, 74°0'0"W 7 4 A, B, J, K
Overall 30 12
Table 3.

Summary of genetic diversity indices for pudu deer of five provinces from southern Chile based on mtDNA Control region (CR) and Cytochrome b (Cyt b) markers. D, Tajima’s D statistic; Fs, Fu’s Fs statistic; H, haplotype diversity ± standard deviation; Hap, number of haplotypes; k, average number of pairwise nucleotide differences; n, sample sizes; Neh, average number of expected haplotypes with 95% confidence interval in parenthesis based on coalescent simulations; S, number of polymorphic sites; Π, nucleotide diversity ± standard deviation. Within a row, means followed by different letters are significantly different. *P < 0.05, **P < 0.01.

Province n Hap S H Π k Neh D Fs
CR
Cautín 5 4 14 0.90000 ± 0.161a 0.01022 ± 0.00208a 6.40000 3.947 (2–4) -0.34706 0.88316
Valdivia 16 14 18 0.98333 ± 0.028a 0.00797 ± 0.00067b* 4.99167 7.783 (4–11) -0.31815 -7.75335**
Osorno 7 6 9 0.95238 ± 0.096a 0.00624 ± 0.00105b* 3.90476 4.223 (2–6) 0.33464 -1.64228
Llanquihue 11 8 12 0.92727 ± 0.066b* 0.00575 ± 0.00119b* 3.60000 5.794 (3–9) -0.52679 -3.93267*
Chiloé 8 4 16 0.64286 ± 0.184b* 0.00668 ± 0.00387b* 4.17857 5.421 (3–8) -1.65898* 1.75610
Overall 47 33 36 0.97600 ± 0.011 0.01126 ± 0.00104 7.04903 16.034 (11–22) -0.50327 -2.13781
Cyt b
Cautín 3 3 2 1.000 ± 0.272a 0.00179 ± 0.00060a 1.333 1.960 (1–3) - -
Valdivia 4 4 5 1.000 ± 0.177a 0.00336 ± 0.00103a 2.500 2.774 (1–4) -0.79684 -1.514
Osorno 7 4 4 0.810 ± 0.130a 0.00192 ± 0.00055a 1.429 3.206 (1–5) -0.59756 -0.780
Llanquihue 9 2 1 0.222 ± 0.166b** 0.00030 ± 0.00022b** 0.222 1.803 (1–4) -1.08823 -0.263
Chiloé 7 4 6 0.810 ± 0.130a 0.00256 ± 0.00102a 1.905 3.718 (2–6) -1.12898 -0.226
Overall 30 12 14 0.798 ± 0.068 0.00315 ± 0.00052 2.338 8.421 (4–13) -1.12937 -4.562
Figure 2. 

Measures of genetic diversity based on CR sequences recorded in this study for pudu deer from southern Chile. A) and B) Genetic diversity observed in Cautín (CAU), Valdivia (VAL), Osorno (OSO), Llanquihue (LLA) and Chiloé Island (CHI) provinces and C) and D) genetic diversity recorded in the Andes range (AND), Central Valley (CEN) and Coastal Range (COA). Columns and bar errors represent the means and standard deviation of the means, respectively. Significant differences between means were calculated according to Welch´s t-test. *P < 0.05, NS = Non significant difference.

Ten new haplotypes for CR were observed, which were named Y26–Y29, Y31–Y32 and Y34–Y36 (Table 2). In addition, two haplotypes were shared among provinces, represented by Y17 among Llanquihue, Valdivia and Osorno provinces and Y05 between Llanquihue and Osorno provinces (Fig. 3a). The remaining haplotypes were unique for each province, composing a ratio of unique haplotypes relative to the total number of haplotypes for each province of 1.0 (4/4) for Cautín, 0.92 (13/14) for Valdivia, 0.66 (4/6) for Osorno, 0.75 (6/8) for Llanquihue and 1.0 (4/4) for Chiloé Island. Therefore, the unique haplotypes consisted of a large proportion of all haplotypes recorded within each province. Haplotype frequencies were heterogeneous for Cautín (chi-square test: χ2 = 12, df = 3, P < 0.001), Chiloé (chi-square test: χ2 = 75, df = 3, P < 0.001) and Llanquihue (chi-square test: χ2 = 37.24, df = 7, P < 0.001) provinces. In fact, Y08 was highly represented in Cautín Province (40%), as were Y13 in Chiloé (62.5%) and Y17 in Llanquihue (27.3%) (Fig. 4a).

Figure 3. 

Venn diagram of haplotypes of the CR sequence found in pudu deer from southern Chile. A) Haplotypes observed in Cautín, Valdivia, Osorno, Llanquihue and Chiloé Island provinces and B) haplotypes recorded in the Andes range, Central Valley and Coastal Range.

Figure 4. 

Haplotype frequencies for the CR sequence found in pudu deer from southern Chile. A) Haplotype frequency observed in Cautín (CAU), Valdivia (VAL), Osorno (OSO), Llanquihue (LLA) and Chiloé Island (CHI)provinces and B) haplotype frequency recorded in the Andes range (AND), Central Valley (CEN) and Coastal Range (COA). Asterisks indicate a heterogeneous frequency based on the chi-square test (P < 0.001). Haplotypes are shown by different color codes as is indicated at the right side of figure.

Twelve Cyt b haplotypes were detected in 30 pudu deer individuals across the five provinces analyzed (Table 2). The coalescent analysis indicated that populations approached the estimated total haplotype number relatively well, since all provinces fell within a 95% confidence interval (Table 3). Fourteen of 743 nucleotide positions were variable (1.88%), including six parsimony informative sites and eight singleton variable sites. Thus, the Cyt b sequence was approximately 70% less variable than the CR sequence. The number of haplotypes per province ranged from 2 to 4, with a mean of 3.4. Haplotype diversity (H) ranged from 0.222 to 1.000 among provinces, and nucleotide diversity (Π) ranged from 0.00030 to 0.00336 (Table 3). Haplotype diversity (H) and nucleotide diversity (Π) were significantly distinct among some provinces, especially among Cautín and Llanquihue in the former (1.000 versus 0.222, Welch´s t-test, P < 0.05) and between Valdivia and Llanquihue in the latter (0.00336 versus 0.00030, Welch´s t-test, P < 0.05).

Four new haplotypes for Cyt b were observed, which were named J–M (Table 2). In addition, only two haplotypes were shared among provinces, with the A haplotype shared among all provinces and the I haplotype between Valdivia and Osorno provinces (Fig. 5a). The remaining haplotypes were unique for each province, composing a ratio of unique haplotypes relative to the total number of haplotypes for each province of 0.66 (2/3) for Cautín, 0.5 (2/4) for Valdivia, 0.5 (2/4) for Osorno, 0.5 (1/2) for Llanquihue and 0.75 (3/4) for Chiloé Island. Therefore, the unique haplotypes involved a large proportion of all haplotypes recorded within each province. Haplotype frequencies were heterogeneous for Chiloé Island (chi-square test: χ2 = 22.49, df = 3, P < 0.001), Llanquihue (chi-square test: χ2 = 60.52, df = 1, P < 0.001) and Osorno (chi-square test: χ2 = 22.49, df = 3, P < 0.001) provinces. In fact, the K haplotype was highly represented in Chiloé Island Province (42.9%), as were the A haplotype in Llanquihue (88.95%) and the I haplotype in Osorno (42.9%) (Fig. 6a).

Figure 5. 

Venn diagram of haplotypes of the Cyt b sequence found in pudu deer from southern Chile. A) Haplotypes observed in Cautín, Valdivia, Osorno, Llanquihue and Chiloé Island provinces and B) haplotypes recorded in the Andes range and Coastal Range.

Figure 6. 

Haplotype frequencies for the Cyt b sequence found in pudu deer from southern Chile. A) Haplotype frequency observed in Cautín (CAU), Valdivia (VAL), Osorno (OSO), Llanquihue (LLA) and Chiloé Island (CHI) provinces and B) haplotype frequency recorded in the Andes range (AND), Central Valley (CEN) and Coastal range (COA). Asterisks indicate a heterogeneous frequency based on the chi-square test (P < 0.001). Haplotypes are shown by different color codes as is indicated at the right side of figure.

Population genetic structure

Global AMOVA results as a weighted average over loci, conducted by partitioning variation among and within populations, revealed that CR haplotype variation was attributed mainly to within-population differences (58.1%) but also to differences among populations (41.9%). In fact, the FST index reached a value of 0.41905 (P < 0.001). Population pairwise FST values ranged from 0.03617 to 0.68651, with seven out of ten pairwise comparisons showing significant differences (P < 0.05) (Table 6). This pattern of genetic differentiation was clearly observed between Cautín and Valdivia and Llanquihue (FST = 0.15385–0.16810) and among Chiloé Island compared to the other four continental provinces (FST = 0.59004–0.68651). The clustering analysis based on allelic frequencies of CR revealed the existence of either two (K = 2, pseudo-F = 29.620) or three (K = 3, pseudo-F = 24.695) groups as the optimal number of clusters. In both simulations, individuals from Cautín, Valdivia, Osorno and Llanquihue provinces overlapped into one or two clusters, in contrast to those of Chiloé Island Province, which mostly formed a distinctive single group (Fig. 7a).

Figure 7. 

Results of the clustering analysis based on the CR sequence for pudu deer from southern Chile. A) Clusters observed for individuals from Cautín (CAU), Valdivia (VAL), Osorno (OSO), Llanquihue (LLA) and Chiloé Island (CHI) provinces and B) clusters recorded for individuals from the Andes range (AND), Central Valley (CEN) and Coastal range (COA). K values correspond to the optimal number of clusters according to the Caliński–Harabasz pseudo-F statistic.

Gene flow among provinces

Pairwise estimate of Nm for CR haploytpe data ranged from 0.22832 to 13.32296. The number of migrants per generation was high among Osorno, Valdivia and Llanquihue provinces (Nm = 5.18197–13.32296), medium between Cautín and Osorno plus Valdivia and Llanquihue (Nm = 2.47437–3.66488), and low for Chiloé Island compared to the other four continental provinces (Nm = 0.2283–0.34740) (Table 6). Thus, patterns in pairwise Nm values were consistent with population structure based on AMOVA and clustering analysis.

Evolutionary neutrality and population expansion

Tajima’s D estimate for CR was negative and significant in one province that corresponds to Chiloé (D = -1.65898, P < 0.05). In contrast, Tajima’s D statistic of other provinces was consistent with population equilibrium (Table 3). Fu’s Fs was negative and significant for Valdivia (Fs = -7.75335, P < 0.001) and Llanquihue (Fs = -3.93267, P < 0.05) provinces, which occurs when an excess of rare haplotypes is present and suggests that either population expansion or genetic hitchhiking has occurred. These results are consistent with deviation from neutrality due to either selection or population expansion in some provinces. Tajima’s D and Fu’s Fs estimations based on Cyt b were negative and not significant for all provinces (Table 3). These results are consistent with population equilibrium.

Analysis at the longitudinal range

Analysis at the longitudinal range based on CR revealed a number of haplotypes per area ranging from 6 to 11 in 32 analyzed individuals (Table 4), with a mean of 7.3. Haplotype diversity (H) ranged from 0.955 to 1.000 among areas and nucleotide diversity (Π) ranged from 0.006513 to 0.008307 (Table 5). Only nucleotide diversity (Π) was significantly different among areas, with the Andes range and Central Valley showing higher diversity than the Coastal Range (Welch´s t-test, P < 0.05) (Fig. 2d). Only four out of 22 haplotypes were shared among the Andes range, Central Valley and Coastal Range. In fact, the Y27 haplotype was shared among the Andes range and Central Valley, as were the Y17 haplotype between the Andes range and Coastal Range, and the Y05 and Y10 haplotypes among the Central Valley and Coastal range (Fig. 3b). Haplotype frequency were heterogeneous only for Coastal range (chi-square test: χ2 = 23.49, df = 10, P < 0.01). In fact, the Y17 and Y32 haplotypes showed high frequencies in this area (21.4 and 14.3%, respectively) (Fig. 4b). In the case of Cyt b, haplotypes per area ranging from 4 to 6 in 22 analyzed individuals were observed (Table 4). Haplotype diversity was significantly higher in the Andes range than in the Coastal Range (0.803 vs. 0.533, Welch´s t-test, P < 0.05) (Table 5). Moreover, the haplotype frequencies were heterogeneous for Andes (chi-square test: χ2 = 58.37, df = 5, P < 0.001) and Coastal (chi-square test: χ2 = 108, df = 3, P < 0.001) ranges. In fact, the A haplotype showed a high frequency in both provinces (41.7 and 70%, respectively) (Fig. 6b). The clustering analysis based on allelic frequencies of CR revealed the existence of either two (K = 2, pseudo-F = 8.802) or three (K = 3, pseudo-F = 8.011) groups as the optimal number of clusters. In both simulations, individuals from the Andes range, Central Valley and Coastal Range overlapped into two or three clusters (Fig. 7b).

Table 4.

Recorded haplotypes for pudu deer of three longitudinal ranges from southern Chile based on mtDNA Control region (CR) and Cytochrome b (Cyt b) markers. Bold indicates new recorded haplotypes.

Longitudinal range Coordinates (Lat., Long.) Sample size (n) No. of haplotypes Haplotypes
CR
Andes range 39°-41°S, 72°7'W 12 9 Y08, Y09, Y17, Y02, Y18, Y15, Y26, Y27, Y28
Coastal range 39°-41°S, 73°30'W 14 11 Y17, Y05, Y10, Y14, Y11, Y03, Y04, Y36, Y32, Y35, Y34
Central Valley 39°-41°S, 72°57'W 6 6 Y27, Y05, Y10, Y07,Y16, Y06
Overall 32 22
Cyt b
Andes range 39°-41°S, 72°7'W 12 6 A, I, C, F, H, L
Coastal range 39°-41°S, 73°30'W 10 4 A, I, E, M
Overall 22 8
Table 5.

Summary of genetic diversity indices for pudu deer of three longitudinal ranges from southern Chile based on mtDNA Control region (CR) and Cytochrome b (Cyt b) markers. D, Tajima’s D statistic; Fs, Fu’s Fs statistic; H, haplotype diversity ± standard deviation; Hap, number of haplotypes; k, average number of pairwise nucleotide differences; n, sample sizes; Neh, average number of expected haplotypes with 95% confidence interval in parenthesis based on coalescent simulations; S, number of polymorphic sites; Π, nucleotide diversity ± standard deviation. Within a row, means followed by different letters are significantly different. *P < 0.05. **P < 0.01.

Longitudinal range n Hap S H Π k Neh D Fs
CR
Andes range 12 9 15 0.955 ± 0.047a 0.007890 ± 0.00101a 4.939 6.413 (3–9) -0.02394 -2.24957
Coastal range 14 11 17 0.956 ± 0.045a 0.006513 ± 0.00110b* 4.077 7.178 (4–11) -0.97905 -6.64752**
Central Valley 6 6 13 1.000 ± 0.096a 0.008307 ± 0.00103a 5.200 4.366 (2–6) -0.52988 -2.08376
Overall 32 22 26 0.968 ± 0.018 0.00741 ± 0.00067 4.639 11.866 (7–17) -0.51095 -3.66028
Cyt b
Andes range 12 6 6 0.803 ± 0.096a** 0.00167 ± 0.00040a 1.242 4.270 (1–7) -1.42890 -2.666
Coastal range 10 4 5 0.533 ± 0.180b** 0.00135 ± 0.00064a 1.000 3.793 (2–6) -1.74110* -0.876
Overall 22 8 9 0.688 ± 0.099 0.00151 ± 0.00039 1.121 6.147 (3–10) -1.82961* -4.230
Table 6.

Pairwise FST values (below the diagonal) and number of migrants per generation (Nm) (above the diagonal) among pudu deer of five provinces from southern Chile based on mtDNA control region. Average F–Statistics over all loci: FST = 0.41905. *P < 0.05.

Cautín Valdivia Osorno Llanquihue Chiloé
Cautín 2.47437 3.66488 2.75000 0.34740
Valdivia 0.16810 * 13.32296 5.18197 0.28967
Osorno 0.12005 0.03617 9.71845 0.27240
Llanquihue 0.15385 * 0.08800 * 0.04893 0.22832
Chiloé 0.59004 * 0.63318 * 0.64733 * 0.68651 *

Global AMOVA results as a weighted average over loci revealed that CR haplotype variation was attributed mostly to within area differences (98.6%) and scarcely to differences between areas (1.3%). In fact, the FST index was not significant and reached a value of only 0.01360 (P > 0.05). Pairwise FST values based on CR haplotype frequencies ranged from -0.01741 to 0.03434, with no significant differences observed in any comparison (P > 0.05). Pairwise estimate of Nm for CR haplotype data ranged from 14.05932 to infinity, revealing the occurrence of an important number of migrants per generation among the Andes range, Central Valley and Coastal Range. Tajima’s D was negative in all areas, but it was not significant (D = from -0.02394 to -0.97905, P > 0.05) (Table 5). Thus, these results were consistent with population equilibrium. However, Fu’s Fs was negative and highly significant for the Coastal Range (Fs = -6.64752, P < 0.001), which suggests that in this area population expansion may have occurred. Tajima’s D estimated based on Cyt b showed a similar result as CR given that the D value was negative and significant for the Coastal Range (D = -1.74110, P < 0.05) (Table 5).

Discussion

Analysis of the CR and Cyt b mitochondrial sequences from pudu deer populations of five provinces from southern Chile revealed different levels of genetic variation, haplotype heterogeneity, a large proportion of unique haplotypes, and strong genetic structuring and clustering of individuals in some provinces. The most likely explanation for our results is that pudu deer from the study area have a low to moderate level of gene flow among populations, which suggests the occurrence of reduced panmixia across the collection range. This genetic pattern is in accordance with our predictions of marked genetic divergence among pudu deer populations, especially those located at the extremes of the distribution or those that inhabit islands, due to the existence of geographic barriers that may have promoted genetic divergence. In fact, individuals from Chiloé Island that are separated from the continental populations by the Chacao channel formed clear single groups in the clustering analysis in concordance with the existence of such a geographic barrier, which should have reduced gene flow with the other continental analyzed populations. In this regard, the separation among continental and island populations of pudu deer, estimated at least 217,500 year ago (Fuentes-Hurtado et al. 2011), would have been key event that has promoted genetic divergence of pudu deer that inhabit Chiloé Island. Pudu deer of other sites, such as those in Cautín Province, also showed a marked divergence with respect to other analyzed provinces (FST = 0.15385–0.59004) and a reduced number of migrants per generation (Nm = 0.34740–3.66488). Since this province is located in the northernmost part of the range of the species distribution, it is likely that geographic barriers, such as rivers, and/or distance constrain gene flow and, therefore, facilitate an increase in genetic divergence. All haplotypes of the CR marker found in Cautín Province were private, a finding that provides further support for a divergent evolutionary process for pudu deer located at this site. Although further analysis is required, genetic differentiation between provinces recorded for pudu deer in our study could be also related to their habitat loss and fragmentation in southern Chile as a consequence of the expansion of human activities (Echeverría et al. 2006, 2007; Miranda et al. 2017). This could be so given that habitat fragmentation contributes to genetic differentiation between populations, since fragmented populations are expected to see more differentiation because migration between them is impaired. Of note is that our previous interpretations should be taken with caution because, as discussed below, the modifications of the genetic structure of pudu deer could also be a consequence of Pleistocene glaciations that affected the landscape in southern Chile in the past. Thus, the modifications of the genetic structure of pudu deer are likely to be a more complex process in which different factors have intervened.

In addition, our results also suggest that some provinces, particularly Valdivia, Llanquihue and Chiloé Island, underwent population expansion given that they showed significant negative Tajima’s D and Fu's Fs values. This expansion process could explain the shared haplotypes among neighboring provinces, such as Valdivia, Osorno and Llanquihue, where the Y05 and Y17 CR haplotypes were highly frequent. In the case of the Cyt b marker, although it did not present any significant population disequilibrium, the occurrence of shared haplotypes for this marker among all provinces analyzed, particularly the A and I haplotypes, suggests that population expansion for pudu deer could be a more extensive process. Of note is that the analysis at the longitudinal range supports that populations located in the Coastal Range would be an important source of individuals for colonization, given that both markers showed significant population disequilibrium in this geographic area. These results are in accordance with the statement that this area was a lowland refuge for vertebrates during the Pleistocene glaciations (Sérsic et al. 2011), where the vertebrates persisted during this geological period to later recolonize the areas where the ice sheet had retreated (Rodríguez-Serrano et al. 2006; Himes et al. 2008; Victoriano et al. 2008). In addition, we also observed that pudu deer of Coastal Range presented lower or similar genetic diversity, measured as haplotype diversity and nucleotide diversity, than those located in the Andes range. It has been observed that vertebrate populations distributed in the Andes range (e.g., 39°S–42°S) tend to exhibit lower genetic diversity, while the opposite is observed for populations distributed in the Coastal Range (Sérsic et al. 2011). This pattern of genetic variation seems to be related to the recolonization process that occurred in the Andes Mountains after the Pleistocene glaciations, where founder effects or selection probably changed the genetic structure of populations. The modifications of the genetic structure as a result of Pleistocene glaciations in the Andes range have been well documented in other vertebrates that inhabit this geographic area, such as the lizard Liolaemus tenuis from 38°S to the southern limit of its distribution (Victoriano et al. 2008), the marsupial Dromiciops gliroides located between 40 and 43°S (Himes et al. 2008) and the lizard Liolaemus pictus (Vera-Escalona et al. 2012). In fact, for L. tenuis, there is a general pattern, whereby populations located in the mesomorphic zone have a higher genetic diversity than those inhabiting in the periglacial and glacial zones. In the case of D. gliroides, southern Clade C, strongly affected by glacial recession, showed lower nucleotide diversity than northern Clade A. Similarly, L. pictus in the northern group (37°S–39°S) was heterogeneous, with well differentiated and deeply divergent internal subgroups, whereas the southern group (39°S–41°S) was more homogeneous. Since Andes range populations of pudu deer analyzed in this work presented high genetic variation, an alternative explanation is required. Among them a mixing origin of this population cannot be ruled out, which implies the occurrence of various recolonization processes with different source populations. This hypothesis is not unlikely because the haplotype composition of CR observed in the Andes range included haplotypes proper to Cautín (Y02, Y08), Valdivia (Y09) and Osorno (Y15, Y26, Y28) provinces. The multiple sources origin hypothesis is not new for vertebrates that inhabit the western Andes range in southern Chile. For example, Vidal et al. (2012) by using the Cyt b marker, found evidence that the lizard Liolemus pictus that inhabits this area may have experienced such a process during Pleistocene glaciations, given their geographically differentiated gene pools. For this lizard, postglacial colonization could have occurred from sources located in ice-free sectors of the Longitudinal Valley adjacent to the Antillanca locality and not from a single western refuge.

A number of studies have indicated that Pleistocene glaciations modified the landscape in southern Chile, which affected the species distribution and created habitat fragmentation of areas that species used as refuges (Lessa et al. 2010; Vera-Escalona et al. 2012; Vidal et al. 2012; Marín et al. 2013). From a population genetic point of view, this glaciation process may have promoted the genetic divergence among populations, with some characteristics of this process being the existence of populations with private haplotypes, haplotype heterogeneity or a reduction in genetic variation as a consequence of a decrease in population effective sizes. Thus, the genetic differentiation between pudu deer from Chilóe Island and surrounding provinces in continental areas fit well with this divergence pattern. In fact, we determined that pudu deer from Chilóe Island do not share CR haplotypes with any continental populations. Moreover, they showed significantly lower haplotype diversity than other analyzed provinces (e.g., with Valdivia Province), which may reflect past processes of population size reduction. Chilóe Island Province also showed significant haplotype heterogeneity, a process that may be related to reduction in population size, since this provokes genetic drift, causing genetic differentiation in populations. In other vertebrates that inhabit this island a similar divergence pattern has been observed. For example, the lizard L. pictus conforms to a well-supported subclade and is composed exclusively of haplotypes proper from Chiloé Island, although the genetic distances between haplotypes from this island and the surrounding mainland populations are low compared to those from the northern clade (Vera-Escalona et al. 2012). The absence of shared haplotypes between populations from the mainland and Chiloé Island, which give rise to divergent haplogroups in L. pictus, might reflect ancient isolation, much older than the LGM (Vidal et al. 2012). Of note is that further niche simulation analyses support that this pattern of divergence in L. pictus may be related mostly to the LGM (Vera-Escalona et al. 2012), given that at that age a greater suitability of the habitat in the extreme north of the island of Chiloé may have existed, which probably acted as a refuge for different vertebrates.

Previous reports indicate that pudu deer from Chiloé Island have marked genetic divergence with respect to continental populations of this species, especially those from the northernmost distribution (Fuentes-Hurtado et al. 2011; Cabello 2019). For example, Fuentes-Hurtado et al. (2011) using the CR marker, found that there were higher genetic distances among continental populations and Chiloé Island populations (2.2%) than between continental pudu deer populations (0.9%). Our clustering analysis results are in accordance with this divergence pattern, since almost all pudu deer individuals from Chiloé Island formed a consistent single cluster when using a clustering algorithm based on a Monte Carlo Markov chain. Moreover, pudu deer from this site showed a marked divergence with respect to those inhabiting the other analyzed provinces (FST = 0.59004–0.68651). Therefore, these results support that they may constitute an evolutionarily significant unit and could be treated as a separate management unit. This finding has conservation implications because pudu deer in Chile is classified as a threatened species according to the Chilean Ministry of Environment, and the identification of an evolutionarily significant unit is of major importance when conservation plans for threatened species are considered (Torres-Florez et al. 2018). Currently, this issue is of public interest because in Chiloé Island Province, pudu deer faces several local threats linked to the expansion of human activities, such as forest loss and fragmentation, predation by domestic dogs and competition with exotic species (Miller et al. 1973; Wemmer et al. 1998; Silva-Rodríguez et al. 2010, 2011; Silva-Rodríguez and Sieving 2012; Jiménez and Ramilo 2013). Pudu deer from the other analyzed provinces, could be also treated as different management units because they presented a large ratio of unique haplotypes (e.g., from 0.66 to 1.0 for CR), which may reflect the existence of local adaptation, a process that may arise in response to localized natural selection pressures. In this regard, we hope that the genetic divergence data presented in this study will help to improve the conservation status of pudu deer in southern Chile.

Conclusion

Genetic and demographic analyses of the mtDNA control region and cytochrome b sequences carried out for pudu deer of five provinces located at different latitudinal sites in southern Chile, indicate the existence of different levels of genetic variation, haplotype heterogeneity, a large proportion of unique haplotypes, strong genetic structuring, clustering of individuals in a province (i.e., Chiloé Island) and a signature of population expansion in some provinces (i.e., Valdivia, Llanquihue and Chiloé Island). These results indicate that pudu deer have a low to moderate level of gene flow, especially among distant continental provinces or between continental provinces and Chiloé Island, which suggests the occurrence of reduced panmixia. This genetic pattern could be related to the existence of geographic barriers and/or paleoclimatic events that may have promoted genetic divergence, such as the Pleistocene glaciations that modified the landscape in southern Chile. Analysis at the longitudinal range within the studied provinces revealed haplotype variation attributed mostly to within areas differences and the existence of a population expansion process in the Coastal mountain range. The last finding supports the statement that this area was a lowland refuge for vertebrates during the Pleistocene glaciations, where the vertebrates persisted during this geological period to later recolonize the areas where the ice sheet retreated. Since pudu deer in Chile is a threatened species, our results have conservation implications because populations of some analyzed geographic areas may constitute significant evolutionary units, which is of major importance when conservation plans are considered.

Acknowledgements

We would like to thank the following people for providing samples of Pudu deer from Osorno Province: Carlos Oyarzún, Museo de Historia Natural de Purranque; Carlos Hernández, Parque Nacional Puyehue; Mario Prussing, Centro de Reproducción del Pudu, Osorno; and Hugo Oyarzo, Sitio Paleontológico de Pilauco, Osorno. This study was supported by grant R25-19 of the Dirección de Investigación of the Universidad de Los Lagos. The mapping support by Alicia Vásquez Parraguez is also appreciated.

References

  • Cabello J (2019) Genética y conservación de tres especies de fauna Chilena amenazadas, distribuidas en la ecorregión bosque Valdiviano: Pudu puda, Lycalopex fulvipes y Rhinoderma darwinii. PhD Thesis, Universidad de Castilla-La Mancha, Ciudal Real, España.
  • CONAF-UACh [Corporación Nacional Forestal-Universidad Austral de Chile] (2014) Monitoreo de cambios, corrección cartográfica y actualización del catastro de recursos vegetacionales nativos de la Región de Los Lagos. Informe final. Laboratorio de Geomática, Instituto de Manejo de Bosques y Sociedad, Universidad Austral de Chile, Valdivia, Chile, 54 pp.
  • Dinerstein E, Olson D, Graham D, Webster A, Primm S, Bookbinder M, Ledec G (1995) A conservation assessment of the terrestrial ecore-gions of Latin America and the Caribbean. Report number 14996. The World Bank, Washington DC. https://doi.org/10.1596/0-8213-3295-3
  • Echeverría C, Coomes D, Salas J, Rey-Benayas JM, Lara A, Newton A (2006) Rapid deforestation and fragmentation of Chilean Temperate Forests. Biological Conservation 130(4): 481–494. https://doi.org/10.1016/j.biocon.2006.01.017
  • Echeverría C, Newton AC, Lara A, Benayas JMR, Coomes DA (2007) Impacts of forest fragmentation on species composition and forest structure in the temperate landscape of southern Chile. Global Ecology and Biogeography 16(4): 426–439. https://doi.org/10.1111/j.1466-8238.2007.00311.x
  • Eldridge WD, MacNamara MM, Pacheco NV (1987) Activity patterns and habitat utilization of pudus (Pudu puda) in south-central Chile. In: Wemmer C (Ed.) Biology and management of the Cervidae. Smithsonian Institution Press, Washington DC, 352–370.
  • Errazuriz A, Cereceda P, González J, González M, Henriquez M, Rioseco R (2000) Manual de geografía de Chile. 3rd edn. Andrés Bello, Santiago de Chile, 443 pp.
  • Excoffier L, Lischer HL (2010) Arlequin suite ver 3.5: A new series of programs to perform population genetics analyses under Linux and Windows. Molecular Ecology Resources 10(3): 564–567. https://doi.org/10.1111/j.1755-0998.2010.02847.x
  • Excoffier L, Smouse PE, Quattro JM (1992) Analysis of molecular variance inferred from metric distances among DNA haplotypes: Application to human mitochondrial DNA restriction data. Genetics 131(2): 479–491. https://doi.org/10.1093/genetics/131.2.479
  • Fuentes-Hurtado M, Marín JC, González-Acuña D, Verdugo C, Vidal F, Vianna JA (2011) Molecular divergence between insular and continental Pudu deer (Pudu puda) populations in the Chilean Patagonia. Studies on Neotropical Fauna and Environment 46(1): 23–33. https://doi.org/10.1080/01650521.2010.537906
  • Gilbert C, Ropiquet A, Hassanin A (2006) Mitochondrial and nuclear phylogenies of Cervidae (Mammalia, Ruminantia): Systematics, morphology, and biogeography. Molecular Phylogenetics and Evolution 40(1): 101–117. https://doi.org/10.1016/j.ympev.2006.02.017
  • Himes CMT, Gallardo MH, Kenagy GJ (2008) Historical biogeography and post-glacial recolonization of South American temperate rain forest by the relictual marsupial Dromiciops gliroides. Journal of Biogeography 35(8): 1415–1424. https://doi.org/10.1111/j.1365-2699.2008.01895.x
  • Hudson RR (1990) Gene genealogies and the coalescent process. In: Futuyma D, Antonovics J (Eds) Oxford Surveys in Evolutionary Biology. Oxford University Press, New York, 1–44. https://books.google.cl/books?id=VuE9jgEACAAJ
  • Hulton NRJ, Purves RS, McCulloch RD, Sugden DE, Bentley MJ (2002) The Last Glacial Maximum and deglaciation in southern South America. Quaternary Science Reviews 21(1–3): 233–241. https://doi.org/10.1016/S0277-3791(01)00103-2
  • Jiménez J (2010) The southern pudu (Pudu puda). In: González S, Barbanti J (Eds) Neotropical cervidology: biology and medicine of Latin American deer. Funep/IUCN, Jaboticabal, Brazil, 140–150.
  • Marín JC, Varas V, Vila AR, López R, Orozco-terWengel P, Corti P (2013) Refugia in Patagonian fjords and the eastern Andes during the Last Glacial Maximum revealed by huemul (Hippocamelus bisulcus) phylogeographical patterns and genetic diversity. Journal of Biogeography 40(12): 2285–2298. https://doi.org/10.1111/jbi.12161
  • Meier D, Merino ML (2007) Distribution and habitat features of southern pudu (Pudu puda Molina, 1782) in Argentina. Mammalian Biology - Zeitschrift fur Saugetierkunde 72: 204–212. https://doi.org/10.1016/j.mambio.2006.08.007
  • Miller S, Rotmann J, Taber RD (1973) Dwindling and endangered ungulates of Chile. Vicugna, Lama, Hippocamelus and Pudu. Transactions of the thirty-eighth North American Wildlife and Natural Resources Conference, March 18, 19, 20, 21, 1973. Washington, Wildlife Management Institute, 55–68.
  • Miranda A, Altamirano A, Cayuela L, Lara A, González M (2017) Native forest loss in the Chilean biodiversity hotspot: Revealing the evidence. Regional Environmental Change 17(1): 285–297. https://doi.org/10.1007/s10113-016-1010-7
  • Myers N, Mittermeier RA, Mittermeier CG, da Fonseca GAB, Kent J (2000) Biodiversity hotspots for conservation priorities. Nature 403(6772): 853–858. https://doi.org/10.1038/35002501
  • Nei M, Kumar S (2000) Molecular evolution and phylogenetics. Oxford University Press, New York, 352 pp.
  • Pavez-Fox M, Estay SA (2016) Correspondence between the habitat of the threatened pudú (Cervidae) and the national protected-area system of Chile. BMC Ecology 16(1): е1. https://doi.org/10.1186/s12898-015-0055-7
  • Rodríguez-Serrano E, Cancino RA, Palma RE (2006) Molecular Phylogeography of Abrothrix olivaceus (Rodentia: Sigmodontinae) in Chile. Journal of Mammalogy 87(5): 971–980. https://doi.org/10.1644/05-MAMM-A-393R2.1
  • Ruxton GD (2006) The unequal variance t-test is an underused alternative to Student’s t-test and the Mann–Whitney U test. Behavioral Ecology 17(4): 688–690. https://doi.org/10.1093/beheco/ark016
  • Sérsic AN, Cosacov A, Cocucci AA, Johson LA, Pozner R, Avila LJ, Sites Jr JW, Morando M (2011) Emerging phylogeographical patterns of plants and terrestrial vertebrates from Patagonia. Biological Journal of the Linnean Society. Linnean Society of London 103(2): 475–494. https://doi.org/10.1111/j.1095-8312.2011.01656.x
  • Silva-Rodríguez EA, Verdugo C, Aleuy OA, Sanderson JG, Ortega-Solís GR, Osorio-Zúñiga F, González-Acuña D (2010) Evaluating mortality sources for the Vulnerable pudu Pudu puda in Chile: Implications for the conservation of a threatened deer. Oryx 44(01): 97–103. https://doi.org/10.1016/j.biocon.2012.03.008
  • Silva-Rodríguez EA, Aleuy OA, Fuentes-Hurtado M, Vianna JA, Vidal F, Jiménez JE (2011) Priorities for the conservation of the pudu (Pudu puda) in southern South America. Animal Production Science 51(4): 375–377. https://doi.org/10.1071/AN10286
  • Sokal R, Rohlf FJ (2009) Introduction to biostatistics. Dover Publications Inc.
  • Torres-Florez JP, Johnson WE, Nery MF, Eizirik E, Oliveira-Miranda MA, Galetti Jr PM (2018) The coming of age of conservation genetics in Latin America: What has been achieved and what needs to be done. Conservation Genetics 19(1): 1–15. https://doi.org/10.1007/s10592-017-1006-y
  • Vanoli T (1967) Beobachtungen an pudus, Mazama pudu (Molina, 1782). Säugetierkundliche Mitteilungen 15: 155–165.
  • Vera-Escalona I, D’Elía G, Gouin N, Fontanella FM, Muñoz-Mendoza C, Sites Jr JW, Victoriano PF (2012) Lizards on ice: Evidence for multiple refugia in Liolaemus pictus (Liolaemidae) during the Last Glacial Maximum in the southern Andean beech forests. PLoS ONE 7(11): e48358. https://doi.org/10.1371/journal.pone.0048358
  • Victoriano PF, Ortiz JC, Benavides E, Adams BJ, Sites Jr JW (2008) Comparative phylogeography of codistributed species of Chilean Liolaemus (Squamata: Tropiduridae) from the central-southern Andean range. Molecular Ecology 17(10): 2397–2416. https://doi.org/10.1111/j.1365-294X.2008.03741.x
  • Villagrán C, Armesto J (2005) Fitogeografía histórica de la Cordillera de la Costa de Chile. In: Smith Ramírez C, Armesto JJ, Valdovinos C (Eds) Historia, biodiversidad y ecología de los bosques costeros de Chile. Editorial Universitaria, Santiago de Chile, 99–116.
  • Villagrán C, Hinojosa L (2005) Esquema biogeográfico de Chile. In: Llorente J, Morrone J (Eds) Regionalización biogeográfica en Iberoamérica y tópicos afines. Ediciones de la UNAM, México DF, México, 551–577.
  • Villagrán C, Abarzúa AM, Armesto JJ (2019) Nuevas evidencias paleobotánicas y filogeográficas de la historia Cuaternaria de los bosques subtropical-templados de la Cordillera de la Costa de Chile. In: Smith-Ramírez C, Squeo FA (Eds) Biodiversidad y Ecología de los Bosques Costeros de Chile. Universidad de Los Lagos, Osorno, Chile, 3–21.
  • Wemmer C, McCarthy A, Blouch R, Moore D (1998) Deer: status survey and conservation action plan. IUCN/SSC Deer Specialist Group, Gland, Switzerland, 106 pp.