###### Abstract

Many natural systems are organized as networks, in which the nodes interact in a time-dependent fashion. The object of our study is to relate connectivity to the temporal behavior of a network in which the nodes are (real or complex) logistic maps, coupled according to a connectivity scheme that obeys certain constrains, but also incorporates random aspects. We investigate in particular the relationship between the system architecture and possible dynamics. In the current paper we focus on establishing the framework, terminology and pertinent questions for low-dimensional networks. A subsequent paper will further address the relationship between hardwiring and dynamics in high-dimensional networks.

For networks of both complex and real node-maps, we define extensions of the Julia and Mandelbrot sets traditionally defined in the context of single map iterations. For three different model networks, we use a combination of analytical and numerical tools to illustrate how the system behavior (measured via topological properties of the Julia sets) changes when perturbing the underlying adjacency graph. We differentiate between the effects on dynamics of different perturbations that directly modulate network connectivity: increasing/decreasing edge weights, and altering edge configuration by adding, deleting or moving edges. We discuss the implications of considering a rigorous extension of Fatou-Julia theory known to apply for iterations of single maps, to iterations of ensembles of maps coupled as nodes in a network.

Real and complex behavior for networks of coupled logistic maps

Anca Rǎdulescu, Ariel Pignatelli

Department of Mathematics, SUNY New Paltz, NY 12561

Department of Mechanical Engineering, SUNY New Paltz, NY 12561

## 1 Introduction

### 1.1 Network architecture as a system parameter

Because many natural systems are organized as networks, in which the nodes (be they cells, individuals, populations or web servers) interact in a time-dependent fashion – the study of networks has been an important focus in recent research. One of the particular points of interest has been the question of how the hardwired structure of a network (its underlying graph) affects its function, for example in the context of optimal information storage or transmission between nodes along time. It has been hypothesized that there are two key conditions for optimal function in such networks: a well-balanced adjacency matrix (the underlying graph should appropriately combine robust features and random edges) and well-balanced connection strengths, driving optimal dynamics in the system. However, only recently has mathematics started to study rigorously (through a combined graph theoretical and dynamic approach) the effects of configuration patterns on the efficiency of network function, by applying graph theoretical measures of segregation (clustering coefficient, motifs, modularity, rich clubs), integration (path length, efficiency) and influence (node degree, centrality). Various studies have been investigating the sensitivity of a system’s temporal behavior to removing/adding nodes or edges at different places in the network structure, and have tried to relate these patterns to applications to natural networks.

Brain functioning is one of the most intensely studied contexts which requires our understanding of the tight inter-connections between system architecture and dynamics. The brain is organized as a “dynamic network,” self-interacting in a time-dependent fashion at multiple spacial and temporal scales, to deliver an optimal range for biological functioning. The way in which these modules are wired together in large networks that control complex cognition and behavior is one of the great scientific challenges of the 21st century, currently being addressed by large-scale research collaborations, such as the Human Connectome Project. Graph theoretical studies of empirical empirical data support certain generic topological properties of brain architecture, such as modularity, small-worldness, the existence of hubs and “rich clubs” [bullmore2009complex, sporns2011non, sporns2002graph].

In order to explain how connectivity patterns may affect the system’s dynamics (e.g., in the context of stability and synchronization in networks of coupled neural populations), and thus the observed behavior, a lot of effort has been thus invested towards formal modeling approaches, using a combination of analytical and numerical methods from nonlinear dynamics and graph theory, in both biophysical models [gray2009stability] and simplified systems [siri2007effects]. These analyses revealed a rich range of potential dynamic regimes and transitions [brunel2000dynamics], shown to depend as much on the coupling parameters of the network as on the arrangement of the excitatory and inhibitory connections [gray2009stability]. The construction of a realistic, data-compatible computational model has been subsequently found to present many difficulties that transcend the existing methods from nonlinear dynamics, and may in fact require: (1) new analysis and book-keeping methods and (2) a new framework that would naturally encompass the rich phenomena intrinsic to these systems – both of which aspects are central to our proposed work.

In a paper with Dr. Verduzco-Flores [radulescu2015nonlinear], one of the authors of this paper first explored the idea of having network connectivity as a bifurcation parameter for the ensemble dynamics in a continuous time system of coupled differential equations. We used configuration dependent phase spaces and our probabilistic extension of bifurcation diagrams in the parameter space to investigate the relationship between classes of system architectures and classes of their possible dynamics, and we observed the robustness of the coupled dynamics to certain changes in the network architecture and its vulnerability to others. As expected, when translating connectivity patterns to network dynamics, the main difficulties were raised by the combination of graph complexity and the system’s intractable dynamic richness.

In order to break down and better understand this dependence, we started to investigate it in simpler theoretical models, where one may more easily identify and pair specific structural patterns to their effects on dynamics. The logistic family is historically perhaps the most-studied family of maps in nonlinear dynamics, whose behavior is by now relatively well understood. Therefore, we started by looking in particular at how dynamic behavior depends on connectivity in networks with simple logistic nodes. This paper focuses on definitions, concepts and observations in low-dimensional networks. Future work will address large networks, and different classes of maps.

Dynamic networks with discrete nodes and the dependence of their behavior on connectivity parameters have been previously described in several contexts over the past two decades. For example, in an early paper, Wang considered a simple neural network of only two excitatory/inhibitory neurons, and analyzed it as a parameterized family of two-dimensional maps, proving existence of period-doubling to chaos and strange attractors in the network [wang1991period]. Masolle, Attay et al. have found that, in networks of delay-coupled logistic maps, synchronization regimes and formation of anti-phase clusters depend on coupling strength [masoller2011complex] and on the edge topology (characterized by the spectrum of the graph Laplacian) [atay2004delays]. Yu has constructed and studied a network wherein the undirected edges symbolize the nodes’ relation of adjacency in an integer sequence obtained from the logistic mapping and the top integral function [yu2013logistic].

In our present work, we focus on investigating, in the context of networked maps, extensions of the Julia and Mandelbrot sets traditionally defined for single map iterations. For three different model networks, we use a combination of analytical and numerical tools to illustrate how the system behavior (measured via topological properties of the Julia sets) changes when perturbing the underlying adjacency graph. We differentiate between the effects on dynamics of different perturbations that directly modulate network connectivity: increasing/decreasing edge weights, and altering edge configuration by adding, deleting or moving edges. We discuss the implications of considering a rigorous extension of Fatou-Julia theory known to apply for iterations of single maps, to iterations of ensembles of maps coupled as nodes in a network.

### 1.2 Networking logistic maps

The logistic map is historically perhaps the best-known family of maps in nonlinear dynamics. Iterations of one single quadratic function have been studied starting in the early 19th century, with the work of Fatou and Julia.

The prisoner set of a map is defined as the set of all points in the complex dynamic plane, whose orbits are bounded. The escape set of a complex map is the set of all points whose orbits are unbounded. The Julia set of is defined as their common boundary . The filled Julia set is the union of prisoner points with their boundary .

For polynomial maps, it has been shown that the connectivity of a map’s Julia set is tightly related to the structure of its critical orbits (i.e., the orbits of the map’s critical points). Due to extensive work spanning almost one century, from Julia [julia1918memoire] and Fatou [fatou1919equations] until recent developments [branner1992iteration, qiu2009proof], we now have the following:

Fatou-Julia Theorem. For a polynomial with at least one critical orbit unbounded, the Julia set is totally disconnected if and only if all the bounded critical orbits are aperiodic.

For a single iterated logistic map[carleson1993complex, devaney2006criterion], the Fatou-Julia Theorem implies that the Julia set is either totally connected, for values of in the Mandelbrot set (i.e., if the orbit of the critical point 0 is bounded), or totally disconnected, for values of outside of the Mandelbrot set (i.e., if the orbit of the critical point 0 is unbounded). In previous work, the authors showed that this dichotomy breaks in the case of random iterations of two maps [pignatelli2016]. In our current work, we focus on extensions for networked logistic maps. Although Julia and Mandelbrot sets have been studied somewhat in connection with coupled systems [isaeva2010phenomena], none of the existing work seems to address the basic problems of how these sets can be defined for networks of maps, how different aspects of the network hardwiring affect the topology of these sets and whether there is any Fatou-Julia type result in this context.

These are some of the questions addressed in this paper, which is organized as follows: In Section 2, we introduce definitions of our network setup, as well as of the extensions of Mandelbrot and Julia sets that we will be studying. In order to illustrate some basic ideas and concepts, we concentrate on three examples of 3-dimensional networks, which differ from each other in edge distribution, and whose connectivity strengths are allow to vary. In Section 3, we focus on the behavior of these 3-dimensional models when we consider the nodes as complex iterated variables. We analyze the similarities and differences between node-wise behavior in each case, and we investigate the topological perturbations in one-dimensional complex slices of the Mandelbrot and Julia sets, as the connectivity changes from one model to the next, through intermediate stages. In Section LABEL:real_maps, we address the same questions for real logistic nodes, with the advantage of being able to visualize the entire network Mandelbrot and Julia sets, as 3-dimensional real objects. In both sections, we conjecture weaker versions of the Fatou-Julia theorem, connecting points in the Mandelbrot set with connectivity properties of the corresponding Julia sets. Finally, in Section LABEL:discussion, we interpret our results both mathematically and in the larger context of network sciences. We also briefly preview future work on high-dimensional networks and on networks with adaptable nodes and edges.

## 2 Our models of networked logistic maps

We consider a set of nodes coupled according to the edges of an oriented graph, with adjacency matrix (on which one may impose additional structural conditions, related to edge density or distribution). In isolation, each node , , functions as a discrete nonlinear map , changing at each iteration as . When coupled as a network with adjacency , each node will also receive contributions through the incoming edges from the adjacent nodes. Throughout this paper, we will consider an additive rule of combining these contributions, for a couple of reasons: first, summing weighted incoming inputs is one simple, yet mathematically nontrivial way to introduce the cross talk between nodes; second, summing weighted inputs inside a nonlinear integrating function is reminiscent of certain mechanisms studied in the natural sciences (such as the integrate and fire neural mechanism studied in our previous work in the context of coupled dynamics). The coupled system will then have the following general form:

where are the weights along the adjacency edges. One may view this system simply as an iteration of an -dimensional map , with (in the case of real-valued nodes), or respectively (in the case of complex-valued nodes). The new and exciting aspect that we are proposing in our work is to study the dependence of the coupled dynamics on the parameters, in particular on the coupling scheme (adjacency matrix) – viewed itself as a system parameter. To fix these ideas, we focused first on defining these questions and proposing hypotheses for the case of quadratic node-dynamics. The logistic family is one of the most studied family of maps in the context of both real and complex dynamics of a single variable. It was also the subject of our previous modeling work on random iterations.

In this paper in particular, we will work with quadratic node-maps, with their traditional parametrization , with and for the complex case and and for the real case. The network variable will be called respectively in the case of complex nodes, and in the case of real nodes. We consider both the particular case of identical quadratic maps (equal values), and the general case of different maps attached to the nodes throughout the network. In both cases, we aim to study the asymptotic behavior of iterated node-wise orbits, as well as of the -dimensional orbits (which we will call multi-orbits). As in the classical theory of Fatou and Julia, we will investigate when orbits escape to infinity or remain bounded, and how much of this information is encoded in the critical multi-orbit of the system.

For the following definitions, fix the network (i.e., fix the adjacency and the edge weights ). To avoid redundancy, we give definitions for the complex case, but they can be formulated similarly for real maps:

###### Definition 2.1.

For a fixed parameter , we call the filled multi-Julia set of the network, the locus of which produce a bounded multi-orbit in . We call the filled uni-Julia set the locus of so that produces a bounded multi-orbit. The multi-Julia set (or the multi-J set) of the network is defined as the boundary in of the filled multi-Julia set. Similarly, one defines the uni-Julia set (or uni-J set) of the network as the boundary in of its filled counterparts.

###### Definition 2.2.

We define the multi-Mandelbrot set (or the multi-M set) of the network the parameter locus of for which the multi-orbit of the critical point is bounded in . We call the equi-Mandelbrot set (or the equi-M set) of the network, the locus of for which the critical multi-orbit is bounded for equi-parameter . We call the th node equi-M set the locus such that the component of the multi-orbit of corresponding to the th node remains bounded in .

We study, using a combination of analytical and numerical methods, how the structure of the Julia and Mandelbrot sets varies under perturbations of the node-wise dynamics (i.e., under changes of the quadratic multi- parameter ) and under perturbations of the coupling scheme (i.e., of the adjacency matrix and of the coupling weights ). In this paper, we start with investigating these questions in small (3-dimensional) networks, with specific adjacency configurations. In a subsequent paper, we will move to investigate how similar phenomena may be quantified and studied analytically and numerically in high-dimensional networks. In both cases, we are interested in particular in observing differences in the effects on dynamics of three different aspects of the network architecture: (1) increasing/decreasing edge weights, (2) increasing/decreasing edge density, (3) altering edge configuration by adding, deleting or moving edges.

While a desired objective would be to obtain general results for all network sizes (since many natural networks are large), we start by studying simple, low dimensional systems. In this study, we focus on simple networks formed of three nodes, connected by different network geometries and edge weights. To fix our ideas, we will follow and illustrate three structures in particular (also see Figure 1): (1) Two input nodes and are self driven by quadratic maps, and the single output node is driven symmetrically by the two input nodes; additionally communicates with via an edge of variable weight , which can take both positive and negative values. We will call this the simple dual model. (2) In addition to the simple dual scheme, the output node is also self-driven, i.e. there is a self-loop on of weight (which can be positive or negative). We will call this the self-drive model. (3) In addition to the self-driven model, there is also feedback from the output node into the node , via a new edge of variable weight . We will call this the feedback model. Unless specified, edges have positive unit weight. Notice that the same effect as negative feed-forward edges from and into can be obtained by changing the sign of , etc. The three connectivity models we chose to study and compare are described by the equations below:

Simple dual model:

Self-drive model:

Feedback model:

For a fixed multi-parameter for example, one can see all three systems as generated by a network map , , , defined as , for any .

We try to classify and understand the effects that coupling changes have on the topology of multi-J and multi-M sets for both complex and real networked maps. We don’t expect all classical topology results on the Julia and Mandelbrot sets for single maps (e.g., Fatou-Julia theorem, or connectivity of the Mandelbrot set) to carry out for networks of coupled maps. However, since the topology of the full sets in is somewhat harder to inspect, we study as a first step their equi-slices and node-wise equi-slices, which are objects in .

We will track and compare in particular the differences between the three models, but also the geometric and topological changes produced on the equi-slices within each one model for different values of the parameters , and . None of these results, however, can be directly extrapolated to similar conclusions on the full sets. To offer some insight into the latter, we study the multi-M and multi-J sets in the context of real maps, for which there objects can be visualized in .

## 3 Complex coupled maps

### 3.1 Equi-Mandelbrot sets

A first intuitive question is when the nodes of the network have similar behavior, and whether if one node-wise orbit is bounded, the others will remain bounded. This relationship is trivial to establish in some cases, such as for example in the simple dual model with independent input nodes (i.e., ). Indeed, in this model, for any fixed , the origin’s orbit in under can be described as:

The projection of the orbit in any of the three components only depends on the previous states of and , and these three sequences are simultaneously bounded in , hence the node-specific equi-Mandelbrot sets are all identical with the traditional Mandelbrot set. Some basic connections between node-wise equi-M sets in each of the three models are stated below. We will prove these incrementally (recall that the dual model is a particular case of self-drive for , and the self-drive is a particular case of feedback model with ).

###### Proposition 3.1.

In the simple dual model, the node-wise equi-M sets for the nodes and are identical subsets of the traditional Mandelbrot set (which is the equi-M set for node ).

###### Proof.

The simple case was already discussed. We will now assume . Suppose the critical orbit for node is bounded by a radius , that is , for all . Hence (omitting the subscript for simplicity):

It follows that:

Hence if the orbit is bounded, then the orbit of is bounded. This applies in particular for the critical orbit, showing that the equi-M set for if a subset of the equi-M set for .

We will next show that, for the simple dual model, corresponding orbits of and are simultaneously bounded. For instance, suppose that an orbit is bounded by . It follows, as shown above, that the corresponding orbit is bounded by some . Then:

Hence the orbit is bounded. The converse is similar, showing that the and equi-M sets are always identical subsets of the equi-M set in the simple dual model. In Figure 2a, we show that these are generally strict subsets, and that a non-symmetric communication can introduce significant differences between the traditional Mandelbrot set of the independent node and the equi-M subsets for and . For example, it is not hard to show that, for (illustrated in Figure 2a), the point belongs to the Mandelbrot set of (the critical orbit has period three), but not to the equi-M set of and .

∎

An additional self-drive applied to the output node changes the balance of inputs to , in the following sense:

###### Proposition 3.2.

In the self-drive model, the node-wise equi-M sets of and remain subsets of the standard Mandelbrot set, but the equi-M set of is strictly contained in the equi-M set of (Figure 2b).