SSR primers used. N, number of alleles; size, allele size range (bp).
The distribution of the endangered species Dioon edule is in populations scattered throughout the Sierra Madre Oriental in San Luis Potosí, Mexico. Its habitat is tropical dry forests at lower elevations and oak forests at higher elevations, mainly disturbed by anthropic activities. We determined and analyzed nine populations’ demographic structure and explored the genetic diversity of five using SSR markers. The population density averaged 2050 individuals ha-1 and have an aggregated distribution pattern. Differences in the plants’ size among localities are due to site quality, based on their adaptation capacity and response to climate and soil traits. Most populations have the highest mortality in the early stages of life, with a low mortality rate for those who survive this bottleneck. Two populations have a relatively constant mortality rate, attributed to disturbance of the habitat. The populations show low genetic diversity and an excess of homozygotes. Their similarity is probably related to the formation of natural corridors favoring connectivity between populations. The deterioration and fragmentation of the habitat have severe effects on the populations’ viability, like reducing gene flow, which has led to inbreeding and genetic drift.
- habitat disturbance
- gene flow
The relevance of the species
This species is distributed in the Sierra Madre Oriental (SMO), within the geographic-cultural region known as La Huasteca, including the south of Tamaulipas, east of San Luis Potosi, Hidalgo, Queretaro, and Veracruz. The elevation goes from sea level to 1525 meters above sea level (m.a.s.l.) [2, 3]. It has stems up to 6 m high in some occasions semi-prostrated, presenting reproductive events with an interval of 10 to 52 years in female individuals and 2.8 to 8.8 in male individuals . Its strobilus is pollinated only by specialist insects (beetles) associated with them .
The species is currently classified as near threatened (NT) by International Union for Conservation of Nature and endangered by the NOM-059-SEMARNAT-2010 in Mexico. This category is mainly because of land-use change to agriculture and livestock, causing habitat fragmentation that adversely affects population dynamics and a possible reduction in gene flow among populations [4, 6].
It is necessary to integrate information from different areas, mainly demography and population genetics, to develop and implement effective conservation strategies  to build general conclusions from the observations derived from each of these areas.
Demographic studies provide the basis for knowledge of populations and allow us to assess populations’ viability, identify factors that reduce the population, and predict the persistence of small, isolated populations and evaluate alternative conservation scenarios . Population genetics allows us to describe the genetic composition of natural populations and predict their changes in response to various evolutionary forces operating on them, allowing us to identify some risks of a genetic nature that affect the persistence of species such as habitat fragmentation and loss of adaptive potential associated with the decline of genetic diversity and inbreeding [8, 9]. Diverse studies related to population genetics and genomics of cycads have been published based on molecular markers such as Aloenzymes, RADseq, chloroplast markers, ISSR [4, 10, 11, 12, 13, 14, 15].
The objective of this research was to study the demographic and genetic status of
2. Materials and methods
2.1 Study area
The Sierra Madre Oriental is fundamentally a group of minor mountain ranges formed by marine sedimentary rocks, mainly limestone, sandstone, and shales. The maximum elevations are up to 3180 m.a.s.l., with small intermontane valleys limited by high mountain ranges with steep slopes (Figure 1). In this region, there are soils of alluvial origin formed in the great plains with sediments mostly coming from limestone and lutites and residual and colluvial origin in the high and lower parts of the mountain ranges constituted by the same type of rock. In general, the soil is medium texture, shallow, and rich in organic matter and nutrients depending on the climate and vegetation . The area’s climate is semi-warm, characterized by average annual temperatures above 18 °C and 600 mm of rainfall, mostly during the summer .
In the state of San Luis Potosi, three important vegetation regions for cycads stand out : Tamasopo with tropical rain forest and oak forest (
The study of demography was conducted in nine populations of
The sampling consists of modifying the line interception method  and used by  in transects of 30 m long and 10 m wide. Rectangular units provide the advantage of evaluating variables by walking in a straight line without moving sideways, starting at one end of the line, and considering only those plants within the plot . Due to the harsh characteristics of the area with rocky slopes of up to 70°, we traced two linear transects 100 m long by 2 m wide, parallel to the slope, locating them where the most significant number of individuals were concentrated, considering that
The different stages of the cycads’ life cycle were defined according to the number of leaves, the number of leaflets and stem height  because it has not yet been possible to define the age of the cycads due to their growth characteristics . Therefore, individuals were classified as follows: seedlings, individuals that present one or more leaves without these forming a crown; juveniles, individuals that present one or more crowns, but with a non-visible stem; adults, individuals that present one or more crowns and have a visible stem. Adults were classified as male, female, or non-reproductive according to the presence of reproductive structures.
Once the cycads’ life cycle stages were defined, each stage’s survival probabilities calculation as the proportion of individuals that survived to enter the next category, relative to the initial number of seeds (lx = nx/n0) . A survival curve describes the survival pattern of individuals in a population over time, plotting the lx on the y-axis, in logarithmic scale, with life cycle stages on the x-axis There are three basic survival curves, the survival curve type I represents a population in which most of the organisms die in the older age categories. The type II curve represents a population in which the mortality rate is almost constant, and type III represents a population in which most individuals die during the first stages of the life cycle, and very few reach the intermediate and late age categories .
The stem’s diameter and the plant’s total height were measured, the total number of leaves was counted, and the length of those on the last crown was measured.
Considering that the species
The data on population density and potential reproductive capacity was analyzed using a one-way Analysis of Variance (ANOVA) to determine the statistical difference among the nine populations studied (localities) with two transects in each (n = 18).
Data on plant size and number of leaves of individuals in reproductive and non-reproductive states were analyzed through a nested ANOVA (random-effects model II) to determine the statistical difference among localities and categories of plants’ reproductive states in all localities and each location. When statistically significant differences between means occurred, a Tukey multiple comparison test was applied (p < 0.05). All analyzes were performed with XLSTAT software (Addinsoft Inc., Paris, France).
2.3 Genetic diversity
The genetic diversity was studied in five populations (Figure 1), covering an altitudinal gradient from 388 to 1050 m.a.s.l. In each population, leaf tissue was collected from 21 individuals (seven from each age category). The DNA extraction was done with the DNeasy Plant minikit® Kit (Quiagen, Valencia, California, USA), following the kit’s protocol and using approximately 100 mg of leaf tissue from which the cuticle was previously detached. A total of seven microsatellite markers designed by  (Table 1) labeled M-13 (−21pb) to obtain fluorescent products [23, 24], and PCR products were separated and genotyped with the LI-COR 4300 DNA analyzer (LI-COR Biosciences, Lincoln, Nebraska, USA) in a 6.5% polyacrylamide gel (LI-COR Biosciences, Lincoln, Nebraska, USA), at a wavelength of 700 nm.
|Locus||Primer sequences (5′–3′)||Repeat||N||Size|
It was calculated the average number of alleles per locus, percentage of polymorphic locus, the expected and observed heterozygosity, the fixation index and it was estimated the deviation of the Hardy–Weinberg Equilibrium using the GenAlex software , for all the samples and grouping the individuals according to their life stage. An Analysis of Molecular Variance (AMOVA) was performed on genetic distances. Wright’s F statistics (Fis, Fst, Fit) was calculated for the whole set of populations and each age category with the GenAlex software , and also the degree of genetic differentiation between pairs of populations (Fst) was estimated for all the samples and grouped by life stage.
The density of the nine locations averaged 2050 individuals ha−1 ± 293.36. Density ranges from 3775 individuals ha−1 in SA to 775 individuals ha−1 in MO. The populations with the highest density of individuals are the localities of SA and GA, and the lowest density is MO, although there are no differences between localities in the number of individuals ha−1 of the populations (ANOVA, df = 8; F = 2.287; P > 0.05).
The stages of the life cycle (Table 2) showed that seedlings were the most abundant, with an average of 49.27% ± 5.00, varying from the lowest in RN to the highest in SA. The proportion of young plants averaged 23.05% ± 2.62, with CH standing out with the lowest proportion and RN with the highest. The proportion of non-reproductive adults was 11.39% ± 3.12, being absent in SA and RN localities, unlike CH, where they represent a third of the population. Reproductive adults represented 16.29% ± 4.48, being scarce in SA and RN, but in the last, they are almost half of the population.
The survival curve based on the states of the life cycle of the populations is an inverted “J” (type III), well defined in most locations (Figure 2F,G,H,I), and scarcely defined in other localities (Figure 2D,E). The localities of SA (Figure 2A), GA (Figure 2B) and AN (Figure 2C) showed a curve like type II in the populations.
The population’s potential reproductive capacity was 250 ha−1 ± 45.26 individuals in the localities with reproductive structures, of which 158 ha−1 ± 35.36 were male and 92 ha−1 ± 17.18 female. The highest number of male individuals with the reproductive structure were GA and RN, and with the lowest number was SA (Figure 3A).
Likewise, the localities with the highest number of female individuals were RN and CH, and those with the lowest number were MO and AN (Figure 3A). The average sex ratio for the locations was 2: 1 (male: female), with the highest ratio of 4: 1 in GA and the lowest of 0.5: 1 in CH (Figure 3B). The ANOVA showed that there are no differences between mean localities in the number of female individuals ha−1 (n = 18; df = 8; F = 1.234; P > 0.05) and male individuals ha-1 (n = 18; df = 8; F = 1.006; P > 0.05) of the populations.
An average of 35.70% ±2.79 of individuals was found in the transect up to 25 m, mainly seedlings and juveniles in the localities. From 25 to 50 m was found on average 34.39% ± 5.60, in the 50 to 75 m 19.89% ± 3.41 and in the 75 to 100 m, 10.02% ± 2.44. The 70% of the individuals were concentrated in the first 50 m from the transect’s starting point and particularly in RN 100%. At sites with a steeper slope, there were plants widely dispersed, while at sites with flatter slope (less than 20%, e.g. SA), there were large groups of individuals aggregated around adult plants or in the shade of oaks (
Considering the plants’ size and the number of leaves as parameters associated with the cycads’ age, the longest and largest plants are in the communities of PO, CH, GA, AN, and SA. If it was impossible to measure the adult plants’ height and diameter because they had a subterranean stem, only the number of leaves was considered, obtaining similar results.
The average height for non-reproductive individuals was 23.4 cm ± 3.4, with a maximum average height of 52.2 cm and the lowest of 14.5 cm. For reproductive plants, the maximum recorded height was 56.3 cm, and the lowest was 9.5 cm, with an average of 28.8 cm ± 8.3. However, the ANOVA showed that there are differences among means in plant height between localities (df = 8; F = 5.428; P < 0.001) and between categories of reproductive status in localities (df = 8; F = 4.992; P < 0.05), but not between categories of reproductive status in each locality (df = 1; F = 2.633; P > 0.05) (Figure 4).
The average stem diameter was 12.2 cm ± 2.5 in non-breeding individuals, with a width of 23.6 cm to 15.0 cm (Figure 4A), and in reproductive individuals, the average was 16.5 cm ± 4.4 with a maximum of 31.2 cm and a minimum of 5.5 cm. The ANOVA showed that there are no differences in the mean diameter between localities (df = 8; F = 1627; P > 0.05) or between categories of reproductive status in each locality (df = 1; F = 2848; P > 0.05), but there are differences between the categories of reproductive status of individuals in the localities (df = 8; F = 4034; P < 0.05) (Figure 4).
The average number of leaves was 17.69 ± 2.18 in all individuals. The average number of leaves was 17.69 ± 2.18 in all individuals. In nonreproductive adults it was 14.3 with a maximum of 21 and a minimum of 3.5 (Figure 4A), for reproductive adults it was 15 with a maximum of 22.5 leaves. The ANOVA showed that there are differences in the mean number of leaves between localities (df = 8; F = 3.166; P < 0.05) and between categories of reproductive status in each locality (df = 1; F = 7.312; P < 0.05) but not between the categories of reproductive status of individuals in the localities (df = 8; F = 1.381; P > 0.05) (Figure 4).
3.2 Genetic diversity
Only three of the first seven proposed were amplified, ED9, TOM5, and CAP5. All the loci were polymorphic, but the population AN have the lowest proportion of polymorphic loci (33.33%) in the seedlings (Table 3). The average number of alleles per locus was 3.66; in the seedling stage, the lowest value corresponds to Los AN (1.33) (Table 3). There were no unique alleles within the adult stage. However, in the juvenile and seedling, four and three unique alleles were found, respectively, and ED9 was the locus that presented the highest number of unique alleles.
|Location||Life cycle stage||N||Na||Polymorphic Loci (%)||Ho||He||F|
The average heterozygosity (Ho) was 0.207. The highest value was in the seedlings in PO with 0.333, and the lowest in the seedlings in PO (Ho = 0.0), and the adult stage in CH (0.048) (Table 3). For the expected Heterozygosity, the average was 0.504 and the highest value was 0.634 in the juvenile category for the PO locality.
In the populations, only two loci resulted in Hardy–Weinberg equilibrium (H-W), ED9 in the populations of RN, PO, and AN, and TOM5 in the population of RN (Table 4). In the seedlings, the ED9 locus resulted in H-W balance in the populations of CH, PO, and SA and TOM5 resulted in balance in RN. The loci CAP5 in PO and SA populations and ED9 in SA, were monomorphic, and CAP5 in PO and SA populations and ED9 in SA. In the juveniles, ED9 resulted in H-W equilibrium in all populations. However, in adults, this same locus resulted in equilibrium in most populations except CH (Table 4).
The AMOVA showed that the categories’ differences are explained by the genetic variation between individuals for all populations and each age category (Seedlings = 50%; juveniles: 59%; adults: 62%: Total: 58%) (Table 4). The inbreeding coefficient (FIS) generally shows a deficiency of heterozygotes in all populations and age categories (Fis ≈ 0.666–0.708) (Table 5).
|Life cycle stage||Source||df||SS||MS||VAR||%|
The fixation index (Fst) indicates a moderate differentiation in all populations, as well as in the age categories (Fst ≈ 0.11–0.115). The populations that share more alleles are CH and PO located in the south of the study area (Fst = 0.057; p = 0.006); those that share fewer alleles are PO and SA also in the south (Fst = 0.265; p = 0.001) (Table 5). The seedlings have a high genetic differentiation between pairs of populations (Fst > 0.26; p < 0.05), while in the juveniles, the differentiation is moderate (Fst > 0.110; p < 0.05), and in adults from moderate to high (Fst > 1.15; p < 0.05).
In the study area the average density was 2103 plants ha−1, lower than the 3200–4633 plants ha−1 of
The absence of statistically significant differences between locations is mainly due to the parameters’ wide variation observed. However, the populations with lower density coincide with other studies [27, 28], where the habitat’s quality has decreased. The low quality of the habitat is due to the recurrent disturbance caused by agriculture, mainly the establishment of sugar cane fields in RN, pastures for livestock in EJ and MO, and the extraction of complete plants or their parts, like seeds and leaf harvest for tamale production [1, 29, 30]. In general, plant density and size are population attributes that indicate habitat quality .
Most populations presented the survival inverted “J” curve (type III) typical of shade-tolerant species  and most populations of
The plant’s reproductive capacity is of great interest in conservation since it significantly impacts management and recovery strategies . It is remarkable the high number of female plants present compared with other studies. The imbalance of females is probably related to climatic factors, soil quality, or the occurrence of fires, which in the study area is very common due to agricultural activities [18, 26, 34, 35, 36, 37, 38]. The plants’ sex ratio average was 2: 1 (male: female) in all locations. However, the ideal ratio for
The differences in height between reproductive and non-reproductive individuals have been attributed to the lowest growth rate of reproductive individuals due to the compression of stem tissues by the increase in weight of stored water and the production of heavy strobili [18, 40]. However, in this study, no significant differences were observed in any location, indicating that this difference does not always apply to all cycad populations.
Considering that there were significant plants’ height differences among populations, it would be an indicator of the site quality. The sites of CH, GA and AB that presented the tallest plants would have the best quality since they present the lowest proportion of organic matter and smooth slopes. Diference in height and diameter of
The plants’ diameter did not present significant differences between localities. It has already been observed in
As highlighted in this study, another determining factor in the development of
The genetic diversity analysis tried to clear the picture of the medium and long-term repercussions of adverse factors identified when analyzing the species’ population demography. However, we consider that the number of populations and molecular markers should be increased since
The populations analyzed show low genetic diversity and an excess of homozygotes in the seedlings, juveniles, and adults. Contrary to that expected for adult plants, which were believed to have become established before the disturbances. Similar behavior was present in a study conducted in Veracruz, mentioning that populations outside the Pleistocene flora and fauna refuges located in Veracruz and Oaxaca probably suffered a bottleneck due to exposure to adverse environmental conditions generating a significant decrease in population size reflected in the low levels of observed heterozygosity detected .
Previous studies suggest that life history issues, such as geographic distribution area, reproduction systems, and pollen and seed dispersal systems, have essential consequences on genetic structure and genetic variability levels in natural populations [4, 6, 10, 45, 46]. These factors are especially relevant in cycads because of their distribution in particular areas. They also present vegetative reproduction forming inside the populations .
In addition to the patches of clonal sprouts, the seedling banks that form around adult plants of
Deviations in H-W were at most of the loci assessed, all tending towards a deficit of heterozygotes. The deviations may be due to the lag in the reproductive season, which can vary from four to 52 years in female plants . This produces a time barrier that prevents individuals’ panmictic crossing, causing groups to form that can function as independent populations.
Five unique alleles were identified with low frequencies in the ED9 and TOM5 locus in seedlings and juveniles. These alleles’ appearance could be attributed to recurrent mutations or gene flow by pollen dispersal at great distances . These alleles may benefit from genetic drift promoted by population size reduction.
AMOVA indicates a high differentiation in seedlings (FST = 0.26), probably due to the reduction in gene flow, which has led to inbreeding within the populations. Simultaneously, in adults and juveniles, the values are very similar and may be maintained by self-incompatibility that creates a dependence on pollinators to facilitate reproduction .
Some populations share some alleles, suggesting that they were connected by natural corridors some time ago, mainly in the vicinity of the Sierra Gorda in the northern of Querétaro and south of San Luis Potosí. However, in the current generations of seedlings, the Fst values indicate that the populations analyzed are becoming isolated. The isolation is probably caused by physical barriers associated with anthropic activities increasing in the area. This situation diminishes the habitat conditions and leads to higher mortality, affecting pollinators and favoring inbreeding and gene drift .
GRM thanks the support of CONACYT (Scholarship for Postgraduate Studies 290685) and the PMPCA-UASLP. The Cycad Society financed the research through the project “The chamal (