Hierarchical dynamics in large assemblies of interacting oscillators

Hierarchical dynamics in large assemblies of interacting oscillators

PHYSICS LETTERS A Physics Letters A 160(1991) 227—232 North-Holland Hierarchical dynamics in large assemblies of interacting oscillators Erik D. Lum...

565KB Sizes 63 Downloads 37 Views

PHYSICS LETTERS A

Physics Letters A 160(1991) 227—232 North-Holland

Hierarchical dynamics in large assemblies of interacting oscillators Erik D. Lumer and Bernardo A. Huberman Department ofApplied Physics, Stanford University, Stanford, CA 94305, USA and Xerox Palo AltoResearch Center, 3333 Coyote HillRoad, Palo Alto, CA 94304, USA Received 9 August 1991; accepted for publication 16 September 1991 Communicated A.R. Bishop

We study a collection of phase-coupled oscillators possessing a hierarchical coupling structure. We establish a necessary condition for the existence ofa phase transition to collective synchrony for finite values ofthe coupling strength in terms of an inequality involving the connectivity between clusters of oscillators, the rate at which coupling strengths decrease with ultrametric distance, and the dispersion of intrinsic frequencies. When the inequality is not satisfied, there is a cascade of discrete transitions to intracluster synchrony as the coupling strength is increased, but no global synchronization is possible in the infinite size limit.

Progress in nonlinear dynamics is leading to an increasing understanding of systems with many degrees of freedom [11. These range from the dynamics of laser arrays [21, and lattices of logistic maps 3], to turbulent fluids [4]. Of particular interest are phenomena of collective synchronization taking place in large assemblies of coupled oscillators [5— 11]. In the biological realm, such regimes characterize situations such as the flickering of swarms of fireflies, and the peristaltic motion of gastrointestinal tracts [5]. Furthermore, recent experimental results in the area of neuroscience lead to the consideration of synchronization in large assemblies of coupled neural oscillators as a possible solution to the so-called binding problem [12—15]. Since most of the models of collective synchronization studied so far considered oscillators interacting via uniform fields, it is of interest to extend them to the case where the coupling strengths vary with distance and are distributed over a hierarchy of values, a situation that seems to be more realistic for heterogeneous systems such as those encountered in the study ofcomplex physical [16,17] and biological structures. This is the purpose of this paper, in which we elucidate the conditions under which collective rythmic activity can exist in a system with a hierarchy of coupling constants. Specifically, we introduce a model of phase-cou-

pled oscillators organized into clusters of coupling strengths, and such that intracluster interactions are weaker than intercluster ones. An ultrametric measure is then introduced that quantifies the interaction distance between any two oscillators in the clustured structure. A similar metric has been used in the past to study the transport properties of hierarchical systems [18]. By varying the shape of the embedding tree of interactions and the functional dependence between coupling strength and ultrametric distance, we are able to consider different classes of interactions in a systematic way. We also show the existence of a phase transition to collectivesynchrony that depends upon the validity of an inequality involving the connectivity between clusters, the rate at which coupling strengths decrease with ultrametric distance, and the dispersion of intrinsic frequencies. By collective synchrony we mean a situation whereby a macroscopic number of oscillators evolve in time with the same instantaneous frequency. These results are supported by extensive numerical simulations. Consider a collection of oscillators with coupling strengths that are weak compared to the forces attracting the oscillators to their limit cycle. Once relaxed to its limit cycle, an oscillator can be described by a single degree offreedom, that is its phase on the cycle. In the limit of weak interactions, a simple

0375-9601/91/S 03,50 © 1991 Elsevier Science Publishers B.V. All rights reserved.

227

Volume 160, number 3

PHYSICS LETTERS A

18 November 1991

model of the dynamics for the phases consists in the following equation [6,8,9,11], dO

ultrametric distance

(1)

where O~represents the phase of the jth oscillator in a population of N elements and ~ its intrinsic frequency. The second term on the right hand side of eq. (1) describes the interactions ofthejth oscillator with some or all ofthe N— 1 others. Accordingly, the coupling between oscillators i and j is proportional to a coupling strength, K,,, and to a periodic function of their phase difference, h (z~9).A convenient assumption which is often made is that of symmetric interactions, which translates in the symmetry of K~ and the antisymmetry of the coupling function for the permutation of the indices i and j. Furthermore, we will allow the intrinsic frequency of the oscillators, w~,to vary from one element to another and we therefore characterize the distribution of its values by a density function f( a). In the above model, the matrix of coupling strength K: {K,~}defines the pattern of interactions among oscillators. Two special and rather extreme connectivity patterns have been thoroughly studied in the past. The first type, a fully interconnected population of infinite size, corresponds to a matrix K whose nondiagonal components are all equal. Such population was shown by Kuramoto [6] to exhibit a phase transition to global synchronism at a critical value K~of the uniform coupling strengths. The second type, a regular lattice of oscillators with nearest neighbor coupling, was studied by Daido [9], who clarified the lattice dimensionality dependence ofthe onset of collective synchronism. In order to study the case of a hierarchy of coupling strengths, it is useful to introduce a visual representation of the structure of interactions between oscillators as done in fig. 1. The oscillators are represented as the leaves of a tree so that an ultrametric distance, l,~,between any pair (i, I) of oscillators can be defined as the level of their nearest common ancestor node. The ancestor node also defines the cluster of oscillators at the bottom of its branching subtree that are within an ultrametric distance 1~,of each other. The coupling strength between i and f is set to be a decreasing function of l~,and is expressed as 228

I~=3

L=2

L=1

spatial distance Fig. 1. Hierarchy ofcoupling strengths. The oscillators are at the leaves of the tree. The ultrametric distance between two oscillators is given by the level of their nearest common ancestor in the tree. For instance, the distance between i and j in fig. 1 is equal to 2. Their ancestor node also defines the cluster of oscillators branching out (in bold in the figure).

Ku = Kd( lii) (2) The profile of coupling strengths, d(1), defines the range ofmacroscopic behaviors exhibited as the control parameter K ofthe dynamics is driven across its critical value(s). The other relevant parameters of the model are the branching ratio b of the embedding tree, its depth L ~‘, and the specific form of the distribution of intrinsic frequencies f( a,). Concerning the latter, Daido [9] showed that the states of collective Oscillations depend crucially on the asymptotic behavior of f( w) in the limit of large deviations from the average th. Without loss of generality, one may write in this case .

f(~~)x lw—WI

a , (3) with 0
~ Notice that N=bL, where N is the size of the population.

Volume 160, number 3

PHYSICS LETTERS A

drops with the corresponding power-law. For instance, a = 1 corresponds to the Lorentzian distribution, often invoked to model the distribution of intrinsic frequencies in a population ofbiological oscillators [9,10]. The statistical properties of aggregates of frequencies indicate that distributions with a finite variance such as the normal law are also correctly modelled by the power-law asymptotic distribution when a is set to 2 *2 Before performing the renormalization-group analysis of the system (1), we must first make explicit our choices of profile functions d(1). We begin by expressing the fact that the total coupling to any oscillator j remains bounded as the population size grows, that is N

~ K~=K (N>> 1)

(4)

.

i= I

This relation, along with eq. (2), constrains the profile function d(1) to drop quickly with ultrametnc distance. It also gives a physical interpretation of the parameter K as a static measure of the total coupling strength between any oscillator and the rest of the population. The left hand side of eq. (4) can be factored into contributions from the (b— 1 )b’—’ ~ cillators at the same ultrametric distance 1 from oscillator j. Using the definition (2), the previousequation then rewrites as

18 November 1991

factor a of the total coupling strength with the oscillators at the new distance. Let us now consider the ensemble of nodes at level 1 in the coupling tree. These nodes divide the whole population in b’~~’ separate clusters so that each one contains b’ oscillators (we assume b’>> 1) and we denote such a cluster as C~.with the index k varying in the range {l, b’~~’}. For each cluster C~,we can define an aggregated phase, in fact that of an imagmary giant cluster oscillator, as Øk~

~ O~. je C,,

Notice that here and in subsequent notations, we make the index 1 implicit. The dynamics of each cluster oscillator is obtained by summing over eqs. (1) for all the oscillators in the cluster and is expressed as døk ~J7= ~ + ~ K1~h(øg 0k + öj,q öt,k), (7) —

i,j,q



where the sum extends over all the clusters C~at level lother than C~and each pair of oscillators (i,J) such that oscillator i belongs to the cluster C~and oscillatorj is in C~.The intrinsic frequency of the cluster oscillator is defined as 1

L

~ (b—l)b’~d(1)=l (L>>l) i

(5)

i

where the sum extends over all the levels in the coupling tree. A typical class of profiles which is found to satisfy the above constraint is called exponential and is defined as d(l) =

1 a—l (b_l)bt_I a’

(a>’l).

(6)

For convenience we will neglect a correcting multiplicative factor in eq. (6), which is due to the finite value of L and whose value converges to unity as L—~cc.To interpret the exponential profile, we nolice that an increment by one of the ultrametric separation corresponds to a reduction by a geometric S2

A qualitative justification of this point is given later in the paper. A mathematical justification is given in ref. [9].

and ôj,q is the deviation of the phase of oscillator i from that of its cluster oscillator, or more explicitly: dz,q

~i

Øq•

In the next step, we rely on a technique used by Daido in ref. [9]. First let the population be of infinite size (i.e. N and L are infinite) and notice the fact that in the limit of n-+~the random variable 1 ‘~

~

~i=1

possesses a stable distribution whose characteristic function is of the form xexp( I zla) provided that a proper choise of the constant Yn ~5 made [19]. We then exploit this property through the following transformations in eq. (7), —

r

tb’~—a)/a

0k = Ok —

(Yb//b’) t.

(8) 229

Volume 160, number 3

PHYSICS LETTERS A

The rescaled dynamics for the cluster oscillator now reads d~. ~ 1,k + —~b

~ K,ih(OqO~k+öiqöjk),

b

(9)

i,j.q

where &I,k comprises a sum over the b’ oscillators in C~.The specific hierarchical structure of coupling between oscillators is such that the coupling strength is uniform for any pair (1,1) of oscillators in which oscillator i belongs to C~and oscillator Ito C~.Furthermore, there are exactly (b 1) bm —1 clusters that connect with 1, C’~,m levels above 1, i.e. at level 1+ m. These facts allow us to decompose the right hand side of eq. (9) into terms of same coupling strength and rewrite it as —

dOk

~

4~.k

(bI

X

+K~ d(1+m)

)b”~

(10)

1~kq~(Oq,~ ~k),

~

where we have made use of relation (2) and an effective coupling function Ii~,,,between clusters is defined as

(

fikq,~ Oqn,



Ok)

~ h (O~q,,,



=b

Ok +

~



~5j,k).

(11)

(

&



I/a)/

a

where c is a constant independent of the level 1. Therefore, if b1 I ~ < a and in the limit of 1—+ ~, the system represented by eqs. (1) converges to the fixed point of the transformation —

=~k, ~k= lim ~b!k. (12) dt Since the random variable £~kpossesses a stable distribution with the characteristic exponent a, the existence of such fixed point indicates the absence of ~

230

global synchrony for any finite value of the coupling constant. Therefore, one can express the necessary conditions for the existence of a phase transition to collective synchronism as ~ l
ieC~,,,~

With the definition of the exponential profile function and the help of some simple arithmetic, it can be shown that the second term on the right hand side of eq. (9) has an upper bound of the form C

18 November 1991

duced by the profile factor a as compared to the coupling between clusters at level 1— 1. Kuramoto [6] has established that for a population of oscillators with infinite range and uniform coupling, the critical value of the coupling constant leading to collective synchrony is inversely proportional to the dispersion of intrinsic frequencies. Therefore, if ~/&~ ~ a, the clusters at level I should also synchronize as soon as K exceeds K~.Since the same argument applies recursively to any two consecutive levels larger than 1, one deduces that a single transition will lead to the bulk synchronization ofall the clusters larger than a certain size. Notice that the condition given by eq. (14) is exactly the one predicted by eq. (13) when a is set to 2. If the effective coupling strength drops faster than the frequency dispersion, i.e. a> \[b, successive threshold values of K that correspond to a unit increment of the maximum ultrametric distance sep-

Volume 160, number 3

PHYSICS LETTERS A

arating two synchronized oscillators are spaced by a multiplicative factor which should converge to a/..fb as one moves up in the coupling tree. Increasing K thus derives the system through an infinite cascade of abrupt transitions to synchrony of increasingly larger clusters. We have simulated a discrete version of the system (1) for a population of 1024 oscillators by approximating the first order derivative d0~/dt with [O~((n+ I )dt)—0~(ndt)]/dt. Eqs. (1) were integrated forward in time for 9000 iterations, with a temporal step dt = 0.05 and a convergence test made throughout the integration. In fig. 2, we show the largest fraction of frequency-locked cells, NS/N, as a function of the coupling constant K. The data are averaged over three successive runs with different mitial conditions. To evaluate NS/N, we computed averages of the instantaneous frequencies d0~/dt,.~,and divided the population in synchronized groups ofoscillators whose members satisfy IJ—fI
(a)

N3

18 November 1991

sition to global synchronization is observed. Specifically, the profile function is exponential with a= 1.2, the distribution of frequencies Gaussian (with variance a= 1.2) and the treebinary (b= 2) with a depth of 10 levels. The coupling function has the form h(A0)=sin(2iti~9)/2~t. In fig. 2b, the profile is exponential with a = 2 and the other parameters are as in fig. 2a. By raising the value ofa, we invalidate the conditions (13) and a cascade of transitions to synchrony, each up to a certain ultrametric distance is observed, takes place as indicated by the heights of the plateaus equal to fractionsof power of 2. Notice the different scale along the K-axis. A discrete transition is recovered in fig. 2c by doubling the branching ratio of the tree (b = 4). The profile factor is still 2 which leads to the equality ofboth sides in relation (13). Finally, in fig. 2d the profile factor is set to its value in fig. 2a (i.e. a= 1.2) but the frequency distribution is now chosen Lorentzian (i.e. corresponding to a= 1). The simulation in this case seems to indicate the absence of global synchrony for any finite value ofK, as predicted by the theory. More defmite conclusions would require simulations involving much larger populations.

±

(b) ~

0.5

1

1.5

2

2.5

3

5

0.5

1

1.5

2

2.5

3

0.5

10

1

15

1.5

20

2

2.5

25

3

Fig. 2. The largest fraction of frequency-locked cells,N,/N, versus coupling constant in a population of 1024 oscillators is shown. T= 9000, dt=0.05 (a) a=l.2, distribution is Gaussian with a=l.2, b=2. (b) a=2 and the other parameters as in (a). (c) a=2, b=4 and other parameters unchanged. (d) a = 1.2, b = 2, distribution is Lorentzian.

231

Volume 160, number 3

PHYSICS LETTERS A

In conclusion, we proposed and studied a model of coupled oscillators with a distribution of intrinsic frequencies and having a hierarchy of coupling strengths. The latter were set to decrease with the ultrametric distance separating the coupled oscillators, As shown, a necessary condition for the existence of a phase transition to collective synchronism at finite values of the coupling strengths is determined by an inequality involving the connectivity between clusters of oscillators, the rateat which coupling strengths decrease with ultrametric distance, and the dispersion of intrinsic frequencies. When the inequality is not satisfied, no global synchrony is possible in a system of infinite size. Instead, increasing the coupling strength drives the system through a. cascade of discrete transitions to synchrony within clusters of increasingly larger dimensions, which reach global synchrony only in the unphysical limit of an infinite coupling constant. These results are quite general. Indeed, the exponential profile of coupling range provides a critical shape to which others might be easily compared. Consider for example a coupling profile which decreases slower than exponentially below a certain distance d and faster beyond, possibly becoming null at finite range. Our analysis indicates that in this case collective synchrony might be reached in clusters of size d or smaller, but with cluster oscillations mutually incoherent. Furthermore, the choice of a regular coupling tree is not restrictive, for our renormalization analysis demonstrates that the interaction between clusters is only a function of their overall coupling strength, and is therefore independent of the detailed structure of interaction among their members. Last but not least, our results have a direct application to the neurological realm. Recent experiments suggest that the brain binds its fragmentary representations of visual events by synchronizing the activities of the stimulated neural oscillators in the visual vortex. If such is the case, our model might provide some basic constraints on the nature of the interactions that must take place between the oscillating cells. We are encouraged in this belief by a number of other recent theoretical studies [14] of

232

18 November 1991

this phenomenon making use of phase models. Furthermore, the distribution of the intrinsic frequencies of the neural oscillators seems to be dependent, at least within a certain range, on that of the values taken locally by the stimulating signal [12]. Since the latter can be of any form, it is important to elucidate the dynamics of the system of oscillators for a wide range of frequency distributions, as was the case with this model. This work was partially supported by ONR contract No. N000 14-82-0699.

References [11 Special issue on spatiotemporal chaos, Prog. Theor. Phys. Suppl. 99 (1989). [2] K. Otsuka, Phys. Rev. Lett. 65 (1990) 329. [3] K. Kaneko, Physica D 37 (1989) 60. [4] P.C. Hohenberg and B.I. Shraiman, Physica D 37 (1989) 109. [5] A.T. Winfree, The geometry of biological time (Springer, Berlin, 1980). [6] Y. Kuramoto, Prog. Theor. Phys. Suppl. 79 (1984) 223. [71G.B. Ermentrout and N. Kopell, Siam J. Math. Anal. 15 (1984) 215. [8] H. Sakaguchi, S. Shinomoto and Y. Kuramoto, Prog. Theor. Phys. 77 (1987) 1005. [9] H. Daido, Phys. Rev. Lett. 61(1988)231. [10] M. Shiino and M. Frankowick, Phys. Lett. A 136 (1989) 103. H. Daido, J. Stat. Phys. 60 (1990) 753. C.M. Gray, P. Konig, A.K. Engel and W. Singer, Nature 338 (1989) 334. Ch. Von der Malsburg, Internal Report 81-2, Department of Neurobiology, Max Planck Institute for Biophysical Chemistry, Gottingen (1981). [14] H. Sompolinsky, D. Golomb and D. Kleinfeld, Proc. NatI. [11] [12] [13]

Acad. Sci. USA 87 (1990) 7200. [151 E.D. Lumer and B.A. Huberman, submitted to Neural Computation (1991). [161 W.P. Keirstead, H.A. Ceccato and B.A. Huberman, J. Stat. Phys. 53 (1988) 733. [17] K. Kaneko, Phys. Rev. Lett. 63(1989) 219. [18] B.A. Huberman and M. Kergsberg, J. Phys. A 18 (1985) 1331. [19] B.V. Gnedenko and A.N. Kolmogorov, Limit distributions for sums of Independent random variables (AddisonWesley, Reading, 1954).