Covariation among multimodal components in the courtship display of the túngara frog

Animals communicate with signals across a wide range of ecological contexts, such as mate attraction, territory defense, danger warnings and parent–offspring interactions (Bradbury and Vehrencamp, 2011). To communicate, signaling individuals must perturb the environment around them. For example, many animals produce sound using internal air to alter the acoustic pressure field around them, a form of communication that has been under selection in many species. The act of producing a sound, however, often includes perturbations that can be perceived via other modalities, such as alterations in reflected light patterns or substrate-borne vibrations (i.e. visual and vibratory cues). In some cases, these additional cues can provide information to receivers, and are sometimes elaborated through selection into additional signals themselves (Halfwerk et al., 2019; Taylor and Ryan, 2013).

Although the multiple components associated with a communication signal may be perceived by receivers using a single sensory system (‘multi-component communication’), receivers often use distinct sensory systems for different components (‘multimodal communication’; Halfwerk et al., 2019; Higham and Hebets, 2013; Partan and Marler, 1999; Partan and Marler, 2005; Taylor and Ryan, 2013). For example, various spiders can detect vibratory, visual and/or chemical cues in courtship displays (Hebets et al., 2013; Uetz and Roberts, 2002; Uetz et al., 2017), bees can detect visual and chemical cues from flowers (Leonard and Masek, 2014), and many primates (including humans) can detect acoustic and visual cues from the face during vocalizations (McGurk and MacDonald, 1976; Micheletta et al., 2013; Partan, 2002; Wilke et al., 2017). These multimodal cues have been proposed to help signalers and/or receivers by improving detection, localization and discrimination, or providing additional information, among other benefits (Mitoyen et al., 2019). Of particular interest has been the courtship displays of birds, and the intricately timed visual and acoustic displays of some bird species that may convey overall motor proficiency to receivers (Byers et al., 2010; Fusani et al., 2014; Hogan and Stoddard, 2018; Miles and Fuxjager, 2018; Ullrich et al., 2016).

The first studies of multimodal communication primarily analyzed behavior using ‘cue-isolation’ experiments, and paid less attention to the natural variation and covariation present in and among cues. Researchers have proposed methods to extend cue-isolation studies and map receiver behavior onto the full range of natural variation among multiple signal components (Halfwerk et al., 2019; Mitoyen et al., 2019; Ronald et al., 2017; Smith and Evans, 2013). Recent studies that utilize such methods have found that comparing across modalities can be essential to determining receiver responses (Halfwerk et al., 2014c; Ręk and Magrath, 2017), and that interactions between different modalities can result in previously unobserved behaviors (Ronald et al., 2017; Stange et al., 2017). Phenotype network and systems approaches have also shown how integrating across components can provide valuable insight into signal form/function relationships (Hebets et al., 2016; Reichert and Höbel, 2018; Wilkins et al., 2015).

Understanding the extent to which different components of a multimodal display covary is critical to understanding the functionality of these displays. But a full understanding requires deciphering how animals respond to multimodal stimuli, which necessitates perceptually focused analyses on how receivers bound the multimodal components into a single percept (Halfwerk et al., 2019; Taylor and Ryan, 2013). As a first step in this pursuit, we measured the natural variation within and covariation among the components of a multimodal communication system and asked how that covariation influences potential information for receivers.

The túngara frog [Engystomops (=Physalaemus) pustulosus] mating display offers an excellent model system to investigate how receivers integrate multimodal components of displays across the full landscape of natural variation contained among those display components. Like many other species of frogs, males produce calls to attract females. Each túngara frog call contains an acoustic signal that consists of a downward sweeping ‘whine’ that is facultatively followed by one or more ‘chucks’. The whine by itself is necessary and sufficient to attract females to males and also for females to discriminate between conspecifics and heterospecifics (Ryan and Rand, 1995). The addition of a chuck increases the male's attractiveness to the female five-fold (Gridi-Papp et al., 2006; Ryan et al., 2019). Females also show a significant, though substantially weaker, preference for lower frequency chucks (Bosch et al., 2002; Ryan, 1980, 1983) and lower frequency whines (Bosch et al., 2000), both of which are produced by larger males. At the mechanistic level, lower frequency chucks elicit greater neural excitation in the female's inner ear (Ryan et al., 1990). At the ultimate level, females typically gain a reproductive advantage choosing larger males as there is a better mechanical fit between the male and the female cloacas during amplexus, which enhances the efficiency of fertilization (Ryan, 1985). Thus, the correlation between the spectral characteristics of the call and male body size is relevant to fitness effects of female mate choice.

As with most frogs, the primary component of the túngara frog's mating display is the acoustic signal. However, additional components can be integrated with the call, resulting in a multimodal display (Grafe et al., 2012; Höbel and Kolodziej, 2013; Reichert and Höbel, 2015; Rosenthal et al., 2004; Starnberger et al., 2014; Taylor et al., 2007, 2008). Specifically, male túngara frogs call by inflating a vocal sac, which results in a highly conspicuous visual cue during each call. This visual cue can affect female perception and decision-making (Rosenthal et al., 2004; Stange et al., 2017; Taylor and Ryan, 2013; Taylor et al., 2008). Furthermore, male túngara frogs call in shallow pools, with their calls creating ripples across the water surface. These ripples have been shown to affect receivers, including the fact that variation in ripple height can influence calling in rival males (Halfwerk et al., 2014b,c). Receivers of the túngara frog call thus potentially have access to at least three distinct components of the sexual display, which are likely to be perceived in distinct sensory modalities: hearing (call), vision (vocal sac) and tactile/vibration (ripple). In this study, we asked how information from these different modalities might allow receivers to predict male body size, and whether the interaction of these cues increases information about body size.

Eavesdroppers abound in this communication system, and potential predators and parasites also have access to these components (Bernal and de Silva, 2015). Indeed, the fringe-lipped bat (Trachops cirrhosus) can hunt male frogs by eavesdropping on their calls (Tuttle and Ryan, 1981), and prefers to attack calls presented with a simultaneous vocal sac inflation cue or ripple cue (Gomes et al., 2016; Halfwerk et al., 2014a,b; Rhebergen et al., 2015). Therefore, this system allows us to ask how multiple receivers with very different goals – searching for a meal versus searching for a mate – respond to the full landscape of multimodal variation within the same stimuli.

To understand the independent and cumulative effects of each of these three display components for receivers, we must first understand their natural variation and covariation. To this end, we have measured these components in the calls of wild túngara frogs, assessed their variation and covariation, and asked to what degree each component can predict body size, a trait that influences male mating success, female mate choice and fecundity (Ryan, 1985).

We collected túngara frogs [Engystomops pustulosus (Cope 1864)] in and around Gamboa, Panamá (9°07.0′N, 79°41.9′W), shortly after sunset. Frogs were tested in a Smithsonian Tropical Research Institute (STRI) laboratory and returned to the collection site within 2 days. To identify any recaptures, we toe-clipped all males that we successfully induced to produce calls in the experimental setup following Guidelines for the Use of Live Amphibians and Reptiles in Field Research (Beaupre et al., 2004). We also measured each frog's mass (in g) and snout–vent length (SVL; in mm). All procedures were approved by the University of Texas at Austin (IACUC: AUP-2019-00067), STRI (IACUC: 2018-0411-2021) and the Ministry of the Environment of Panamá (MiAmbiente: SE/A-40-19).

We developed an apparatus to simultaneously measure three focal features of the túngara frog display: call amplitude, vocal sac inflation size and ripple height (Fig. 1A). We placed frogs within a circular pool of water (diameter 42.5 cm) inside a ‘ripple cage’ developed to be transparent to ripples propagated from the calling frog (Halfwerk et al., 2016). The water depth was maintained at 2.2 cm from the surface to the bottom of the ripple cage for each frog. Túngara frog males typically begin a call bout with very low amplitude calls and then quickly ramp up their calling to a relatively stable plateau of consistent calls (Fig. 1B). To capture inter- and intra-individual variation in call parameters for correlation analyses (see below), we measured the first 15, middle 15 and final 15 calls produced by each male (n=54 males, n=2430 calls). To obtain average peak calling measurements for each individual frog (to compare with body size), we took the mean for each feature from the final set of calls. However, when frogs ceased calling on their own, they often produced a very weak final call, so we excluded the last call for this measure of peak calling, averaging across the penultimate 14 calls.

We recorded calls using an omnidirectional AKG microphone (model C417) placed 19.6 cm above the ripple cage (Fig. 1C). While the placement of the microphone and circular shape of the cage ensured that the distance from the frog to the microphone was quite consistent, there could at most be a 1 cm difference in distance, which would cause an ∼0.43 dB decrease in measured sound level. To calibrate call amplitude, we broadcast a tone from a speaker and measured the sound pressure level (SPL) using a GenRad SPL meter just above the microphone. The microphone was connected to a one-channel Behringer preamp, and logged into a laptop computer using a National Instruments data acquisition device (DAQ) and LabVIEW 2019. We digitized the acoustic signal and the ripple height vibrometer signal (see below) simultaneously at a sample rate of 25.6 kHz.

As noted above, each túngara frog call consists of a downward sweeping whine and is sometimes followed by one or more chucks. Calls that that are followed by chucks are referred to as ‘complex calls’. Because all calls contain a whine, we took the maximum amplitude of the whine only. To calculate the amplitude metric, we first squared the recorded pressure signals. A Gaussian window 5 ms in length with a time constant of 2 ms was applied and the area under the window was calculated. This area has the units of Pa² s and is often referred to as a windowed sound exposure. The window was shifted by one sample and the calculation was repeated for the length of each whine. We conducted this smoothing procedure in order to ensure that artifacts in the recordings could not have large impacts on our measurements of peak amplitude. The amplitude was converted to dB re. 20 µPa calibrated by comparison to the tone played back immediately after each frog's recording session.

We also measured the dominant frequency (DF) of each whine, defined as the peak frequency of a fast Fourier transformation of the entire whine recording. Because the DF was unreliable in the quietest whines, we limited this analysis to the middle and final 15 calls analyzed. The DF of a whine is nearly always in the fundamental of the downward sweeping harmonic; however, three measurements were extreme outliers because the DF was in the second harmonic. We excluded these three whines from analysis (n=1617 total DF measurements for analysis).

Ripple height was measured using a Polytec PDV 100 portable laser Doppler vibrometer simultaneously sampled along with the acoustic signal, as described above. We aimed the vibrometer at a small piece of reflective tape adhered to a piece of aluminium foil to provide flotation positioned 12 cm from the edge of the ripple cage and 8 cm from the edge of the pool. We attached this tape to a weight using nylon thread to maintain the tape's position. Using the velocity signal from the vibrometer and the distance from the vibrometer to the water surface, we calculated the maximum ripple height associated with each call. For this, we took a 1.5 s window from the vibrometer signal following the onset of each whine (as annotated from the airborne audio signal; Fig. 1D). We took the integral of this signal (to convert mm s⁻¹ to mm) and measured the maximum peak-to-peak amplitude (i.e. the difference between the greatest positive and negative displacement of the water surface following each call). We also calculated and plotted the energy spectral density of the displacement for a series of example ripple measurements (Fig. 1E).

Ripple amplitude declines with propagation; therefore, we measured the distance from the frog to the reflective tape (using the overhead camera stills, the distances ranged from 12 to 15.4 cm). We then included distance (d; Fig. 1A) to the vibrometer as a covariate in the analyses comparing the ripple height with whine amplitude and vocal sac size. Distance did not significantly correlate with ripple height in any of the models where it was included, nor did distance affect the significance of other covariates.

Any call that was produced right before, during or after the frog kicked its legs or repositioned itself within the cage was excluded from the analysis. These excluded values were typically outliers because the ripples created by swimming were much larger than those created simply by calling. This resulted in n=2134 ripple measurements. We observed the video (see below) before and after each ripple measurement to determine which ripple data points met the inclusion criterion.

We used an overhead camera to measure the maximum volume of each vocal sac inflation. We used a Baesler camera (model: acA1920-150um) and captured video using the Pylon viewer program. We captured 25 frames s⁻¹ with 39 ms per exposure.

The shape of the vocal sac is roughly a spheroid. Thus, the volume is estimated as: ⁠, where a is the semi-minor axis (as viewed from the top camera) and b is the semi-major axis. In this regard, it is assumed that the height is equal to the minor axis from above (Dudley and Rand, 1991) (Fig. 1D).

Calibration videos were taken with each frog tested. A transparent ruler was floated on the water's surface just next to the ripple cage. Using frames from this calibration, we used ImageJ to convert pixels to millimeters, and then we manually drew lines on the frames containing the peak vocal sac inflation to calculate each of a and b as described above. Two naïve researchers performed these video analyses, and data from 135 images analyzed by both researchers suggested reasonable inter-rater reliability (linear correlation: r=0.85, P<0.0001).

We also confirmed the calibration by measuring each frog's SVL from a still frame of the frog floating and compared with the SVL manually measured in the field. These two measures of SVL were highly correlated (r=0.92, P<0.0001) and we used those measured in the field for analysis.

We induced males to call by broadcasting a field recording of a large túngara frog chorus shortly after placing them individually into the recording apparatus (Fig. 1A). This apparatus was contained within a noise-attenuating walk-in chamber (Acoustic Systems, ETS-Lindgren, Austin, TX, USA; 2.7×1.8 m and 2 m high). The low amplitude of this chorus and the microphone placement directly above the focal frog were such that one could easily extract the spectral detail of the focal male's calls. We generally allowed males that were waiting to be tested to interact with any females with whom they were collected. We tested males found in amplexus before those that remained single because we observed that frogs recently in amplexus were much more likely to begin calling in the experimental setup. Each frog was given at least 10 min to begin calling before moving on to the next individual. The experiment was terminated if or when a frog produced at least 50 complex calls in rapid succession.

To analyze the covariation among the call amplitude, vocal sac size and ripple height, we ran linear mixed-effects models in R with the ‘lme4’ package (https://CRAN.R-project.org/package=lme4) on measurements from the first 15 calls produced by each male, the last 15 calls and the middle 15 calls, to capture intra-individual variation. Specifically, we used one measure as the independent variable, a second as the dependent variable and frog ID as a random factor. We assessed the significance of each model with type III ANOVA F-tests with Satterthwaite's method using the ‘lmerTest’ package (https://CRAN.R-project.org/package=lmerTest). We also calculated and visualized simple linear correlations between pairs of features within each individual frog using MATLAB. Furthermore, we report the marginal correlation coefficients (representing variance explained by the fixed effect) for each correlation using the ‘MuMIn’ package (https://CRAN.R-project.org/package=MuMIn).

To compare these three call features to body size, we focused on each individual's ‘peak calling’, measured as the average value for each feature across the final 15 calls, with the very last call excluded (see above). We compared these peak calling values with our two measures of body size (SVL and mass) using simple linear correlations. We also used model selection and dimensionality reduction techniques to compare unimodal, bimodal and trimodal models to predict body size. We conducted model selection in R using the ‘MuMIn’ package (https://CRAN.R-project.org/package=MuMIn) with all possible models (including all possible interaction terms) and report those with a delta <6. Because the call features are correlated among themselves, we conducted a principal component analysis (PCA) with z-scores and included a model with the first component (PC1) as the sole predictor of body size in the model selection process. In this way, we could assess whether a model with all three call features included in a composite score was a better predictor of body size than a unimodal model. Because we measured body size in two ways, we ran all of these analyses for each measurement separately. Data are available through figshare at: https://doi.org/10.6084/m9.figshare.14488971.

In this study, we measured three main components from 45 calls produced by each of 54 male túngara frogs (n=2430 calls). Namely, we compared the peak amplitudes across three modalities: auditory (peak whine amplitude), visual (peak vocal sac size) and seismic/tactile (peak ripple height). To maximize within-individual variation, we analyzed the first 15 calls produced by each male, the last 15 calls and the middle 15 calls. Fig. 1B depicts all 45 calls produced by an example male in this study. We also measured DF, a secondary auditory component that is known to affect some receivers (albeit with less robust effects compared with amplitude). Table 1 contains the data demonstrating the natural variation we observed across all males. Furthermore, we compared the average coefficient of variation (CV) of each measurement within each frog with the overall CV using data from all frogs. We observed substantial variation in each feature within each frog, with the average CV in each frog quite similar to the overall CV.

Next, we asked what natural covariation existed between pairs of these call components. The data for ripple height and vocal sac size were log transformed. We found significant and positive correlations between each pairwise comparison of components using mixed-effects models with frog ID as a random effect (whine amplitude and vocal sac size: r=0.789, F_1,2405.7=6865.9, P<0.0001; whine amplitude and ripple height: r=0.092, F_1,2128.3=13.6, P=0.0002; vocal sac size and ripple height: r=0.116, F_1,2118.1=20.3, P<0.0001; Fig. 2A–C). Among these analyses, the relationship between whine amplitude and vocal sac size was clearly the strongest. We also analyzed these correlations within each frog and plotted each frog's slope (Fig. 2D–F) and Pearson's correlation coefficient r (Fig. 2G–I). For the relationship between vocal sac size and whine amplitude, every frog had a positive slope and the mean r-value was 0.89, which was substantially higher than the r-values for the two comparisons involving ripple height.

We also conducted a PCA to represent the range of variation across these three modalities (Fig. 3). We found that PC1 explained 61.5% of the variation in the data, with coefficients of 0.68, 0.68 and 0.28 for whine amplitude, vocal sac size and ripple height, respectively. PC2 explained 31.1% of the variation in the data, and had coefficients of −0.16, −0.23 and 0.96 for whine amplitude, vocal sac size and ripple height, respectively.

Finally, we asked to what degree these three call components could predict the body size of the calling frog. We measured body size in two ways – SVL and mass – and these two measures were correlated (n=54, r=0.7430, P<0.0001). We compared these sizes with each frog's ‘peak calling’, for which we took his final 15 calls, excluded the very last call (as some frogs produced a very quiet final call) and calculated the mean. For SVL, we detected a significant relationship for vocal sac size (n=54, r=0.379, P=0.0047), but not whine amplitude (n=54, r=0.136, P=0.3285) nor ripple height (n=54, r=–0.014, P=0.9177; Fig. 4A–C). For mass, we found stronger correlations with all three call components, including significant relationships for whine amplitude (n=54, r=0.531, P<0.0001) and vocal sac size (r=0.573, P<0.0001), and a non-significant relationship for ripple height (n=54, r=0.184, P=0.1826; Fig. 4D–F).

To assess whether multimodal cues could provide additional information, we compared unimodal, bimodal and trimodal models as predictors of body size. Because we found significant correlations among call features, we also ran a PCA that included the peak calling values of all three features across every frog, and conducted models with PC1 as the sole predictor of body size. Thereafter, we used model selection techniques with every possible unimodal, bimodal and trimodal model (as well as every possible combination of interaction terms) included alongside the model with PC1. For SVL, we found that a unimodal model with just the vocal sac size fit best (Table 2) and was similar in explanatory power to bimodal or trimodal models. However, for mass, we found that bimodal models including both sac size and whine amplitude, as well as the PC1 and a trimodal model, all fit the data better than the unimodal model with just the vocal sac size (Table 2).

Finally, we analyzed how DF compared with the three primary measures of amplitude (Fig. S1). We observed significant negative correlations between whine amplitude and DF (r=–0.143, F_1,1197.6=12.8, P=0.0004; Fig. 2A) and vocal sac size and DF (r=–0.119, F_1,1390.7=11.4, P=0.0007). We observed a non-significant positive trend for the relationship between ripple height and DF (r=0.035, F_1,1461.3=3.2, P=0.0750). On a frog-by-frog basis, we measured the mean DF for each frog's peak calling (Fig. S1). We observed no relationship between SVL and DF (n=54, r=0.005, P=0.9689), but we did observe a significant negative relationship between mass and DF (n=54, r=–0.310, P=0.0223). When we ran model selection using all four features (and a PC1 on all four features), we observed very similar results. Sac size alone was still the best model for predicting SVL, and multimodal models were still better for predicting mass than unimodal models.

When animals communicate, they often create disturbances in the environment that receivers could perceive using multiple modalities. Some of these disturbances have been favored by selection and are thus evolved signals, such as the túngara frog's acoustic call. In contrast, other disturbances are merely incidental artifacts or cues, such as the inflating vocal sac and the water-borne ripples generated by vocal sac movements, which are unlikely to have been under selection directly. Regardless, receivers may still garner information from all of these components. To understand the range of information potentially available to the receiver, it is thus important to know to what degree these different components covary.

Here, we found that each of three call components, whine amplitude, vocal sac size and ripple height, were correlated with each other across individual frogs. However, we also found that both within and across frogs, the relationship between whine amplitude and vocal sac size was much stronger than either of the relationships with ripple height, indicating that vocal sac size relates to sound propagation. Interestingly, the weak relationship between vocal sac size and ripple height suggests that the ripple is not just a byproduct of the vocal sac inflation, but may also result from other body movements during the production of the call (for example, movement of the arms and/or legs, and movement of the body wall owing to inflation and deflation of the lungs). If so, it is possible that calling males might alter the heights of their ripples independent of changes to vocal sac size and whine amplitude.

We also asked whether any of these three components could predict body size. Although whine amplitude significantly correlated with mass, we found that the best predictor of body size, when measured as either the SVL or mass, was the vocal sac size during peak calling. These results contrast with the fact that the acoustic signal has been shown to be the most important component in preference tests for both túngara frog females and fringe-lipped bats (Gomes et al., 2016; Rosenthal et al., 2004; Taylor et al., 2008). Although it is logical that predatory bats should prefer larger prey to smaller ones, it is yet unknown whether the small variation we see in male túngara frog size is salient to hunting bats; we do know that differences in the chuck's dominant frequency, which is correlated to male body size and which influences female choice, has no influence on the bats' attraction to calls (M. J. Ryan and M. D. Tuttle, unpublished data). Future work assessing how receivers respond to variation in the vocal sac cue will be useful in determining whether the information contained in this cue is used by receivers. Indeed, it is also possible that variation in sensory acuity or sensitivity among receivers affects successful choices for larger males (Cummings et al., 2008; Ronald et al., 2018; Ryan et al., 1992).

We also asked whether access to multimodal components could provide additional information about male body size. For male mass, we found that bimodal and trimodal models had better fit than any unimodal model. Even when we accounted for the covariation among the trimodal components with a PCA, we still found that the composite trimodal dataset fit better than any unimodal model. That being said, how receivers perceptually bind these components together, which is unlikely to directly resemble a PCA, remains unknown. We also note that these components will be important to receivers beyond their potential information on body size. For instance, how these components together influence saliency and the ability to localize signalers will be important to test in future studies.

Túngara frog females gain an advantage in reproductive success by choosing larger males. Previous studies have shown that females prefer lower frequency chucks as well as lower frequency whines, both of which are correlated with male size (Bosch et al., 2000; Ryan, 1980, 1983; Wilczynski et al., 1995). Preferences for these lower frequency calls are favored by selection as they result in females mating with larger males, who fertilize more of the female's eggs (Ryan, 1980, 1985). Interestingly, the inner ear tuning that seems to contribute to the preference for lower frequency chucks did not evolve in túngara frogs but is an ancestral trait that existed long before these types of frogs evolved chucks (Wilczynski et al., 2001). The results from the present study highlight that preferences for features other than frequency could be used to predict body size, and that, even when accounting for dominant frequency, multimodal integration can better predict mass.

Finally, there are other measurements within the túngara frog call that could affect receivers. Here, we primarily analyzed only one measurement (peak amplitude) per component. Within each of these components, however, there are other measurements that may be important. For instance, we know that call frequency can affect female preference in túngara frogs (Bosch et al., 2000; Wilczynski et al., 1995), and we found that dominant frequency of the whine was correlated with other call components and body mass. How call complexity interacts with these features remains to be explored (Ryan, 1980, 1983). Much less is known about what features within the vocal sac inflation and ripple propagation could be important, and future work could ask whether receivers respond to differences in aspects such as vocal sac inflation rate or ripple propagation speed. Furthermore, additional modalities can affect receivers; hormones released into the water, for example, affect calling in rival males (Still et al., 2019).

Overall, these data demonstrate how multiple components of a multimodal communication system relate to one another. We also demonstrate that the integration of multimodal components can provide additional information about body size, an important trait for reproductive success. We hope this descriptive dataset will provide the foundation for future studies addressing how this natural variation and covariation is integrated to inform perception and decision-making in both mating and foraging contexts.

We thank Emily G. Roberts for help with video analysis and the members of the Ryan lab for feedback on the project. Gregg Cohen and the Smithsonian Tropical Research Institute provided equipment and logistical support.