Introduction
In fisheries and marine biology, appraisal of reproductive potential in broadcast-spawning marine invertebrates such as most cnidarians, molluscs and echinoderms conventionally considers various gonad-related traits at the population level. These include sex ratio, frequency of the gonad developmental stages, gonad size increase, and the intensity and period of the reproductive season (‘breeding seasonality’, ‘reproductive event’). Such data give essential information about the reproductive seasonality of these taxa, which can often be of economic significance, especially if they are fishery resources.
Since these organisms are sedentary, reproductive success is characteristically an increasing function of gamete release via broadcast spawning, and such organisms typically have enormous gonads (Parker et al., Reference Parker, Ramm, Lehtonen and Henshaw2018). Their reproductive potential (‘reproductive activity’, ‘reproductive effort’, ‘energy invested in reproduction’ or ‘reproductive output’) is often measured simply through gonad weight expressed as a percentage of total body weight, i.e. ‘gonadosomatic index’ or ‘GSI’, where the highest GSI value corresponds in most cases to maximum (peak) ripe stage (i.e. maximum gonad investment, ‘MGI’) and determines the reproductive season. However, there are some exceptions; for example, in sea urchins, the gonads first increase in mass due to the accumulation of nutritive phagocytes, which later transfer resources to the gametes (Unuma & Walker, Reference Unuma and Walker2009; Walker et al., Reference Walker, Lesser and Unuma2013; Benítez-Villalobos et al., Reference Benítez-Villalobos, Avila-Poveda, Diaz-Martinez and Bravo-Ruiz2015; Kalachev & Yurchenko, Reference Kalachev and Yurchenko2017; Eckelbarger & Hodgson, Reference Eckelbarger and Hodgson2021). In both sexes, a cautious interpretation of the highest GSI value is, therefore, necessary: the gonads increase greatly in mass due to accumulation of nutritive reserves, or development of mature gametes, or both effects simultaneously, according to the cellular dynamics of reserve exchange between the nutritive and gametogenic cells (Harrington et al., Reference Harrington, Walker and Lesser2007; Hernandez et al., Reference Hernandez, Vázquez, Torruco and Rahman2020; Eckelbarger & Hodgson, Reference Eckelbarger and Hodgson2021).
Invertebrate reproduction studies generally attempt to understand the duration of the reproductive season at the population level in response to specific factors of each sampling location. Since external sexual dimorphism is often not evident, intrinsic sex differences often remain unstudied. However, the gonad-size increase is a highly fitness-related trait in both sexes of many marine invertebrates with broadcast spawning (Mercier & Hamel, Reference Mercier and Hamel2009; Abadia-Chanona et al., Reference Abadia-Chanona, Avila-Poveda, Arellano-Martinez, Ceballos-Vazquez, Benítez-Villalobos, Parker, Rodriguez-Dominguez and Garcia-Ibañez2018; Parker et al., Reference Parker, Ramm, Lehtonen and Henshaw2018; Avila-Poveda et al., Reference Avila-Poveda, Abadia-Chanona, Alvarez-Garcia and Arellano-Martinez2021). Broadcast-spawning marine invertebrates regularly exhibit a maximum GSI that is approximately equal in both sexes under high levels of sperm competition or sperm limitation (Levitan, Reference Levitan, Westneat and Fox2010; Parker et al., Reference Parker, Ramm, Lehtonen and Henshaw2018; Parker, Reference Parker2020). However, each sex can exhibit different life-history strategies when males and females have different optimal gonadal investments (Charnov, Reference Charnov1982; Levitan, Reference Levitan2005; Parker, Reference Parker2014; Parker & Pizzari, Reference Parker, Pizzari and Hoquet2015; Parker et al., Reference Parker, Ramm, Lehtonen and Henshaw2018). In echinoderms of the genus Isostichopus and Holothuria, there is generally a female-biased GSI asymmetry, though male-biased GSI does occur (see table 1 and figure 1 in Parker et al., Reference Parker, Ramm, Lehtonen and Henshaw2018).
Sampling locations of H. inornata in the southern Mexican Pacific. CC, Caleta de Campos, Michoacán (18°04′N 102°45′W); PM, Puerto Madero, Chiapas (14°42′N 92°24′W).

Some examples of population variation in holothurian maximum GSI with latitude

a Both sexes combined; – without data.
In the build-up to spawning, gonads typically increase greatly in mass as the gametes are generated, develop and mature physiologically, with MGI occurring just before gamete release. Spawning physiology is synchronous in many broadcast spawning invertebrates, an adaptation related to increasing fertilization success (McHugh & Rouse, Reference McCune, Grace and Urban1998), quantifiable through the population frequencies of the gonad developmental stages (GDSs), particularly the ripe and spawning stages (i.e. the reproductive phase or reproductive season). The intensity and period of MGI (representing a mass investment) reflect the gonad developmental stages (representing physiological investment) during the reproductive season. Variation in intensity (GSI %) and period (months) of the MGI peaks determine the reproductive season, setting up cycles that can be annual, semi-annual or continuous (Brewin et al., Reference Brewin, Lamare, Keogh and Mladenov2000; Oyarzun et al., Reference Oyarzun, Toro, Garces-Vargas, Alvarado, Guiñez, Jaramillo, Briones and Campos2018; Machado et al., Reference Machado, Moura, Pereira, Vasconcelos and Gaspar2019; figure 1 in Avila-Poveda et al., Reference Avila-Poveda, Abadia-Chanona, Alvarez-Garcia and Arellano-Martinez2021) among populations of the same species.
The present investigation considers the sea cucumber Holothuria (Halodeima) inornata (Echinodermata: Holothuroidea). This ectothermic species constitutes a potential fishery in Mexico, despite not being referenced in law. As a ubiquitous member of neritic benthic communities in the Eastern Pacific Ocean, including the Mexican tropical Pacific, it represents an ideal broadcast spawner model for examining how MGI is modulated by sex. We tested four aspects of population-level reproductive seasonality separately for each sex by investigating: (i) the intensity and duration of GSI variations, (ii) the frequency of specimens in each GDS through time, (iii) sexual asymmetry in GSI by examining the natural logarithm of the ratio of male MGI to female MGI (i.e. RGES: see Abadia-Chanona et al., Reference Abadia-Chanona, Avila-Poveda, Arellano-Martinez, Ceballos-Vazquez, Benítez-Villalobos, Parker, Rodriguez-Dominguez and Garcia-Ibañez2018; Parker et al., Reference Parker, Ramm, Lehtonen and Henshaw2018), (iv) sexual asymmetry in gonads development, measured through the natural logarithm of the ratio of the frequency of each of the GDSs per sex. The adult sex ratio (ASR) was also estimated. ‘Adult’ refers to an individual physiologically capable of producing offspring. ASR is here equivalent to the operational sex ratio, the ratio of sexually active males to females, i.e. individuals exhibiting the ripe/spawning GDS (see Szekely et al., Reference Szekely, Weissing and Komdeur2014). Additionally, we evaluated the possible association of GSI with the local and regional seawater temperature in terms of the sea surface temperature (SST) and the Oceanic Niño Index (ONI).
Materials and methods
Sampling areas, specimen collection and measurements
Two coastal locations in the southern Mexican Pacific were sampled, distant from each other by four latitudinal degrees (Figure 1): Caleta de Campos, Michoacán (CC; 18°04′N 102°45′W) and Puerto Madero, Chiapas (PM; 14°42′N 92°24′W). Both locations were sampled during two consecutive years (May 2011–June 2013). Although the sampled months were not the same in both populations, they are consecutive, and for a given locality, these sampled months have similar thermal deviations from the El Niño-Southern Oscillation (Oceanic Niño Index) ranging from −1.1 (mildly cold) to 0.3 (regular climatic event: NOAA, 2017); so we anticipate that their thermal effects on gonadal development will have been consistent across each of the two populations, and hence comparable concerning regional temperature.
A mean of 20 H. inornata individuals of adult status (54 ≤ TBW ≤ 390 g in total body weight) were randomly collected monthly at each location by scuba diving from a depth interval of 4–10 m. Unfortunately, due to adverse sea conditions and/or logistical difficulties, some monthly sampling events were prevented between January and April 2012 in Caleta de Campos.
Individuals were anaesthetized by submersion in seawater with 5% KCL added for at least 10 min at room temperature before dissection (see Lincoln & Sheals, Reference Lincoln and Sheals1979). The entire gonadal tubules (here categorized simply as gonads) were removed by opening the ventral side of each specimen, and the body, viscera and gonads were allowed to drain until damp-dried (around 5 min) before being weighed (Benítez-Villalobos et al., Reference Benítez-Villalobos, Avila-Poveda and Gutierrez-Mendez2013). The total body weight (TBW, including viscera and gonad) and the gonad weight (GW) were recorded using an analytical balance (± 0.001 g). Later, the gonads were fixed in 10% neutral formalin-saline solution for 2 weeks and later preserved in 70% ethanol for histology (Humason, Reference Humason1962).
Assignment of the microscopic gonad developmental stages (GDS)
Histological procedure and microscopic examination previously used with H. fuscocinerea were followed, in which gonad developmental stages (GDS) for each sex were assigned to one of the following five GDS (Benítez-Villalobos et al., Reference Benítez-Villalobos, Avila-Poveda and Gutierrez-Mendez2013): 1 – gametogenesis, 2 – ripe, 3 – spawning, 4 – post-spawning and 5 – gonadal recovery or ‘sexually undifferentiated’. These GDS take account of three criteria: development of sex cells, the wall thickness of gonadal tubules, and regenerative or resorptive processes of gonads (for criteria to establish GDS, see Benítez-Villalobos et al., Reference Benítez-Villalobos, Avila-Poveda and Gutierrez-Mendez2013). The 5-gonadal recovery or ‘sexually undifferentiated’ stage is characterized by gonads with a macroscopically crystalline appearance and very short tubules (about 6 mm in length). At the same time, microscopically, they comprised only a large amount of connective tissue containing germinal cells, often termed gonadal-basis (Benítez-Villalobos et al., Reference Benítez-Villalobos, Avila-Poveda and Gutierrez-Mendez2013).
Mass investment: intensity and period of maximum GSI values (MGIs)
The intensity of the reproductive season was examined by the traditional gonadosomatic index (GSI = [GW/(TBW – GW)] × 100), where GW is gonad weight and TBW is total body weight. Peak GSI values obtained indicate that most individuals have gonads full of mature gametes, i.e. maximum gonad investment (MGI) or maximum ripe stage.
The period of the reproductive season refers to the timing over which animals reproduce, characterized by a high proportion of individuals in ripe and spawning stages (i.e. the reproductive phase) (see figure 1 in Avila-Poveda et al. Reference Avila-Poveda, Abadia-Chanona, Alvarez-Garcia and Arellano-Martinez2021). In GSI terms, it refers to the group of consecutive months showing no significant monthly differences in the highest value of GSI.
Physiological investment: development and period of the GDSs
Gonad development was measured as the percentage of specimens in each GDS throughout the sampling period, indicating the ‘entire’ reproductive cycle. Consequently, in broadcast-spawning invertebrates, an entire reproductive cycle refers to the time (months) in which the full succession of GDSs is completed by the population, therefore, in that period of time all the GDSs appear simultaneously for each sex (e.g. Giese, Reference Giese1959). The entire reproductive cycle was summarized as two phases to more easily visualize trends for males and females: (1) the non-reproductive phase with the GDS-4 (post-spawning), GDS-5 (recovery) and GDS-1 (gametogenesis) amalgamated, and (2) the reproductive phase with the GDS-2 (ripe) and GDS-3 (spawning) amalgamated (see figure 1 in Avila-Poveda et al., Reference Avila-Poveda, Abadia-Chanona, Alvarez-Garcia and Arellano-Martinez2021). Therefore, in broadcast-spawning marine invertebrates, a reproductive phase (or reproductive season) forms part of an entire reproductive cycle, but it is also setting up its own cycle that can be annual, semi-annual or continuous.
Sampling periods in which more than 50% of specimens occurred in the reproductive phase, i.e. joint duration of GDS-2 (ripe) and GDS-3 (spawning) stages were considered to be a reproductive event (the reproductive season).
Sexual asymmetry in MGIs: mass investment
While overall mean GSI values can give a good indication of the reproduction season at the population level, much can be gained by comparing the gonad expenditure of the two sexes in broadcast-spawning marine invertebrates; sexual asymmetry in MGI can be related to the selective forces of sperm competition and sperm limitation (Abadia-Chanona et al., Reference Abadia-Chanona, Avila-Poveda, Arellano-Martinez, Ceballos-Vazquez, Benítez-Villalobos, Parker, Rodriguez-Dominguez and Garcia-Ibañez2018; Parker et al., Reference Parker, Ramm, Lehtonen and Henshaw2018). Since both sexes may show large variations in GSI with locality (see Table 1), a ratio of male:female GSI measures is more informative. The relative gonad expenditure of the sexes (RGES) was calculated for H. inornata as the natural log of male GSI/female GSI according to Abadia-Chanona et al. (Reference Abadia-Chanona, Avila-Poveda, Arellano-Martinez, Ceballos-Vazquez, Benítez-Villalobos, Parker, Rodriguez-Dominguez and Garcia-Ibañez2018) and Parker et al. (Reference Parker, Ramm, Lehtonen and Henshaw2018). This index gives a symmetric estimate of sexual asymmetry in GSI around zero for equal expenditure, with three possible outcomes (Parker et al., Reference Parker, Ramm, Lehtonen and Henshaw2018): M > F (male-biased asymmetry, RGES > 0), F = M (equal GSI; RGES ≈ 0), and F > M (female-biased asymmetry, RGES < 0). The standard deviation (σ) of the ratio RGES was estimated according to the method of Karlen et al. (Reference Karlen, McNair, Perseguers, Mazza and Mermod2007) and Wolter (Reference Wolter2007), and ln(σ) was used to correspond with ln in RGES (see Abadia-Chanona et al., Reference Abadia-Chanona, Avila-Poveda, Arellano-Martinez, Ceballos-Vazquez, Benítez-Villalobos, Parker, Rodriguez-Dominguez and Garcia-Ibañez2018).
Sexual asymmetry in GDSs: physiological investment
Reproductive cycles in marine invertebrates generally show variations in frequencies of each microscopic GDSi between sexes, particularly ripe and spawning stages during each MGI peak. Therefore, we used an index that illustrates how the GDS data translates into a relative gonadal physiological investment of the sexes in reproduction (RGPIS) to estimate sexual asymmetry in the frequency of GDS, particularly the ripe (GDS-2) and spawning stages (GDS-3), i.e. the reproductive phase. This index was calculated as the natural log-ratio between sexes of the frequency of each microscopic GDSi (FGDSi) through time, i.e.:
where subscripts M and F refer to males and females, respectively. This index again gives three possible outcomes: M > F (male-biased asymmetry, RGPIS > 0), F = M (equal, RGPIS ≈ 0) and F > M (female-biased asymmetry, RGPIS < 0).
Adult sex ratio ‘ASR’
The sex ratio of individuals that could be determined as male or female in each population was calculated monthly and annually. To assess monthly deviations from expected equality (1:1 ratio), the χ2 goodness of fit test adjusted by Yates’ correction (χ2c) was used (Yates, Reference Yates1934; Zar, Reference Zar2010).
Environmental temperature
The main regulator of the physiological aspects of reproduction in ectotherm animals is the temperature of the surrounding seawater. Taking this into account, for both geographic locations, we used the sea surface temperature (SST) and the Oceanic Niño Index (ONI) as proxy variables of the thermal conditions experienced by surface-benthic sea cucumbers to assess associations between temperature (°C) and reproductive intensity (i.e. GSI%). The available dataset was gathered from NOAA (2017, 2018). ONI index gives a symmetric estimate of warm and cold periods for the El Niño–Southern Oscillation features consistent with coupled ocean-atmosphere phenomenon periods. A regular period is defined as ONI index values between 0.4 and −0.4, a cold period exhibits values of −0.5 downward and a warm period shows values of 0.5 upward, according to the description in NOAA (2017).
Data analysis
GSI (proportion data) were normalized by arcsine square root transformation (McCune et al., Reference McCune, Grace and Urban2002) before testing with ANOVA or factorial ANOVA to evaluate differences in GSI according to sex, months, geographic locations and temperature (SST and ONI), followed by Fisher's Least Significant Difference (LSD) test post hoc for multiple comparisons. Statistically significant differences were accepted at P < 0.05 (Sokal & Rohlf, Reference Rubilar and Crespi-Abril1995). Statistical analyses were carried out using Statistica® and SigmaPlot ®. Data are reported as mean ± standard deviation (SD).
Results
Population characteristics: size and temperature effects on reproductive intensity
The size (total length ‘TL’) frequency distribution of the sea cucumbers ranged between 13 and 36.5 cm at CC, Michoacán (N = 115, median = 19.3) and between 7 and 28.3 cm at PM, Chiapas (N = 123, median = 17). There are significant differences in the TL median values between the two populations (Mann–Whitney U tests, U = 4810, P = <0.001).
In both geographic locations, the reproductive intensity (i.e. GSI%) was weakly inversely correlated with TL (CC, r = −0.19, P = 0.04; PM, r = −0.21: P = 0.02) (Figure 2A). Changes in reproductive intensity were not underpinned by changes in local temperature (CC, r = 0.12, P = 0.19; PM, r = 0.10, P = 0.24) (Figure 2B). However, reproductive intensity did show a strong positive correlation with regional temperature at CC (r = 0.69, P = 0.00), but not at PM (r = 0.05, P = 0.53) (Figure 2C).
Reproductive intensity (GSI %) in relation to size (A), local temperature by SST (B) and regional temperature by ONI index (C) for the populations of H. inornata in Caleta de Campos, Michoacán (18°04′N), and Puerto Madero, Chiapas (14°42′N). In the ONI index, the regular period lies between 0.4 and −0.4, cold period has index values of −0.5 downward, and warm period has values of 0.5 upward (NOAA, 2017). Each scatterplot shows a linear fit and regression bands with a confidence level of 0.95.

Mass investment: intensity and period of maximum GSI values (MGIs)
At Caleta de Campos (18°04′N), the population GSI values ranged between 0.01 and 3.46. The effect of sex on GSI is highly significant (F 1,101 = 10.09, P = 0.0019) and summed throughout the entire sampling period, mean male GSI was higher than the mean female GSI (Figure 3A, right panel). Only one MGI for both sexes occurred in June (males 1.46 ± 1.03, females 0.66 ± 0.33), with male MGI significantly higher (F 1, 20 = 5.4509, P = 0.0301). The sampling months with the highest male GSI were significantly different (P < 0.05), confirming that June is the only month of peak reproduction for males. Contrastingly, in females, the highest GSI sample months showed no significant difference between June and July (P = 0.23) and June and August (P = 0.48), suggesting three months of peak reproduction by females (Figure 3A); i.e. a reproductive season of 1 month in males to 3 months in females.
Intensity (%) and period (months) of the reproductive season for the populations of H. inornata in (A) Caleta de Campos, Michoacán (18°04′N), and (B) Puerto Madero, Chiapas (14°42′N). The dark grey areas represent the principal maximum gonad investment (MGI). The light grey areas represent the reproductive season in each site according to the rule of no significant differences between maximum GSIs (P > 0.05, see text). Means ± SD are shown for females and males (for clarity, population means are omitted). The right panels show the overall effect of sex on GSI (combined data for all samples at each location).

At Puerto Madero (14°42′N), the population GSI values ranged between 0.12 and 13.49. In this site, the effect of sex on GSI summed throughout the entire sampling period was not significant (F 1,108 = 1.11, P = 0.2937), though females exhibited a slightly higher mean GSI than males (Figure 3B, right panel). July was the main MGI month for both sexes (females 4.79 ± 4.40, males 3.91 ± 5.29), and their values were not significantly different (F 1, 14 = 0.3126, P = 0.5848). November was a secondary MGI month for females (2.30 ± 0.26). The highest female GSIs showed no significant difference between July and May (P = 0.12), July and June (P = 0.54), July and August (P = 0.24), July and September (P = 0.17), July and October (P = 0.06), and July and November (P = 0.39), suggesting seven months as the reproductive season for females. In males, the highest GSIs showed no significant difference across all months sampled (P > 0.05) except January to March (P = 0.03), suggesting an extended reproductive activity throughout the year (Figure 3B); i.e. a reproductive season of 9 months in males and 7 months in females.
In both locations, the GSI values fluctuate similarly through time. However, MGI was ~2.5 (in males) and 7 (in females) times higher at PM, Chiapas than at CC, Michoacán (Figure 3). At the population level, the mean GSI at PM was significantly higher than the mean GSI at CC (F 1,211 = 99.25, P = 0.0000) (Figure 3, right panel), suggesting a highly significant effect of the geographic location on the intensity and duration of the MGI.
Physiological investment: frequencies and period of GDSs
The main histological features of each gonad development stage in males and females of H. inornata are illustrated in Figure 4, and the respective monthly frequencies of the non-reproductive phase (amalgamated GDS-4, GDS-5 and GDS-1) and the reproductive phase (amalgamated GDS-2 and GDS-3) from CC and PM populations are presented in Figure 5.
Histological sections of the testes (left panel) and the ovaries (right panel) of Holothuria inornata. (A) gametogenesis, (B) mature, (C) spawning, (D) post-spawning.

Frequency (%) and period (months) of the non-reproductive phase (amalgamated GDS-4, GDS-5 and GDS-1) and the reproductive phase (amalgamated GDS-2 and GDS-3) of H. inornata females and males. (A) Caleta de Campos, Michoacán (18°04′N), and (B) Puerto Madero, Chiapas (14°42′N). The dashed line indicates 50% frequency. Sampling periods in which more than 50% of all specimens occurred in the reproductive phase indicated the number of reproductive seasons.

Figure 5 shows how the transition from the non-reproductive to reproductive phase varies between the sexes and the geographic locations, apparently resulting in more than one reproductive season (RS). At CC, Michoacán, reproductive phase (GDS-2 plus GDS-3) frequency larger than 50% suggested two reproductive seasons of different duration for the female population (first, June–September; second, during November). In contrast, reproductive phase-frequency in the male population was not larger than 33%, although the data also suggested two reproductive seasons (first, during June; second, during three months: August–October) (Figure 5A). At PM, Chiapas, reproductive phase-frequency larger than 50% suggested three reproductive seasons for the female population (first, during 3 months: May–July; second, during September; third, during November–December). For the male population, reproductive phase-frequency suggested two reproductive seasons (first, during 3 months: July–September; second, during 3 months: November–January) (Figure 5B).
That the number of reproductive seasons in the female population appears to be larger than those in the male population, is likely to be an artefact of sample size. In the CC population, the number of reproductive seasons appears to be limited to just one MGI, although its duration (months) increases, i.e. a single, annual, extended reproduction. In the PM population, the number of reproductive seasons increases (more than one MGI) although the duration of each MGI decreases, i.e. more than one discontinuous reproduction per year.
Sexual asymmetry in the MGIs: mass investment
Seasonal shifts in RGES are illustrated in Figure 6. The annual mean RGES showed slight female-biased asymmetry (−0.12) in the PM population, while it showed a male-biased asymmetry (0.26) in the CC population. The most significant biases for both locations were: female-biased asymmetry (−3.00) during November in the CC, Michoacán and male-biased asymmetry (2.12) during December in the PM, Chiapas (Figure 6).
Variation in sexual asymmetry in GSI (relative gonad expenditure RGES, with standard deviation, see text) of Holothuria inornata at each location. (A) Caleta de Campos, Michoacán (18°04′N). (B) Puerto Madero, Chiapas (14°42′N). The RGES variation includes all GDS (I to V) for each sex. F > M refers to female-biased asymmetry in RGES (circular symbols) and M > F to male bias (triangular symbols). F = M refers to GSI equality (i.e. RGES ≈ 0) (square symbols).

Sexual asymmetry in the GDSs: physiological investment
Figure 7 illustrates how the data in Figure 6 translate into our measure of the relative gonadal physiological investment of the sexes in reproduction (RGPIS) through time. Small sample sizes preclude firm conclusions, but sexual equality (mean RGPIS = 0) appeared more recurrent in the PM population (three times) than in the CC population (once). These equalities were for the spawning stage (July and November) and post-spawning stage (October) in PM, Chiapas, but just for the gametogenesis stage (May) in the CC, Michoacán.
Variation in sexual asymmetry in GDSs (relative gonadal physiological investment of the sexes in reproduction, RGPIS) for each GDS of Holothuria inornata at each location, as explained in the text. (A) Caleta de Campos, Michoacán (18°04′N). (B) Puerto Madero, Chiapas (14°42′N). The arrows in September for both locations show RGPIS almost equidistant from zero between GDS. F > M refers to female-biased asymmetry in RGPIS and M > F to male bias. F = M refers to the equality of the physiological investment of the GDSs (i.e. RGPIS ≈ 0).

Particularly for both geographic locations during September, some RGPIS values were almost equidistant from the zero value between GDS. Comparatively, the RGPIS showed several values with male-biased asymmetry (>0) in the PM population, while it showed few values with male-biased asymmetry in the CC population (Figure 7).
Adult sex ratio
Both H. inornata populations showed an annual sex ratio (m/f) that did not differ significantly from 1:1, although the CC population was closer to being male-biased (CC 1.7:1, χc2 = 3.44, P = 0.06; PM 1.3:1, χc2 = 1.38, P = 0.24). The monthly sex ratio also did not differ significantly from 1:1, except during June in the CC population and January in the PM population (Table 2, Figure 8). Individuals with recovering gonads that are sexually undifferentiated were much more evident in the CC population (47% of the annual total population) than the PM population (11% of the annual total population) (Table 3, Figure 8).
Adult sex ratio of H. inornata for (A) Caleta de Campos, Michoacán (18°04′N), and (B) Puerto Madero, Chiapas (14°42′N).

Male to female adult sex ratio (excluding sexually undifferentiated individuals) of Holothuria inornata by location, month, and the total annual sample

Notes: χ2 goodness of fit test using the Yates' correction for continuity (χ2c) monthly and annual. Nm, number of males; Nf, number of females; Nt, total number of males and females; %m, percentage of males; %f, percentage of females; m:f, male to female ratio.
* Deviation from the expected ratio (1:1) indicates a significant difference at P < 0.05. Other P-values in bold and without a superscript ‘*’ show closeness to P = 0.05. Sea cucumbers with recovering gonads recorded as sexually undifferentiated were not included in this analysis. The grey shading represents the period of the reproductive season respectively for each place.
Numbers of males, females and sexually undifferentiated sea cucumbers of Holothuria inornata by location, month, and the total sample

Notes: Nm, number of males; Nf, number of females; Ns-U, number of sea cucumbers with recovering gonads or sexually undifferentiated; %m, percentage of males; %f, percentage of females; %s-U, percentage of sea cucumbers with recovering gonads or sexually undifferentiated. The grey shadow represents the period of the reproductive season respectively for each place.
Discussion
Our examination of maximum gonad investment (MGI) in mass (intensity and period) and developmental frequency through time demonstrates that the broadcast spawning sea cucumber Holothuria inornata exhibits different MGI responses for males and females. Sex exerts a strong influence on both the gametogenic cycle and its duration in the reproductive season. From life-history theory, we expect these traits to reflect optimum reproductive strategies (Stearns, Reference Stearns1992).
The influence of temperature on reproductive season flexibility
Changes in the intensity and duration of the reproductive season (MGI) of H. inornata were not underpinned by changes in local temperature, but rather by regional temperature, suggesting that though reproduction is temperature-sensitive, H. inornata has the ability to cope with local fluctuations. Temperature variation (Figure 2) at Puerto Madero appeared to be within the regular period of ONI, and had little effect on GSI. During the sampling period at Caleta de Campos, it was noticeably colder and within the cold period of ONI, and appeared to reduce GSI intensity. GSI intensity decreased as the temperature dropped. However, this was compensated by extension of the reproductive season, demonstrating possible phenotypic reproductive flexibility to temperature changes, particularly to regional thermal events, as has been suggested for a wide range of ectotherm invertebrates of the intertidal and subtidal zone (Manríquez et al., Reference Manríquez, Guiñez, Olivares, Clarke and Castilla2018; Avila-Poveda et al., Reference Avila-Poveda, Abadia-Chanona, Alvarez-Garcia and Arellano-Martinez2021; Flores et al., Reference Flores, Wiff, Ahumada, Queirolo and Apablaza2021).
Mass investment: MGI and its sexual asymmetry
The relationship between gonad weight and body weight (i.e. gonadosomatic index, GSI) is the most generalized and widely used method for evaluating reproductive season, particularly in broadcast spawning invertebrates (Giese, Reference Giese1959). Typically, the highest GSI value corresponds to a peak of the ripe stage (i.e. maximum gonad investment, MGI). We have here used a recently proposed GSI sexual asymmetry index that indicates the relative cost of gonad investment of the sexes, and which can be related to selection due to the level of sperm competition and/or sperm limitation (Abadia-Chanona et al., Reference Abadia-Chanona, Avila-Poveda, Arellano-Martinez, Ceballos-Vazquez, Benítez-Villalobos, Parker, Rodriguez-Dominguez and Garcia-Ibañez2018; Parker et al., Reference Parker, Ramm, Lehtonen and Henshaw2018).
Intensity of reproductive season
The PM population showed ~4 times larger MGI than the CC population, suggesting that (for both sexes; Figure 3) MGI is influenced by site-specific factors. The fact that both higher MGI values and an extended reproductive season are attained at PM, Chiapas, suggests that resources are more abundant and available for more extended periods of the year than at CC, Michoacán. Sea cucumbers are benthic, surface-deposit feeders responsible for significant removal of the particulate organic carbon loading on the seabed. Field and laboratory experiments have shown a correlation between the density of holothurian species, geographic location and substrate composition for achieving optimum reproduction (Chao et al., Reference Chao, Chen and Alexander1995; Gianasi et al., Reference Gianasi, Hamel and Mercier2019). For male broadcast spawners, GSI intensity is predicted to depend on the number of males simultaneously releasing gametes during a spawning event (the mating density: Levitan, Reference Levitan2005) and on the operational sex ratio ‘OSR’ (the average ratio of fertilizable females to sexually active males at any given time: Emlen & Oring, Reference Emlen and Oring1977; Szekely et al., Reference Szekely, Weissing and Komdeur2014). These two factors together determine sperm densities and the levels of sperm competition and/or sperm limitation (Levitan, Reference Levitan, Westneat and Fox2010; Parker et al., Reference Parker, Ramm, Lehtonen and Henshaw2018; Evans & Lymbery, Reference Evans and Lymbery2020).
Duration of the reproductive season
In both H. inornata populations, the main period of the reproductive season occurred once a year, with the same month of ripe peak (MGI) being June. Nevertheless, there were differences in duration of the reproductive season (i.e. those months without a significant GSI difference to the MGI), which lasted 3 months in the CC population (June–August), while in the PM population, there was a second MGI and an additional four more high GSI months, generating a 7-month duration (May–November) (Figure 3). Although these MGI analyses were based on only two populations distant from each other by four latitudinal degrees along a relatively narrow latitudinal gradient, we speculate that in H. inornata, the duration of the reproductive season tends to increase as latitude decreases, changing from annual to almost continuous semi-annual. These differences could, however, arise from site-specific resource-level factors unrelated to latitude. Similar patterns have been observed in other broadcast spawning marine invertebrates (Benítez-Villalobos et al., Reference Benítez-Villalobos, Avila-Poveda, Diaz-Martinez and Bravo-Ruiz2015; Oyarzun et al., Reference Oyarzun, Toro, Garces-Vargas, Alvarado, Guiñez, Jaramillo, Briones and Campos2018; Avila-Poveda et al., Reference Avila-Poveda, Abadia-Chanona, Alvarez-Garcia and Arellano-Martinez2021). In ectothermic organisms, this pattern typically reflects both latitude and temperature (which changes with latitude; Giese, Reference Giese1959).
GSI sex asymmetry
Observed data on sedentary broadcast spawners show that approximate GSI symmetry is most typical. However, when significant differences are observed, female GSI usually exceeds male GSI (as in our holothurian PM population). Male GSI higher than female GSI is less common overall, except possibly in tropical molluscs (see Abadia-Chanona et al., Reference Abadia-Chanona, Avila-Poveda, Arellano-Martinez, Ceballos-Vazquez, Benítez-Villalobos, Parker, Rodriguez-Dominguez and Garcia-Ibañez2018; table 1 and figure 1 in Parker et al., Reference Parker, Ramm, Lehtonen and Henshaw2018) and echinoderms (Engstrom, Reference Engstrom1980; Guzman et al., Reference Guzman, Guevara and Hernandez2003), including our holothurian CC population. Sea cucumbers of the genus Holothuria generally show female-biased GSI asymmetry, though male-biased GSI does occur (see table 1 and figure 1 in Parker et al., Reference Parker, Ramm, Lehtonen and Henshaw2018).
Sperm competition and sperm limitation are two selective forces that can generate sexual asymmetry in MGI and affect gonad expenditure (Levitan, Reference Levitan, Westneat and Fox2010; Parker et al., Reference Parker, Ramm, Lehtonen and Henshaw2018; Evans & Lymbery, Reference Evans and Lymbery2020; Parker, Reference Parker2020). In the absence of sperm limitation, increasing sperm competition favours increasing male GSI towards that shown by females, so that under very high sperm competition, GSIs become equal (RGES = 0). Similarly, in the absence of sperm competition, increasing sperm limitation has a similar effect: male GSI should be lower than female GSI when virtually all available eggs are fertilized by just one male, increasing towards female GSI only when very few available eggs are fertilized, and all by the same male. Both selective forces appear to generate either female-biased or equal GSI (RGES ≤ 0). Because high resource levels lead to higher GSI values (and possibly higher population densities), our holothurian PM population might be expected to experience significant sperm competition, leading to slightly female-biased RGES values, as we observed. However, the male-biased GSI asymmetry (RGES > 0) found in the CC population is more difficult to explain.
Two possibilities have been proposed for male-biased GSI (Parker et al., Reference Parker, Ramm, Lehtonen and Henshaw2018). Most plausibly, higher male GSI can be generated when the cost (per unit mass) of producing ova is higher than the cost of producing sperm. Different life-history strategies can arise because males and females have different optimal gonadal investments (Charnov, Reference Charnov1982; Levitan, Reference Levitan2005; Parker et al., Reference Parker, Ramm, Lehtonen and Henshaw2018), and a second possibility is that under certain conditions, sperm competition can cause males to mature earlier than females (Parker et al., Reference Parker, Ramm, Lehtonen and Henshaw2018). Higher mean GSI for males can arise if immature females are included in samples (Parker et al., Reference Parker, Ramm, Lehtonen and Henshaw2018). Because all individuals in our samples are adult that had reached sexual maturity and are physiologically capable of producing offspring, this second possibility can be ruled out. However, there are reasons to suppose that egg production can be costlier than sperm production (Crean & Marshall, Reference Crean and Marshall2008; Parker et al., Reference Parker, Ramm, Lehtonen and Henshaw2018); if this effect is particularly acute under scarcer resources, as likely experienced in the CC population, this may explain the divergence towards male-biased GSI in that study site.
Both holothurian populations showed (non-significant) overall trends towards male-biased sex ratios, and more so in the CC population (1.7:1) than in the PM population (1.3:1). If the survival costs of producing each unit mass of ova exceeds that for sperm, and if restricted resources at the CC population exacerbate this effect, this could generate both the significantly male-biased GSI asymmetry and the trend for male-biased sex ratio. The model of Parker et al. (Reference Parker, Ramm, Lehtonen and Henshaw2018) derives evolutionarily stable strategies (ESSs; Maynard Smith, Reference Maynard Smith1982) for gonad expenditures,
$x_f^\ast $ and
$x_m^\ast $ for females and males respectively (equations 18 and 19, and figure 4 in Parker et al., Reference Parker, Ramm, Lehtonen and Henshaw2018), when gonad investment (GSI values expressed as proportions) trades off against survival by increasing adult mortality rates at an accelerating rate. Male-biased GSI (
$x_m^\ast > x_f^\ast $) can occur if ova are physiologically more expensive than sperm, and sperm competition is high. Under high sperm competition, equations 14 and 15 in the model of Parker et al. become equal, and the ESS mortality rate for sex i,
$p_i^\ast ( x_i^\ast ) $, using the model for figure 3 of Parker et al. (Reference Parker, Ramm, Lehtonen and Henshaw2018) becomes
where p 0 is the random mortality rate and c i is a cost constant for sex i that increases the rate of acceleration of mortality due to gonad expenditure (this solution for
$p_i^\ast $ is obtained by substituting equation 18 into equation 16 in Parker et al., Reference Parker, Ramm, Lehtonen and Henshaw2018). Thus if ova are more costly than sperm, i.e. (c f > c m), this implies a higher ESS mortality rate for females, leading to a male-biased sex ratio (figure 4 in Parker et al., Reference Parker, Ramm, Lehtonen and Henshaw2018), conforming with the observed trend for the male-biased sex ratio and the male-biased GSI asymmetry in the CC, Michoacán population.
Physiological investment: GDSs and their sexual asymmetry
Another criterion for assessing reproductive season is observing ripe gametes by histological examination (i.e. through gonad developmental stages, GDS). Some confusion occurs over the terms ‘reproductive cycle’ and ‘reproductive season’: ‘reproductive cycle’ often relates to events over an entire year, which may contain more than one ‘reproductive season’, each involving the change from unripe to ripe gametes, quantified by population (or sex) frequency of GDSs (see Giese, Reference Giese1959 and figure 1 in Avila-Poveda et al., Reference Avila-Poveda, Abadia-Chanona, Alvarez-Garcia and Arellano-Martinez2021). As expected, the reproductive phase of the gonads (i.e. ripe and spawning stages greater than 50%) matches MGI in both PM and CC populations.
In many broadcast spawning invertebrates, ripe and spawning physiology is synchronous between sexes, an adaptation related to increasing fertilization success (McHugh & Rouse, Reference McCune, Grace and Urban1998). Since broadcast sperm survive between hours to at most a few days in seawater, we expect selection to match sperm release rates closely with ovum release rates (see Discussion in Abadia-Chanona et al., Reference Abadia-Chanona, Avila-Poveda, Arellano-Martinez, Ceballos-Vazquez, Benítez-Villalobos, Parker, Rodriguez-Dominguez and Garcia-Ibañez2018). Inevitably, when sample sizes are reduced due to subdivision into male and female GDSs, any firm conclusions about the sexual synchrony of GDSs are impossible. However, some tendency for synchrony related to sex ratio seemed apparent in H. inornata. For example, in the June MGI period, males outnumbered females (m/f ratios: 8.0:1 and 2.3:1, respectively for CC and PM populations), while female ripe stage frequency outnumbered male ripe stage frequency (f to m: 100%:12% and 50%:28%, respectively for CC and PM populations). Whether the high male-male competition in June favoured some males delaying spawning until competition reduces remains to be determined.
Specific holothurian populations have exhibited reproductive peaks of different intensities with variable ripe stage frequency. For example, in an H. grisea population, during the first peak, the ovaries were heavier and riper than testes (female MGI was 44–48% with a frequency of 87–100% in ripe stage and 0% in spawning stage, while male MGI was 37–39% with a frequency of 71–87% in ripe stage and 0–13% in spawning stage). However, at the next peak, the trend disappeared (male MGI 45–50% and female MGI 40–45%, and both sexes with a frequency of 100% in ripe and 0% in spawning stages) (Leite-Castro et al., Reference Leite-Castro, Souza-Junior, Salmito-Vanderley, Ferreira-Nunes, Hamel and Mercier2016). A theoretical model would derive the ESS distribution of spawning period by individuals in a population, such that the reproductive success of all individuals of sex is equal (‘ideal free’ fertilizations through time), and would take into account sex ratio and the constraints on the frequency of gamete release by each sex (see Discussion in Abadia-Chanona et al., Reference Abadia-Chanona, Avila-Poveda, Arellano-Martinez, Ceballos-Vazquez, Benítez-Villalobos, Parker, Rodriguez-Dominguez and Garcia-Ibañez2018).
Overall, our results suggest that, across narrow latitudinal gradients within the tropical zone, our two natural populations of Holothuria (Halodeima) inornata exhibit adaptive variation in the intensity, period and number of the maximum gonad investment (MGI) peaks per annual cycle for each sex associated with ripe stage frequency and adult sex ratio. However, we do not rule out the possibility that extrinsic aspects, in addition to temperature, such as variation in food supply or specific site differences (e.g. substrate composition) can also exert an influence; this question remains to be investigated.
Acknowledgements
Both H. Cancino-Guzman and E. Ramos-Ramirez received scholarships and the results presented here are part of their respective bachelor's thesis projects. Histological procedures and analyses were carried out at the histology laboratory of the Universidad del Mar (UMAR). OH Avila-Poveda is commissioned as CONACYT Research Fellow/UAS-FACIMAR (project No. 2137) to the Academic group ‘Manejo de Recursos Pesqueros UAS-CA-132; UAS-FACIMAR’.
Financial support
This research was funded by the Programa de Mejoramiento del Profesorado (ProMeP) of the Secretaria de Educacion Publica (SEP) [grant number 2IR0707] awarded to F. Benítez-Villalobos, and partially by the INE-UMSNH project [agreement number INE/PS-031/2011] shared to F. Benítez-Villalobos by Francisco Alonso Solis-Marin from ICMyL, UNAM.
Conflict of interest
The authors declare that they have no known competing financial interests or personal relationships that could have influenced the work reported in this paper.
Ethical standards
The specific research collection permit (No. DGOPA.02934.250308.0646) was acquired with the ProMeP project. All sampling procedures and manipulations of individuals was according to welfare and ethically responsible research with sea cucumbers (Rubilar & Crespi-Abril, Reference Rubilar and Crespi-Abril2017; Crespi-Abril & Rubilar, Reference Crespi-Abril and Rubilar2018), and with the ethical recommendations for the humanitarian killing of animals as established under Mexican law (NOM-033-SAG/ZOO-2014 derogating to NOM-033-ZOO-1995).






