We used ecological niche models (ENM) to identify regions in Cuba most likely to be colonized by L. fragilis. The ENM was first built using records of the native distribution of the species and bioclimatic variables, before being projected onto the Cuban archipelago to predict the potential distribution. The ENM was generated using Maximum entropy model (MaxEnt v.3.3.3, Phillips et al., 2006), that predicts habitat suitability as a function of environmental variables and species occurrence data. It relies on presence-background data (i.e., randomly selected absences from areas that have been accessible to the species; Phillips and Dudík, 2008), and is particularly suitable to explore the likely distributions of invasive species (Jiménez-Valverde et al., 2011; Cordier et al., 2020).

We obtained 551 georeferenced records describing the native distribution of L. fragilis in mainland America from Medina et al. (2020). To reduce clusters of localities that might create bias in environmental space, we used the R package spThin (Aiello-Lammens et al., 2015) to spatially filter occurrence records with a minimum distance of at least 5km. After filtering, 462 occurrence records remained for modeling distributions (Fig. 1a). To build ENM, we choose nine climatic variables: annual mean temperature (bio 1), mean diurnal range (bio 2), temperature seasonality (bio 4), maximum temperature of warmest month (bio 5), minimum temperature of coldest month (bio 6), annual precipitation (bio 12), precipitation seasonality (bio 15), precipitation of wettest quarter (bio 16), and precipitation of driest quarter (bio 17). Such bioclimatic variables were selected, as they are critical for amphibian survival and reproduction (Wells, 2007). Indeed, due to complex life cycle with aquatic and terrestrial life stages, shell-less eggs, permeable and exposed skin, and ectothermic condition, amphibians are very sensitive to changes in environmental quality (Hillman et al., 2008). All variables were extracted from Worldclim database 2.0, based upon weather conditions recorded between 1970 and 2000, with a grid cell resolution of 30 arc seconds (Fick and Hijmans, 2017; http://www.worldclim.org). Percent contribution of each variable was then calculated as detailed in Phillips et al. (2006). In order to reduce model complexity and multicollinearity, only six variables (bio 2, 4, 6, 12, 15 and 17) were retained for further analysis, based on their high contribution and low degree of correlation between themselves (|Pearson correlation|≤0.75). We restricted background sampling to a 100km radial buffer zone around each occurrence record. We considered this extent appropriate for background selection, as it does not include large regions that the species does not inhabit due to dispersal limitations and/or biotic interactions (Barve et al., 2011). We extracted 10 000 background random points within these buffers using ArcGIS v. 10.2 (ESRI, Redlands, California, USA).

Fig. 1.

Ecological niche model for Leptodactylus fragilis based on 462 occurrence records (orange points) from its native range, green shading indicates the suitable climatic area using “10th percentile training presence” threshold (a). Potential distribution of L. fragilis onto the Cuban archipelago (green area) based on ENM and suitable land cover, the orange spot indicates the known localities of L. fragilis; the records of presence of Peltophryne empusa (yellow triangles) from Alonso Bosch (2011) and Rivalta et al. (2014) (b). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

(0.86MB).

We tuned the value of the regularization multiplier and determined the optimal selection of feature classes using the R package ENMeval (Muscarella et al., 2014). Regularization values ranged from 0.5 to 3, in increments of 0.5. Five settings of feature classes were evaluated: linear (L), linear and quadratic (LQ), hinge (H), LQH, and LQH plus product (P). The complexity of each model was evaluated using Akaike Information Criterion corrected for small sample sizes (AICc) (Warren and Seifert, 2011). We also inspected the predictive power (omission rate E=10%) and the overfitting by calculating AUC difference, which is the difference between AUC calculated on training localities and AUC calculated on evaluation localities (Warren and Seifert, 2011). To evaluate the selected model, we used the area under the curve (AUC) of the partial ROC (receiver operating characteristic; pROC), which allows differential weighting of omission and commission errors, and more accurately assesses the quality of the model (Peterson et al., 2008). Using program NicheA (Qiao et al., 2016), pROC was used as a measure of model performance, with pROC ranging from 0 to 2, where 2 indicates a perfect prediction and values <1 indicate predictions no better than random (Peterson et al., 2008).

We used clamping which treats variables outside the training range as if they were at the limit of the training range (Elith et al., 2010). The similarity between calibration areas and Cuban archipelago (novel area) was assessed using Multivariate Environmental Similarity Surface (MESS) implemented within Maxent. This procedure allows the identification of areas with high extrapolation for which obtained results should be interpreted with caution (Elith et al., 2010). The final model was run in MaxEnt using 50 bootstrapping replicates and the combination of regularization multiplier and feature classes with lowest AICc using all thinning occurrence records; this model was projected onto the Cuban archipelago. To obtain presence – absence map from the logistic values of environmental suitability (continuous probability from 0 to 1), we used the 10th percentile threshold, which is poorly sensitive to extreme environmental values (Radosavljevic and Anderson, 2014), reduces commission errors, and produces more conservative distributional maps (Liu et al., 2005). Finally, we applied a habitat filtering approach, removing areas with unsuitable land cover types from the climatic niche model projected to Cuban archipelago; this procedure acts to reduce commission errors inherent to the ENM based exclusively on climatic data (Rondinini et al., 2006). After extracting cover types from GlobCover project v. 2.3 (E.S.A., 2010), we reclassified the map with all categories suitable for L. fragilis (e.g. croplands, grasslands and shrublands) being grouped in a single class.

Recordings were obtained from six individuals at Sandino (22.075837N, −84.219685W, 9.5masl), Pinar del Río province, using a digital audio field recorder (Sony PCM-M10) with a directional microphone (Senheiser MEE/K6). Advertisement calls were obtained at a 44,100Hz sampling frequency, with a 16 bit resolution. The analysis of calls was performed with the Raven 1.3 software (Bioacoustics Research Program, Cornell Laboratory of Ornithology, 2012), using Hanning window, FFT size 2048, and overlap of 95%. From the oscillograms, call duration, interval between calls (error 0.001s), and number of pulses in each call were recorded, while dominant frequency to the nearest 0.02kHz at the peak of maximum amplitude were obtained from the power spectrum.

In order to assess the possible overlap in mean values of temporal (call duration, pulses per call and pulse rate) and spectral (dominant frequency) features of advertisement calls between L. fragilis, and native species, we examined the acoustic features of the advertisement calls of more than 60 Cuban species of frogs and toads (Díaz and Cádiz, 2006, 2008; Alonso Bosch, 2011). Special attention was devoted to species potentially occurring in sympatry with the invasive species in open areas, particularly those who vocalize and breed in temporary and permanent water bodies (Alonso Bosch et al., 2007; Díaz and Cádiz, 2008).

Pre-recorded signals used as playbacks were broadcasted using a loudspeaker placed in the natural environment in order to measure propagation efficiency, similar to previous experiments with other anuran amphibians (Kuczynski et al., 2010; Llusia et al., 2013; Bleach et al., 2015; Penna and Moreno-Gómez, 2015; Schwartz et al., 2016). Propagation experiments were carried out during the wet season, between 18:00 and 20:00h. The audio file for the experiment was designed using the software Adobe Audition CC 6.0 (Compilation 732, 64 bits). The 3-min long playback (44,100Hz and 16 bit) included sequences of pure tones of 1kHz, and a series of three different stimuli: consecutive calls of one male, vocal interaction between two males, and a chorus of L. fragilis (Appendix S1 file). Sound pressure level (SPL) was standardized for all the stimuli at 0.5m from the source.

The playback stimuli were broadcast with a self-powered loudspeaker (Pignose No.7-100) connected to a laptop computer and placed at positions typically occupied by calling males on the ground. The broadcast signals were recorded with a directional microphone (Sennheiser MEE/K6) pointing at the source of the playback, placed successively at distances of 0,5; 1; 2; 4; 8m from the loudspeaker, and a digital audio field recorder (Sony PCM-M10). At the same positions, a sound level meter (RadioShack, Error=2dB) was placed to obtain dB values of the signal at each distance. Relative humidity and air temperature were measured at each distance using a digital thermo-hygrometer (Winflex, Error=±1°C; ±5%RH).

We selected six localities to perform sound propagation experiments, considering the potential distribution of L. fragilis, the description of its advertisement call, and the potential overlapping with the structure of advertisement calls of native syntopic species. In order to compare the excess attenuation and degradation of advertisement calls of L. fragilis, sound propagation experiments were conducted in forest (PMf) and open (PMo) habitats of the ancient Botanical Garden of Havana (23.1001N, −82.3985W, 32.3masl.). The first habitat is a perturbed mesophyllous semideciduous forest with herbaceous stratum bears grasses, ferns and abundant lianas, while the second is located in an open area on the border of the forest. Two other localities in the main island were included in the study: Sandino, Pinar del Río, province (22.075837N, −84.219685W, 9.5masl.) and near to Campo Florido town in La Habana province (23.128944N, −82.157216W, 31.1masl.). Sandino (SA), one of the two previously known localities of L. fragilis in Cuba, corresponds to the shore of a temporary small pond surrounded by grass pastures and bushes, while Campo Florido (CF) is a very irregular terrain near a temporary stream surrounded by grass pastures and several small bushes. Additionally, two open habitat localities were selected from Isla de la Juventud: surroundings of El Colony (21.633332N, −82.983118W, 6.7masl.) and outside of Nueva Gerona city (21.870908N, −82.79291W, 22.5masl.). The study area in El Colony (CO) is a temporal flooded depression with abundant thin grass forming like a mattress, while Nueva Gerona (NG) is a grass pasture area, with a sandy soil substrate.

In order to assess the possible degradation of the signals, sequences of advertisement calls of one male and the interaction between two males for all distances (0,5; 1; 2; 4; 8m) were analyzed. Call duration (CD) was measured at zero amplitude level on the oscillogram (error=0.1ms), while the number of pulses per call (NP) was counted. The dominant frequency (DF) was measured to the nearest 0.02kHz, at the peak of maximum amplitude in the power spectrum. All of these features were measured on the advertisement calls of one male for 10 alternate calls, starting on the second call of the sequence. For the sequence of the interaction between two males, the same procedure was followed, using the call of the closest male. Considering that signal attenuation is a major effect of sound transmission that may interrupt effective communication if signals are attenuated below auditory thresholds of recipients, attenuation effects were calculated with the equation: spherical transmission loss (dB)=20log[dB Sound Pressure Level (SPL) at far distance (m)/dB SPL at 0.5 (m)] (Llusia et al., 2013). Excess attenuation was calculated by subtracting values of spherical spreading from the actual transmission loss (differences between SPL at 0.5m from the loudspeaker and the corresponding higher distances) measured with the sound level meter for 1, 2, 4 and 8m relative to measurement at 0.5m.

We obtained values of intensities for all studied localities and three broadcasted stimuli only when testing at a distance of 0.5–1m. At 0.5–2m values were obtained only for the three open localities for three stimuli, whereas the sensitivity of the equipment was too low to perform experiments at distances of 0.5–4m and 0.5–8m. Increasing the sensitivity of the sound level meter to its maximum value at the farthest distances from the loudspeaker had no effect on detectability.

We calculated mean and standard deviation for each parameter included in our analyses. For the analysis of variation of three acoustic variables along five distances among localities, Inverse Gaussian response distribution was selected by the lowest values of Akaike Information Criterion (AIC) after SEVERITY procedure, and reciprocal link function, while the analysis of excess attenuation along five distances, Gamma distribution and log link function were selected. We used a two-way Analysis of Variance, and Type III sums of squares, in a generalized linear mixed model with PROC GLIMMIX (Littell et al., 2006). We considered distance, locality and their interaction as fixed-effect factors. If there was an overall significant effect, we conducted Tukey–Kramer adjusted multiple comparison tests between all levels. For all analyses, the best-fit model was selected using AICC (Burnham and Anderson, 1998). All statistical analyses were conducted using SAS software version 9.3 (SAS Institute Inc., Cary, NC, USA), with significance level at α=0.05.

The ecological niche model of L. fragilis (Fig. 1a) with lowest values of AICc (delta AICc=0) resulted from regularization multiplier of 2.5 and LQHP feature class. The mean AUCdiff (0.025±0.0006 SD) and mean 10% training omission rate (0.12±0.0003 SD) were indicative of good model performance and low overfitting. pROC value was equal to 1.49 and significantly different from 1 (P<0.001), indicating significant predictive ability of model. The analyses of variable contributions revealed that minimum temperature of coldest month, temperature seasonality and mean diurnal range had the highest explanatory power (Appendix S2 File). Geographic projections of ENM showed climatic suitability for L. fragilis throughout the Cuban archipelago (Fig. 1b). Clamping and multivariate environmental similarity surfaces (MESS) showed that no areas had variables outside of the training data range. The ENM filtered with suitable land cover indicated that 66,930km2 are highly susceptible to invasion by L. fragilis, of which 5680km2 are included within the Cuban system of protected areas (Fig. 1b).

Males of L. fragilis emit an amplitude-modulated advertisement call with harmonic structure. The advertisement call consists of a single pulsed note with 19.4±1.5 pulses (16–22), and two or three indistinct pulses added at the end. Its dominant frequency ranged between 1.38–2.07kHz (Mean±SD, 1.76±0.14kHz). Acoustic features of advertisement calls of L. fragilis, and those of eight native species whose vocal activity and mating take place more frequently in lentic water bodies in open areas, are provided in Table 1. The advertisement calls of Peltophryne empusa appear to be most similar in temporal and spectral features to L. fragilis relatively to the other species of toads we analyzed (Table 1). According to the previously obtained ecological niche model of L. fragilis in Cuba, this introduced species could be widely sympatric with the endemic P. empusa (Fig. 1b). In the hypothetical case where both species coexist in the same locality, share the microhabitat for vocalizing, and their calling activities coincide in time and space, the physical features of their advertisement calls could largely overlap (Table 1).

Although we did not observe any severe degradation in the signals of L. fragilis, significant effects of distance and locality (P<0.001) were observed (Table 2, Appendix S3 File). Call duration was the most stable acoustic feature, when considering either the signal of one single male or the interaction between two males (Table 2, Appendices S4 and S5). Locality had no significant influence on the three acoustic features of one single male call (Table 2, Appendix S4). In terms of call duration of one single male, Nueva Gerona (NG) was the most different locality (with less variation in relation to the stimuli), whereas this variable for the interaction between males, was very similar among localities (Appendix S5). The number of pulses in both cases (one single male and interaction) followed the same pattern: all localities were significantly different from Nueva Gerona and Sandino (SA), whereas the latter two did not differ between themselves. Finally, the comparative analysis of the dominant frequency among localities did not reveal a clear pattern. However, Campo Florido (CF) and both areas in the ancient Botanical Garden of Havana [PMo (Open) and PMf (Forest)] appeared to be the most distinct localities, with the highest variation for both the advertisement calls of one single male and the interaction between males (Appendices S4 and S5).

The comparative analysis of L. fragilis signal degradation along five distances to the loudspeaker showed significant differences for the three studied acoustic variables, for both the advertisement calls of one single male and the interaction between males (Table 2). Noticeable variations were detected in terms of call duration and number of pulses when examining both stimuli (one male call and interaction between males) at the farthest distances (Appendices S4 and S5).

There was a significant interaction between distance to the loudspeaker and locality (Loc×Dist) on signal degradation (Table 2, Appendices S3–S5). The most relevant pattern was observed for the dominant frequency, with more structurally complex localities (higher number of obstacles and denser vegetation) being significantly different from more open ones at all distances (0.5, 1, 2, 4, and 8m to the loudspeaker). A similar pattern was detected both in the analysis of the signal degradation of the one single male call and in the interaction between males (Appendix S3). Differences in dominant frequency among open localities appeared only at the farthest distances. No particular pattern was observed for the remaining acoustic variables analyzed (Appendices S3–S5). We also checked for the effect of the interaction between distance and locality, when establishing comparisons between localities, considering different distances (Appendix S3).

Measurements of playback experiments between 0.5 and 1m from the loudspeaker in the forest of the ancient Botanical Garden of Havana (PMf) showed increased attenuation (P<0.001) of advertisement calls of one single male, interaction between two males, and chorus, compared to open areas (Fig. 2, Appendix S6). The excess of attenuation was very similar between Colony (CO), Campo Florido (CF), the open area in the ancient Botanical Garden of Havana (PMo) and Nueva Gerona (NG) at the first evaluated distance (Fig. 2, Appendix S6). However, the excess of attenuation in these open areas was also significantly different in relationship to the obtained values in Sandino (SA) at the same distance from the loudspeaker. This locality showed the lowest excess of attenuation for one single male call when performing tests at the 0.5–1m distance, while the excess of attenuation for both the interaction and chorus was minimal in Nueva Gerona (Fig. 2, Appendix S6). The excess of attenuation at the 0.5–2m distance differed among the open localities where the three stimuli were perceived by the instrument (CF, NG, SA), with Nueva Gerona showing the lowest values of attenuation for the three types of stimuli (Appendix S6).

Fig. 2.

Excess attenuation of sequences of the advertisement calls of one single male, an interaction between two males, and a chorus of Leptodactylus fragilis between 0.5 and 1.0m from the loudspeaker. Filled rhombus: advertisement call of one single male; Open rhombus: interaction between two males; Asterisks: chorus. CF: Campo Florido, CO: Colony, NG: outside Nueva Gerona, SA: Sandino, PMf: Forest in Parque Metropolitano of Havana city, Pmo: Open area in Parque Metropolitano of Havana city. For each locality mean and 95% confidence intervals are shown. Gray scale bar represents schematically the density of vegetation and obstacles for each locality. The results of the post hoc Tukey-Krammer test are also showed in the inserted table.

(0.11MB).