Interaction Mechanisms Quantified from Dynamical Features of Frog Choruses

TL;DR

This study combines mathematical modeling and empirical data to quantify interaction mechanisms in frog choruses, revealing how male frogs allocate attention based on behavioral strategies to optimize mate attraction.

Contribution

It introduces a Bayesian-identified phase oscillator model for frog interactions and links attention levels to behavioral parameters, providing new insights into acoustic communication strategies.

Findings

Identified a model that reproduces frog chorus dynamics

Discovered a negative correlation between attention and inter-frog distance

Indicated frogs selectively pay attention to less attractive males

Abstract

Interaction mechanism in the acoustic communication of actual animals is investigated by combining mathematical modeling and empirical data. Here we use a deterministic mathematical model (a phase oscillator model) to describe the interaction mechanism underlying the choruses of male Japanese tree frogs (Hyla japonica) in which the male frogs attempt to avoid call overlaps with each other due to acoustic communication. The mathematical model with a general interaction term is identified by a Bayesian approach from multiple audio recordings on the choruses of three male frogs. The identified model qualitatively reproduces the stationary and dynamical features of the empirical data, supporting the validity of the model identification. In addition, we quantify the magnitude of attention paid among the male frogs from the identified model, and then analyze the relationship between the attention and behavioral parameters by using a statistical model. The analysis demonstrates the biologically valid relationship about the negative correlation between the attention and inter-frog distance, and also indicates the existence of a behavioral strategy that the male frogs selectively pay attention towards a less attractive male frog so as to utilize the advantage of their attractiveness for effective mate attraction.

Equations

$${\varphi }_n(t_i)=2\pi \frac{t_i-T_{n,k}}{T_{n,k+1}-T_{n,k}}.$$

$$\frac{d{\varphi }_n(t)}{dt}$$

$$\displaystyle L_{n}=\prod_{i}\mathcal{N}\Big{[}\frac{\phi_{n}(t_{i+1})-\phi_{n}(t_{i})}{t_{i+1}-t_{i}}\Big{|}\omega_{n}+\sum_{m=1,m\neq n}^{3}\Gamma_{nm}(\psi_{nm}(t_{i})),\frac{\sigma_{n}^{2}}{t_{i+1}-t_{i}}\Big{]},$$

$$\Gamma_{nm}(\psi_{nm})=\sum_{k=1}^{M}\big{(}a_{nm}^{(k)}\cos k\psi_{nm}+b_{nm}^{(k)}\sin k\psi_{nm}\big{)}.$$

$$p_{\mathrm{post}}({\mathbf{c}}_n)\propto L_n({\mathbf{c}}_i)p_{\mathrm{prior}}({\mathbf{c}}_n).$$

$$G_{nm}({\psi }_{nm})\equiv \frac{d({\varphi }_n-{\varphi }_m)}{dt}={\omega }_n-{\omega }_m+{\Gamma }_{nm}({\psi }_{nm})-{\Gamma }_{mn}(-{\psi }_{nm}),$$

$$\frac{d{\varphi }_1}{dt}={\omega }_1+{\Gamma }_{12}({\varphi }_1(t)-{\varphi }_2(t))+{\xi }_1(t),$$

$$\frac{d{\varphi }_2}{dt}={\omega }_2+{\xi }_2(t).$$

$$\frac{d{\psi }_{12}}{dt}={\omega }_1-{\omega }_2+{\Gamma }_{12}({\psi }_{12})+{\xi }_1(t)+{\xi }_2(t).$$

$$\frac{\partial f({\psi }_{12},t)}{\partial t}=-\frac{d{\Gamma }_{12}}{d{\psi }_{12}}({\psi }_{12})f-[{\omega }_1-{\omega }_2+{\Gamma }_{12}({\psi }_{12})]\frac{\partial f}{\partial {\psi }_{12}}+\frac{{\sigma }_1^2+{\sigma }_2^2}{2}\frac{{\partial }^{2}f}{\partial {\psi }_{12}^2},$$

$$D_{\mathrm{KL}}(\hat{f}||u)={\int }_0^{2\pi }\hat{f}({\psi }_{nm})\log \frac{\hat{f}({\psi }_{nm})}{u({\psi }_{nm})}d{\psi }_{nm}.$$

$$\log \alpha ={\beta }_0+{\beta }_{\mathrm{dis}}{X}_{\mathrm{dis}}+{\beta }_{\mathrm{int}}{X}_{\mathrm{int}}+{\beta }_{\mathrm{prob}}{X}_{\mathrm{prob}}+{\xi }_{\mathrm{frog}}+{\xi }_{\mathrm{date}},$$

$$Y_{\mathrm{att}}\sim \mathrm{Gamma}(\alpha \beta ,\beta ).$$

Full text

**Interaction Mechanisms Quantified from Dynamical Features of Frog Choruses **

Kaiichiro Ota1,2, Ikkyu Aihara3,∗, Toshio Aoyagi2,4

**1 Cybozu, Inc., Tokyo, Japan

2 JST CREST, Tokyo, Japan

3 Graduate School of Systems and Information Engineering, University of Tsukuba, Tsukuba, Japan

4 Graduate School of Informatics, Kyoto University, Kyoto, Japan

$\ast$ E-mail: aihara@cs.tsukuba.ac.jp

 

Key Words: Acoustic communication, Nonlinear dynamics, Bayesian approach, Synchronization, Japanese tree frogs**

Abstract

Interaction mechanism in the acoustic communication of actual animals is investigated by combining mathematical modeling and empirical data. Here we use a deterministic mathematical model (a phase oscillator model) to describe the interaction mechanism underlying the choruses of male Japanese tree frogs (Hyla japonica) in which the male frogs attempt to avoid call overlaps with each other due to acoustic communication. The mathematical model with a general interaction term is identified by a Bayesian approach from multiple audio recordings on the choruses of three male frogs. The identified model qualitatively reproduces the stationary and dynamical features of the empirical data, supporting the validity of the model identification. In addition, we quantify the magnitude of attention paid among the male frogs from the identified model, and then analyze the relationship between the attention and behavioral parameters by using a statistical model. The analysis demonstrates the biologically valid relationship about the negative correlation between the attention and inter-frog distance, and also indicates the existence of a behavioral strategy that the male frogs selectively pay attention towards a less attractive male frog so as to utilize the advantage of their attractiveness for effective mate attraction.

1 Introduction

Animals show various types of behavior in the form of a swarm. For instance, fish and birds construct a robust and flexible flock [1]. To maintain the flock, they need to synchronize their velocity and direction with each other. On the other hand, various animals (e.g., mammals, birds, anurans, and insects) aggregate in a certain area and utilize acoustic signals for communication and also mate identification [2, 3, 4, 5]. Experimental studies demonstrate that these animals tend to alternate their acoustic signals with each other [6, 7, 4, 5]. Because such alternating behavior reduces the acoustic interference of their signals, they are likely to effectively communicate or make themselves more conspicuous within a swarm. Thus, synchronization and alternation are abundant in an animal swarm, and can work as an important indicator that determines the quality of their behavior.

To synchronize or alternate behavior, animals must recognize a specific target in a swarm. Such selective attention is reported in various systems. For example, humans pay attention to one of talking people in noisy environment like a party [8]; fish and birds attend their neighbor in a flock [9, 10]; bats pay attention to specific targets during a prey capture [11]; male insects and frogs call alternately with their neighbor [4, 5]. To understand the roles of the selective attention, it is essential for us to quantify their interaction mechanisms and compare it with the behavior of an individual animal.

This study aims to quantify interaction mechanisms inherent in the choruses of male Japanese tree frogs (Hyla japonica) by identifying a deterministic mathematical model (a phase oscillator model) from empirical data. Japanese tree frogs are observed widely in Japan, and breed mainly at a paddy field from April to July [12]. The male frogs form a lek, and produce successive calls to attract conspecific females. Indoor experiments with multiple male Japanese tree frogs demonstrate that they avoid call overlaps with each other in the forms of anti-phase synchronization of two frogs, tri-phase synchronization of three frogs, and clustered anti-phase synchronization of three frogs [13, 14] (see Figure 1 for examples of three frogs). Moreover, field observations reveal two-cluster synchronization with a larger number of the male frogs in which each pair of neighbors tend to call alternately in their natural habitat [15, 16]. Because the temporal overlap of acoustic signals generally mask the information included in each call, these alternating behavior would be important for the male frogs to effectively advertise themselves towards a conspecific female [17, 13, 14, 15, 16].

This paper is organized as follows: (1) a deterministic mathematical model (a phase oscillator model [18]) is identified for each pair of male Japanese tree frogs from our empirical data according to a Bayesian approach, (2) the attention paid among the male frogs is quantified on the basis of the identified model, and (3) the relationship between the quantified attention and behavioral parameters is analyzed by using a statistical model.

2 Results

To investigate interaction mechanisms inherent in the acoustic communication of actual animals, we analyze empirical data of male Japanese tree frogs that were obtained from our previous experiments and data analysis [14, 19] (see Sec. 4.1 for details). In each experiment, we randomly captured three male frogs at a field site, and then placed them along a straight line at the interval of $50$ cm in a stationary room. Spontaneous calling behavior of the male frogs was recorded by three microphones that were placed in the vicinity of each frog. The audio data were analyzed by the method of independent component analysis, and were separated into call signals of respective frogs. Here we use four datasets of the audio data in which all the frogs stably called more than 1400 times in four hours, allowing us to precisely estimate a phase oscillator model by utilizing the large sample size of call timing.

Figure 2 explains how we identify a phase oscillator model from the empirical data. A phase oscillator model is a well-known deterministic mathematical model derived from simple assumptions about periodicity and interaction [18] (see Sec.4.2 for the details of a phase oscillator model) that are essential for the acoustic communication of animals. Various experimental and theoretical studies have shown that a phase oscillator model can qualitatively reproduce the synchronization phenomena in actual systems [20, 21, 22, 23] including alternating chorus patterns of male Japanese tree frogs [13, 14, 15, 19]. In this study, we first calculate a phase $\phi_{n}(t_{i})$ ($n=1$, $2$, $3$) with discrete time $t_{i}$ for the $n$th frog from the separated audio data according to Eq. (1). Then, a Bayesian method [24] is applied to the time series data of $\phi_{n}(t_{i})$ so as to estimate the parameters of a phase oscillator model (i.e., an interaction term $\Gamma_{nm}(\phi_{n}(t_{i})-\phi_{m}(t_{i}))$, a natural frequency $\omega_{n}$, and a noise intensity $\sigma_{n}$). Figure 3 shows a representative result of the model identification that is obtained from a single dataset with three frogs. Green bold lines represent the interaction terms with the $95\%$ confidence interval.

To confirm the validity of the model identification, we perform numerical simulation by using the identified model and compare it with empirical data. Note that we set only one parameter $\sigma_{n}$ to be slightly larger than its estimated value throughout the following analysis because this parameter was very likely to be estimated at a smaller value because of the intermittency of the chorus over a long time scale (see Sec. 4.1 for details) that is out of the scope of the phase oscillator model. Figure 4A represents the scatter plot of phase differences that are obtained from numerical simulation of the identified model. By contrast, Figure 4B shows the scatter plot of phase differences that are directly calculated from our empirical data. Here we plot the phase differences only when one of the phases hits [math], which is consistent with our method for calculating a phase difference from discrete call timing [14, 16]. The comparison of Fig. 4A and B demonstrates that the identified model can qualitatively reproduce the experimental result, supporting the validity of our model identification. In addition, our empirical data of Fig. 4B shows a complex transition among clustered anti-phase synchronization and tri-phase synchronization, which is consistent with our previous study [14]. To investigate the mechanism of such a transition, we further analyze the relationship between a phase difference $\phi_{n}-\phi_{m}$ and its differential $d(\phi_{n}-\phi_{m})/dt$ by using the identified model (see Sec.4.4). Figure 5 demonstrates that $d(\phi_{1}-\phi_{3})/dt$ and $d(\phi_{2}-\phi_{3})/dt$ have critical states while $d(\phi_{1}-\phi_{2})/dt$ has equilibrium states. Because the transitions among multiple states can be driven by added noise under the existence of critical states in various dynamical systems, this property about critical states is likely to be the origin of the transitions observed in the empirical data of the frog choruses.

Next, we quantify the magnitude of the attention paid among the male frogs based on the result of the model identification. The stationary distribution of each phase difference is first calculated by numerically solving the Fokker-Planck equation of the identified model (see Sec. 4.5). Figure 6 demonstrates that some distributions have a sharp peak. For instance, the distribution of $\phi_{2}-\phi_{1}$ has a sharp peak around $\pi$, meaning that the $2$nd frog attempted to call alternately with the $1$st frog. Then, we calculate the Kullback-Leibler divergence of the calculated distribution from uniform distribution for each pair of the male frogs (see Sec. 4.5), which corresponds to the quantification of the attention. The right panel of Fig. 6 shows the result of the quantification in which line width represents the magnitude of the attention. Consequently, an asymmetric structure is observed in the attention paid among the male frogs. For example, the $1$st frog pays strong attention to the $2$nd frog while the $2$nd frog pays just weak attention to the $1$st frog.

To further examine the validity of the model identification, we analyze the relationship between the attention and behavioral parameters by using a generalized linear mixed model (see Sec. 4.6). In this analysis, the magnitude of the attention is treated as a response variable, and three behavioral parameters (i.e., an inter-frog distance, the probability for taking a chorus leader (leader probability), and an inter-call interval) are treated as explanatory variables of fixed factors. These variables are calculated from four datasets with $12$ frogs in which all the frogs stably produced calls. Figure 7 shows the posterior mean and the 95$\%$ confidence interval of the coefficients of respective fixed factors. This analysis demonstrates that (1) the inter-frog distance has a negative effect on the attention, meaning that male frogs pay more attention to their neighbor, and also indicates that (2) the leader probability and the inter-call interval have negative and positive effects on the attention, respectively.

3 Discussions

In this study, we quantify interaction mechanisms in the chorus of male Japanese tree frogs by identifying a phase oscillator model from empirical data. The identified model qualitatively reproduces the stationary and dynamical features of the frog choruses, which supports the validity of our model identification. Then, the magnitude of the attention paid among the male frogs is quantified on the basis of the identified model. To our knowledge, this is the first study that shows the evidence of selective attention inherent in an animal chorus by combining empirical data with a deterministic mathematical model.

The relationship between the attention and the behavioral parameters (Figure 7) gives perspectives on the choruses of male Japanese tree frogs.

• •

Inter-frog distance: Our analysis using a statistical model demonstrates the negative relationship between the attention and the inter-frog distance, which means that male Japanese tree frogs pay more attention to their neighbor. This is consistent with previous studies reporting that a neighboring pair of males alternate their calls in various species of frogs and insects [4, 5, 16, 25, 26]. Because sounds attenuate depending on distance, the calls of a neighboring pair should arrive at a higher intensity than the calls of a distant pair. We speculate that such a sound attenuation depending on distance is the origin of the negative relationship between the attention and the inter-frog distance.

• •

Inter-call interval and leader probability: In various species of frogs and insects, females prefer a conspecific male that leads their chorus [27, 28] and produces calls at a higher repetition rate [29]. Our analysis indicates that male Japanese tree frogs pay more attention to a male that rarely leads their chorus and produces calls at a lower repetition rate. These results suggest that the male frogs pay more attention to a less attractive male and then call alternately with him. Given that alternating chorus patterns can reduce the acoustic interference of their calls [17, 13], this feature would be important for male frogs to effectively advertise themselves towards conspecific females by utilizing the advantage of their attractiveness over a neighboring male.

Thus, our analysis is likely to show the important features in frog choruses, which is valid in the contexts of the acoustic communication of male frogs as well as the strategy for mating. However, it should be noted that the 95$\%$ confidence intervals of the coefficients of the inter-call interval and leader probability do not include but very close to [math] while the 95$\%$ confidence interval of the coefficient of the inter-frog distance is obviously beyond [math] (Figure 7). This suggests that the effects of the leader probability and inter-call interval on the attention are marginal compared to the effect of the inter-frog distance.

The present methodology is widely applicable to the analysis on various types of communication in animals because of the following reasons: (1) a phase oscillator model is derived from the simple assumptions about periodicity and interaction both of which are ubiquitous in the communication of animals relying on various signals (e.g., sounds [4, 5], lights [30, 31], visual display [32, 33], and electric fields [34]), (2) our methodology only requires the timing of signal emissions so as to identify the mathematical model, (3) our methodology allows us to separately identify the unidirectional interaction term of a phase oscillator model, and (4) the identified model qualitatively reproduces the stationary and dynamical features of empirical data. An important point is that the interaction term of the phase oscillator model varies depending on the value of the phase difference that corresponds to the change in the inter-signal interval among individual animals. This study demonstrates that such a dynamical property of the phase oscillator model can precisely capture not only the stationary distribution of the phase difference but also the dynamical feature of the transition among multiple synchronization states. We believe that this property of a dynamical model is advantageous compared to traditional methods (e.g., the calculation of the histogram of the inter-signal interval) when studying selective attention in animal communication that shows complicated deterministic dynamics.

4 Materials and Methods

4.1 Estimation of phase dynamics from empirical data

In this study, we use the empirical data of male Japanese tree frogs that were obtained from our previous experiments [14]. In each experiment, we randomly captured three individuals of the male frogs that were calling in a paddy field of Kyoto University (35°01’ 57.2”N, 135°47’00.0”E). Then, we put each frog in a small mesh cage, and placed three cages along a straight line at the interval of $50$cm in a stationary room. After sunset, the spontaneous calling behavior of the male frogs was recorded in the stationary room by three microphones that were placed in the vicinity of each cage. After the recording, the male frogs were released at the paddy field where they were captured. The experiment was carried out $44$ times in total between 2008 and 2009. In most trials of the experiments, some of three frogs rarely or did not call at all. Because this study aims to estimate the parameters of a mathematical model due to a statistical approach and the precision of the parameter estimation generally depends on a sample size of the empirical data, we choose four datasets in which all the three frogs stably produced calls more than 1400 times in four hours; the four datasets were obtained from the experiments performed on $26$th May, $16$th June, $17$th June in 2008, and $29$th May in 2009, respectively. All the experiments and collection of Japanese tree frogs were carried out within the facility of Kyoto University in accordance with the guideline approved by the Animal Experimental Committee of Kyoto University. Then, the audio data of the four datasets were analyzed by the method of independent component analysis, and were separated into call signals of respective frogs [14]. In addition, we carefully checked the quality of the sound-source separation and excluded just three choruses from 143 choruses in which the separation did not work well [19]. From the separated audio data, we estimated the call timing of respective frogs according to the method of Ref.[14].

For each dataset, the call timing of three male frogs is described as $T_{n,k}$ that represents the timing of the $k$th call vocalized by the $n$th frog ($n=1$, $2$ or $3$). By using the call timing $T_{n,k}$, we calculate a phase $\phi_{n}(t_{i})$ for the $n$th frog at discrete time $t_{i}$ as follows [35, 14, 15]:

[TABLE]

This calculation is performed for the set of $T_{n,k}$ and $T_{n,k+1}$ that satisfy the conditions of $T_{n,k}\leq t_{i}<T_{n,k+1}$ and $T_{n,k+1}-T_{n,k}<0.9$ sec. The first condition (i.e., $T_{n,k}\leq t_{i}<T_{n,k+1}$) is necessary for restricting the phase in the range of $0\leq\phi_{n}(t_{i})<2\pi$. The second condition (i.e., $T_{n,k+1}-T_{n,k}<0.9$ sec) is assumed because of the property of the frog chorus in which male Japanese tree frogs intermittently start and stop their periodic calling behavior over a long time scale. Namely, each frog periodically calls almost at the interval of $0.3$ sec for several tens of seconds, stays silent for several minutes, and then repeats this cycle [19]. Because such an intermittency over a long time scale is out of the scope of the phase oscillator model, we assume the second condition and then calculate the phase $\phi_{n}(t_{i})$ only during the periodic calling behavior. Consequently, the phase $\phi_{n}(t_{i})$ increases from [math] to $2\pi$, and then reset to be [math] when the $n$th frog produces a call during the periodic calling behavior.

4.2 Phase oscillator model

To reproduce the stationary and dynamical features in the choruses of male Japanese tree frogs, we use a phase oscillator model [18] with additive noise as follows:

[TABLE]

where $\phi_{n}(t)\in[0,2\pi)$ ($n=1$, $2$, …, $N$) is the phase of the $n$th frog. $\omega_{n}$ is a positive parameter that describes the intrinsic angular velocity of the $n$th frog. We then assume that the $n$th frog produces a call at $\phi_{n}(t)=0$, which is consistent with the definition of Eq. (1). Consequently, $2\pi/\omega_{n}$ gives the intrinsic inter-call interval of the $n$th frog [13, 14, 15, 19]. $\Gamma_{nm}(\phi_{n}(t)-\phi_{m}(t))$ is the unidirectional interaction term of the phase oscillator model, and is defined as a $2\pi$-periodic function of $\phi_{n}(t)-\phi_{m}(t)$ [18]. In the context of acoustic communication amomg male frogs, this term represents how the $n$th frog controls his call timing by responding to the calls of the $m$th frog [13, 14, 15]. $\xi_{n}(t)$ is the term of additive noise. We assume that this term is given by independent white Gaussian noise with magnitude $\sigma_{n}$ that satisfies the relationship $<\xi_{n}(t)\xi_{n}(s)>=\sigma_{n}\delta(t-s)$. Note that we utilize this term $\xi_{n}(t)$ so as to better fit the model parameters from empirical data that contain noisy component.

4.3 Identification of a phase oscillator model

Given the time series data of the phase $\phi_{n}(t_{i})$ with discrete time $t_{i}$ that was obtained from our experiments using actual frogs (see Sec. 4.1), we estimate the unknown parameters of the phase oscillator model according to the Bayesian method of Ref. [24]. In the method, we first define the following likelihood function from the time series data of the phase $\phi_{n}(t_{i})$:

[TABLE]

where $\mathcal{N}(x|m,s^{2})$ represents a Gaussian function with mean $m$ and variance $s^{2}$, and $\psi_{nm}$ denotes the phase difference $\phi_{n}-\phi_{m}$. Because each interaction term is defined as a $2\pi$-periodic function of the phase difference $\psi_{nm}$ [18], it can be expanded into a Fourier series as follows:

[TABLE]

We apply a standard model comparison method to determine the maximum order $M$ in Eq.(4). That is, we choose $M$ such that a marginal likelihood function $\mathcal{L}_{M}$ is maximized. Technically, we computed and compared $\mathcal{L}_{M}$ for $M=1,2,\dots,30$. For our data, the maximum order $M$ was chosen between 1 and 9.

Then, the parameters to be estimated are $\omega_{n},a_{nm}^{(1)},\dots,a_{nm}^{(M)},b_{nm}^{(1)},\dots,b_{nm}^{(M)}$ and $\sigma_{n}$, which are denoted by a shorthand notation $\bm{c}_{n}$. To estimate $\bm{c}_{n}$ in a Bayesian framework, we use a reasonable conjugate prior distribution $p_{\mathrm{prior}}(\bm{c}_{n})$, which is a Gaussian-inverse-gamma function. Bayes’ theorem then gives the posterior parameter distribution as follows:

[TABLE]

Because the prior distribution $p_{\mathrm{prior}}(\bm{c}_{n})$ is conjugate to the likelihood $L_{n}(\bm{c}_{i})$, we can easily calculate the posterior distribution $p_{\mathrm{post}}(\bm{c}_{n})$ (see Ref.[24] and its supplemental material for detail).

4.4 Examination of critical states

Empirical data on the choruses of male Japanese tree frogs demonstrate the complicated transition among multiple synchronization states (Fig.4B) [14]. Here we examine the origin of this transition on the basis of the estimated model. Equation (2) without a noise term yields the time evolution of the phase difference $\psi_{nm}\equiv\phi_{n}-\phi_{m}$ as follows:

[TABLE]

where $\Gamma_{nm}(\psi_{nm})$ is the maximum a posteriori (MAP) estimation of the interaction term that was obtained from the analysis of Sec. 4.3. This function $G_{nm}(\psi_{nm})$ quantifies how the phase difference changes in time under the mutual interaction between the $n$th and $m$th frogs without the effect of another frog. If $G_{nm}(\psi_{nm})$ has a zero-crossing (a stable equilibrium) at $\psi_{nm}=\psi_{nm}^{*}$ , the phase difference approaches $\psi_{nm}^{*}$ and then stays near the value forever. Such stationary dynamics typically come from strong interaction. By contrast, if $G_{nm}(\psi_{nm})$ has no zero-crossing but is very close to zero at a certain point $\psi_{nm}=\psi_{nm}^{*}$ (a critical state), the phase difference stays near $\psi_{nm}^{*}$ for a certain amount of time and eventually departs the point. Such critical dynamics typically come from moderate interaction.

4.5 Quantification of selective attention

To assess selective attention among male frogs, we evaluate the stochastic feature of the estimated model by taking the effect of noise into consideration. For example, to evaluate the attention from the $1$st frog to the $2$nd frog, we analyze the following model only with the unidirectional interaction term $\Gamma_{12}(\phi_{1}(t)-\phi_{2}(t))$:

[TABLE]

These equations describe the situation that the $1$st frog exposed to the noise term $\xi_{1}(t)$ pays attention to the $2$nd frog according to the interaction term $\Gamma_{12}(\phi_{1}(t)-\phi_{2}(t))$ (Eq. (7)) but the $2$nd frog exposed to the noise term $\xi_{2}(t)$ does not pay any attention to the $1$st frog (Eq. (8)). Hence, this is a concise mathematical model capturing the attention from the $1$st frog to the $2$nd frog under the effect of the noise $\xi_{1}(t)$ and $\xi_{2}(t)$. Subtracting Eq. (8) from Eq. (7) yields the time evolution of the phase difference $\psi_{12}\equiv\phi_{1}(t)-\phi_{2}(t)$ as follows:

[TABLE]

Then, we calculate the stationary distribution of the phase difference $\psi_{12}$ by numerically solving the Fokker-Planck equation of Eq.(9) that gives the time evolution of the distribution of the phase difference (i.e., $f(\psi_{12},t)$),

[TABLE]

until it converges. Subsequently, we can calculate the stationary distribution of the phase difference for all the pairs of male frogs, and then describe it as $\hat{f}(\psi_{nm})$.

The stationary distribution of a phase difference ($\hat{f}(\psi_{nm})$) allows us to quantify the degree of attention. If there is no attention from the $n$th frog to the $m$th frog, the distribution of $\psi_{nm}$ is almost uniform; on the other hand, if the $n$th frog calls synchronously with the $m$th frog as the result of paying strong attention to him, the distribution has a sharp peak. In this study, we quantify such a sharpness of the distribution $\hat{f}(\psi_{nm})$ by using the Kullback-Leibler divergence from uniform distribution $u(\psi_{nm})\equiv\frac{1}{2\pi}$ as follows:

[TABLE]

Consequently, a nearly uniform distribution is represented by $D_{\mathrm{KL}}(\hat{f}||u)\sim 0$ while a sharply peaked distribution is represented by $D_{\mathrm{KL}}(\hat{f}||u)$ that is much larger than [math].

4.6 Relationship between selective attention and behavioral parameters

To further examine the validity of the model identification, we analyze the relationship between selective attention and behavioral parameters. Here we focus on the following behavioral parameters of male Japanese tree frogs: (1) an inter-frog distance, (2) an inter-call interval, and (3) leader-follower relationship. An inter-frog distance represents the distance between each pair of male frogs that was measured in our experiments [14] (see Sec. 4.1). An inter-call interval is calculated as $\delta T_{n,k}=T_{n,k+1}-T_{n,k}$ using the sequences of the call timing $T_{n,k}$ only when the condition $T_{n,k+1}-T_{n,k}\leq 0.6$ sec is satisfied. The leader-follower relationship is determined according to the following definition: the leader, the $1$st follower, and the $2$nd follower are defined as the frogs that start calling first, second, and third within the same chorus, respectively. Note that male Japanese tree frogs start calling with low-intensity sound, making it difficult for us to automatically determine the leader-follower relationship. Hence, we manually determined the leader-follower relationship by carefully looking at all the separated audio data and calculated the probability for taking a chorus reader for each frog.

Next, we analyze the relationship between the magnitude of attention and the three behavioral parameters by using a generalized linear mixed model (GLMM). GLMM is a well-known statistical model that is used in various research areas so as to analyze the effects of multiple explanatory variables on a response variable [36, 37]. Here, we treat the magnitude of the attention paid from the $n$th frog to the $m$th frog (see Eq.(11)) as a response variable $Y_{\mathrm{att}}$. We then treat the three behavioral parameters (i.e., the inter-frog distance between the $n$th and $m$th frogs, the inter-call interval of the $m$th frog, and the leader probability of the $m$th frog) as explanatory variables of fixed factors, and describe them as $X_{\mathrm{dis}}$, $X_{\mathrm{int}}$ and $X_{\mathrm{prob}}$, respectively. In addition, we treat frog index and experimental date as explanatory variables of random factors ($\xi_{\mathrm{frog}}$ and $\xi_{\mathrm{date}}$) because these factors are difficult to be quantified but very likely to affect the calling behavior of male frogs. Here, we assume that the random factors $\xi_{\mathrm{frog}}$ and $\xi_{\mathrm{date}}$ follow a normal distribution of zero mean value with standard deviations of $\sigma_{\mathrm{frog}}$ and $\sigma_{\mathrm{date}}$, respectively. By using these variables, we construct the following GLMM:

[TABLE]

Here, the response variable $Y_{\mathrm{att}}$ is always positive because of its definition (see Eq.(11)). To reproduce this feature, we assume that (1) the response variable follows a gamma distribution that always takes a positive value, and (2) the parameter $\alpha$ of the gamma distribution (this parameter gives the mean value of the distribution) is linked to the explanatory variables with a log function; this framework is consistent with a traditional gamma regression with multiple explanatory variables. We confirmed that there is no multicollinearity among any pairs of the explanatory variables of fixed factors (the absolute value of Pearson’s correlation coefficient is less than $0.36$). Posterior distributions of all the unknown parameters (i.e., $\beta_{0}$, $\beta_{\mathrm{dis}}$, $\beta_{\mathrm{int}}$, $\beta_{\mathrm{prob}}$, $\beta$, $\sigma_{\mathrm{date}}$, and $\sigma_{\mathrm{site}}$) were estimated from MCMC samples generated by R Statistical Software (ver. 3.4.2) and Stan (ver. 2.17.2). Note that we normalized the explanatory variables of fixed factors from [math] to $1$ prior to the calculation of the MCMC samples, and confirmed the convergence of the MCMC samples by using $\hat{R}$ with a threshold of $1.01$ [38]. The posterior mean and the 95$\%$ confidence interval of $\beta_{0}$, $\beta_{\mathrm{dis}}$, $\beta_{\mathrm{int}}$, and $\beta_{\mathrm{prob}}$ are shown in Figure 7.

Acknowledgement

We thank K. Aihara, H.G. Okuno, R. Takeda, T. Mizumoto, H. Awano, K. Itoyama, and Y. Bando for their valuable comments on this study.

Author Contributions

K.O., I.A. and T.A. designed the research; K.O. and I.A. analyzed the data; K.O. performed simulation; K.O., I.A. and T.A. wrote the paper.

Competing Interests

We declare we have no competing interests.

Funding

This study was partially supported by JSPS Grant-in-Aid for Challenging Exploratory Research (No. 16K12396) and Grant-in-Aid for Young Scientists (No. 18K18005) to I.A.