Study area
We performed this study in San Nicolas Tecoaco village (20° 2′ 38.2ʺ N, 98° 35′ 16ʺ W), Hidalgo State, central Mexico, from March 2014 to October 2014. This location has an annual average temperature of 16 °C and an average altitude of 2600 m above sea level. The type of vegetation occurring in this area is classified as a xerophilous shrubland42.
Study species
Opuntia robusta (Cactaceae) is an endemic plant found in Meridional Altiplano, México43, which exists in the following three sexual forms: hermaphrodite, dioecious (male and female), and trioecious44. In a parallel study, Sandoval-Molina45 found that the most common herbivores of this plant were leaf-footed bugs, Chelinidea sp., Narnia sp. (Hemiptera: Coreidae), the cactus long-horned beetle, Moneilema sp. (Coleoptera: Cerambycidae), and mining insects. Before 2017, this population was considered to be gynodioecious; thus, we did not collect samples from male individuals in this study. In 2018, fewer than 15 male individuals were reportedly present in a population of more than 800, and most of these were hermaphrodites (Supplementary Information).
Determination of plant sex
White empty anthers, short style, and well-developed lobular stigma characterised female flowers, while a relatively longer style compared to that of the female and functional anthers characterised hermaphrodite individuals44.
Comparison of tissue cost between female and hermaphrodite individuals
In March 2017, we undertook a census in San Nicolas Tecoaco, to identify the number of female and hermaphrodite plants with cladode and flower sprouts from the set of plants studied in the previous years. We selected 1–2 m tall plants, located 5–10 m apart for sampling. Finally, we randomly selected 19 plants (eleven female and eight hermaphrodite individuals) bearing flower buds and young cladodes on different branches for analysis and tagged the cladodes and flower sprouts using a permanent marker. We marked the flower sprouts on the adjacent side of their parental cladode surface.
Between March 2017 and June 2017, we obtained sufficient data to estimate the relative growth rates of the species, in order to explore possible differences in the energy costs of cladodes and flower buds between the two sexual forms of O. robusta. We measured the length, width, and thickness of each cladode and flower bud twice during the study, once at the beginning, and once at the end of the study. Additionally, we also measured the lengths of the flowers from the base to the beginning of the sepals. Since the flower buds were spherical, we considered the thickness to be equal to the width. Subsequently, we calculated the flower volume immediately after the emergence of cladodes and flower buds, and the final volume after anthesis. We estimated the initial and final volumes (Vx) of the cladodes using the formula Vx = ((a/2))/((b/2)π)c, and those of the flowers using the formula Vx = 4/3πa2b. Here, x represents the time of measurement (initial or final), a and b represent the major and minor axes of the ellipsis, while c represents the cladode thickness. We measured all estimators to the nearest 1.0 mm and represented values in centimetres. We estimated the relative growth rate (RGR) using the formula proposed by Hunt46: RGR = (lnVf – lnVi)/(t2 – t1). Here, Vf represents the final volume [cm3], Vi represents the initial volume [cm3], t1 represents the initial time [day], and t2 represents the final time [day].
We compared relative growth rate data using a generalized linear model (GLM) with gamma error distribution in the R software, using the log link function47. The explanatory variables included sex, type of structure, and their interactions. We performed partial regression using the ggeffects package in R48.
We obtained meteorological variables, including total precipitation [mm], maximum temperature [°C], minimum temperature [°C], mean temperature [°C], global radiation [W(m2)−1], relative humidity [%], reference evapotranspiration [mm], and potential evapotranspiration [mm] for the Singuilucan municipality from March–October 2014, from the official Mexican Government weather station database of the Instituto Nacional de Investigaciones Forestales Agrícolas y Pecuarias49. We summed up the data for the per-day total precipitation, and that for the reference and potential evapotranspiration, from the beginning of each month through the sampling day. In the months (March, April, and May) or days when values from the meteorological database were underestimated, we averaged the values for the closest preceding and following days. If we lacked the data for more than one day and the data for such days could not be acquired, we considered a repetition of the averaged value for the days for which we lacked data, between the existing days. For July, we considered the values for the previous day (11/07/14), since we lacked the data for the days on which sampling was performed and the subsequent days. For the additive variables (total precipitation, reference, and potential evapotranspiration), we summed up data for 30 days, excluding data for one day, for the 31-day period.
To determine the effects of the environmental variables on the concentration and presence/absence of secondary metabolites, we used R to formulate a structural equation model (SEM) in piecewiseSEM47,50. For concentrations, we fitted linear mixed-effects models using the nlme package51 and used the plant ID as a random factor. To evaluate the presence or absence of substances, we fitted generalized linear models with binomial error distributions and logits as the link functions. The concentration and presence/absence of 4-HBA, CGA, and QUE were dependent variables, and total precipitation, average temperature, global radiation, relative humidity, and potential evapotranspiration were explanatory variables. We analysed the sexes separately, and the substance concentration variables were log + 1 transformed. We assessed the goodness-of-fit using the Fisher function in the piecewiseSEM package50, where a larger p-value implies better data adjustment to the model. We conducted a visualisation of the SEM models using Biorender52, flaticon53, and CorelDRAW54.
We estimated fruit traits (biomass [g], volume [cm3], and tissue density [g × cm−3]) and the number of fruits eaten by fructivores and compared them between the sexual forms using data reported by Janczur et al.18. The former comparison enabled the assessment of the possible differences in reproduction per fruit biomass between the sexual forms. The latter comparison enabled the assessment of the differences in preference for fruits eaten by animals in relation to the different sexual forms, and thus, the mechanisms by which this may increase the probability of seed dispersal. Higher zoochory of one sexual form may occur not only because of differences in fruit biomass density [g × cm−3], but also because of differences in the volatile substance content between the sexual forms.
To test the effects of sexual form on fruit traits, we used generalized linear models in R. To analyse the number of fruits eaten, we used the negative binomial error distribution and log link function, and the Gaussian error distribution and identity link function for the other fruit traits47,55,56. We performed all post-hoc contrasts for fruit traits using the emmeans package47,57, and generated plots using the ggplot R package47,58. We compared the average number of fruits produced by the two sexual forms using the Kruskal–Wallis test.
Comparison of secondary metabolite occurrence/concentration between female and hermaphrodite individuals
We obtained plant samples for secondary metabolite analysis using 100 m long Canfield lines, which were parallel to the contours of the hill and located 60 m from each other, and selected plants that were located near the lines and were 10 m apart for analysis. We randomly assigned each plant to one of the eight groups established herein, with three female plants and twelve hermaphrodite plants. The uneven number of individuals of each sex was attributable to the low proportion of females in the population. We tagged examined cladodes on their surface using a permanent marker.
We used a stainless-steel punch (Ø = 0.5 cm) to remove two samples of vegetative tissue from cladodes belonging to the same order of each plant. We perforated the mid-section of the arc delimited by the border of the upper quarters of the cladodes, approximately 1 cm away from the edge. We placed samples in labelled Ziploc bags, stored them in a cooler containing ice, and then transported them to the laboratory in a portable refrigerator at − 20 °C. The samples were stored in the laboratory at − 40 °C until extraction.
We performed homogenisation of approximately 1 g of the sample containing the cuticle in 35 mL of 100% methanol in an ultrasonic 6 L bath for 30 min at room temperature (21 °C). We filtered the methanol extracts, placed them in amber bottles, and stored the bottles at − 20 °C until further analysis59. We determined the types and concentrations of secondary metabolites in these tissues using high-performance liquid chromatography (HPLC), in accordance with the procedure described by Janczur and González Camarena59, using the following: Waters 717 liquid chromatograph with autosampler, Waters 2487 HPLC Absorbance UV–Vis Detector, Waters 1525 Binary HPLC Pump, Waters control module with SAT/IN Bus (Waters, Milford, MA, USA), Symmetry HPLC C18 column (particle size 5 µm, length 250 mm, internal Ø = 4.6 cm; Waters, Milford, MA, USA). We filtered the extracts using a 0.45 µm pore size nylon-membrane filter. The mobile phase consisted of 0.1% v/v acetic acid (A) together with 100% acetonitrile (B). For the mobile phase A, we dissolved 1 mL of glacial acetic acid with HPLC water, until the volume was 1 L. For the mobile phase B, we used 100% acetonitrile. We filtered both mobile phases using a 0.45 µm nylon membrane. We degasified them with an ultrasonic bath for 30 min. We set the column temperature at 25 °C, used the 254 nm UV detector, and established the flow of the mobile phase, injection volume, and run time as 0.2–0.8 mL/min, at 8 µL, and 35 min, respectively. To wash the piston seals, we used MeOH : H2O (60 : 50). To generate the calibration curves, we used standards for salicylic acid (SA), 4-hydroxybenzoic acid (4-HBA), chlorogenic acid (CGA), and quercetin (QUE) (Sigma-Aldrich). We generated the following calibration curves: yi = 1109.4xi + 481.67, yi = 296.01xi + 133.74, yi = 551.41xi + 263.64, and yi = 919.96xi + 201.64; here, yi represents the area below the absorbance curve, xi represents the concentration of the secondary metabolite, and i = 1, 2, 3, and 4 for 4-HBA, CGA, QUE, and SA, respectively. SA was not present in any of the samples tested (Table S1 online31).
We used a logistic regression model to test the effect of the sexual form, month of study, cladode age category, cladode size, the number of cladodes above a given cladode, and the cladode order above the soil level, on the probability of detecting secondary metabolites in the cladodes. Since the latter data were ordinal, the sexual form and month were considered as discrete variables and treated the other traits as continuous variables60. We applied the generalized linear mixed model (GLMM) with a logit link function [ln(P/(1-P)], where P indicated the probability of detecting a given metabolite, binomial response distribution, maximum likelihood estimation technique, Newton–Raphson optimisation algorithm, and Person Chi-Square/df fit criterion. We used the GLIMMIX procedure in SAS statistical software61 (Methods S1).
We used generalized linear models (GLMs) in R47 to determine the relationship between cladode length, width, thickness, months, age, cladode order from the soil, and cladodes above a given cladode, and the concentrations of the different secondary metabolites. Since many concentrations were null, we analysed only the positive concentrations (Methods S1).
Comparison of damage between female and hermaphrodite individuals
We used the same plants as those used for relative growth rate analysis. We analysed the extent of damage caused by herbivorous insects on both sexes of O. robusta from March–June 2017. We selected two branches, one with flowers and the other with cladodes, from each plant. We estimated two types of damage caused by herbivores using image analysis, to determine the total percentage of tissue removed and other types of damage, such as scars or necrosis. We acquired photographs of one randomly selected face of each structure, using a Nikon D3200 with an AF-S DX NIKKOR 18–55 mm f/3.5–5.6G VR lens (Nikon Corporation, Tokyo, Japan) mounted on a tripod, using a 1-cm piece of millimetre paper as a reference for size. We analysed all images using ImageJ62 to estimate the total proportion of damaged areas.
We analysed data on herbivore damage and other damages using a GLM procedure with the Gaussian error distribution and identity link function47 in R. The response variables were the logit transformed proportion of damage (ln[P/(1-P)]), where P represents the proportion of tissue damaged. In our statistical models, the transformation improved the distribution of residuals. The explanatory variables were sex, type of structure, and their interactions. We performed partial regressions using the ggeffects package in R48.
Comparison of the occurrence/concentrations of secondary metabolites between younger and older vegetative tissues
We named the oldest cladodes (closest to the soil) as ‘first-order cladodes,’ those growing on the oldest cladodes as ‘second-order cladodes’ etc. We selected each plant branch with the largest number of cladodes. We measured the length, width, and thickness of each cladode. We sampled vegetative tissues from plants belonging to each of the eight groups; the first group on the 10th March, the second group on the 12th April and so on, through the 10th May, 14th June, 12th July, 10th August, 13th September, and 11th October 2014. We measured the length and width of each cladode to the nearest 0.5 cm, using a measuring tape, and their thickness to the nearest 0.01 mm, using a calliper. We conducted the latter measurement in the apical part of the cladodes in the case of apical cladodes, or at the point of ramification of the daughter cladode when it grew on its apex.
During eight years of observations prior to the commencement of this study, we observed that the age of the cladodes in the studied zone could be estimated by examining the following colour patterns of their spines: 1—yellowish, 2—yellow, white base, 3—white-yellowish, 4—white, 5—greyish, 6—black, with ‘1’ being the youngest, and ‘6’ being the oldest. We assigned each cladode to one of the classes. We used the HPLC procedure described by Janczur and González Camarena59 to determine the concentrations of different secondary metabolites in the plant tissues.
To test whether different estimators of cladode age were parallel (to test whether younger cladodes were mostly apical, and thus bore fewer cladodes above), we examined the relationship between the cladode order from the soil or cladode number above a given cladode and cladode age, using ordinary least squares regression (OLS). We used a numerical algorithm applied to the SMATR software for R63. We included a test for the determination of the effects of cladode age estimators on the SMSs occurrence/concentration in the same GLM models, as described in the previous section.
Trade-off between investments in defence, growth, and reproduction
We tested the relationship between cladode length and cladode order or cladode age to determine whether cladode size was parallel to cladode age. We performed OLS analysis and slope comparison between sexual forms using the Wald test (WT—test statistic) and tested the significance of differences between the intercepts. We used a numerical algorithm applied in the SMATR software63. To estimate the relative investment in growth and reproduction, we counted the number of flower and cladode buds on parental cladodes of the same plants used in the study performed by Sandoval and Janczur (Dataset online29). We used generalized linear models in R, with a negative binomial error distribution and log link function47,55,56, to test the effects of sexual form on the average number of flower and cladode buds. Significant differences between the number of flowers and cladodes for certain sexual forms implies a higher relative reproductive investment.
We used the same method of quantification for the standardized major axis and GLM models for intersexual comparisons, as described in the previous Sect. 59. For example, larger relative allocations for reproduction and secondary metabolites together with lower allocation to growth in one sexual form, compared to lower allocations for reproduction and secondary metabolites, and higher allocations for growth in the other sexual form imply that the production of secondary metabolites does not compete with either growth or reproduction; rather, growth competes with reproduction, and allocation to the production of secondary metabolites is an outcome of the gain in terms of fitness from such an allocation.
Effects of the existence of trade-offs between different secondary metabolites on the predictions of the plant defence hypothesis
We used ordinary least squares regression (OLS), coefficient of determination, and t-tests to determine the existence of possible trade-offs in the proportion of cladodes harbouring different secondary metabolites. We performed the t-test to determine the significance of correlation between cladode order and cladode age64.
Ethics statement
This research did not involve any human or animal measurements. We obtained permission from the head of the Singuilucan municipality, State of Hidalgo, Mexico, to conduct research activities at the selected sites of the municipality. The owners of the lands permitted us to conduct the study and were informed of the permission granted by the municipality. MKJ obtained a permit (09,448/14) from the Ministry of Environment and Natural Resources of the United States of Mexico (SEMARNAT), which stated that no permission is necessary to conduct field studies on plants belonging to the genus Opuntia. The study site was not considered to be a protected area65, and O. robusta was not considered to be an endangered species66. During this study, we did not affect or involve any endangered species. As we did not sample all plants, we did not deposit specimens in a public herbarium. No plant was killed or severely damaged as a result of our research activity; the plant material used for this study was sampled at a limited scale, and therefore, the sampling presented with negligible effects on the functions of the broader ecosystem. All the methods were carried out in accordance to relevant guidelines and regulations.
Source: Ecology - nature.com