Finding Quasi-Optimal Network Topologies for Information Transmission in Active Networks

This work clarifies the relation between network circuit (topology) and behaviour (information transmission and synchronization) in active networks, e.g. neural networks. As an application, we show how one can find network topologies that are able to transmit a large amount of information, possess a large number of communication channels, and are robust under large variations of the network coupling configuration. This theoretical approach is general and does not depend on the particular dynamic of the elements forming the network, since the network topology can be determined by finding a Laplacian matrix (the matrix that describes the connections and the coupling strengths among the elements) whose eigenvalues satisfy some special conditions. To illustrate our ideas and theoretical approaches, we use neural networks of electrically connected chaotic Hindmarsh-Rose neurons.


Introduction
Given an arbitrary time dependent stimulus that externally excites an active network formed by systems that have some intrinsic dynamics (e.g. neurons and oscillators), how much information from such stimulus can be realized by measuring the time evolution of one of the elements of the network ? Determining how and how much information flows along anatomical brain paths is an important requirement for the understanding of how animals perceive their environment, learn and behave [1,2,3].
Even though the approaches of Ref. [1,2,3,4,5,6] have brought considerable understanding on how and how much information from a stimulus is transmitted in a neural network, the relation between network circuits (topology) and information transmission in a neural as well as an active network is still awaiting a more quantitative description [7]. And that is the main thrust of the present manuscript, namely, to present a quantitative way to relate network topology with information in active networks. Since information might not always be easy to be measured or quantified in experiments, we endeavour to clarify the relation between information and synchronization, a phenomenon which is often not only possible to observe but also relatively easy to characterize.
We initially proceed along the same line as in Refs. [8,9], and study the information transfer in autonomous systems. However, instead of treating the information transfer between dynamical systems components, we treat the transfer of information per unit time exchanged between two elements in an autonomous chaotic active network. Thus, we neglect the complex relation between external stimulus and the network and show how to calculate an upper bound value for the mutual information rate (MIR) exchanged between two elements (a communication channel) in an autonomous network. Ultimately, we discuss how to extend this formula to non-chaotic networks suffering the influence of a time-dependent stimulus.
Most of this work is directed to ensure the plausibility and validity of the proposed formula for the upper bound of MIR (Sec. Results) and also to study its applications in order to clarify the relation among network topology, information, and synchronization. We do not rely only on results provided by this formula, but we also calculate the MIR by the methods in Refs. [10,11] and by symbolic encoding the trajectory of the elements forming the network and then measuring the mutual information provided by this discrete sequence of symbols.
To illustrate the power of the proposed formula, we applied it to study the exchange of information in networks of coupled chaotic maps (Sec. Methods) and in Hindmarsh-Rose neural networks bidirectionally electrically coupled (Sec. Results). Our formula can be used to a larger class of active networks than the ones here considered. As the networks formed by elements coupled both electrically and chemically (see Ref. [12]). Still, the studied network topologies are much simpler than the ones found in the brain [13,14]. Nevertheless, we do believe our approaches can be used to better understand how information is transfered in more realistic networks as the scale-free networks [15], the small-world networks [16], or power-law networks [17].
The analyses are carried out using quantities that we believe to be relevant to the treatment of information transmission in active networks: a communication channel, the channel capacity, and the network capacity (see definitions in Sec. Methods).
A communication channel represents a pathway through which information is exchanged. In this work, a communication channel is considered to be formed by a pair of elements. One element represents a transmitter and the other a receiver, where the information about the transmitter can be measured.
The channel capacity is defined in terms of the proposed upper bound for the MIR. It measures the local maximal rate of information that two elements in a given network are able to exchange, a point-to-point measure of information exchange. As we shall see, there are two network configurations for which the value of the upper bound can be considered to be maximal with respect to the coupling strength.
The network capacity is the maximum of the KS-entropy, for many possible network configurations with a given number of elements. It gives the amount of independent information that can be simultaneously transmitted within the whole network, and naturally bounds the value of the MIR in the channels, which concerns only the transmission of information between two elements.
While the channel capacity is bounded and does not depend on the number of elements forming the network, the network capacity depends on the number of elements forming the network.
As a direct application of the formula for the upper bound value of the MIR, we show that an active network can operate with a large amount of MIR and KS-entropy and at the same time it is robustly resistant to alterations in the coupling strengths, if the eigenvalues of the Laplacian matrix satisfy some specified conditions (Sec. Results). The Laplacian matrix describes the connections among the elements of the network.
The conditions on the eigenvalues depend on whether the network is constructed in order to possess communication channels that are either self-excitable or non-self-excitable (see definition in Sec. Methods). Active networks that possess non-selfexcitable channels (formed by oscillators as the Rössler, or the Chua's circuit) have channels that achieve their capacity whenever their elements are in complete synchrony. Therefore, if a large amount of information is desired to be transmitted point-to-point in a non-self-excitable network, easily synchronizable networks are required. On the other hand, networks that possess self-excitable channels (as the ones formed by neurons), achieve simultaneously its channel and network capacities when there is at least one unstable mode of oscillation (time-scale) that is out of synchrony.
While non-self-excitable channels permit the exchanging of a moderate amount of information in a reliable fashion, due to the low level of desynchronization in the channel, self-excitable channels permit the exchange of surprisingly large amounts of information, not necessary reliable, due to the higher level of desynchronization in the channel.
We do not intend to find the best network topology among all possible ones. But rather, we aim at finding classes of network topologies that can not only transmit large amounts of information but are also robust under alterations in the coupling strengths. We arrive at two relevant eigenvalues conditions which provide networks that satisfy all these requirements. Either the network has elements that remain completely desynchronous for large variations of the coupling strength, forming the self-excitable channels, or the network has elements almost completely synchronous, forming the non-selfexcitable channels. In fact, the studied network, a network formed by electrically connected Hindmarsh-Rose neurons [18], can have simultaneously self-excitable and non-self-excitable channels.
Self-excitable networks, namely those that have a majority number of self-excitable channels, have the topology of a perturbed star, i.e., they are composed of a central neuron connected to most of the other outer neurons, and some outer neurons sparsely connected among themselves. The networks that have non-self-excitable channels have the topology of a perturbed fully connected network, i.e., a network whose elements are almost all-to-all connected. The self-excitable network has thus a topology which can be considered to be a model for mini-columnar structure of the mammalian neocortex [19].
In order to find quasi-optimal network topologies, we have used (Sec. Results) a Monte Carlo evolution technique [20], assuming equal bidirectional coupling strengths. This evolving technique simulates the rewiring of a neuron network that maximizes or minimizes some cost function, in this case a cost function which produces quasi-optimal networks to transmit information.
Finally, we discuss how to extend these results to networks formed by elements that are non-chaotic (Sec. Results), and to non-autonomous networks, that are being perturbed by some time-dependent stimuli (Sec. Results).

Upper bound for the Mutual Information Rate (MIR) in an Active Network
In a recent publication [10], we have argued that the mutual information rate (MIR) between two elements in an active chaotic network, namely, the amount of information per unit time that can be realized in one element, k, by measuring another element, l, regarded as I C , is given by the sum of the conditional Lyapunov exponents associated with the synchronization manifold (regarded as l I ) minus the positive conditional Lyapunov exponents associated with the transversal manifold (regarded as l H ). So, As shown in [11], if one has N = 2 coupled chaotic systems, which produce at most two positive Lyapunov exponents l 1 , l 2 with l 1 .l 2 , then l I = l 1 and l H = l 2 . Denote the trajectory of the element k in the network by x k . For larger number of elements, N, the approaches proposed in [10] remain valid whenever the coordinate transformation X klI = x k +x l (which defines the synchronization manifold) and X klH = x k 2x l (which defines the transversal manifold) successfully separates the two systems k and l from the whole network. Such a situation arises in networks of chaotic maps of the interval connected by a diffusively (also known as electrically or linear) all-to-all topology, where every element is connected to all the other elements. These approaches were also shown to be approximately valid for chaotic networks of oscillators connected by a diffusively all-to-all topology. The purpose of the present work is to extend these approaches and ideas to active networks with arbitrary topologies.
Consider an active network formed by N equal elements, x i (i = 1,…,N), where every D-dimensional element has a different set of initial conditions, i.e., x 1 ?x 2 ?…?xN. The network is described by where G ij is the ij element of the coupling matrix. Since we choose P j G ij~0 in order for a synchronization manifold to exist by the subspace g = x 1 = x 2 = x 3 = … = x N , we can call this matrix the Laplacian matrix. The way small perturbations propagate in the network [21] is described by the i (i = 1,…,N) variational equations of Eqs. (1), namely writing x i = g+dx i and expanding Eq.
Making x i = j, which can be easily numerically done by setting the elements with equal initial conditions and taking H(x j ) = x j , Eq. (2) can be made block diagonal resulting in where c i are the eigenvalues (positive defined) of the Laplacian matrix ordered such that c i+1 $c i . Note that c 1 = 0. Notice that the network dynamics is described by Eq. (1), which assumes that every element has different initial conditions and therefore different trajectories (except when the elements are completely synchronized). On the other hand, Eq. (3) that provides the conditional exponents considers that all the initial conditions are equal. While Eq. gives the upper bound for the amount of information between modes of oscillation, for some simple network geometries, as the ones studied here, we can relate the amount of information exchanged between two vibrational modes to the amount of information between two elements of the network, and therefore, Eq. (4) can be used to calculate an upper bound for the MIR exchanged between pairs of elements in the network. For larger and complex networks, this association is non-trivial, and we rely on the reasonable argument that a pair of elements in an active network cannot transmit more information than some of the i21 values of I i{1 P . The inequality in Eq. (4) can be interpreted in the following way. The right hand side of Eq. (4) calculates the amount of information that one could transmit if the whole network were completely synchronous with the state j, which is only true when complete synchronization takes place and when all the nodes have equal dynamics. Typically, we expect that the elements of the network will not be completely synchronous to j and in realistic networks, the nodes will not be equal. Thus, the amount of information provided by the right part of Eq. (4) overestimates the exact MIR which, due to desynchronization in the network, should be smaller than the calculated one. Equation (5) allows one to calculate the MIR between oscillation modes of larger networks with arbitrary topology rescaling the MIR curve (I 1 P vs. s) obtained from two coupled elements. Denoting s * (N = 2) as the strength value for which the curve for l 2 reaches a relevant value, say, its maximum value, then the coupling strength for which this same maximum is reached for l i in a network composed by N elements is given by where c i (N) represents the ith largest eigenvalue of the N-elements network. If the network has an all-to-all topology, thus, s * (N = 2) represents the strength value for which the curve of I 1 P reaches a relevant value, and s * (N) the strength value that this same value for I i P is reached. Notice that symmetries in the connecting network topology leads to the presence of degenerate eigenvalues ( = equal eigenvalues) in the Laplacian matrix, which means that there are less independent channels of communication along which information flows. Calling Q the number of degenerate eigenvalues of the Laplacian matrix, Eq. (4) will provide N2Q different values.
As the coupling strength s is varied, the quantities that measure information change correspondingly. For practical reasons, it is important that we can link the way these quantities (see Sec. Finally, there are a few more relevant coupling strengths, which characterize each communication channel. First, s i min , for which l i equals the value of l 1 , with i$2. For svs i min , the communication channel c i21 (whose upper rate of information transmission depends on the two oscillation modes j 1 and j i ) behaves in a self-excitable way, i.e., l 1 ,l i . For s §s i min , l 1 $l i . Secondly, s i* indicates the coupling strength at which I i{1 P is maximal. Thirdly, s i CS indicates the coupling strength for which the communication channel c i21 becomes ''stable'', i.e., l i ,0. At s = s i* the selfexcitable channel capacity of the channel c i21 is reached and at s~s i CS , the non-self-excitable channel capacity is reached. Finally, s C is the coupling for which the network capacity is reached, and then, when the KS-entropy of the network is maximal.

The MIR in networks of coupled Hindmarsh-Rose neurons
We investigate how information is transmitted in self-excitable networks composed of N bidirectionally coupled Hindmarsh-Rose neurons [18]: The parameter r modulates the slow dynamics and is set equal to 0.005, such that each neuron is chaotic. The index i?j assumes values within the set [1,…,N]. S k represents the subsystem formed by the variables (x k , y k , z k ) and S l represents the subsystem formed by the variables (x l , y l , z l ), where k = [1,…,N21] and l = [k+1,…,N]. The Laplacian matrix is symmetric, so G ji = G ij , and sG ji is the strength of the electrical coupling between the neurons, and we take for I i the value I i = 3.25.
In order to simulate the neuron network and to calculate the Lyapunov exponents through Eq. (2), we use the initial conditions x = 21.3078+g, y = 27.3218+g, and z = 3.3530+g, where g is an uniform random number within [0,0.02]. To calculate the conditional Lyapunov exponents, we use the equal initial conditions, x = 21.3078, y = 27.3218, and z = 3.3530.
All-to-all coupling. Here, we analyze the case where N neurons are fully connected to every other neuron. The Laplacian matrix has N eigenvalues, c 1 = 0, and N21 degenerate ones c i = N, i = 2,…,N. Every pair of neurons exchange an equal amount of MIR. Although, there are N6(N21)/2 pairs of neurons, there is actually only one independent channel of communication, i.e., a perturbation applied at some point of the network should be equally propagated to all other points in the network. In Fig. 1(A), we show the MIR, I C , calculated using the approaches in Refs. [10,11], I P , calculated using the right hand-side of Eq. (4), and I S , calculated encoding the trajectory between pair of neurons, and the Kolmogorov-Sinai entropy, H KS , for a network composed by N = 2 neurons. In (B), we show these same quantities for a network formed by N = 4 neurons.
While for s>0 and s$s CS , we have that I C >I P >I S , for s>s 2* (when the self-excitable channel capacity is reached) it is clear that I P should be an upper bound for the MIR, since not only I P .I C but also I P .I S . Notice the good agreement between I C and I S , except for s%s 2 min , when I S .H KS , which violates Eq. (11). The star symbol indicates the value of the coupling, s BPS , for which burst phase synchronization (BPS) appears while the spikes are highly desynchronous. The appearance of BPS coincides with the moment where all the quantifiers for the MIR are large, and close to a coupling strength, s C , for which the network capacity is reached (when H KS is maximal).
At this point, the network is sufficiently desynchronous to generate a large amount of entropy, which implies a large l i , for i$2. This is an ideal configuration for the maximization of the MIR. There exists phase synchrony in the subspace of the slow time-scale z variables (which is responsible for the bursting-spiking behavior), but there is no synchrony in the (x,y) subspace. This supports the binding hypothesis, a fundamental concept of neurobiology [19] which sustains that neural networks coding the same feature or object are functionally bounded. It also simultaneously supports the works of [22], which show that desynchronization seems to play an important role in the perception of objects as well. Whenever l 2 approaches zero, at s = s CS , there is a drastic reduction in the value of H KS as well as I P , since the network is in complete synchronization (CS), when all the variables of one neuron equals the variables of the other neurons.
Therefore, for coupling strengths larger than the one indicated by the star symbol, and smaller than the one where CS takes place, there is still one time-scale, the fast time-scale, which is out of synchrony.
For sws 2 min , the only independent communication channel is of the non-self-excitable type. That means l i #l 1 (i$2), and as the coupling strength increases, H KS decreases and I P increases.
Note that the curve for I P shown in Fig. 1(B) can be obtained by rescaling the curve shown in Fig. 1(A), applying Eq. (5).
Star coupling. We consider N = 4. There is a central neuron, denoted by S 1 , bidirectionally connected to the other three (S k , k = 2,3,4), but none of the others are connected among themselves. The eigenvalues of the Laplacian matrix are c 1 = 0, c 2,3 = 1, To treat general types of networks, it is useful to define two quantities related to the excitability of the communication channels. The here called non-self-excitable (NSE) robustness parameter of the channel c i21 (i$2) as Ds i NSE~s CS {s i min and the self-excitable (SE) robustness parameter for the communication channel c i21 as Ds i SE~s i min (i$2). It is also useful to define a quantity that measures the distance between the eigenvalues, the normalized spectral distance (NED) between the two eigenvalues.
Having a large NED between the ith largest and the first largest eigenvalues (c i 2c 2 )/N, results in a non-self-excitable channel, c i21 , with a large NSE robustness parameter that implies that the channel preserves its NSE character under large alterations of the coupling strength. On the other hand, having a large NED between the largest and the ith largest eigenvalues (c N 2c i )/N, results in a self-excitable channel, c i21 , with a large self-excitable robustness parameter that implies that the channel preserves its SE character under large alterations of the coupling strength.
So, for the star topology network, not only the NED between c N and c N21 is large but also between c N and c N22 , and therefore, Ds N{1 SE and Ds N{2 SE are large. This provides a network whose channels c 1 and c 2 have a large MIR for a large coupling strength alteration. Note that if c N21 is far away from c N that implies that c N22 is also far away from c N . Thus, a reasonable spectral distance between c N21 and c N is a ''biological requirement'' for the proper function of the network, since even for larger coupling strengths there will be at least one oscillation mode which is desynchronous, a configuration that enables perturbation (meaning external stimuli) to be propagated within the network [23].
The largest eigenvalue is related to an oscillation mode where all the outer neurons are in synchrony with each other but desynchronous with the central neuron. So, here it is clear the association between |l 1 2l 4 | and the MIR between the central neuron with an outer neuron, since l 1 represents the amount of information of the synchronous trajectories among all the neurons, while l 4 is the amount of information of the desynchronous trajectories between the central neuron and any outer neuron. The other eigenvalues (c 2 ,c 3 ) represent directions transverse to the synchronization manifold in which the outer neurons become desynchronous with the central neuron in waves wrapping commensurately around the central neuron [21]. Thus, l 2 and l 3 are related to the error in the transmission between two outer neurons, k and l, with k,l?1. Notice that this network topology has two independent channels of communication.
Note that the MIR between S 1 and an outer neuron (upper bound represented by I 3 represented by I S (k, l), in Fig. (2), for small coupling (for when the channel c 3 is self-excitable, and s §s 4 min ). Similar to what happens to nearest-neighbour networks, the self-excitable and the non-selfexcitable channel capacities of the channel associated with the transmission of information between closer elements (the channel c 3 ) are achieved for a smaller value of the coupling strength than the one necessary to make the channels associated with the transmission of information between more distant elements (the channel c 1 ) to achieve its two channel capacities. That property permits this network, for s%s 4 min , to transmits simultaneously reliable information using the channel c 3 and with a higher rate using the channel c 1 .
Notice, in Fig. 2, that s 2 1 %s 4 min %s BPS %s C . So, when the channel capacity of the channel c 1 is reached, also H KS of the network is maximal, and the network operates with its capacity.
Another point that we want to emphasize in this network is that while a large NED between c N and c N21 provides a network whose channel c 1 is self-excitable and can transmit information at a large rate for a large coupling strength interval, a large NED between c 3 and c 2 leads to a non-self-excitable channel c 3 even for small values of the coupling amplitudes, and it remains non-selfexcitable for a large variation of the coupling strength. Thus, while a large NED between the second and the first largest eigenvalues leads to a network whose channels are predominantly of the selfexcitable types, a large NED between the second largest and the third largest eigenvalues provide a network whose communication channels are predominantly of the non-self-excitable types.

Eigenvalues conditions
Finding network topologies and coupling strengths in order to have a network that operates in a desired fashion is not a trivial task (see Sec. Methods). An ideal way to proceed would be to evolve the network topology in order to achieve some desired behaviour. In this paper, we are interested in maximizing simultaneously I P , the KS-entropy, and the average AEI P ae, for a large range of the coupling strength, characteristics of a quasioptimal network. However, evolving a network in order to find a quasi-optimal one would require the calculation of the MIR in every communication channel and H KS for every evolution step. For a typical evolution, which requires 10 6 evolution steps, such an approach is impractical.
Based on our previous discussions, however, a quasi-optimal network topology can be realized by only selecting an appropriate set of eigenvalues which have some specific NED. Evolving a network by the method of Sec. Methods using a cost function which is a function of only the eigenvalues of the Laplacian matrix is a practical and physible task.
The present section is dedicated to describe the derivation of this cost function.
We can think of two most relevant sets of eigenvalues which create quasi-optimal networks, and they are represented in Fig. 3. Either it is desired eigenvalues that produce a network predominantly self-excitable [SE, in Fig. 3] or predominantly non-self-excitable [NSE, in Fig. 3].
In a network whose communication channels are predominantly self-excitable, it is required that the NED (c N 2c N21 )/N is maximal and (c N21 )/N minimal. Therefore, we want a network for which the cost function is maximal. A network whose eigenvalues maximize B 1 has self-excitable channels for a large variation of the coupling strength. As a consequence, AEI P ae as well as H KS is large for s[ts N min ,s 2 min s. In a network whose communication channels are predominantly non-self-excitable, it is required that the NED (c 3 2c 2 )/N is maximal and (c 2 )/N minimal. Therefore, we want a network for which the cost function is maximal. A network whose eigenvalues maximize the condition in Eq. (8) have non-self-excitable channels for a large variation of the coupling strength. As a consequence, AEI P ae is large for s[ts N min ,s 3 min s, which is a small coupling range, but since there is still one oscillation mode that is unstable (the mode j 2 ), H KS is still large for a large range of the coupling strength svs 2 min À Á . Most of the channels will transmit information in a reliable way, since the error in the transmission, provided by l i (i$2), of most of the channels will be zero, once l i ,0.
Since degenerate eigenvalues produce networks with less vibrational modes and therefore less independent channels of communication, we assume in the following the absence of such degenerate eigenvalues. In addition, we assume that there is a finite distance between eigenvalues so that the network becomes robust under rewiring, and therefore, perturbing G ij will not easily create degenerate eigenvalues.
A network that is completely synchronous and has no unstable modes does not provide an appropriate environment for the transmission of information about an external stimulus, because they prevent the propagation of perturbations. Networks that can be easily completely synchronized (for small coupling strengths) requires the minimization of c N 2c 2 , or in terms of the eigenratio, the minimization of c N /c 2 . We are not interested in such a case. To construct network topologies that are good for complete synchronization, see Refs. [21,24,25,26].

Quasi-optimal topologies for information transmission
Before explaining how we obtain quasi-optimal network topologies for information transmission, it is important to discuss the type of topology expected to be found by maximizing either B 1 , in Eq. (7) or B 2 , in Eq. (8). Notice that Laplacians whose eigenvalues maximize B 1 are a perturbed version of the star topology, and the ones that maximize B 2 are a perturbed version of the all-to-all topology. In addition, in order to have a network that presents many independent modes of oscillations it is required that the Laplacian matrix presents as much as possible, a large number of non-degenerate eigenvalues. That can be arranged by rewiring (perturbing) networks possessing either the star or the nearest-neighbour topology, breaking the symmetry.
In order to calculate a Laplacian from a quasi-optimal network, we propose an approach described in Sec. Methods, based on the reconstruction of the network by evolving techniques, simulating the process responsible for the growing or rewiring of real biological networks, a process which tries to maximize or minimize some cost function.
In order to better understand how a network evolves (grows) in accordance with the maximization of the cost functions in Eqs. (7) and (8), we first find the network configurations with a small number of elements. To be specific, we choose N = 8 elements. To show that indeed the calculated network topologies produce active networks that operate as desired, we calculate the average upper bound value of the MIR [Eq. (10)] for neural networks described by Eqs. (6) with the topology obtained by the evolution technique, and compare with other network topologies. Figure 4 shows AEI P ae, the average channel capacity, calculated for networks composed of 8 elements, using one of the many topologies obtained by evolving the network maximizing B 1 (circles, denoted in Fig. by ''evolving  1''), all-to-all topology (squares), star topology (diamonds), nearestneighbor (upper triangle), and maximizing B 2 (down triangle, denoted in Fig. by ''evolving 2''). The star points to the value of s 2 min , when c 1 , the most unstable communication channel (a selfexcitable channel), becomes non-self-excitable.
As desired the evolving network 1 has a large upper bound for the MIR (as measured by AEI P ae) for a large range of the coupling strength, since the network has predominantly self-excitable channels. The channel c 1 has a large robustness parameter Ds 2 SE , i.e., it is a self-excitable channel for svs 2 min , where s 2 min~2 :0. In contrast to the other topologies, in the star, nearestneighbour, and all-to-all topologies, Ds 2 SE is smaller and Ds 2 NSE is larger. Even though most of the channels in the evolving 2 topology are of the non-self-excitable type, AEI P ae remains large even for higher values of the coupling strength. That is due to the channel c 1 which turns into a self-excitable channel only for s.2.
The KS-entropies of the 5 active networks whose AEI P ae are shown in Fig. 4 are shown in Fig. 5. Typically, the network capacities are reached for roughly the same coupling strength for which the maximum of AEI P ae, is reached. In between the coupling strength for which the network capacities and the maximal of AEI P ae are reached, l 3 becomes negative. At this point, also BPS appears in the slow time-scale, suggesting that this phenomena is the behavioral signature of a network that is able to transmit not only large  amounts of information between pairs of elements (high MIR) but also overall within the network (high H KS ).
Note however, that since the evolving networks have a small number of elements, the cost function cannot reach higher values and therefore, the networks are not as quasi-optimal as they can be. For that reason, we proceed now to evolve larger networks, with N = 32.
Maximization of the cost function B 1 leads to the network connectivity shown in Fig. 6(A) and maximization of the cost function B 2 leads to the network connectivity shown in Fig. 6(B). In (A), the network has the topology of a perturbed star, a neuron connected to all the other outer neurons, thus a hub, and each outer neuron is sparsely connected to other outer neurons. The arrow points to the hub. In (B),the network has the topology of a perturbed all-to-all network, where elements are almost all-to-all connected. Note that there is one element, the neuron S 32 , which is only connected to one neuron, the S 1 . This isolated neuron is responsible to produce the large spectral gap between the eigenvalues c 3 and c 2 .
AEI P ae for the network topology represented in Fig. 6(A) is shown in Fig. 7 as circles, and AEI P ae for the network topology represented in Fig. 6(B) is shown in Fig. 7 as squares. We see that the star topology, whose connectivity is represented in 6(A), has larger AEI P ae for a larger coupling strength than the topology whose connectivity is represented in 6(B). Other relevant parameters of the network whose topology is represented in 6(A) are s 2 min~0 :8468, s 3 min~0 :8249, s N min~0 :0278, s CS = 0.9762 and for the topology represented in 6(B) are s 2 min~0 :8512, s 3 min~0 :042, s N min~0 :031, and s CS = 0.9761.
It is worth to comment that the neocortex is being simulated in the Blue Brain project, by roughly creating a large network composed of many small networks possessing the star topology. By doing that, one tries to recreate the way minicolumnar structures [19] are connected to minicolumnar structures of the neocortex [27]. Each minicolumn can be idealized as formed by a pyramidal neuron (the hub) connected to its interneurons, the outer neurons in the star topology, which are responsible for the connections among this minicolumn (small network) to others minicolumn. So, the used topology to simulate minicolumns is an good topology in what concerns the transmission of information.

Active networks formed by non-chaotic elements
The purpose of the present work is to describe how information is transmitted via an active media, a network formed by dynamical systems. There are three possible asymptotic stable behaviours for an autonomous dynamical system: chaotic, periodic, or quasiperiodic. A quasi-periodic behaviour can be usually replaced by either a chaotic or a periodic one, by an arbitrary perturbation. For that reason, we neglect such a state and focus the attention on active channels that are either chaotic or periodic.
Equation (4) is defined for positive exponents. However, such an equation can also be used to calculate an upper bound for the rate of mutual information in systems that also possess negative Lyapunov exponents. Consider first a one-dimensional contracting system being perturbed by a random stimulus. Further consider that the stimulus changes the intrinsic dynamics of this system. This mimics the process under which an active element adapts to the presence of a stimulus.
Suppose the stimulus, h n , can be described by a discrete binary random source with equal probabilities of generating '0' or '1'. Whenever h n = 0, the system presents the dynamics x n+1 = x n /2, otherwise x n+1 = (1+x n )/2. It is easy to see that the only Lyapunov exponent of this mapping, l 1 , which is equal to the conditional exponent, l 1 , is negative. Negative exponents do not contribute to the production of information. From Eq. (4) one would arrive at I P = 0. However, all the information about the stimulus is contained in the trajectory. If one measures the trajectory x n , one knows exactly what the stimulus was, either a '0' or a '1'. The amount of information contained in the stimulus is log(2) per iteration which equals the absolute value of the Lyapunov exponent, |l 1 |. In fact, it is easy to show that I C = I P = |l 1 | = |l 1 | = log(2), or if we use the interpretation of [28], I C = I P = l, where l = |l 1 | is the positive Lyapunov exponent of the time-inverse chaotic trajectory, x n+m , x n+m21 , …, x 0 , which equals the rate of information production of the random source. So, in this type of active communication channel, one would consider in Eq. (4) the positive Lyapunov exponents of the time-inverse trajectory, or the absolute value for the negative Lyapunov exponent.
Another example was given in [11]. In this reference we have shown that a chaotic stimulus perturbing an active system with a space contracting dynamics (a negative Lyapunov exponent) might produce a fractal set. We assume that one wants to obtain information about the stimulus by observing the fractal set. The rate of information retrieved about the stimulus on this fractal set equals the rate of information produced by the fractal set. This amount is given by D 1 |l|, where D 1 is the information dimension of the fractal set and |l| the absolute value of the negative Lyapunov exponent. In fact, D 1 |l| is also the rate of information produced by the stimulus. So, if an active system has a space contracting dynamics, the channel capacity equals the rate of information produced by the stimulus. In other words, the amount of information that the system allows to be transmitted equals the amount of information produced by the chaotic stimulus.

The role of a time-dependent stimulus in an active network
The most general way of modelling the action of an arbitrary stimulus perturbing an active network is by stimulating it using uncorrelated white noise. Let us assume that we have a large network with all the channels operating in non-self-excitable fashion. We also assume that all the transversal eigenmodes of oscillations except one are stable, and therefore do not suffer the influence of the noise. Let us also assume that the noise is acting only on one structurally stable ( = far from bifurcation points) element, S k . To calculate the upper bound of the MIR between the element S k and another element S l in the network, we assume that the action of the noise does not alter the value of l 1 . Then, the noise on the element S k is propagated along the vibrational mode associated with the one unstable transversal direction, whose conditional exponent is l 2 . As a consequence, the action of the noise might only increase l 2 , while not affecting the negativeness of all the other exponents (l m , m.2), associated with stable transversal modes of oscillation. That means that the channels responsible for transmitting large amounts of information (associated with l m , with m large) will not be affected. So, for such types of noises, Eq. (4) of the autonomous network is an upper bound for the non-autonomous network. Consider now a situation where the noise acts equally on all the elements of an active network. The mapping was proposed as a way to understand such a case. In this mapping, we consider r$0 and x n , y n M[0,1], which can be accomplished by applying the mod(1) operation. Note that the term r x 2 n zy 2 n À Á that enters equally in all the maps has statistical properties of an uniformly distributed random noise. Calculating I P for r = 0 (the noise-free map) we arrive at I P >2s, for small s, while the true MIR I C >2(s2r). These results are confirmed by exact numerical calculation of the Lyapunov exponents of Eq. (9) as well as the calculation of the conditional exponents of the variational equations. So, this example suggests that Eq. (4) calculated for an autonomous non-perturbed network gives the upper bound for the mutual information rate in a nonautonomous network.

Discussion
We have shown how to relate in an active network the rate of information that can be transmitted from one point to another, regarded as mutual information rate (MIR), the synchronization level among elements, and the connecting topology of the network. By active network, we mean a network formed by elements that have some intrinsic dynamics and can be described by classical dynamical systems, such as chaotic oscillators, neurons, phase oscillators, and so on.
Our main concern is to suggest how to construct a quasioptimal network. A network that simultaneously transmits information at a large rate, is robust under couplings alterations, and further, it possesses a large number of independent channels of communication, pathways along which information travels.
We find that there is not the best topology but many that can be classified in two classes. Self-excitable [maximizing Eq. (7)] or nonself-excitable [maximizing Eq. (8)] (see definition of self-excitability in Sec. Methods). Self-excitable networks have communication channels that transmit information in a higher rate for a large range of the coupling strength. Most of the oscillation modes in these networks are unstable, and therefore, information is mainly propagated in a desynchronous environment. Non-self-excitable networks have communication channels that transmit information in a higher rate for a small range of the coupling strength, however, they have channels that transmit reliable information in a moderate rate for large range of coupling strengths. Most of the oscillation modes in these networks are stable, and therefore, information is mainly propagated in a synchronous environment, a highly reliable environment for information transmission.
One of the main results of our work, the Eq. (4), which relates synchronization, topology and information in active networks, can only be used in networks composed of nodes that have equal dynamics. We have reasons to believe that if the nodes have non equal dynamics, Eq. (4) provides an upper bound for the value of the mutual information rate that modes in the network exchange. That was shown in Ref. [11] for two linear coupled maps. Another reason is given in the following. When the nodes are not completely synchronous, networks of nodes with equal dynamics but randomly coupled (as the networks in [12]) in Ref. [12]), are good models of networks with nodes that have different dynamics. We have found that these random networks with nodes electrically connected usually become more non-self-excitable than the networks with nodes being connected with equal bidirectional couplings. As a consequence, both the network capacity and the channel capacities become smaller. It remains still to be verified if that is so for networks whose nodes are connected with chemical synapses. As shown in Ref. [12], chemical couplings make the network to become highly excited. As a consequence, it might be that as the nodes are made non-equal, the network gains a selfexcitable character, resulting in an increase of the information capacities. In such a case, Eq. (4) would provide a lower bound for the mutual information rate of networks with nodes that have non equal dynamics. If brain-networks somehow grow in order to maximize the amount of information transmission, simultaneously remaining very robust under coupling alterations, the minimal topology that small neural networks must have should be similar to the one in Fig. 6(A), i.e., a network with a star topology, presenting a central element, a hub, very well connected to other outer elements, which are sparsely connected.

Self-excitability
In Ref. [11] self-excitability was defined in the following way. An active network formed by N elements, is said to be selfexcitable if H KS (N, s).H KS (N, s = 0), which means that the KSentropy of the network increases as the coupling strength is increased. Thus, for non self-excitable systems, an increase in the coupling strength among the elements forming the network leads to a decrease in the KS-entropy of the network.
Here, we adopt also a more flexible definition, in terms of the properties of each communication channel. We define that a communication channel c i behaves in a self-excitable fashion if l i .l 1 . It behaves in a non-self-excitable fashion if l i #l 1 .

Mutual Information Rate (MIR), channel capacity, and network capacity
In this work, the rate with which information is exchanged between two elements of the network is calculated by different ways. Using the approaches of Refs. [10,11], we can have an estimate of the real value of the MIR, and we refer to this estimate as I C . Whenever we use Eq. (4) to calculate the upper bound for the MIR, we will refer to it as I P . Finally, whenever we calculate the MIR through the symbolic encoding of the trajectory, we refer to it as I S .
We define the channel capacity of a communication channel formed by two oscillation modes depending on whether the channel behaves in a self-excitable fashion or not. So, for the studied network, every communication channel possess two channel capacities, the self-excitable capacity and the non-selfexcitable one. A channel c i operates with its self-excitable capacity when I i P is maximal, what happens at the parameter s (i+1)* . It operates with its non-self-excitable capacity when l i+1 = 0.
We also define the channel capacity in an average sense. In that case, the averaged channel capacity is given by the maximal value of the average value The network capacity of a network composed of N elements, C N (N), is defined to be the maximum value of the Kolmogorov-Sinai (KS) entropy, H KS , of the network. For chaotic networks, the KSentropy, as shown by Pesin [29], is the sum of all the positive Lyapunov exponents. Notice that if I denotes the MIR then As shown in Ref. [11] and from the many examples treated here, C N (N)/N, and so, the network capacity grows linearly with the number of elements in an active network.
Understanding Eq. (4): Positiveness of the MIR for selfexcitable channels in the (non-linear) HR network To show that indeed I i P should be positive in case of a selfexcitable channel in the HR network, one can imagine that in Eq. (1) the coupling strength is arbitrarily small and that N = 2. At this situation, the Lyapunov exponent spectra obtained from Eq. (2) are a first-order perturbative version of the conditional exponents, and they appear organized by their strengths. One arrives at l 1 >l 2 and l 2 >l 1 , which means that the largest Lyapunov exponent equals the transversal conditional exponent and the second largest Lyapunov exponent equals the conditional exponent associated with the synchronous manifold. Using similar arguments to the ones in Refs. [10,11,30], we have that the MIR is given by the largest Lyapunov exponent minus the second largest, and therefore, I C = l 1 2l 2 , which can be put in terms of conditional exponents as I P #l 2 2l 1 , or as represented in Eq. (4), I P #|l 1 2l 2 |.

Understanding Eq. (4): The inequality in Eq. (4)
To explain the reason of the inequality in Eq. (4), consider the following two coupled maps: x nz1~2 x n {rx 2 n z2ss y n {x n ð Þ , y nz1~2 y n {ry 2 n z2ss x n {y n ð Þ , with s = 1 and x n , y n M[0,1]. For this mapping, the MIR can be written in terms of the Lyapunov Exponents [11,31]. For two coupled systems, the MIR can be exactly calculated by I C = l 1 2l 2 , since l I = l 1 and l H = l 2 , assuming that both l 1 and l 2 are positive. Calculating the conditional exponents numerically, we can show that I P $I C , and thus I P is an upper bound for the MIR. For more details on this inequality, see [12] Evolutionary construction of a network In our simulations, we have evolved networks of equal bidirectional couplings [32]. That means that the Laplacian in Eq. (1) is a symmetric matrix of dimension N with integer entries {0,1} for the off diagonal elements, and the diagonal elements equal to { P j G ij , with i?j. Finding the network topologies which maximize B in Eq. (7) is impractical even for moderately large N. Figuring out by ''brute force'' which Laplacian produces the desired eigenvalue spectra would require the inspection of a number of 2 N N{1 ð Þ =2

N!
configurations. To overcome this difficulty, Ref. [20] proposed an evolutionary procedure in order to reconstruct the network in order to maximize some cost function. Their procedure has two main steps regarded as mutation and selection. The mutation steps correspond to a random modification of the pattern of connections. The selection steps consist in accepting or rejecting the mutated network, in accordance with the criterion of maximization of the cost function B, in Eq. (7).
We consider a random initial network configuration, with N elements, which produce an initial Laplacian G 0 , whose eigenvalues produce a value B 0 for the cost function. We take at random one element of this network and delete all links connected to it. In the following, we choose randomly a new degree k to this element and connect this element (in a bidirectional way) to k other elements randomly chosen. This procedure generates a new network that possesses the Laplacian G9, whose eigenvalues produce a value B9. To decide if this mutation is accepted or not, we calculate De = B92B 0 . If De.0, the new network whose Laplacian is G9 is accepted. If, on the other hand, De,0, we still accept the new mutation, but with a probability p(De) = exp(2De/ T). If a mutation is accepted then the network whose Laplacian is G 0 is replaced by the network whose Laplacian is G9.
The parameter T is a kind of ''temperature'' which controls the level of noise responsible for the mutations. It controls whether the evolution process converges or not. Usually, for high temperatures one expects the evolution never to converge, since new mutations that maximizes B are often not accepted. In our simulations, we have used T>0.0005.
These steps are applied iteratively up to the point when |De| = 0 for about 10,000 steps, being that we consider an evolution time of the order of 1,000,000 steps. That means that the evolution process has converged after the elapse of some time to an equilibrium state. If for more than one network topology |De| = 0 for about 10,000 steps, we choose the network that has the larger B value. This constraint avoids the task of finding the best network topology. However, we consider that a reasonably low number of mutations would recreate what usually happens in real networks.