Collective dynamics in two populations of noisy oscillators with asymmetric interactions
We study two intertwined globally coupled networks of noisy Kuramoto phase oscillators that have the same natural frequency, but differ in their perception of the mean field and their contribution to it. Such a give-and-take mechanism is given by asymmetric in- and out-coupling strengths which can be both positive and negative. We uncover in this minimal network of networks intriguing patterns of discordance, where the ensemble splits into two clusters separated by a constant phase lag. If it differs from , then traveling wave solutions emerge. We observe a second route to traveling waves via traditional one-cluster states. Bistability is found between the various collective states. Analytical results and bifurcation diagrams are derived with a reduced system.
pacs:05.40.-a, 05.45.Xt, 87.10.Ca
From neurons in the brain to cells in the heart, from electrons in superconductors to planets in the universe, collective oscillations are ubiquitous Strogatz (2003); ?. Clearly, the underlying mutual synchronization critically relies on the presence of interactions. If individuals are allowed to interact non-uniformly, the collective behavior remains particularly elusive. Beyond incoherence and synchronization, what other collective states can emerge, and under which conditions? Motivated by Daido’s seminal work on “oscillator glasses” Daido (1992), Hong and Strogatz tackled this question recently in a series of papers Hong and Strogatz (2012, 2011a); ?. They investigated two mutually globally coupled populations of Kuramoto phase oscillators that differed in their coupling strengths. Two scenarios of mixed attractive and repulsive interactions were distinguished. In the first case some oscillators’ phases repel the phases of all the others, while the remaining attract all the other phases Hong and Strogatz (2012). This situation resembles neural networks with excitatory and inhibitory connections Wilson and Cowan (1972); ?. In the second case some oscillators tend to align with the mean field, while others oppose it favoring an antiphase alignment Hong and Strogatz (2011a); ?. This version is analogous to sociophysical models of opinion formation Miguel et al. (2005). Surprisingly, only the second scenario led to enriched dynamics beyond the traditional order-disorder transition.
In this paper, we unify both coupling scenarios. We unveil the existence of two routes to traveling waves, which in parameter space are surrounded either by states of diametral two-cluster synchronization or by one-cluster partially synchronous states. Differently than in Refs. Hong and Strogatz (2012, 2011a); ? we consider temporal fluctuations acting on the frequencies. Our main contribution is the derivation of bifurcation diagrams for all possible collective states.
One of the most prominent models describing phenomena of mutual synchronization is due to Kuramoto Kuramoto (1984). It describes how the phases of coupled oscillators evolve in time, and is applicable to systems of nearly identical, weakly coupled limit-cycle oscillators. We consider a stochastic version with twofold disordered coupling strengths:
where . Using the notion of give-and-take as a metaphor, oscillator contributes to the mean field with weight and at the same time incorporates the mean activity with weight into its own dynamics. Accordingly, we call in- and out-coupling strength, respectively. Grouping together oscillators with the same coupling strengths, the number of different pairs coincides with the number of subpopulations.
We consider two equally-sized subpopulations denoted by “” and “.” Hence, the oscillators are distinguished by a pair of coupling strengths, or , and all of those can be positive or negative. We choose the parametrization
and are average in- and out-coupling strengths, while and give corresponding mismatches. If or , then half of the couplings are positive (attractive) and half negative (repulsive). In such cases we speak of mixed interactions. Note that (2) leads to point symmetries, because changing or yields the same situations. In the following, all oscillators have the same constant natural frequency . Therefore, by virtue of the rotational symmetry, we can set without loss of generality. Time-dependent disorder is included as Gaussian white noise, The angular brackets denote averages over different realizations of the noise and the single non-negative parameter denotes the noise intensity. The noise terms can be regarded as an aggregation of various stochastic processes Anishchenko et al. (2007). In Refs. Hong and Strogatz (2012, 2011a); ? it is found for the deterministic case that mixed out-couplings alone do not enable more than partial synchronization, whereas mixed in-couplings yield traveling waves reached through diametrically synchronized states. We intermingle both types of mixing and explore in particular whether traveling waves persist in the presence of noise . Without loss of generality all subsequent results are obtained with , but for illustration we keep in the derivations.
“Discordant synchronization” is used here as an umbrella term for situations where the ensemble splits into two partially synchronized clusters. This will include traveling waves and -states, the latter being the extreme form of discordance with two oscillator populations anti-aligned to each other.
We investigate the thermodynamic limit , where propagation of molecular chaos Kac (1956) allows us to describe each population by a one-oscillator probability density . Normalization requires . For given coupling strengths and , denotes the fraction of oscillators with phase between and at time . The densities are governed by the nonlinear Fokker-Planck equations Sakaguchi (1988); ?:
The global mean-field amplitude and phase follow from a superposition:
Note that subpopulations of different sizes can be treated simply by rescaling . The local mean-field variables obey
The level of synchrony in the two subpopulations is measured separately by , whereas for the global measure we take the classical Kuramoto order parameter which differs from Eq. (4) by taking out . The order parameter values lie between zero (incoherence) and (complete synchronization). The variables and stand for the corresponding mean phases. Of special interest is the phase lag , i.e. the difference in the mean phases of the two populations, .
In the sequel, variables without a dot or other indicated time dependence refer to the long-time limit. Let us outline the four qualitatively different self-organized states observed here after some transient dynamics. (i) In the incoherent state the whole population of oscillators rotates asynchronously, . (ii) The classical partially synchronized state has zero phase lag, . We use “zero-lag sync” as a shortcut to denote this state. (iii) The -state describes a partially synchronized state, where the two subpopulations are anti-aligned to each other, . (iv) In the traveling wave state the whole population is also partially synchronized, but oscillates with a frequency different from the frequency of single oscillators. This spontaneous change in rhythm is induced by a phase lag that is neither zero nor , . We calculate the wave speed as
where represents a long-time average Hong and Strogatz (2011a); ?.
In Fig. 1 example snapshots from simulations of oscillators are shown ( are equidistantly chosen for visualization). Figure 1(a) displays a -state with , , , and . Figure 1(b) shows a traveling wave state with , , , and . It is equally possible that the wave runs in the other direction, depending on initial phases and realization of the noise. Note that perfect synchrony, , cannot be achieved with finite coupling strengths, if an infinitesimal amount of noise is present.
In order to analytically investigate the collective dynamics that are governed by (3)-(5), we approximate the phase distributions in the two populations by time-dependent Gaussians. This well-known method is motivated by numerical observations Kurrer and Schulten (1995); ?; Sonnenschein and Schimansky-Geier (2013). It has been generalized recently to encompass heterogeneities in couplings Sonnenschein and Schimansky-Geier (2013); Sonnenschein et al. (2013); ?. Extending those derivations to the present case, we obtain the following three-dimensional system of ODE’s:
All the four aforementioned collective states are fixed points of (7) with . Two types of fixed point solutions have to be distinguished, because there are two possibilities that holds:
Equation (8) describes zero-lag and -states, whereas Eq. (9) underlies traveling waves. Intermediate phase lags cause spontaneous drifts, because according to Eqs. (7) and (9) the common frequency of the traveling waves obeys
Two more equations are obtained from imposing in (7):
With (8), (9) and (11) we have three coupled equations for three unknowns, and . No stationary solution with can be found, but the singularities and turn out to have a special meaning. In particular, if one of the in-coupling strengths vanishes, the corresponding population remains incoherent. Numerical continuation around this point shows that in order to avoid a negative local order parameter, , which is unphysical, the whole population transfers to a -state. In our parametrization this first critical condition can be written as
One can show that in general the incoherent state, , loses linear stability, if the noise intensity falls below a certain value (cf. the Appendix). This happens at
Finally, as outlined below, the intersection given by and (13) coincides with the origin of bistability. Note that the aforementioned conditions are exact. Figures 2(a) and 2(b) depict bifurcation diagrams in the planes spanned by the coupling mismatches and the average coupling strengths , respectively. Solid lines are given by the critical conditions (12) and (13). Dashed lines are obtained on the basis of the reduced system (7) with the help of MATCONT Dhooge et al. (2003). We detect branch and limit points, since at all lines one eigenvalue vanishes, except at the ones given by (12), because those do not correspond to real bifurcations, but delineate two analogous partially synchronous states: zero-lag and -states. We additionally test all these findings by numerically calculating the eigenvalues of the Jacobian of (7) with the fixed points given by (8)-(11). When the two lines given by (12) and (13) intersect, the boundaries (12) cease to exist. We emphasize two distinct routes to TW states; TW1 is surrounded by -states, TW2 by classical zero-lag sync states. Delimiting lines approach each other, see insets for enlarged areas. We further find bistability between incoherence and zero-lag or -states, see panels (a) and (b) in Fig. 2, respectively. The bistable areas are circumscribed by two lines that intersect at the points given by and (13). The location of this intersection determines the type of bistability. Interestingly, traveling and non-traveling wave states can coexist in small parameter regions. We show this in Fig. 2(c). In particular it is observed that traveling waves can coexist with complete incoherence, as well as with -states and zero-lag partially synchronous states.
Figure 2 suggests certain conditions for the various collective states. In order to observe -states, mixed attractive-repulsive in-couplings have to be included. Traveling waves surrounded by -states are possible, if one includes mixed in-couplings without mixed out-couplings. In contrast, traveling waves surrounded by zero-lag sync states exist, if there are mixed out-couplings and a non-zero mismatch without mixing in the in-couplings. Bistability between incoherence and zero-lag sync requires lack of mixed in-couplings, but the presence of mixed out-couplings. Finally, bistability between incoherence and -states is possible by combining mixed in-couplings with vanishing mixing in the out-couplings. These conditions appear to supplement consistently the observations in Refs. Hong and Strogatz (2012, 2011a); ?. In particular, we verified numerically that the traveling waves surrounded by zero-lag sync can also be observed in the setting studied in Ref. Hong and Strogatz (2012), if a small mismatch in the in-couplings is present, as suggested by our bifurcation diagram in Fig. 2(a).
It is worth asking how crucial asymmetric interactions are for the discordant synchronization patterns discussed here. As it is easily seen from Eq. (2), the interactions are symmetric if in- and out-coupling strengths balance each other such that the equation holds. This condition can be projected onto straight cuts through the parameter space. The demarcations would not cross the traveling wave areas [see Fig. 2(a) and 2(b)]; instead they would divide the bifurcation diagrams into parameter regions that contain both routes to traveling waves, i.e. TW1 and TW2. In other words, asymmetric interactions are needed to get traveling waves. Interestingly however, those straight cuts would go through the -state regimes. Thus, asymmetry in the interactions is not a necessary ingredient to observe -states. This conclusion is not evident from the works presented in Hong and Strogatz (2012, 2011a); ?.
Additional insights can be gained by numerically solving the three coupled equations (8)-(11) to get and . The only subtlety is that one has to factor in the bifurcation values previously obtained in order to correctly choose between (8) and (9). Solutions are shown in Fig. 3, and compared with the results from numerical simulations. For the latter, initial phases are randomly chosen from the uniform distribution . For each value a long-time average is taken over with integration time step . Upper panels depict the order parameters, while in the lower panels the corresponding phase lags are shown. The colored regions match those in Fig. 2 and discriminate the different collective states. In Fig. 3(b) one can see that at the transition from incoherence to -state the suborder parameters abruptly jump from zero to high values in a hysteretic manner. For very long time averages it is expected that the hysteresis is washed out due to noise-induced jumps between the two stable steady states. We do not report this here. In Fig. 3(c), around , the abrupt change in the phase lag and the non-vanishing order parameters signal extended stability of traveling waves, as discussed for Fig. 2(c). Such bistable dynamics appears to be a promising topic for future studies. In general, the abrupt changes and the local minima in the order parameters as a function of the average in-coupling strength , as presented in Fig. 3, are of vital interest on their own, see e.g. Refs. Pazó (2005); ? and Omel’chenko and Wolfrum (2012), respectively.
In Fig. 4 we compare the common frequency obtained from the reduced system, Eq. (10), with numerical simulations. It serves as an alternative measure to the wave speed (6), which is calculated from the individual instantaneous frequencies that do not exist in the analytical treatment. Therefore for the wave speed no comparison with theory is being made. In Fig. 4 zero-lag synchronous and -states become unstable in the regime of traveling waves. One can observe that both measures, the common frequency and the wave speed, highlight equally well the onset of traveling waves. As mentioned before, the waves emerge in frequency pairs, meaning that they can travel in both directions, depending on realization of random numbers.
In Fig. 5 results of numerical simulations are superimposed on smaller theoretical bifurcation diagrams. In Figs. 5(a)-5(c) local and global order parameters are depicted, while in Fig. 5(d) the wave speed is plotted. For each of the data points in the grid, the equations of motion (1) are integrated with oscillators and observables are then averaged over time, . Different initial conditions are chosen in Figs. 5(e)-5(g) in order to find bistability in numerical simulations ( data points there). Specifically, in Figs. 5(e)-5(g) the area circumscribed by dashed and solid lines shows -states: the local order parameters attain large values, but due to the anti-phase alignment given by the phase lag , the total order is small. In Figs. 5(a)-5(c) this area is filled with incoherence. Hence, numerical simulations agree again very well with theoretical results. We remark that bistability between incoherence and zero-lag synchronous states, as predicted by the theory, can be analogously found by varying the initial conditions (not shown).
We explored the rich dynamics that emerge from asymmetric in- and out-coupling strengths among two mutually globally coupled oscillator populations. As an illustrative example we considered identical noisy Kuramoto phase oscillators with non-uniform and mixed attractive-repulsive interactions. We observed that the two populations can partially synchronize in-phase and with a constant phase lag to each other. We referred to the latter as “discordant synchronization.” The phase lags induced spontaneous drifts. As a result, traveling waves were formed in which the whole population oscillated with a different frequency than the individual units. However, in the state of maximal discordance, where the two partially synchronized populations are anti-aligned to each other, the spontaneous drifts disappeared. Correspondingly, we revealed two distinct routes to traveling waves, one through diametral two-cluster states, the other one through classical one-cluster states. Since the latter are ubiquitously investigated in the literature, we expect the second route to traveling waves to be more prevalent in general. Appropriate experimental setups made up of two constituents are conceivable, realized e.g. with laser systems Wünsche et al. (2005); ?, metronomes Martens et al. (2013) or chemical Belousov-Zhabotinsky oscillators Tinsley et al. (2012).
With the help of a Gaussian approximation in the infinite system size limit we derived a three-dimensional system of coupled ODE’s. This reduced system allowed a thorough bifurcation analysis and further analytical treatment, in excellent agreement with numerical simulations of a large but finite number of oscillators. We found that physically relevant singularities constitute a significant part of the bifurcation scenario. We further showed which collective states can coexist. Our results help to understand the emergence of discrepancies between individual and collective rhythms, as is observed e.g. in neuronal networks Wilson and Cowan (1972); ?. Specifically, if the connection strengths were capable to slowly vary in time, one could expect temporal patterns reminiscent of the high-frequency rhythmic events observed in hippocampal networks Gloveli et al. (2005); ?. Attractive (positive) and repulsive (negative) couplings are often associated with excitatory and inhibitory connections among neurons. This is reasonable, since positive couplings tend to increase synchrony, which is also the case for excitatory connections in the brain. In contrast, negative and inhibitory connections have in common that they tend to decrease synchrony Cumin and Unsworth (2007). Exceptions however exist, see e.g. Ref. Van Vreeswijk et al. (1994). The individual out-coupling strengths considered here are particularly suitable to emulate the role of excitatory or inhibitory neurons. This is pointed out in Ref. Hong and Strogatz (2012) by referring to Dale’s principle, according to which a neuron releases the same set of neurotransmitters at all its synapses Eccles (1976). Such a comparison would become even more applicable by including an excitation threshold into the system Kurrer and Schulten (1995); ?; Sonnenschein et al. (2013); ?; Hong (2014). Whether the combination of mixed attractive-repulsive interactions on the level of in- and out-coupling strengths is experimentally relevant, remains an interesting topic for the future.
Furthermore, future work should approach real networks of networks by considering multiple populations Kivelä et al. (2014); ?, see the Appendix for first steps. It would also be interesting to extend the present framework towards inertia Ji et al. (2014); ? and imposed phase shifts Iatsenko et al. (2014); ?; Maistrenko et al. (2014); ?. Finally, in small oscillator populations additional peculiarities can be expected Maistrenko et al. (2014); ?.
Acknowledgements.Thanks to S. Milster for comments on the manuscript. B.S. and L.SG. acknowledge funding from the Bundesministerium für Bildung und Forschung (BMBF) (BCCN II A3, grant 31401211). T.K.DM.P. acknowledges FAPESP (grant 2012/22160-7) and IRTG 1740. F.A.R. acknowledges CNPq (grant 305940/2010-4), FAPESP (grant 2013/26416-9) and IRTG 1740. J.K. acknowledges IRTG 1740 (DFG and FAPESP). *
Appendix A Arbitrary number of populations
Here we discuss the transition from incoherence to partial synchrony in an arbitrary number of interacting populations of arbitrary sizes. To this end, oscillators with the same pair of in- and out-coupling strengths are again grouped into one population. Then the nonlinear Fokker-Planck equation for the one-oscillator probability density reads
The global mean-field amplitude and phase are given by
The averages take into account all in- and out-coupling strengths, and , via the corresponding joint probability distribution . In Eq. (15) the averages are further taken over the local mean-field variables, By a linear stability analysis of the incoherent state, , one can exactly derive the critical noise intensity for the synchronization transition:
This result follows from adapting the derivation Sakaguchi (1988); ? to the present case. Above the whole ensemble is incoherent, below the incoherent state loses stability and partial synchrony is observed. Now we give an alternative derivation on the basis of the Gaussian approximation. It recovers the exact condition (16) and yields an instructive generalization of the equations in the main part of the manuscript. Inserting the Gaussian approximation for all phase distributions into Eqs. (14)-(15), we obtain a set of coupled differential equations for the local mean-field amplitudes and mean phases :
Let us consider small perturbations of the incoherent state, . The perturbations may give rise to zero-lag synchronous or -states, . Accordingly, we can separate the network of networks into two groups, and , which contain subpopulations with coupling strength pairs that lead to the same mean phases . At the same time, the mean phases of the subpopulations in the groups and differ by . Linearizing around the perturbations in the two groups separately, we obtain from Eqs. (17):
Here, the two separated networks of networks are labeled by the indices and . According to Eq. 15 the two perturbations can be put together as , where is some constant coming from an arbitrary global mean phase. As a result we obtain
which leads to the critical condition (16).
- S. Strogatz, Sync: The Emerging Science of Spontaneous Order (Hyperion, New York, 2003).
- A. Pikovsky, M. Rosenblum, and J. Kurths, Synchronization: A Universal Concept in Nonlinear Sciences (Cambridge University Press, Cambridge, 2003).
- H. Daido, Phys. Rev. Lett. 68, 1073 (1992).
- H. Hong and S. H. Strogatz, Phys. Rev. E 85, 056210 (2012).
- H. Hong and S. H. Strogatz, Phys. Rev. Lett. 106, 054102 (2011a).
- H. Hong and S. H. Strogatz, Phys. Rev. E 84, 046202 (2011b).
- H. R. Wilson and J. D. Cowan, Biophys. J. 12, 1 (1972).
- N. Brunel, J. Comput. Neurosci. 8, 183 (2000).
- M. S. Miguel, V. M. Eguíluz, R. Toral, and K. Klemm, Comput. Sci. Eng. 7, 67 (2005).
- Y. Kuramoto, Chemical Oscillations, Waves, and Turbulence (Springer-Verlag, Berlin, 1984).
- V. S. Anishchenko, V. Astakhov, A. Neiman, T. Vadivasova, and L. Schimansky-Geier, Nonlinear Dynamics of Chaotic and Stochastic Systems (Springer-Verlag, Berlin, 2007).
- M. Kac, in Proc. Third Berkeley Symposium on Mathematics, Statistics, and Probability, Vol. 3 (University of California Press, 1956) p. 171.
- H. Sakaguchi, Prog. Theor. Phys. 79, 39 (1988).
- S. H. Strogatz and R. E. Mirollo, J. Stat. Phys. 63, 613 (1991).
- C. Kurrer and K. Schulten, Phys. Rev. E 51, 6213 (1995).
- M. A. Zaks, A. B. Neiman, S. Feistel, and L. Schimansky-Geier, Phys. Rev. E 68, 066206 (2003).
- B. Sonnenschein and L. Schimansky-Geier, Phys. Rev. E 88, 052111 (2013).
- B. Sonnenschein, M. A. Zaks, A. B. Neiman, and L. Schimansky-Geier, Eur. Phys. J. Special Topics 222, 2517 (2013).
- B. Sonnenschein, T. K. D. Peron, F. A. Rodrigues, J. Kurths, and L. Schimansky-Geier, Eur. Phys. J. B 87, 182 (2014).
- A. Dhooge, W. Govaerts, and Y. A. Kuznetsov, ACM Trans. Math. Softw. 29, 141 (2003).
- D. Pazó, Phys. Rev. E 72, 046211 (2005).
- J. Gómez-Gardeñes, S. Gómez, A. Arenas, and Y. Moreno, Phys. Rev. Lett. 106, 128701 (2011).
- O. E. Omel’chenko and M. Wolfrum, Phys. Rev. Lett. 109, 164101 (2012).
- H.-J. Wünsche, S. Bauer, J. Kreissl, O. Ushakov, N. Korneyev, F. Henneberger, E. Wille, H. Erzgräber, M. Peil, W. Elsäßer, and I. Fischer, Phys. Rev. Lett. 94, 163901 (2005).
- J. Zamora-Munt, C. Masoller, J. Garcia-Ojalvo, and R. Roy, Phys. Rev. Lett. 105, 264101 (2010).
- E. A. Martens, S. Thutupalli, A. Fourrière, and O. Hallatschek, Proc. Natl. Acad. Sci. USA 110, 10563 (2013).
- M. R. Tinsley, S. Nkomo, and K. Showalter, Nature Phys. 8, 662 (2012).
- T. Gloveli, T. Dugladze, S. Saha, H. Monyer, U. Heinemann, R. D. Traub, M. A. Whittington, and E. H. Buhl, J. Physiol. 562, 131 (2005).
- K. T. Hofer, Á. Kandrács, I. Ulbert, I. Pál, C. Szabó, L. Héja, and L. Wittner, Hippocampus 25, 169 (2015).
- D. Cumin and C. Unsworth, Physica D 226, 181 (2007).
- C. Van Vreeswijk, L. Abbott, and G. B. Ermentrout, J. Comput. Neurosci. 1, 313 (1994).
- J. Eccles, Notes Rec. R. Soc. Lond. 30, 219 (1976).
- H. Hong, Phys. Rev. E 89, 062924 (2014).
- M. Kivelä, A. Arenas, M. Barthelemy, J. P. Gleeson, Y. Moreno, and M. A. Porter, J. Complex Netw. 2, 203 (2014).
- S. Boccaletti, G. Bianconi, R. Criado, C. Del Genio, J. Gómez-Gardeñes, M. Romance, I. Sendina-Nadal, Z. Wang, and M. Zanin, Phys. Rep. 544, 1 (2014).
- P. Ji, T. K. D. Peron, F. A. Rodrigues, and J. Kurths, Phys. Rev. E 90, 062810 (2014).
- M. Komarov, S. Gupta, and A. Pikovsky, Eur. Phys. Lett. 106, 40003 (2014).
- D. Iatsenko, P. V. E. McClintock, and A. Stefanovska, Nat. Commun. 5, 4118 (2014).
- V. Vlasov, E. E. N. Macau, and A. Pikovsky, Chaos 24, 023120 (2014).
- Y. Maistrenko, B. Penkovsky, and M. Rosenblum, Phys. Rev. E 89, 060901 (2014).
- O. Burylko, Y. Kazanovich, and R. Borisyuk, Phys. Rev. E 90, 022911 (2014).