Hierarchy and Polysynchrony in an adaptive network

We describe a simple adaptive network of coupled chaotic maps. The network reaches a stationary state (frozen topology) for all values of the coupling parameter, although the dynamics of the maps at the nodes of the network can be non-trivial. The structure of the network shows interesting hierarchical properties and in certain parameter regions the dynamics is polysynchronous: nodes can be divided in differently synchronized classes but contrary to cluster synchronization, nodes in the same class need not be connected to each other. These complicated synchrony patterns have been conjectured to play roles in systems biology and circuits. The adaptive system we study describes ways whereby this behaviour can evolve from undifferentiated nodes.

of the paper is as follows. In section II the adaptive network introduced in [23] is defined. The dynamics of the network topology has a stochastic element driven by a homophilic principle, so nodes in similar states 'like' to be connected together. At each time step the network topology can change according to a set of rules, thus changing the inputs to the dynamics at the nodes. We refer to this process as rewiring. In section III we explore the dynamics of the network numerically. We show that the network reaches a frozen state where the rewiring stops. The transient times to the frozen state are evaluated as a function of the network size and the coupling strength. The different dynamical regimes are also described in this section in terms of the synchronization of the nodes. We provide several examples of polysynchronous networks and study the probability of finding polysynchrony as a function of the coupling strength. The numerical observations show that the synchronization effect is very strong, and the dynamics at different nodes can become indistinguishable at machine accuracy. This effect, which we believe is interesting in that it reflects what any finite measurement could discern, means that some of the final topological states observed are extremely unlikely from a mathematical point of view. In section IV we show formally that a closely related network rule that eliminates these mathematically unlikely states must lead to a stationary topology. In section V we summarize and discuss the main results of this work and their potential applications. The detailed stability analysis of the fully synchronous and the polysynchronous states is given in two appendices.
Many accounts of adaptive networks concentrate on the increased complexity of the evolving network topology (to scale-free networks for example). In contrast, the systems described in this paper evolve towards a stationary network topology with some striking features such as a strong hierarchical structure and polysynchronous dynamics at the nodes. Our models therefore point the way to rather different application areas: the evolution to networks with relatively simple structure having dynamics correlated in different nodes that are not directly connected by the network lends itself to interpretations in terms of functional differentiation of initially equivalent units, where the differentiated systems are distributed across the network rather than clustered. This and other possible applications in biological and social systems is commented on further in the final section of this paper.

II. THE MODEL
The model consists of a directed network of N nodes where the dynamics of the ith node (i = 1, . . . , N ) are given by We choose f to be the fully-chaotic logistic map f (x) = 4x(1 − x) and A n is the adjacency matrix of the network at time step n, so A ij n = k if there are k directed edges from j to i. In the figures we represent the directed edges by an edge with an arrow indicating the direction of the flow of information. Thus the head of the arrow is the node that receives the input and the tail of the arrow is the node that influences the node at the head, i.e. if A ij n = 0 there will be a directed edge (an arrow) from node j to node i. Each node is assigned the same fixed number m of incoming links so N j=1 A ij n = m, and we choose m = N − 1 throughout this paper. The input degree of the nodes is therefore fixed. This is particularly important for the interpretation of the examples we show throughout the paper where we have avoided labelling the weights of the connections; they always sum to m. Moreover, we will not allow a link from a node to itself so A ii n = 0 for all n. At each iteration the ith node is influenced by the dynamics of those nodes to which it is connected by an incoming arrow. We will call these nodes the neighbours of node i. Due to the condition imposed by (2), a node can have at most m neighbours.
As indicated in the introduction, the network topology changes according to a homophilic principle. At each time step the node dynamics evolves according to (1). The values of the map f at each node is compared to the values of f at its neighbours and then a 'bad' set of neighbours is identified. These are those with f values far from that at the node they influence. The connections to the bad node are then changed at random to nodes that are not bad, then the process repeats. More precisely, the nodes rewire their links through the following mechanism. At each iteration n we compute the distance matrix D ij and calculate from it the mean distance of a node to all its neighbours where a i n is the unweighted number of neighbours of node i at time step n, i.e. the sum over j of sign(A ij n−1 ). We have choosen the rewiring to be homophily-driven, so nodes prefer to be connected to nodes being in a similar state. Therefore, we identify the bad neighbours B i n of each node i at iteration n Thus a neighbour j is considered bad if its distance D ij n to the node is larger than the average distance of the neighbourhood D i n . The good neighbours of node i are then given by Once the good and bad neighbours have been identified node i will break the links coming from B i n and randomly rewire them to nodes in G i n . Let b i n be the number of bad connections, i.e. the sum of the connections to i from bad neighbours: Now choose b i n elements of G i n at random and suppose that r ik n is the number of times node k is chosen. The adjacency matrix at the next time step is It is worth noting that G i n contains two sets of nodes: those that were neighbours of i at time n − 1 and which were not bad according to the criterion (5) at time n − 1, and those that were not neighbours of i at time n − 1. This means that at each time step with B i n non-empty, connections from outside the set of previous neighbours becomes possible, and also that there is no memory of whether a node has been bad in the past.
In all the cases described here the initial connectivity is the symmetric all-to-all connectivity where each node in the network is connected to all the possible m = N − 1 neighbours and A ii 0 = 0.

A. Asymptotic network topology
The first main observation is that, contrary to other models of adaptive networks of chaotic maps [5,6], in this model the network reaches a frozen state where the rewiring stops for all values of ε ∈ [0, 1]. The existence of the frozen state is partly explained by the rewiring mechanism chosen (see section IV for further explanation and mathematical proof). If, for instance, a node i receives all its incoming links from one single neighbour k at some iteration n ′ , then D i n = D ik n and B i n = ∅ for all n > n ′ . Therefore i will remain locked to this neighbour and there will be no further change to this part of the network topology. The duration of the transient to the frozen state appears to increase exponentially with the system size N (Fig. 1) and depends on the value of ε (Fig. 2). The dependence of the transient length on the coupling constant is a sign of the influence of the dynamics in the rewiring and freezing processes. The exponential increase of the transient time with the system size is similar to that described in [32] for the case of a coupled map lattice with diffusive coupling although the definition of the transient is different. In the lattice case the topology is fixed and the transient is defined as the time it takes to reach a certain attractor.
In both Fig. 1 and Fig. 2 there seems to be a marked difference between parameters ε in the interval [0.25, 0.75] and parameters outside this interval. The transient times appear significantly shorter for parameters inside this central interval, and as we shall see (although this is, of course, not an explanation) the dynamics of the nodes for the stationary network is different in these two cases too. In Fig. 3 we show six examples of final topologies of a network of N = 10 nodes for different values of ε. The most clear feature of these network examples is the strong hierarchical structure. This model does not allow a tree structure as a final topology since all nodes have input links by definition and therefore the network will have at least one cycle. However, the structure is very close to the hierarchy of a tree structure if we consider strongly connected components of the network as roots. (We say a set of nodes is strongly connected if there is a path in the graph following the directed edges or arrows between any two nodes.) Inspired by the definitions of 'trophic level' and 'trophic height' introduced in [33] for the study of food webs, we can define the 'level' of a node as the minimum (directed) path length from the root to the node and the 'height' of a node as the average distance over all possible directed paths from the root to the node. We say a network is strongly hierarchical if level and height coincide for all the nodes in the network. We can see that following this definition all the topologies shown in Fig. 3 are strongly hierarchical.
The observation of these topologies also allows us to deduce some dynamical properties of the network. In examples a. and b. (ε = 0.1) all nodes are locked to one single neighbour. In the remaining the examples there are nodes with inputs coming from two different neighbours. As we shall see in the following sections, this is due to synchronization phenomena in the strongly connected components. When a node i has only two neighbours j, k and these are synchronized (x j n = x k n , for all n), then D ij n = D ik n = D i n for all n and the node remains locked to its neighbourhood. Mathematically this is highly unlikely since usually the orbits synchronize only eventually and are therefore never exactly the same. However, we find such states due to machine-precision effects in the numerical computations. In section IV, where we prove the freezing of the network, a extra rule for the rewiring mechanism is added to avoid such situations.
In Fig. 4 we show the probability that a node belongs to a strongly connected component of a certain size as a function of ε. As could be already appreciated in Fig. 3, the most common strongly connected components are pairs and triplets except in the region of large ε (ε > 0.9) where bigger strongly connected components are possible. This and the rest of the variations of the probabilities with ε can be better understood by studying the synchronization dynamics. Therefore, we will come back to this figure in the next section.

B. Dynamics
We study now the dynamics of the nodes in the final network as a function of the coupling constant. Our main interest is to see if any nodes in the network synchronize and how they do it.
To measure the synchronization of nodes we find it useful to define the matrix where θ(x) is the Heaviside step function, n ′ is a long transient that we allow in order to be sure the network has frozen and the dynamics have stabilized and δ is a small quantity that we introduce to properly detect eventually synchronous dynamics (numerically, machineprecision effects do the work). The element β ij is equal to zero if the trajectories of nodes i and j are fully synchronized (x i n = x j n ) during τ iterations after the transient and is equal to one otherwise.
We can now define a measure of full synchronization of the network as This measures the percentage of synchronized pairs of nodes (connected or not) over the total number of pairs. If α = 1 all nodes in the network are synchronized in the same trajectory while if α = 0 no two nodes in the network are synchronized. Since our network can split in several disconnected components and each connected component could be fully synchronized in a different trajectory, we introduce a second quantity to take this into account and measure the synchronization only between pairs of connected nodes. We can thus define the connected component synchronization as where C is the set of pairs of connected nodes and |C| is the cardinality of this set.
We can see the values of α and α C as a function of ε in Fig. 5. To explain the different regimes in this figure it is very useful to study first the dynamics of the most common small strongly connected components such as the completely connected pair, the triplet with transposition symmetry and the 3-cycle shown in Fig. 6. Since these act as roots from which the rest of the network takes their inputs, the dynamics of these components is what determines the behaviour of the rest of the nodes. In appendix A we detail the calculations. Here we will only report the results that are of interest for the discussion. The Lyapunov exponent of the fully-chaotic logistic map (r = 4) is λ = ln 2. Substituting this in (A3) we find that the synchronous chaotic state of the completely connected pair is stable in the interval 0.25 < ε < 0.75. Similarly, the triplet with transposition symmetry (A6) has a stable synchronous state if 0.5 < ε < 0.75. On the other hand, for the fully chaotic logistic map (r = 4) used here the 3-cycle (A8) has no stable synchronous state. Another important fact is that a node locked to a synchronized set of nodes (all following an orbit of the uncoupled logistic map), as in (A12), will synchronize to them if ε > 0.5. This results makes the interpretation of Fig. 5 much more straightforward. The change of regime at ε = 0.25 is explained by the strongly connected pairs becoming synchronized. Also, in Fig. 4 we can see that the probability of finding pairs in the final network greatly increases. At ε = 0.5 the synchronized state of the triplet with transposition symmetry becomes stable. This is likely to be the cause of the increase, at ε = 0.5, of the probability of being in a strongly connected component of size n = 3. Moreover, a node locked to a synchronized pair or triplet will become synchronized with it and due to the hierarchical structure of the networks, this opens the possibility for the whole network to synchronize in the same orbit. When ε = 0.75 the pair and the triplet synchronized states both lose stability. However, Fig. 5 shows that in both 0.12 ε 0.2 and ε > 0.75 there is a considerable amount of synchronized nodes in the final networks even though none of the most common strongly connected components has a stable synchronous state. This is partly caused by the phenomena of polysynchrony that we explain in detail in the next section.

C. Polysynchrony
In most studies of synchronization on networks, if two or more nodes synchronize then they are connected directly in the network, and the synchronized states form clusters. The term polysynchrony [24,25] was introduced to describe a more general form of synchronization on networks for which the synchronized states are not necessarily directly connected within the network. Examples, and further analysis of general conditions for the existence of such states can be found in [24][25][26][27][28].
We illustrate now the phenomenon of polysynchrony in our model with several examples from the simplest case of fixed point dynamics to more involved examples of quasiperiodic and chaotic polysynchronous dynamics. In Fig. 7, for ε = 0.85, we find that each synchrony class has a fixed point as the final attractor. These fixed points correspond to the fixed point dynamics of the completely connected pair since the root of the network in this example is composed of three completely connected pairs. In fact, for most of the examples of polysynchrony provided the quotient system of the network, obtained by identifying synchronized nodes, reduces to a completely connected pair [23]. Therefore, the available dynamics are those of the completely connected pair (see Fig. 8). For ε = 0.18, in Fig. 9 we find a network with polysynchronous period-2 dynamics. The dynamics is divided into two synchrony classes following the same period-2 orbit in antiphase. Since for a given ε the dynamics of all pairs is the same, the nodes of a pair are synchronous with the corresponding nodes of the other pair.
In Fig. 10 (for ε = 0.861), the network has divided into two separate clusters. In one of them the dynamics of the nodes is quasiperiodic while in the other it is periodic with period 3. Both clusters have a triplet with transposition symmetry as a root. As in the previous examples no two synchronized nodes are connected and all nodes with equivalent inputs are synchronized. Although the two clusters have different dynamics, their quotient systems are completely connected pairs and therefore, both the period-3 and the quasiperiodic orbit are attractors of the completely connected pair (see Fig. 8) when ε = 0.861.   11 shows a slightly different example of polysynchrony. As in Fig. 10, the network has split into two clusters. The roots of the clusters are completely connected pairs. The dynamics of the nodes in the pairs is periodic with period 2 as in Fig. 9. The two nodes inside the pair follow the same orbit but they are out of phase. The dynamics of the rest of the nodes in the network is periodic with period 4. Interestingly, in this example we can see how nodes with the same input (such as nodes 1 and 4) do not necessarily synchronize. Nodes filled with the same pattern are synchronous. In this case the nodes in the completely connected pair follow period-2 orbits as in Fig. 9 and the rest of the nodes follow different period-4 orbits.
This particular case is more involved because we find here an instance of spatial route to chaos in an open flow similar to that described in [34][35][36][37]. In Fig. 12 we observe the prototypical open flow system, consisting of a chain of unidirectionally coupled maps. In this case the system is closed on one side by a completely connected pair. The period-2 orbit of the pair for ε = 0.14 is fed into the chain as a fixed boundary condition and we observe a spatial period-doubling bifurcation. What we observe in Fig. 11 is merely the beginning of this route-to-chaos. to an open flow system with a completely connected pair at one end that forces the system with a period-2 orbit. We can see how the dynamics of the nodes follow a spatial route-to-chaos along the chain. In this example ε = 0.14. Fig. 13 shows an example of polysynchronous network for ε = 0.78 where the dynamics of the nodes is chaotic. It is important to note that the synchronized trajectories do not correspond to trajectories of the uncoupled logistic map. If this were the case, nodes 7 and 9 would synchronize to nodes 2 and 3 since ε > 0.5 (see appendix A 4).
Finally, in Fig. 14 we show the probability of finding polysynchrony in the network as a function of the coupling constant ε. As we have already detected in the study of synchrony in Fig. 5, there are two intervals of The regions of polysynchrony in Fig. 14 can be better understood by studying the stability of the different polysynchronous dynamics. In Appendix B we study the stability of the simplest (and more common) polysyn-chronous states in the triplet with transposition symmetry which is one of the smallest structures that can show polysynchrony. As explained in the appendix, the obtained results are independent of the value of m ′ (see Fig. 6B.) and therefore, are also valid for the case of a completely connected pair with a third node unidirectionally coupled to one of the nodes in the pair (m ′ = m or m ′ = 0). The polysynchrony in the interval 0.12 ε 0.20 is mainly explained by the stability of the period-2 polysynchronous state. This state is stable in the range 0.140375 ε 0.193814. Although this state occupies most of the interval, quasiperiodic polysynchronous states and periodic polysynchronous states of higher period can be numerically witnessed (Fig. 8). As we have also determined analytically, chaotic polysynchrony with coupled fully-chaotic logistic maps (λ = ln 2) is only possible for ε > 0.75. This fact, together with the stability of the fixed-point polysynchronous state in the range 0.806186 ε 0.86 accounts for most of the polysynchrony found for ε > 0.75. However, as before, other periodic and quasiperiodic polysynchronous states can be found in that region.

IV. ANALYSIS OF THE FROZEN STATE
As noted earlier in the discussion of Fig. 3 some numerical simulations lead to stationary (frozen) networks in which one or more nodes has more than one input. Mathematically this can only happen if there is complete synchronization between the input nodes, not simply that the values approach each other. This is highly unlikely, and the fact that we find these configurations in examples reflects the speed of convergence to synchronization and the finite precision of the simulations. To avoid this possibility, and to make our proofs simpler, we use a slightly modified rewiring rule in this section. Equations (5)-(8) specify how to change the network structure if there exists j with D ij n > D i n . We now add an additional rule to resolve the ambiguity that arises if node i has more than one neighbour and D ij n = D i n for all of these neighbours. The additional rule, which we will refer to as rule (R), states: (R)If B i n = ∅ and i has more than one neighbour then choose a neighbour k i at random and set Thus if all the neighbours are 'good', which would have led to no rewiring in the previous rule, we choose one of these at random and rewire all inputs to node i from this choice. Of course, once this is done there can be no further rewiring (as the node k is compared only to itself) and so the connection to node i is from a single node k i . In terms of the polysynchronous states observed numerically, this rule would lead to further evolution in the network topology of the polysynchronous networks shown in the previous section but polysynchrony would still occur (although in simpler networks).
Note that since the dynamics has a stochastic component it is not surprising that the freezing theorem will also be probabilistic: we will prove that the probability that the network has not frozen by time n, i.e. the dynamics of the network (but not necessarily the dynamics on the nodes) is stationary from time n onwards, tends to zero as n tends to infinity.
Before writing down the detailed calculations we will describe the strategy of the proof. We begin by considering a slightly modified system; one that is realized with non-zero probability in the dynamics described above. The finite probability system analyzed here is a subsystem of the general case in which at each time step, every connection that is rewired is rewired to the good node to which it already has the most connections (or one of these at random if there are two or more such nodes). At each time step, either a node has only one neighbour, and there can be no rewiring, or the number of connections to the most connected node increases by at least one. Since each node has m inputs, there will be one neighbour within m time steps. At each time step this happens in the real system with a probability that is bounded below by a fixed non-zero p. Hence for any finite T there is a finite probability (greater than p T ) that this revised rule will be used at each of the next T time steps and hence, as n goes to infinity, the probability of freezing goes to one. Returning to the original system this modified system occurs for m times steps with a finite probability, and hence the probability that this modified rule is applied is non-zero and the probability that the original system does not freeze tends to zero as time goes to infinity. We will now provide the details.
The modified system is specified as follows. At each time step n the mean distance D i n is calculated which determines the good set and bad set, G i n and B i n for each i ∈ {1, . . . , N } as explained in section II, equations (3)- (6). Now, the rewiring condition (5) implies that at least one of the nodes that is 'good' for a given node at time n − 1 is also 'good' for that node at time n, so it is always possible to choose a node k(i, n) ∈ G i n such that where the maximum is over j ∈ G i n−1 ∩ G i n , and if the maximum is attained by more than one node, one of these is chosen at random. Then if b i n is the valency of the bad nodes as defined in (7) then If b i n = 0 then by definition A ik(i,n) n > A ik(i,n−1) n−1 and so (since they are bounded by m) after at most m iterations for each i there exists k, and r ≤ m such that A ik r = m, b i r = 0 and A ik s = m for all s > r. In other words the network has frozen. This rule could be the outcome of the original rules when the bad set is non-empty if all but one of the numbers r ij n were zero and so the remaining r ij n equals b i n and this final connection is to a particular chosen node (that with the largest current connectivity to i). If there are s bad nodes and N − s good nodes, then for a given i the probability of picking the 'right' good node is 1/(N − s) and so the probability of rewiring all the bad nodes to this node is 1/(N − s) b i n . Now, b i n ≤ m and N − s ≤ N so the probability of making this choice is greater than (1/N ) m . This is true for each of the N nodes labelled by i and so the probability of the original system behaving in this way in one time step is greater than (1/N ) mN . Now consider using the modified rule (13) together with the additional rule (R). Then at each time step either the number of connections of the most connected node to i increases by at least one, or there is only one node connected to i and so there can be no further changes to the connections to i. Since there are a total of m connections to each node, this latter state must be achieved within m time steps of this modified system, after which it is frozen (and it is frozen whichever rules are used after this stage). Now, the probability of applying this modified rule to the original system for m consecutive time steps is just (1/N ) m 2 N , so if we return to the original system with our additional rule (R), time can be divided up into segments of length m, and so in time rm there are r independent opportunities to use the modified rule that leads to freezing, each with probability greater than (1/N ) m 2 N , where the extra factor of m in the exponent reflects the fact that the modified rule is applied at most m times to obtain the frozen state. So the probability that the system does not freeze in time rm is less than which obviously tends to zero as r → ∞, completing the proof.
The estimate of the probability could be improved considerably, for example by considering overlapping time intervals, but we are only interested in whether the probability of this not happening tends to zero, and for this the argument above suffices and has the virtue of simplicity.
Note that each node of the frozen network topology has precisely one neighbour. From this it is easy to show that each connected component of the network has one and only one strongly connected component (a cycle) and then trees based on the elements of the cycle. This means that the network eventually has precisely the hierarchical structure of [33] when the cycle is considered as the root of the network.

V. DISCUSSION
In this paper we have studied the dynamics of a simple adaptive network model as a function of the coupling parameter. We have rigorously proved that the network reaches a frozen state where the rewiring stops. We have shown that the final topologies are usually hierarchical and that polysynchronous dynamics appear in the frozen networks for certain parameter values. The hierarchical structure of the networks facilitates the appearance of polysynchrony as a stable attractor of the dynamics by making it easier to establish a balanced equivalence relation on the nodes. The stability study of different polysynchronous states explains the concrete coupling parameter ranges for which polysynchrony can be observed.
Unlike many adaptive network studies, the system described here evolves from a totally connected initial state to a much more constrained final topology. This simplifying structure could be relevant to the formation of functional groups in social interactions of biological systems. The dynamics on the network also has rich features; so far as we are aware this is the first network which can evolve naturally to a polysynchronous state. Such states could describe a form of functional evolution where a uniform population separates into different functional units described by different synchrony classes. The novel feature of polysynchrony is that these groups do not have to separate spatially as in the standard clusters, which are directly connected within the network. From this point of view, a fast time adaptive network of the type described here (to establish differentiated populations amongst a uniform set of initial nodes) followed by a slow differentiation process to lock in the differences created by the different synchrony classes could be a model for processes that require a mixed heterogeneous population from an initially homogenous set. Our models bear some resemblance to models in population dynamics (metapopulations, see [38]) and so this may be another area where polysynchrony might arise.

Appendix A: Synchronization dynamics of strongly connected components
To study the stability of the synchronized state of the different strongly connected components we will follow the approach exposed in [39]. Some of these results are well documented in the literature and are shown here for the sake of completeness.
In all the cases we study the coupling is linear and can be written, in general, as where L is the coupling linear operator. The synchronous state exists if the operator L has an eigenvalue σ 1 = 1 corresponding to the eigenvector e 1 = (1, 1, . . . , 1). Since our coupling is dissipative, the rest of the eigenvalues of the coupling operator are in modulus less than one. The stability of the synchronized state is then given by the condition where λ ⊥ is the transverse lyapunov exponent, λ is the lyapunov exponent of the uncoupled map and σ 2 is the second largest eigenvalue of the coupling operator.

Dynamics of the completely connected pair
The completely connected pair (Fig. 6A.) forms a system of two coupled logistic maps . (A3) This system has been thoroughly studied as a model of population dynamics in [40][41][42].

Dynamics of the triplet with transposition symmetry
The linear operator in the case of the triplet with transposition symmetry (Fig. 6B.) is with eigenvalues corresponding to the eigenvectors e 1 = (1, 1, 1), e 2 = ( m ′ −m m ′ , 0, 1), e 3 = (1, −1, 1). We should note that which eigenvalue has the second largest modulus depends on the value of ε and the stability condition has to be evaluated for both σ 2 and σ 3 . It is an easy calculation to deduce that the synchronized chaotic state will be stable in the range

Dynamics of the 3-cycle
The linear operator of the 3-cycle (Fig. 6C.)reads and has eigenvalues σ 1 = 1, Thus, the stability condition of the synchronous state reduces to Solving this for ε provides us with the condition where Therefore, the stability region for the synchronous state of the 3-cycle is an interval centered around ε = 0.5 of a width depending on the lyapunov exponent λ of the map. When λ = ln 2, B vanishes and the synchronous state becomes unstable for all ε.

Dynamics of the unidirectional coupling
Apart from the dynamics of the strongly connected components, it is necessary to study the case where a node is influenced by a single neighbour following an orbit of the uncoupled logistic map or, equivalently, by a fully synchronized neighbourhood. In both cases the dynamics is given by where x n is the variable of the node being influenced and y n the trajectory of the synchronized neighbourhood. Note that it is implied in the equation that the input is a trajectory of the uncoupled map (y n+1 = f (y n )). If this were not the case we could not perform this analysis. The eigenvalues of the linear operator are σ 1 = 1 (e 1 = (1, 1)) and σ 2 = 1 − ε (e 2 = (1, 0)). Therefore, the influenced node will synchronize to its input if It is very interesting to note that the eigenvalues are independent of m ′ and therefore our conclusions are also valid for a completely connected pair with an outgoing link to a third node (as in the 3-node subsystem of Fig. 3.e.).

Period-2 polysynchronous state
As described in [40], the period-2 orbit of the completely connected pair has as its elements This period-2 dynamics is stable in the range 0.13925 ε 0.193814. Similarly to the fixed-point case, the stability of the period-2 polysychronous state can be studied as the stability of a period-2 orbit of a three dimensional system. Therefore, we should observe the eigenvalues of the jacobian matrix J 2 (z 1 , z 2 ) = J(z 1 , z 2 , z 1 ) · J(z 2 , z 1 , z 2 ). (B1) The eigenvalues of J 2 have been graphed in Fig. 16 as a function of the coupling strength. From this figure we see that the period-2 polysynchronous state is stable when 0.140375 ε 0.193814.
Interestingly, and contrary to the fixed-point case, there is an interval of coupling strengths (0.13925 ε 0.140375) where the period-2 orbit is stable in the pair but the period-2 polysynchronous state of the triplet is unstable. It is in this interval of values where the spatial period-doubling phenomena appears (see Fig. 12 and description in the text).

Chaotic polysynchronous state
To study the stability of the chaotic polysynchronous states we start by making the following change of variable In these new variables, the dynamics is given by The polysynchronous chaotic state corresponds to the case W n = 0 with U n and V n following a non-synchronous chaotic orbit that we denote U * ,V * . By expanding the equation for W n+1 around (U * , V * , 0) we obtain W n+1 ≈ −2(1 − ε)U * n W n .
Assuming ergodicity, the transverse Lyapunov exponent (transverse to the surface of R 3 where the polysynchronous orbit lies) can be written as [39] λ ⊤ = ln | − 2(1 − ε)| + λ * , where λ * is the average Lyapunov exponent of the orbit U * n and has the upper bound ln 2, that would correspond to the complete synchronization of U n and V n (or equivalently, x i and x j ). We need λ ⊤ < 0 for the polysynchronous chaotic state to be stable. This condition provides us with a relation between the Lyapunov exponent of the chaotic orbit of U * and the minimum coupling strength necessary for polysynchrony to be stable It is easy to see that ε min is always greater than 0.5 for λ * ≤ ln 2. This is in agreement with the numerical experiments, which do not witness chaotic polysynchrony for ε < 0.5.