Steering complex networks toward desired dynamics.

We consider networks of dynamical units that evolve in time according to different laws, and are coupled to each other in highly irregular ways. Studying how to steer the dynamics of such systems towards a desired evolution is of great practical interest in many areas of science, as well as providing insight into the interplay between network structure and dynamical behavior. We propose a pinning protocol for imposing specific dynamic evolutions compatible with the equations of motion on a networked system. The method does not impose any restrictions on the local dynamics, which may vary from node to node, nor on the interactions between nodes, which may adopt in principle any nonlinear mathematical form and be represented by weighted, directed or undirected links. We first explore our method on small synthetic networks of chaotic oscillators, which allows us to unveil a correlation between the ordered sequence of pinned nodes and their topological influence in the network. We then consider a 12-species trophic web network, which is a model of a mammalian food web. By pinning a relatively small number of species, one can make the system abandon its spontaneous evolution from its (typically uncontrolled) initial state towards a target dynamics, or periodically control it so as to make the populations evolve within stipulated bounds. The relevance of these findings for environment management and conservation is discussed.

Introduction

Controlling the dynamics of ensembles of units networking via irregular topologies is one of the foremost challenges of modern science, and, in fact, the literature of the last two decades abounds with proposals for network control. In some of the earliest contributions to the field, a pinning method based on applying linear feedback injections to some nodes of a network with the objective of stabilizing a given global fixed point was explored[1,2]. Pinning controllability was further studied in Ref.[3] as a way to synchronize to a given, time-dependent, network evolution. Similar approaches that expand or modify these initial efforts were developed in more recent contributions, see e.g.[4]. Later on, the introduction of multi-layer network representations[5] opened up new avenues, such as the study of complex-network targetability[6], based on considering an identical copy of the graph undergoing a desirable evolution, and gradually creating unidirectional actions from nodes of the copy to the corresponding nodes in the original network, until the latter becomes fully synchronized with the former. These and other related works follow the master stability function approach[7] in assuming that dynamical units are identical, and that their coupling function at each link is the same, in order to derive analytical criteria for controllability.

A different approach was proposed in Ref.[8], where conditions based on classical control and graph theories were given for the identification of the minimal set of nodes that, if forced to follow a prescribed time evolution, suffice to drive the entire network to the target dynamics. This is applicable to graphs whose dynamics is unknown, and for directed and weighted connectivities, when weights may possibly be unknown too. This framework was also used to investigate network properties and their connection to structural controllability[9-11]. However, all these results come at the price of introducing drastic restrictions, as the study focuses on scalar state variables governed by linear equations of motion. In a recent contribution, moreover, control frameworks that are purely based on network topological properties (and completely ignore dynamical considerations), such as this one or the one proposed in Ref.[12], have been shown to fail in Boolean networks and models of biochemical regulation[13]. Later developments along these lines include Ref.[14], which develops a perturbation approach to optimize the structural controllability of a complex network, and Ref.[15], which generalizes the approach to a wider set of topologies via spectral techniques. In one way or another, all such methods rely on rather restrictive sets of assumptions that are not always fulfilled in applications. This does not just mean that, at some level, there is always some degree of approximation: it may just be the case that different assumptions lead to radically different results. For instance, according to Ref.[8] the denser and more homogeneous a network is, the fewer nodes are needed to control its dynamics (though the conclusion appears to be different in later refinements aimed at an efficient choice of driver nodes[16]), whereas the method proposed in Ref.[6] comes to diametrically opposed conclusions. Other network control schemes that, strictly speaking, do not belong to any of the categories above have also appeared recently (see e.g. Refs.[17,18]).

The most common assumption found in the network control literature is that the dynamical units are identical (which greatly simplifies both analytic and numerical treatments). Depending on the problem at hand, this may or may not be a drastic simplification: in physics it is sometimes a sensible approach, while it is clearly not in other disciplines, e.g. in the study of ecological systems. On the other hand, the assumption that the dynamics can be captured by linear ordinary differential equations is certainly not realistic for most applications, as it effectively bans limit-cycle or chaotic oscillations. Moreover, assuming linearly interacting units constitute also a severe limitation, as in most circumstances systems interact non-linearly: many-body gravitational and electrostatic problems in physics include, for instance, forces that are inversely proportional to the square distance between the interacting bodies, and a full treatment of solid-state and molecular systems frequently requires the incorporation of anharmonicities. In other areas of science nonlinear interactions are also the norm: in the modelling of ecosystems predator-prey couplings or competition for resources among species take the form of products of different populations, or more elaborated functional forms, see e.g. Ref.[19]. While in some cases a linear (first-order) approximation might be justified, in some others it may even be not possible at all, as the coupling functions might not be analytic (as in models of neurons, whose action potential is fired when the membrane potential reaches a threshold). Lastly, another common (and quite drastic) assumption is that of having identical coupling functions represented by either directed (unidirectional case) or undirected (bidirectional case) networks, whereas many systems (particularly in biology and social sciences) display in fact mixed couplings implying a combination of bidirectional and unidirectional links, with strengths and even functional forms that may vary from one link to another.

In this work, we introduce a technique for pinning control of networks that does not rely on any of these assumptions and is thus of wide applicability. The basic mechanism, previously introduced in a considerably more restrictive setting[6], consists in establishing unidirectional pinning actions from a copy of the networked system (in practice it may be an experimental recording, or just a simulation of the dynamics) to the system on which one wants to impose the dynamics of the copy. In the jargon of multilayer networks, this is an inter-layer synchronization problem[20,21]: while individual nodes on a layer (the original network) may not be synchronized to each other, each of them is synchronized to its counterpart on the other layer (the copy). By considering synthetic mixed networks of nonlinearly coupled chaotic oscillators, we first derive some general results on the correlations between the nodes that need to be pinned and their topological properties. In essence, we find that those nodes that are influential on the dynamics of many other nodes but are simultaneously less influenceable by the rest of the network are by far the most efficient in setting inter-layer synchronization already with a small number of actions. This analysis also serves to illustrate the method in a relatively simple setting, yet including several features that violate the assumptions used in the above-discussed references.

We then illustrate the applicability of our method to real-world networks by steering the dynamics of a trophic web containing 12 species toward a desired evolution. This allows us to obtain information on which are the appropriate species to target, i.e. which species are keystone in the environment, as well as the best strategies to impose given dynamics on an ecosystem. We discuss how these results can be used as a basis for adaptive management of ecosystems. Such a method can be effective to foster the implementation of adaptive ecosystem management as requested by the application of Malawi principles of the Convention for Biological Diversity, http://www.uni-kiel.de/ecology/users/fmueller/salzau2006/studentpages/Malawi_Principles/index.html. From a formal point of view, this is a challenging networked system to control: its units (the species) are governed by different nonlinear equations, they are nonlinearly coupled via different coupling functions, and the pattern of connections is highly asymmetric and irregular (including different functional forms). This implies a strong departure from the set of assumptions used in all previous methodologies. After almost two decades of intense activity in the field, it is fair to say that none of the previous methods, as far as we are aware, can be applied to such a control problem despite its great environmental interest.

Description of the method and application to networks of chaotic oscillators

We consider a two-layer network. One layer is the slave layer, which corresponds to the original network over which one wants to impose the desired dynamics (i.e. a given evolution compatible with the equations of motion). The other layer is the master layer, which is identical to the slave layer, but starts from a different initial condition (i.e. the one generating the specific desired dynamics towards which the state of the slave layer is to be steered), and evolves autonomously. In applications, the master layer may just be an experimental recording or a simulation of the original system—as long as it can be coupled to the slave network its physical nature is irrelevant. Our control method consists then in establishing directed inter-layer links from nodes in the master layer to their counterparts in the slave layer. Once they are established, these links remain in place as more nodes are connected in sequential control steps. At each step the selected node is the one whose pinning causes the most rapid approach towards inter-layer synchronization (i.e. the imposition of the evolution followed by the master layer on the slave layer). While the two layers have to be identical, the nodes (i.e. the dynamical units) and links (the coupling structure connecting the dynamical systems) on each layer can be completely different, as we will see below. This is thus a generalization of the method proposed in Ref.[6].

We illustrate our method by applying it to networks of identical chaotic oscillators, and leave the applicability to more challenging real-world systems to the next section. Specifically, we consider networks of nodes whose topology is that of a mixed random graph, i.e. containing both bidirectional and unidirectional links. These graphs are realizations of the configuration model[22] with the in-degree (i.e. the number of links pointing to a given node) and the out-degree (i.e. the number of links emanating from a given node) uniformly distributed in . Each node evolves autonomously in time as a chaotic Rössler oscillator, which we simply denote as , where and , with parameters , and . Nodes are coupled quadratically via their z variables, a nonlinear coupling form that was previously considered in Ref.[23].

Before the first control step is applied (prior to the creation of the first inter-layer connection) both master and slave layers evolve spontaneously as followswhere if there is a directed link from node j to node i, and is zero otherwise (for bidirectional links ). As we do not consider self-links, the diagonal terms vanish, i.e. , and the in-degree of node i is . The graph can thus be alternatively represented by the Laplacian matrix . The vector field governs the dynamics of node i, which would evolve autonomously (if uncoupled from its neighbors) simply as , and the parameter is the intra-layer coupling strength

When the control procedure starts, each node i in the master network keeps evolving according to the dynamics in Eq. (1), . In the slave layer dynamics, however, one has to consider an additional term which accounts for the inter-layer coupling from the master layer (without loss of generality, we here take a linear coupling through the y variable). One has

Here is a binary variable that is one if there is a link coupling node i in the master layer to node i in the slave layer (i.e. if the targeting procedure includes a pinning action from master to slave at node i) and is zero otherwise. The parameter is the inter-layer coupling strength. We emphasize that the coupling that is linear (in fact, diffusive) is the externally-imposed inter-layer coupling, which does not restrict in any way the form of the (intra-layer) couplings between the nodes of the system under study. Such diffusive inter-layer coupling is chosen as it is the simplest form that makes the inter-layer synchronization manifold into an invariant set of the dynamics (for a detailed mathematical treatment of invariant sets and related concepts, see e.g. Ref.[24]).

The results of applying our method to the network of Rössler chaotic oscillators are shown in Fig. 1. Two observables are employed to characterize the inter-layer synchronization between master and slave as more and more inter-layer links are established in successive targeting steps. One is the maximum Lyapunov exponent, , computed from the dynamics of the slave network linearized around that of the master network as in Ref.[6]. For a review of the theory and numerical computation of Lyapunov spectra, see e.g. Refs.[25,26]. The other observable is the synchronization error, which is the time average of the Euclidean distance in phase space (N is the number of nodes, m is the phase space dimensionality of the node dynamics—in our case , ) between the full state of the master layer and that of the slave layer, . In practice, T is finite, but orders of magnitude larger than the characteristic timescales of oscillation (thus the numerical convergence to an asymptotic value is guaranteed). In the top panel of Fig. 1, we show the maximum Lyapunov exponent as a function of the targeting step, which is seen to progressively decrease as more and more nodes are pinned. Analogous results in terms of the synchronization error are reported in the inset, which shows how the synchronization error becomes zero when becomes negative.

Figure 1

Controlling the dynamics of a mixed network with uniform and distributions comprising nonlinearly-coupled Rössler oscillators with intra-layer coupling and inter-layer coupling . (Top) Maximum Lyapunov exponent (main panel) and synchronization error (inset) as functions of the targeting step. (Bottom) Influence index of the node that is pinned at each targeting step. The curves are averages of 20 different network realizations. A 4th-order Runge-Kutta method with a step of 0.01 time units has been employed for the numerical integration of the systems of ordinary differential equations corresponding to each layer.

The maximum Lyapunov exponent is also used to identify the node to be targeted at each step: of all the nodes that remain unconnected to their counterparts in the other layer, the one that, when a master-slave connection is established, leads to the largest decrease in is targeted next. An exploration of possible correlations between the resulting targeting sequence (i.e. the ordered list of nodes that are targeted at successive steps) and local topological properties yields a remarkable correspondence between the targeting sequence position and the ranking of nodes in terms of their influence index , as shown in the lower panel of Fig. 1. This index is large when a node has a privileged position for influencing other nodes, while receiving very little influence from the rest of the network. No such correlations are observed for connectivity indices that are insensitive to the directionality of connections, such as , while correlations only based on or give considerably poorer results that those shown in the figure. Other measures of connection directionality that we have inspected, such as , show weaker correlations with the targeting sequence than the influence index does. While these results are based on networks with uniform distributions of and , which have been chosen precisely because a large variety of possible degree values is desirable, a strong correlation between the influence-index ranking and the targeting ranking is also observed for Barabási-Albert scale-free networks[27] and Erdös-Rényi random graph[28] topologies, as shown in Sect. A of the Supplementary Information.

This correlation is most clearly seen for small values of the intra-layer coupling strength, such as the value considered in Fig. 1. For larger values of , which make inter-layer synchronization possible with a very small number of steps, the correlation is less strong, while no obvious correlation between the targeting sequence and local topological properties are found for very large , see again Sect. A of the Supplementary Information. This might be related to the enhanced contribution of next-nearest neighbors and other relative distant nodes as the coupling strength is increased. Despite its limited range of validity, this correlation is nontheless remarkable, as it is very robust, and quite different from the situation observed in undirected networks, where the topological observable correlating with the targeting sequence is the degree[6]. On the other hand, there is an intriguing parallel between the correlation reported in the lower panel of Fig. 1 and the fact that, in undirected networks, nodes with a higher dynamic vulnerability are those with less influence from the rest of the network, followed by those that have the strongest ability to influence the rest[29]. In fact, both aspects of a node position are combined in the influence index in the case of directed or mixed networks.

Controlling ecological networks

We next apply our method to a model of a trophic web involving 12 species. This model describes the dynamics of a generic trophic web, including several categories of consumers such as top predators (, and ), mesopredators ( and ), several large herbivores (from to ), small herbivores ( and ) and also intermediate omnivourous consumers ( and ) which, in the real world, may also rely on predation and scavenging[30] (cf. Sect. B of the Supplementary Information for a full description of the model). The model represents a simplified food web inspired by holarctic ecosystems (see Ref.[30], and references therein). Controlling such a trophic web by means of only pinning a limited number of species, and/or implementing desired control policies for specific populations, are tasks of great societal relevance. As a matter of fact, there are many situations where wildlife agencies aim to control populations in order to reduce crop riding, depredation, as well as to control transmissible disease, to reduce extinction risks, or mitigate conflicts among stakeholders (e.g. conservationists, farmers, hunters).

The trophic web is viewed as a network where the species are the different nodes, and the links stand for the interactions among them. From the point of view of network control, this is a challenging model: each node (species) evolves autonomously following different population dynamics, the links (inter-species interactions) are also species-dependent and vary widely both in number, character (some are directed, some are undirected, and they are assigned different strengths) and in the mathematical form of the couplings, which are usually nonlinear. While the details of model are described in Sect. B of the Supplementary Information, we here briefly summarize its salient qualitative features. Each of the 12 species is described by a scalar that measures the population density at a given time. The coupling between species is given by nonlinear predator-prey response functions and competition-for-resource terms, which are proportional to products of the populations of the competing species. Moreover, there are logistic growth terms for each of the herbivores. A key feature of this model is the periodic nature of masting, which represents the quasi periodic production of forest fruits, such as acorns. Here masting acts as a forcing agent on the growth rate of one of the populations. The forcing makes the dynamics chaotic, with a (numerically calculated) maximum Lyapunov exponent . A representative sample of the highly irregular oscillations of the populations is shown in Sect. B of the Supplementary Information.

In order to apply the pinning procedure we need to construct a copy of the trophic web from which to establish unidirectional links to the original web. Below we clarify how this can be practically achieved in real ecosystems by monitoring populations along time, but for the time being we assume this to be a feasible task. We then apply pinning actions sequentially until, when a sufficiently high number of pinning actions have been established, the slave layer (the ecosystem of interest) follows the dynamics of the master system. As in the previous section, the key information is contained in the sequence of pinned nodes, as this reveals which are the species whose population one must preserve or modify in order to maintain a desired dynamics or disrupt an undesired one. In actual management, wildlife agencies are often requested to purchase action plans for removing or reintroducing individuals, to increase recruitment or reduce natural mortality by supplementary feeding, or to modify to some extent the natural dynamics of the system. This makes sense if the action provides long-lasting results, meaning that the ecosystem would attain a new equilibrium.

In the top panel of Fig. 2, we show the maximum Lyapunov exponent as a function of the targeting step, with labels in the horizontal axis indicating the targeted species at each step (see Sect. B of the Supplementary Information for a detailed description of these codes). It appears that here it is enough to target two species of herbivores, such as deer ( and ), a small herbivore, such as a species of hare (), and a mesopredator, such as jackals or foxes (), to control the network, to which it is probably necessary to add a control on one large predator, such as the wolf or the coguar if allowed by laws, . Analogous results in terms of the synchronization error are reported in the bottom panel of Fig. 2, which shows how the synchronization error becomes zero when becomes sufficiently negative. After pinning four species, is only midly negative and, given the relatively high value of the synchronization error, there appear to be regions in phase space where trajectories diverge despite the fact that the phase space averaging given by the Lyapunov exponent shows that the trend is overall converging. After the fifth pinning (i.e. the large predator), both measures unmistakably show that inter-layer synchronization has been achieved. If the sequence of species to be pinned is instead chosen at random, approximately twice as many species have to be pinned to achieve inter-layer synchronization for the same parameter values. The inter-layer parameter is here chosen to be , as considerably larger values lead to instabilities in the dynamics. In this case the intra-layer coupling is not a free parameter that one can modify at will as in the network of Rössler oscillators of the previous section—in fact, it is determined by the different parameters of the trophic web model and varies from link to link (see Sect. B of the Supplementary Information for more information). The length of the sequence of the species needed to achieve synchronization is largely dependent, but the ranking is robust across a range of .

Figure 2

Controlling the dynamics of a trophic web comprising 12 species with inter-layer coupling . (Top) Maximum Lyapunov exponent (main panel) as a function of the targeting step (codes indicating the species targeted at each step are described in Sect. B of the Supplementary Information). (Bottom) Synchronization error as a function of the targeting step. A 4th-order Runge–Kutta method with a step of 0.01 time units has been employed for the numerical integration of the systems of 12 ordinary differential equations corresponding to each layer.

In fact, the possibility of imposing on a system a given dynamics compatible with the equations of motion from an initial time onwards may not be always realistic. If the master layer were an exact physical replica of the original system (as could approximately be achieved in networks of nonlinear oscillating circuits or other technological systems) or a faithful simulation of its dynamics, one could hope to achieve that. In most cases, however, one can only expect to obtain (finite) recorded segments of the activity of the system in the form of a time series of some of its observables. Fortunately such a recording periodically repeated may suffice to maintain the system evolution within certain desired region of phase space. This is certainly true in the case of a trophic web, so we next illustrate with our model how the pinning method can be based on a short segment of recorded activity.

To do so, we simulate our trophic web system over a time window of tens of thousands of units and record the species populations every 5 time units. From this time series (i.e. our recording of “observational” data), we choose a time window of 325 time units where the populations oscillate relatively regularly within certain bounds that are of course species-dependent. Assuming these are the bounds that, for instance, on one side allow the conservation of a given species, but on the other reduce the amount of economic damage to crops, we take this to be our desired dynamics. Our master layer is simply this segment of recorded activity periodically repeated (i.e. when we come to the last sample, we start again from the first one, and so on), which we impose on the system by pinning a sufficiently high number of species (we choose this number to be 5, following the results in Fig. 2, which also determine the species chosen for the pinning actions). The results are shown in Fig. 3 for two species of the master layer (see black discs), which is just the periodically repeated time window of the recording, and the corresponding species in the slave layer (see red dotted line), which is the actual ecosystem in our model. We see how the pinning rapidly brings the slave system into the desired dynamical regime, despite the fact that it has started from an initial condition which is quite far from it. What we illustrate here for just two species for the sake of brevity, is similarly observed for the remaning ones. While the periodic repetition of the recorded segment introduces some discontinuities in the dynamics, the slave network does not take long (relatively speaking) to follow the master dynamics. In fact the length of the period where the trajectories are visibly different at each start of the cycle is related to the (inverse of the) Lyapunov exponent displayed in Fig. 2, and in general is expected to become smaller as more species are pinned (at least before the exponent saturates, as happens in the results shown in that figure eventually).

Figure 3

Evolution of two species of the trophic web (red dotted line) pinned to a periodically repeated segment of recorded population dynamics (black discs). (Top) Population of one species of deer as a function of time in the recording and under pinning of 5 species using the recording as master dynamics. (Bottom) Population of another species of deer under the same conditions.

To conclude this section, we should mention that the pinning strategy is also expected to work when not all of the nodes of the network are included in the time series (for instance, when only a subset of the interacting populations are tracked). See the relevant discussion in Ref.[6], which is also applicable in the present setting.

Discussion

We have proposed a quite general procedure for controlling the dynamics of complex networks of nonlinear dynamical units that are coupled in a nonlinear fashion, and possibly through mathematically different coupling functions, according to network schemes that may include unidirectional, bidirectional and weighted links. The method, which was first proposed in a much more restrictive setting[6], is based on the establishment of unidirectional inter-layer couplings (pinning actions) that reach the nodes of the original networked system from their counterparts in the other layer, which is an identical copy of the system, with the aim of imposing the dynamics of the copy on the network under consideration. While invoking a two-layer structure might sound quite remote from any practical application, in fact the copy can be simply a set of experimental data characterizing the state of the system across time, which may well be finite (in fact quite short) and even contain information on only some of the nodes.

We first illustrate the method on a network of nonlinearly-coupled (chaotic) Rössler oscillators coupled through a network topology that includes both directed and undirected links. The sequence of pinned nodes that is found to bring the system closer to the desired dynamics at each step shows a remarkable correlation with a ranking of the network nodes in terms of their influence index . This index measures the ability of a node to influence other nodes (as given by the, possibly weighted, out-degree ) normalized by the influence other nodes have on it (as given by the in-degree ).

We then move on to a study of a trophic web model containing 12 species inspired by European and North-American ecosystems. This case is much more challenging, as different nodes evolve according to different dynamic rules, and are coupled via nonlinear mathematical functions that depend on the pair of species involved. The method is shown to be perfectly applicable on such systems, and has the potential to yield valuable information on which species are key in maintaining a given ecosystem dynamics. Moreover, we illustrate the method by using as desired dynamics a short segment of activity where the populations evolve within stipulated bounds. Studying correlations between targeting sequences and topological properties of the kind observed in the network of Rössler oscillators will require further work in the case of this type of systems, where distinct links represent interactions that are mathematically completely different, and therefore it is hard to give precise operational meanings to measures such as the influence index.

In conclusion, we have presented a versatile method for steering networks toward desired dynamics. This method has shown to be valuable for unveiling correlations between node controllability and topological properties, which provide theoretical insight into the interplay of structure and function in complex systems. Most importantly, the method is of practical value for the control of systems that do not satisfy the highly idealized requirements of network control methods in the literature, systems whose dynamics may not be even fully understood nor amenable to realistic theoretical modelling. In fact, the main challenge facing the application of our method to environmental management is that the dynamics of ecosystems are always imperfectly known and in many cases scarcely documented. Thus it is necessary to joint this theoretical approach with adaptive management. Adaptive management is a method of learning by doing. Initially the model used as master shall be very rough, but with subsequent refinements based on management actions and monitoring, the method would improve becoming more and more effective. Our results open a very new avenue to apply adaptive management to nature conservation in the framework of the Convention of Biological Diversity as summarised by the Malawi Principles for the management of whole ecosystems.

Supplementary Information.