Genetic structure of coexisting wild and managed agave populations: implications for the evolution of plants under domestication.

Domestication is a continuous evolutionary process guided by humans. This process leads to divergence in characteristics such as behaviour, morphology or genetics, between wild and managed populations. Agaves have been important resources for Mesoamerican peoples since prehistory. Some species are domesticated and others vary in degree of domestication. Agave inaequidens Koch is used in central Mexico to produce mescal, and a management gradient from gathered wild and silvicultural populations, as well as cultivated plantations, has been documented. Significant morphological differences were reported among wild and managed populations, and a high phenotypic variation in cultivated populations composed of plants from different populations. We evaluated levels of genetic diversity and structure associated with management, hypothesizing that high morphological variation would be accompanied by high genetic diversity in populations with high gene flow and low genetic structure among managed and unmanaged populations. Wild, silvicultural and cultivated populations were studied, collecting tissue of 19-30 plants per population. Through 10 nuclear microsatellite loci, we compared population genetic parameters. We analysed partition of variation associated with management categories to estimate gene flow among populations. Agave inaequidens exhibits high levels of genetic diversity (He = 0.707) and moderate genetic structure (FST = 0.112). No differences were found in levels of genetic diversity among wild (He = 0.704), silviculturally managed (He = 0.733) and cultivated (He = 0.698) populations. Bayesian analysis indicated that five genetic clusters best fit the data, with genetic groups corresponding to habitats where populations grow rather than to management. Migration rates ranged from zero between two populations to markedly high among others (M = 0.73-35.25). Natural mechanisms of gene flow and the dynamic management of agave propagules among populations favour gene flow and the maintenance of high levels of variation within all populations. The slight differentiation associated with management indicates that domestication is in an incipient stage. Published by Oxford University Press on behalf of the Annals of Botany Company.

Introduction

Domestication is a gradual and continuous process through which plants undergo phenotypic and genetic changes, mainly resulting from artificial selection favouring organisms with features desirable, necessary or interesting to humans (Darwin 1859; Casas ; Pickersgill 2007). Since Darwin's time (Darwin 1859, 1868), domestication has been a useful model for studying a variety of questions in evolutionary biology (Gepts 2004; McKey ). Understanding the relationship between domesticated and wild individuals provides valuable opportunities for studying evolutionary processes in nature, especially in areas where the ancestor and descending organisms coexist and are available for comparison (Olsen and Wendel 2013). The most common pattern of domestication involves a reduction of genetic variation in domesticated populations compared with wild ones, mainly due to artificial selection (Sauer 1972; Doebley 1989; Olsen and Wendel 2013). However, some studies in plants have documented that such a trend may vary according to the life history of the plant species analysed, as well as the context of human cultures involved in use and selection of those species (Casas , 2007; Vargas-Ponce ; Parra ; Aguirre-Dugua ).

Mesoamerica is a main centre of diversity of useful plant species (Caballero ). It is also the centre of domestication of nearly 200 economically important plant species worldwide, such as maize, beans, chilli peppers, cocoa, cotton, prickly pears (Vavilov 1951; Harlan 1975) and among them several species of agaves (Gentry 1982; Colunga-GarcíaMarín ). In addition, numerous native plant species (>800) are at various stages of domestication in traditional farming systems, such as home gardens and other agroforestry systems, most of them still unstudied (Colunga-GarcíaMarín and Zizumbo-Villarreal 1993; Casas , 2007; Caballero ; Blancas ). In these systems, it is possible to find coexisting wild, semi-domesticated or incipiently domesticated, and fully domesticated plants belonging to a single species (Casas ). In Mesoamerica, people practice in situ management of wild plants in forests, forest patches and agroecosystems, which includes several types of interaction: letting stand, encouraging growth and special care and protection of more favourable plants. These interactions in some cases involve artificial selection, resulting in documentable domestication processes (Colunga-GarcíaMarín and Zizumbo-Villarreal 1993; Casas , 2007).

Agaves belong to the family Asparagaceae, which comprises 9 genera and ∼300 species (Eguiarte ; Rocha ). The genus Agave is relatively recent, about 10 million years (my) old, with >160 species, ∼75 % of them occurring in Mexico, where it has shown a high adaptive radiation in diverse ecosystems (Eguiarte ). Among the most meaningful adaptations of the genus, several authors (Rocha ; Eguiarte and Souza 2007) have highlighted the morphological and physiological defences against predators, and substantial diversity of reproductive mechanisms and pollinators. Agaves have high chromosome number and large genome size, nearly 64.3 % of species being polyploid (Goldblatt 1980; Palomino ), showing different ploidy levels. For instance, for Agave cupreata Trel. et Berg., Palomino reported diploids, tetraploids, pentaploids and hexaploids. For A. tequilana Weber, Bennett , Robert and Simpson recorded that genome size varies from 2940 to 77 458 million base pairs.

In a review considering 22 species of Agave, Eguiarte estimated that, on average, genetic diversity of this genus is high and genetic differentiation is moderate compared with other outcrossing monocots (Hamrick and Godt 1996).

About 102 taxa of agaves have been reported as being used in Mexico (Colunga-GarcíaMarín ; Torres ), most notably as human food. Before maize was adopted as the main Mesoamerican crop, agaves were a principal source of carbohydrates for peoples of western Mexico and an area that is now the southwest USA, who consumed the stems, leaf bases and flower stalk cooked in stone ovens (Callen 1965; Smith 1986; Hodgson 2001; Zizumbo-Villarreal ).

At present, flower buds, flower stalk and leaf bases are edible, whereas the sap, called aguamiel, is collected for drinking fresh or for preparing the fermented drink called pulque. In addition, various distilled beverages (generically called mescal, including the appellations of origin Tequila, Bacanora and Raicilla) are prepared with fermented cooked stems. Other important uses are as a source of fibre, and as fodder. The activity that has gained the highest economic importance in recent decades is the production of mescal, spirit drinks produced from the central corm of the plant. For the production of this beverage, the use of 53 taxa has been recorded in Mexico (Colunga-GarcíaMarín ; Torres ), most of them extracted from forests. The selection and harvesting of agave corms are carried out just before sexual reproduction takes place, which poses a high risk for numerous populations of several species (Delgado-Lemus ).

Only four Mexican Agave taxa and two of the USA referred to as crops have been studied in terms of the consequences of divergence in genetic diversity caused by domestication and management practices. Among the species studied are A. fourcroydes Lem., which was domesticated by the Pre-Columbian Maya of the Yucatán Peninsula for fibre, and their putative ancestor A. angustifolia Haw. (Gentry 1982). Using isozyme markers, Colunga-GarcíaMarín found that A. fourcroydes exhibits lower levels of genetic diversity than wild A. angustifolia populations, probably due to the predominant asexual reproduction of a clone selected for commercial monoculture plantations since the beginning of the 20th century, in addition to the disappearance of the traditional management. In southern areas of the state of Jalisco, the complex of taxa related to A. angustifolia has been used to produce mescal and tequila. Currently, it is possible to find wild varieties, traditionally managed varieties (of the species A. angustifolia and A. rhodacantha Trel.) and one variety that is predominately found in commercial monoculture plantations (A. tequilana var. azul). Using inter simple sequence repeat markers, Vargas-Ponce found that traditional landraces of A. angustifolia have genetic diversity levels similar to those recorded in wild populations. In contrast, genetic diversity in A. tequilana was markedly low. These trends can be explained as a result of the differences between traditional management and commercial monocultures, and the predominant asexual method used to propagate agave crops. In southeast Arizona, A. parryi Engelm. and A. p. var. huachucensis (Baker) Little were cultivated since prehistory, as sources of food and fibre. A number of relict populations from ancient cultivated areas remain in modern landscapes. Parker , 2014) evaluated the genetic divergence between wild and anthropogenic populations (relicts of ancient managed populations) through isozymes and microsatellites. In this case, the genetic diversity in the anthropogenic populations was lower than that of wild populations, a trend that can be explained due to dominance of asexual reproduction in anthropogenic populations.

In the cases mentioned, the general trend observed is that cultivated populations exhibit lower levels of diversity than wild populations. The most frequent explanation for this pattern is that genetic diversity in wild populations has been modelled for millions of years, whereas managed populations most commonly involve a fraction of such variation in time periods much shorter. However, in the case of the complex A. angustifolia in southern Jalisco, it is clear that traditional management allows the occurrence of high genetic diversity levels in the crop lands, similar or even higher than in the wild, through the constant let standing or encouraging of plants from wild populations already growing in the crop land, or introducing them to these areas. A similar pattern was also reported for several columnar cacti species (Casas ; Parra , 2010).

Agave inaequidens Koch is distributed mainly in pine and pine-oak forests of the Trans-Mexican Volcanic Belt. Historically, this species was consumed as food and its sap extracted for drinking fresh or fermented (Gentry 1982). At present, it is used mainly for producing mescal and fibre. In Michoacán, this species is found in a gradient of management intensity with populations occurring in wild habitats as part of natural forests, but also under silvicultural or in situ management, through which people leave some individuals standing when the forest is cleared, and deliberately propagate agaves in the cleared areas in order to increase their population density (Torres ). In addition, some people cultivate this agave away from its natural habitat (ex situ cultivation), by transplanting saplings taken from their wild populations and others produced in seedbeds or nurseries, but this practice is relatively new, no more than 30 years old (Figueredo ; Torres ).

Based on previous ethnobotanical and morphological studies, we (Figueredo ) documented the morphological variation and differential use of variants in populations of A. inaequidens. We also documented that cultivated populations are composed of individuals from different wild populations, propagated by both seeds and vegetative propagules (Figueredo ; Torres ). This background leads us to hypothesize that the high morphological variation observed in cultivated and managed populations would be accompanied by high genetic diversity, similar to that found in wild populations of this species. Furthermore, given the dynamics of formation of crops of this species, and the relatively short time this activity has been carried out, we expected to record high levels of genetic diversity in both managed in situ and cultivated populations. We also expected to document high levels of gene flow among managed and unmanaged populations and low genetic structure among populations.

Our study aimed to (i) assess levels of genetic diversity in populations of A. inaequidens within a gradient from low to intensive management; (ii) evaluate the levels of genetic structure and divergence between wild, silvicultural managed and cultivated populations of A. inaequidens and (iii) estimate gene flow among populations under different management practices, in order to identify possible sources of the gene pools.

Methods

Study species

Agave inaequidens is recognized with the common name of ‘maguey alto’, which makes reference to its large size, or ‘maguey bruto’, which makes reference to the caustic property of its tissue, because of the presence of saponins and other secondary metabolites that can cause dermatitis (Gentry 1982). This species is endemic to Mexico, naturally growing throughout the Trans-Mexican Volcanic Belt (Gentry 1982). Its distinctive characteristic is teeth dimorphism (successively one shorter and one larger) along the leaves' edges. It is monocarpic, with yellow flowers, protandrous and xenogamous, with pollination conducted by bats and also diurnal animals (León-Jacinto 2013). Fruits are capsules producing shiny black seeds dispersed by the wind (Gentry 1982; León-Jacinto 2013).

This species has >30 different documented uses (Torres ). In the past, this agave was more important as food, since its inflorescences and stems were cooked in underground ovens and consumed. Its sap was extracted for drinking fresh or fermented (Gentry 1982). Nowadays, those uses are rather rare, and its main use is mescal production and, in some areas, the extraction of its fibres for manufacturing cords (Valenzuela-Zapata ).

Study area

We sampled a total of 16 populations of A. inaequidens: 6 wild, 7 cultivated and the only 3 silvicultural in situ managed populations accessible in the state of Michoacán (Fig. 1). Wild populations were growing in pine-oak, oak and subtropical forests. Seven cultivated populations are found growing in orchards together with other species of agaves and fruit trees, where part of the forest was removed for cultivation. Two silvicultural managed populations found in the region grew in secondary forests (natural forest cover has decreased) and one in pastureland (Table 1). We collected tissue of young healthy leaves of 19–30 individuals per population. Tissue samples were kept in silica gel until extraction of total DNA.

Figure 1.

Locations of the 16 populations of A. inaequidens studied in Michoacán state. Triangle, wild; square, managed; inverted triangle, cultivated populations.

Table 1.

Sampling sites studied for A. inaequidens in Michoacán State. The localities were classified according to their category in wild, cultivated and managed. Acronyms of localities and sampling size are shown in parentheses.

Category	Population name and municipality	Acronym (n)	Elevation (m)	Vegetation
Wild	Piedra del Indio, Morelia	SPIE (30)	2414	Pine-oak forest dominated by species of Pinus and Quercus
	Cuanajo, Pátzcuaro	SCUA (30)	2541	Oak forest dominated by Q. castanea, Q. candicans, Q. laeta, Q. crassipes and Q. rugosa
	Pino Real, Queréndaro	SPR (30)	2323	Pine-oak forest dominated by species of Pinus and Quercus
	Icuacato, Quiroga	SICU (30)	2475	Oak forest dominated by Q. castanea, Q. candicans, Q. laeta, Q. crassipes and Q. rugosa
	Pino Bonito, Queréndaro	SPB (30)	2328	Pine-oak forest dominated by species of Pinus and Quercus
	La Manga Manseña, Sahuayo	SSAH2 (19)	1913	Subtropical scrub dominated by Euphorbia tanquehuete, Bursera fagaroides, Eysenhardtia polystachya and Ipomea murucoides. Approximately 100 plants
	Ejido de Parras, Queréndaro	ROC1 (30)	2127	Orchard of 2 ha with fruit trees and other Agaves species (A. angustifolia, A. cupreata and A. americana). Approximately 7000 plants about 4 years of being transplanted
Cultivated	Los Alamos, Queréndaro	ROC2 (30)	2691	Orchard 5 ha with fruit trees. Around 5000 plants about 8 years of being transplanted
	El Salto, Queréndaro	ROC3 (30)	2383	Orchard with few fruit trees, species of Pinus sp., 4 ha, with pepper crops. Approximately 4000 plants about 11 years after being transplanted
	Traspatio, Queréndaro	CTC (30)	2059	Home gardens, plants about 8 years after being transplanted
	La Huertilla, Queréndaro	CLH (30)	2078	Orchard of ¼ ha growing only A. inaequidens, with about 9 years after being transplanted
	Barranca del añil, Sahuayo	CSAH1 (30)	1909	Orchard 5 ha with fruit trees. Approximately 50 plants about 3 years of being transplanted
	Lindero Don Tarsicio, Sahuayo	CSAH2 (20)	1710	Plants arranged as living fences, and mature individuals are used to extract ‘aguamiel’
Managed	Aguacatillos, Indaparapeo	M1 (30)	2508	Secondary vegetation with few oaks. The management in this area is seedlings transplants
				When they are grouped (i.e. aggregate spatial distribution), they are more dispersed
	Salecillo, Indaparapeo	M2 (30)	2486	Same as referred to above
	La Paja, Sahuayo	SSAH1 (30)	1901	Pastureland dominated by Mimosa sp. Some individuals aligned as living fence

DNA extraction, amplification, marker screening and data quality

Total DNA was extracted following the cetyltrimethylammonium bromide protocol for plants (Lefort and Douglas 1999). We used nuclear microsatellites as genetic markers. Eleven different oligonucleotides or primers were assayed, which were originally developed for the studies of A. parryi, A. parryi var. huachucensis and A. palmeri (Parker , 2014; Lindsay ).

Polymerase chain reaction (PCR) was conducted using the QIAGEN (www.qiagen.com) Multiplex kit at a final volume of 5 µL, including 2.5 µL of Master Mix (with HotStartaq DNA polymerase, Multiplex PCR buffer, 3 mM MgCl2 and deoxynucleotides), 1 mM of each primer, 1.5 µL sterilized water and 0.5–1.0 µL of 50–100 ng µL−1 DNA. For all microsatellite loci, amplifications were carried out with a EscoSwiftMaxPro thermal-cycler, using the following protocol: initial heat activation for 15 min at 95 °C, and then 35 cycles of denaturing at 95 °C for 1 min, annealing at 60 °C for 1 min 30 s for all microsatellite loci and extension at 72 °C for 1 min. We included a final step of extension at 72 °C for 7 min. PCR products were mixed with formamide and the size standard Gene Scan LIZ-500 (Applied Biosystems) and denatured for 5 min at 95 °C. PCR products were analysed in a 3130xl Genetic Analyser (Applied Biosystems) sequencer. The resulting electropherograms were analysed using the program Peak Scanner (Applied Biosystems).

Possible genotyping errors due to the presence of null alleles, large alleles or stuttering were identified using the MicroChecker 2.2.3 software (Van Oosterhout ) with 1000 bootstrap simulations and a confidence interval (CI) of 95 %. Deviations from Hardy–Weinberg Equilibrium (HWE) were examined for all loci in each population using the exact test with Arlequin version 3.11 (Excoffier ). Deviations from linkage equilibrium (LE) were estimated using Genepop on the Web program with Fisher's method for each pair of loci (Raymond and Rousset 1995).

Genetic diversity

Number of alleles per locus (Na), effective number of alleles (Ne) per locus, expected heterozygosity (He), unbiased expected heterozygosity (uHe) and observed heterozygosity (Ho) were estimated using GenAlEx (Peakall and Smouse 2006). Kruskal–Wallis one-way analysis of variance range tests were performed to identify significant differences in levels of diversity according to the management type for each of the estimated parameters.

Genetic structure

Following Weir (1996), the fixation index (FST) was calculated with the program FreeNA using the excluding null alleles method assuming null alleles, with 10 000 bootstrap repetitions (Chapuis and Estoup 2007). The inbreeding coefficient (FIS) was calculated correcting for null alleles with the INEst program (Chybicki and Burczyk 2009) using the Bayesian model IIM assuming inbreeding. Every run consisted of 10 000 burn-in and 50 000 periods of Markov Chain Monte Carlo simulations (MCMC). The genetic distances (DC) of Cavalli-Sforza and Edwards (1967) were estimated for each pair of populations using the including null alleles correction described in Chapuis and Estoup (2007). From a matrix of Nei's genetic distances (D), we constructed a dendrogram through the unweighted pair group method with arithmetic mean (UPGMA) method with 1000 bootstrap replicates of the original matrix with TFPGA 1.3 (Miller 1997) and the MEGA program (Tamura ). Isolation by distance and Mantel tests were examined with isolations by distance webb service (Jensen ) using Nei's genetic distances.

Genetic differentiation and genetic flow

The program STRUCTURE version 2.3.4 (Pritchard ; Falush ) was used to perform the Bayesian clustering, in which individuals are probabilistically assigned to one of K predefined groups to identify the optimal number of genetic groups (Evanno ). The optimum group number (K) was determined varying K from 1 to 17, with 10 runs for each K value, in order to determine the maximum value of the posterior probability [ln P(K)]. Every run consisted of 5.0 × 104 burn-in and 106 periods of MCMC repetitions after the burn-in. We used the admixture model with correlated allelic frequencies without prior information. The number of subpopulations was additionally estimated based on the approach of Evanno using the software STRUCTURE Harvester (Earl 2011). In order to align the cluster membership coefficients of the 10 structure runs and to graphically display the results, we used the programs CLUMPP version 1.1.2 (Jakobsson and Rosenberg 2007) and Distruct version 1.1 (Rosenberg 2004). Analyses of molecular variance (AMOVAs) were used to test for genetic differences firstly among the populations under the three types of management and then according to the genetic grouping resulting from the Bayesian clustering. For these tests, we used the stepwise mutation models with the program Arlequin version 3.11 (Excoffier ).

Migration rate (M = m/µ, where m is the migration rate per generation and µ is the mutation rate) paired in both directions, and Theta (Θ = 4Neµ, where Ne is the effective size of population) was estimated through MIGRATE-N (Beerli and Felsenstein 1999, 2001), based on maximum likelihood using the Brownian method and a constant mutation rate (µ). From the values of M and Θ, we estimated the gene flow or number of migrants per population (Nm). The effective population size (Ne) per population was estimated through an average mutation rate for microsatellites (5 × 10−4, Schlötterer 2000; Selkoe and Toonen 2006).

Results

In total, we successfully amplified 10 microsatellite loci. No evidence for large-allele dropout was detected in any locus using the software Micro-Checker. Null alleles seemed to be present in 14 populations for loci 6 and 7 (APARLC20, APARLC21, Lindsay ), whereas populations CSAH1 and M1 had no evidence of null alleles. Significant deviations from HWE were recorded associated with heterozygote excess in the wild populations SPIE, SCUA, SICU and SPB; the cultivated populations CROC1, CROC2, SAHC1 and LH and in the silvicultural managed populations M1, M2 and SSAH1. Heterozygote deficiency was recorded in two wild populations (SPR and SSAH2) and three cultivated populations (CSAH2, CROC3 and TC). A global test of LE indicated that genotypes at one pair of loci (APAR2-12 and APARLC28, Lindsay ) are not independent (P ≤ 0.05), but identical results were found in the rest of the analysis when either of these loci was removed.

The Na ranged from 4.0 to 8.6, with the highest value found in the wild population SPIE and the lowest value recorded in the wild population SSAH2 and the cultivated population CSAH2 (Table 2). The Ne ranged from 2.37 to 5.04, with the highest value found in the wild population SPB. The values of Ho ranged from 0.2 to 0.87, with the highest value recorded for the cultivated population CSAH1. He ranged from 0.449 to 0.754, with the highest value observed in the cultivated population CLH. For uHe, estimates ranged from 0.461 to 0.774; in this case, the managed population M1 had the highest values. An interesting point to emphasize is that the cultivated population CSAH2 exhibited the lowest values in all the diversity parameters estimated (Table 2). Differences in mean values of Na, Ne, Ho, He and uHe among wild, cultivated and silvicultural managed populations were not statistically significant (P > 0.05).

Table 2.

Summary of genetic diversity estimates (means ± SE) at the population and management levels for A. inaequidens based on 10 microsatellite loci. Na, mean number of alleles per locus; Ne, mean effective number of alleles per locus; Ho, mean observed heterozygosity; He, mean expected heterozygosity; uHe, mean unbiased expected heterozygosity.

Population	Na	Ne	Ho	He	uHe
SPIE	8.600 ± 1.118	4.750 ± 0.771	0.778 ± 0.082	0.734 ± 0.042	0.747 ± 0.043
SCUA	8.300 ± 1.309	4.594 ± 0.729	0.764 ± 0.066	0.724 ± 0.044	0.738 ± 0.045
SPR	8.000 ± 1.211	4.601 ± 0.727	0.678 ± 0.096	0.731 ± 0.041	0.745 ± 0.042
SICU	8.300 ± 1.055	4.454 ± 0.537	0.742 ± 0.103	0.734 ± 0.039	0.747 ± 0.040
SPB	7.700 ± 1.230	5.041 ± 1.042	0.756 ± 0.067	0.742 ± 0.037	0.756 ± 0.037
SSAH2	4.700 ± 0.731	2.937 ± 0.458	0.311 ± 0.095	0.556 ± 0.077	0.572 ± 0.079
Mean wild population	7.600 ± 0.593	4.396 ± 0.303	0.672 ± 0.074	0.704 ± 0.029	0.717 ± 0.029
CROC1	7.800 ± 0.786	4.540 ± 0.541	0.783 ± 0.062	0.753 ± 0.025	0.767 ± 0.026
CROC3	7.500 ± 0.898	4.393 ± 0.492	0.692 ± 0.091	0.741 ± 0.032	0.754 ± 0.033
CROC2	7.400 ± 1.087	4.474 ± 0.631	0.755 ± 0.071	0.725 ± 0.044	0.739 ± 0.045
CTC	7.400 ± 0.968	4.176 ± 0.455	0.730 ± 0.082	0.733 ± 0.030	0.754 ± 0.029
CSAH1	7.300 ± 0.955	4.184 ± 0.539	0.872 ± 0.044	0.730 ± 0.029	0.743 ± 0.029
CSAH2	4.000 ± 0.683	2.372 ± 0.442	0.200 ± 0.064	0.449 ± 0.086	0.461 ± 0.088
CLH	8.500 ± 0.654	4.434 ± 0.476	0.817 ± 0.069	0.754 ± 0.023	0.769 ± 0.023
Mean cultivated population	7.129 ± 0.544	4.082 ± 0.290	0.693 ± 0.085	0.698 ± 0.041	0.712 ± 0.042
M2	8.100 ± 0.924	4.444 ± 0.469	0.802 ± 0.052	0.750 ± 0.029	0.767 ± 0.029
M1	8.000 ± 0.856	4.408 ± 0.405	0.784 ± 0.057	0.752 ± 0.026	0.774 ± 0.027
SSAH1	7.200 ± 0.646	3.682 ± 0.438	0.828 ± 0.060	0.697 ± 0.031	0.709 ± 0.031
Mean in situ managed population	7.767 ± 0.285	4.178 ± 0.248	0.805 ± 0.012	0.733 ± 0.017	0.750 ± 0.020

Genetic structure and gene flow

The global FST was 0.112 (CI: 0.0783–0.151), indicating moderate genetic structure among populations. The mean FIS across loci was significantly different from zero (FIS = 0.054; 95 % CI: 0.041–0.069), indicating inbreeding.

The UPGMA dendrogram based on genetic distances among populations (Fig. 2) did not cluster populations according to their management type, but several clusters were partially concordant with the geographic location of several populations. Populations SSAH2 and CSAH2, located in the western part of the sampling area, were separated in a group, different from the remaining populations. These populations were also separated in the dendrogram from the SSAH1 and CSAH1. The next group included the four populations of region of Sahuayo in the western part of the study area, separated from the remaining populations. Another group included population CTC and the two silvicultural managed populations M1 and M2. The remaining populations were grouped neither according to their geographic distances nor according to their management category. The regression coefficient of the linear regression analysis of FST/(1 − FST) on ln of geographic distance (Fig. 3) was positive (β = 0.015) and explained 11 % of the variation in FST/(1 − FST); the association between pairwise in geographic distances and pairwise FST/(1 − FST) values was significant (Mantel test, Z = 6065,23, r = 0.62, P ≤ 0.0010).

Figure 2.

UPGMA cluster based on genetic distances (D) estimated among 16 populations of A. inaequidens. Numbers are Nei's genetic distances, and numbers in parentheses are the results of one thousand bootstrapping random replicates.

Figure 3.

Differentiation among A. inaequidens populations. Multilocus estimates of pairwise differentiation (FST/[1 − FST]) are plotted against the natural logarithm distances (in kilometres).

The highest value of ΔK revealed that the most probable number of genetic groups was K = 5 (Fig. 4A). Figure 4B shows the proportion of ancestry of each population and individual plants based on these groups. Three wild populations and one cultivated population (SPIE, SCUA, SICU and CROC1, respectively) had a higher proportion of the orange genetic group. The green genetic group comprised two wild and two cultivated populations (SPR, SPB, and CROC3, CROC2, respectively). The brown genetic group (SSAH1 and CSAH1) and the pink genetic group (SSAH2 and CSAH2) are markedly different, because their members had a very low proportion of ancestry from other genetic groups. The remaining populations were included in the grey genetic group (CLH, CTC, M2 and M1) (Fig. 5). Overall, the grouping pattern under this analysis was very similar to that found through the UPGMA, in which the groups did not correspond to the categories of management but rather are associated with geographic location and other ecological aspects like life-history characteristics.

Figure 4.

(A) Estimated number of genetic groups (K) derived from clustering analysis using STRUCTURE. Delta K was calculated using the method described by Evanno . (B) Genetic 741 clusters obtained with five groups (K = 5). Each individual plant is represented by one vertical line with K segments coloured proportionally to their belonging to a genetic cluster.

Figure 5.

Pie charts showing proportion of ancestry assigned to individuals of each population by Bayesian clustering analysis using STRUCTURE 2.3.4.

The AMOVA carried out to test whether populations had hierarchical genetic structure according to the management type categories was not significant (RST = 0; P ≥ 0.05). When AMOVA was conducted using the STRUCTURE-defined populations, results indicated that most of the variation occurs among groups (59.32 %; P ≤ 0.001, Table 3), followed by the variation within populations (36 %, Table 3). Total differentiation was high and significant among groups (ΦCT = 0.5926; P ≤ 0.001) and within populations (ΦST = 0.6322; P ≤ 0.001), whereas differentiation among populations (ΦSC = 0.0972; P ≤ 0.001) was moderate but significant.

Table 3.

Analysis of molecular variance using RST estimates for five genetic groups obtained with STRUCTURE for populations of A. inaequidens in Michoacán State. *P ≤ 0.001.

Source of variation	Sum of squares	Components of variance	Percentage of variation	Φstatistics
Among genetic groups	1 393 355.745	2140.946	59.259	ΦCT = 0.592*
Among populations within groups	94 475.027	143.053	3.959	ΦSC = 0.097*
Within populations	1 069 513.519	1328.846	36.781	ΦST = 0.632*
Total	2 557 344.291	3612.845

The migration rate (M) among pairs of populations ranged from 0.703 to 35.25 migrants per generation. The effective population size ranged from 36 in CSAH2 to 709 individuals in the managed population M1 [see ]. The gene flow or number of migrants (Nm) was below one per generation in all cases considering CSAH2 as the receiver population, and the highest value was 29 migrants from managed populations M2 to M1. In general, the wild populations SICU and SPIE and the cultivated population CROC1 had the highest number of migrants to other populations [see ].

Discussion

Among all Agave species so far examined for genetic diversity, A. inaequidens has the highest genetic variation (wild populations averaged He = 0.703), showing estimates comparable with those recorded for the wild populations of A. potatorum (He = 0.549, Félix-Valdez ), wild populations of A. parryi (He= 0.621, Parker ; He = 0.659, Lindsay ) and wild populations of A. palmeri (He = 0.550, Lindsay ), A. utahensis (He = 0.490; Byers ) and A. utahensis subsp. kaibabensis (He = 0.408; Byers ), which have also been studied by using microsatellite loci.

The dynamic movement of propagules among populations of A. inaequidens involved in silvicultural management and cultivation could be one of the causes determining high levels of genetic diversity in these populations. Overall high levels of genetic diversity and the moderate levels of genetic structure among populations of A. inaequidens are consistent with several life-history traits of this species: predominant sexual reproduction, the outcrossing breeding system (León-Jacinto 2013) and involvement of long-distance pollinators (bats and birds) (León-Jacinto 2013).

Although morphological studies have reported significant phenotypic differences between wild and cultivated populations of A. inaequidens (Figueredo ), this pattern was not consistent with the results of our genetic structure analysis using microsatellite loci. No genetic divergence was detected among populations according to management practice, which indicates that management and artificial selection documented to occur on A. inaequidens by other studies (Figueredo ; Torres ) have not affected genetic diversity within the species. Artificial selection could be a recent process in this species or not strong enough to have caused genetic differentiation between wild populations and populations exposed to management practices, because human-induced and naturally occurring gene flow may rapidly counteract the effects of artificial selection. Also, it could be that there may be phenotypic plasticity in the genus since phenotypic differences documented previously could be largely influenced by local environments, and therefore, artificial selection would be acting on features that have a low heritability. Or, finally, the neutral markers used are not efficient to identify artificial selection on the phenotypic features evaluated.

Genetic diversity in A. inaequidens

Despite its relatively recent origin (approximately 10 my old), the genus Agave has a high adaptive radiation and an active process of speciation, with a diversification rate (0.32 species per my), only comparable with the rosetophyllous silversword (0.56 species per my) (Baldwin and Sanderson 1998; Eguiarte and Souza 2007). Often, the highest levels of genetic diversity have been found in the centre of origin, and levels decrease as long as the distance increases until the extremes of the area of distribution (Eckert ; Figueredo and Nassar 2011; Parra ). The main area of diversification of the genus Agave is central Mexico, where A. inaequidens occurs. This fact could help to explain the high levels of genetic diversity of this species.

Polyploidy is a source of genetic variation in plants. Nearly 70 % of Agave species are polyploid. We do not discard the possibility that A. inaequidens populations have different ploidy levels, similar to what has been recorded in A. cupreata (Palomino ) and A. parryi (Parker ). However, this supposition is yet to be confirmed through cytogenetic and flow cytometric studies. Generally, polyploids show the giant phenotype, which includes large size plants. It is possible that in A.inaequidens, polyploidy has an influence on one of the target features of artificial selection in favour of more productive agaves, originally for more edible matter, more aguamiel and pulque production, and more recently for more efficient mescal production. However, the microsatellites tested showed patterns typical of diploid plants.

Other aspects that can contribute to the high levels of genetic diversity and substantial gene flow are life-history traits (Hamrick ; Hamrick and Godt 1996). Agave inaequidens is self-incompatible, with a xenogamous pollination system, and bat pollinated. All these traits can promote gene flow and high levels of genetic variability within populations (Horner ; Rocha ; Eguiarte and Souza 2007; León-Jacinto 2013; Parra ).

Although, in general, we found high levels of genetic diversity in all studied populations, we identified two populations with comparably lower levels of genetic variation: SSAH2 and CSAH2. Possibly, the differences are associated with anthropogenic effects. The wild population SSAH2 has been under an intensive extraction regime for mescal production for >40 years, suffering significant population decline. In addition to agave extraction, the population has also been affected by livestock raising, because animals may consume saplings and trample seedlings, thus affecting recruitment of new plants (Torres ), as also documented for A. potatorum (Delgado-Lemus ; Torres ).

Population CSAH2 is composed of individual agaves that form a live fence, and plants of the fence are used for extracting sap. The process of sap extraction enhances asexual reproduction, forming axillary suckers, and apparently, the new plants forming the fence are clones maintained for several generations (Torres ). In other species of agave (A. tequilana, A. fourcroydes and A. parryi), it has been similarly reported that the low levels of genetic diversity are due to asexual propagation (Colunga-GarcíaMarín ; Vargas-Ponce ; Parker , 2014).

Structure and gene flow among studied localities

We recorded a significant, although slight, pattern of isolation by distance for the set of localities sampled; however, we noticed that this pattern seems mainly driven by distances from the most distant localities. In addition, Bayesian grouping appears to be explained by geographic distance and the type of ecosystem where agaves grow. These results are supported by the AMOVA, indicating that 60 % of the variation occurs among these five genetic groups. Populations from Sahuayo grow in subtropical scrub forest at the western area of Michoacán. These populations are the most genetically distinctive, and there is substantial gene flow among them. These populations could be adapted to a drier and warmer environment. Morphological differences between these populations and populations occurring in oak and pine-oak forests were reported by Figueredo .

High levels of historical gene flow were detected among most populations examined, independently of their management type. As referred to above, bat pollination makes possible a long-distance exchange of genes among populations, since Leptonycteris yerbabuenae can potentially fly nearly 100 km in one single night (Horner ). But it is also relevant to indicate that the agave handlers are moving seeds and vegetative propagules from forest to forest, from forests to plantations and among plantations themselves, thus promoting gene flow among populations. The populations with fewer incoming migrants per generation were located in the Sahuayo region, which can be explained by the relatively long distance from other populations.

Implications for management and conservation

An important aspect associated with management of agaves is that cultivation is conducted mainly in agroforestry systems where sexual and asexual propagules and even plantlets transplanted from different wild and managed sources are maintained together. This management system determines that these managed populations represent important reservoirs of genetic variation. Some handlers sow seeds and then transplant saplings to the forest (Torres ), which increases the probability of establishment of new agave plants in forests with this management regime. In addition, cultivation of agave together with other crops and other wild plants allows high plant heterogeneity, which can significantly contribute to control pests, an aspect that has been documented for other Agave species (Zizumbo-Villarreal ) and crops (Root 1973). Populations of the western area that grow in drier ecosystems are more vulnerable, particularly because of low population density, continual plant extraction and absence of plans for population recovering. The forms of management developed in the eastern areas should be considered for adoption in the western areas. Exchange of management experiences is possible and necessary. This study identified some effects of the management regimes carried out in the eastern zone of the study area. It is still necessary to evaluate the population genetic patterns of A. inaequidens in other regions of Mexico. Particularly interesting would be those areas where the species is used mainly for the production of fibre, which apparently is an activity that has decimated the agave populations, and where no reforestation or cultivation planning has occurred (Valenzuela-Zapata ).

Conclusions

Agave inaequidens is a species with high levels of genetic diversity along the entire management gradient. Population genetic structure in relation to management practices was moderate, because management is relatively recent, and artificial selection can be rapidly counteracted by high gene flow.

Although we recorded high levels of genetic diversity in general in this species, some populations appear to be at risk, because of low population sizes and levels of genetic diversity, in addition to the continuous removal of individuals for the production of mescal, without a strategy for reforestation or farming.

Sources of Funding

This research was funded by Posgrado en Ciencias Biológicas, UNAM and the Consejo Nacional de Ciencia y Tecnología CONACYT, the Programa de Apoyo a Proyectos de Investigación e Innovación Tecnológica (PAPIIT IN205111-3 and IN209214) and CONACYT (CB-2013-01-221800).

Contributions by the Authors

C.J.F., main author, involved in the study design, conducting of field and laboratory work, literature review and general data collection, systematization and analysis, wrote the first draft and concluded the final version of this paper. A.C., main coordinator-supervisor of the research project, contributed with original data and the designing of all the researches providing the information for the current analysis; participated in systematization and analysis of data and reviewed several drafts of the manuscript and concluded the final version this paper. P.C.-G., J.M.N. and A.G.-R. contributed to designing and following progress of the research and fieldwork and data analyses. All authors contributed to writing, reading and approved the final manuscript.

Conflict of Interest Statement

None declared.

Supporting Information

The following additional information is available in the online version of this article –

Appendix S1. Migration rates (M = m/µ) paired in both directions and parameter Theta (Θ = 4Neµ) for 16 populations of A. inaequidens. + = receiving population. Top row of the table shows donor populations, and first column to the left shows populations that receive migrants.

Appendix S2. Effective population size (Ne) and the number of migrants (Nm) for 16 populations of Agave inaequidens. + = receiving population. In the table, the donor population is shown in the top of the table and the population that receives is shown in the left of the table.