Feigenbaum Graphs: A Complex Network Perspective of Chaos

The recently formulated theory of horizontal visibility graphs transforms time series into graphs and allows the possibility of studying dynamical systems through the characterization of their associated networks. This method leads to a natural graph-theoretical description of nonlinear systems with qualities in the spirit of symbolic dynamics. We support our claim via the case study of the period-doubling and band-splitting attractor cascades that characterize unimodal maps. We provide a universal analytical description of this classic scenario in terms of the horizontal visibility graphs associated with the dynamics within the attractors, that we call Feigenbaum graphs, independent of map nonlinearity or other particulars. We derive exact results for their degree distribution and related quantities, recast them in the context of the renormalization group and find that its fixed points coincide with those of network entropy optimization. Furthermore, we show that the network entropy mimics the Lyapunov exponent of the map independently of its sign, hinting at a Pesin-like relation equally valid out of chaos.


Introduction
We expose a remarkable relationship between nonlinear dynamical systems and complex networks by means of the horizontal visibility (HV) algorithm [1][2][3] that transforms time series into graphs. In low-dimensional dissipative systems chaotic motion develops out of regular motion in a small number of ways or routes, and amongst which the period-doubling bifurcation cascade or Feigenbaum scenario is perhaps the better known and most famous mechanism [4,5]. This route to chaos appears an infinite number of times amongst the family of attractors spawned by unimodal maps within the so-called periodic windows that interrupt stretches of chaotic attractors. In the opposite direction, a route out of chaos accompanies each period-doubling cascade by a chaotic band-splitting cascade, and their shared bifurcation accumulation points form transitions between order and chaos that are known to possess universal properties [4][5][6]. Lowdimensional maps have been extensively studied from a purely theoretical perspective, but systems with many degrees of freedom used to study diverse problems in physics, biology, chemistry, engineering, and social science, are known to display lowdimensional dynamics [7].
The horizontal visibility (HV) algorithm converts the information stored in a time series into a network, setting the nature of the dynamical system into a different context that requires complex network tools [8][9][10][11][12] to extract its properties. This approach belongs to an emerging corpus of methods that map series to networks (see for instance [2,[13][14][15][16][17] or a recent review [18]). Relevant information can be obtained through the family of visibility methods, including the characterization of fractal behavior [19] or the discrimination between random and chaotic series [1,20], and it finds increasing applications in separate fields, from geophysics [21], to finance [22] or physiology [23]. Here we offer a distinct view of the Feigenbaum scenario through the specific HV formalism, and provide a complete set of graphs, which we call Feigenbaum graphs, that encode the dynamics of all stationary trajectories of unimodal maps. We first characterize their topology via the order-of-visit and self-affinity properties of the maps. Additionally, a matching renormalization group (RG) procedure leads, via its flows, to or from network fixedpoints to a comprehensive view of the entire family of attractors. Furthermore, the optimization of the entropy obtained from the degree distribution coincides with the RG fixed points and reproduces the essential features of the map's Lyapunov exponent independently of its sign. A general observation is that the HV algorithm extracts only universal elements of the dynamics, free of the peculiarities of the individual unimodal map, but also of universality classes characterized by the degree of nonlinearity. Therefore all the results presented in this work, while referring to the specific Logistic map for illustrative reasons apply to any unimodal map.

Model: Feigenbaum graphs
The HV graph [1] associated with a given time series fx i g i~1,:::,N of N real data is constructed as follows: First, a node i is assigned to each datum x i , and then two nodes i and j are connected if the corresponding data fulfill the criterion x i ,x j wx n for all n such that ivnvj. Let us now focus on the Logistic map [4] defined by the quadratic difference equation x tz1~f (x t )~mx t (1{x t ) where x t [½0,1 and the control parameter m[½0,4. According to the HV algorithm, a time series generated by the Logistic map for a specific value of m (after an initial transient of approach to the attractor) is converted into a Feigenbaum graph (see figure 1). Notice that this is a well-defined subclass of HV graphs where consecutive nodes of degree k~2, that is, consecutive data with the same value, do not appear, what is actually the case for series extracted from maps (besides the trivial case of a constant series). While for a period T there are in principle several possible periodic orbits, and therefore the set of associated Feigenbaum graphs is degenerate, it can be proved of all these Feigenbaum graphs fulfill k k(T)~4(1{ 1 2T ) and d d (T)~1 3T respectively, yielding a linear relation d d( k k)~(4{ k k)=6 that is corroborated in the inset of figure 1. Observe that aperiodic series (T??) reach the upper bound mean degree k k~4.

Results
A deep-seated feature of the period-doubling cascade is that the order in which the positions of a periodic attractor are visited is universal [24], the same for all unimodal maps. This ordering turns out to be a decisive property in the derivation of the structure of the Feigenbaum graphs. See figure 2 where we plot the graphs for a family of attractors of increasing period T~2 n , that is, for increasing values of mvm ? . This basic pattern also leads to the expression for their associated degree distributions, , k~2,4,6,:::,2n, and zero for k odd or kw2(nz1). At the accumulation point m ? the period diverges (n??) and the distribution is exponential for all even values of the degree, , k~2,4,6,:::, ð2Þ and zero for k odd. Observe that these relations are independent of the order of the map's nonlinearity: the HV algorithm sifts out every detail of the dynamics except for the basic storyline. We turn next to the period-doubling bifurcation cascade of chaotic bands that takes place as m decreases from m~4 towards m ? . For the largest value of the control parameter, at m~4, the attractor is fully chaotic and occupies the entire interval ½0,1 (see figure 1). This is the first chaotic band n~0 at its maximum amplitude. As m decreases in value within m ? vmv4 bandnarrowing and successive band-splittings [4][5][6]24] occur. In general, after n reverse bifurcations the phase space is partitioned in 2 n disconnected chaotic bands, which are self-affine copies of the first chaotic band [25]. The values of m at which the bands split are called Misiurewicz points [24], and their location converges to the accumulation point m ? for n??. Significantly, while in the chaotic zone orbits are aperiodic, for reasons of continuity they visit each of the 2 n chaotic bands in the same order as positions are visited in the attractors of period T~2 n [24]. In figure 3 we have plotted the Feigenbaum graphs generated through chaotic time series at different values of m that correspond to an increasing number of reverse bifurcations. Since chaotic bands do not overlap, one can derive the following degree distribution for a Feigenbaum graph after n chaotic-band reverse bifurcations by Figure 2. Periodic Feigenbaum graphs for mvm ? . The sequence of graphs associated to periodic attractors with increasing period T~2 n undergoing a period-doubling cascade. The pattern that occurs for increasing values of the period is related to the universal ordering with which an orbit visits the points of the attractor. Observe that the hierarchical self-similarity of these graphs requires that the graph for n{1 is a subgraph of that for n. doi:10.1371/journal.pone.0022411.g002 using only the universal order of visits , k~2,4,6,:::,2n, and zero for k~3,5,7,:::,2nz1. We note that this time the degree distribution retains some dependence on the specific value of m, concretely, for those nodes with degree k §2(nz1), all of which belong to the top chaotic band (labelled with red links in figure 3).
The HV algorithm filters out chaotic motion within all bands except for that taking place in the top band whose contribution decreases as n?? and appears coarse-grained in the cumulative distribution P m (n,k §2(nz1)). As would be expected, at the accumulation point m ? we recover the exponential degree distribution (equation 2), i.e. lim n?? P m (n,k)~P(?,k). Before proceeding to interpret these findings via the consideration of renormalization group (RG) arguments, we recall that the Feigenbaum tree shows a rich self-affine structure: for mwm ? periodic windows of initial period m undergo successive perioddoubling bifurcations with new accumulation points m ? (m) that appear interwoven with chaotic attractors. These cascades are selfaffine copies of the fundamental one. The process of reverse bifurcations also evidences this self-affine structure, such that each accumulation point is the limit of a chaotic-band reverse bifurcation cascade. Accordingly, we label G(m,n) the Feigenbaum graph associated with a periodic series of period T~m : 2 n , that is, a graph obtained from an attractor within window of initial period m after n period-doubling bifurcations. In the same fashion, G m (n,m) is associated with a chaotic attractor composed by m : 2 n bands (that is, after n chaotic band reverse bifurcations of m initial chaotic bands). Therefore, graphs depicted in figures 2 and 3 correspond to G(1,n) and G m (1,n) respectively and for the first accumulation In order to recast previous findings in the context of the renormalization group, let us define an RG operation R on a graph as the coarse-graining of every couple of adjacent nodes where one of them has degree k~2 into a block node that inherits the links of the previous two nodes (see figure 4.a). This is a real-space RG transformation on the Feigenbaum graph [26], dissimilar from recently suggested box-covering complex network renormalization schemes [27,28,29]. This scheme turns out to be equivalent for mvm ? to the construction of an HV graph from the composed map f (2) instead of the original f , in correspondence to the original Feigenbaum renormalization procedure [30,6]. We first note that RfG(1,n)g~G(1,n{1), thus, an iteration of this process yields an RG flow that converges to the (1st) trivial fixed point R (n) fG(1,n)g~G(1,0):G 0~R fG 0 g. This is the stable fixed point of the RG flow Vmvm ? . We note that there is only one relevant variable in our RG scheme, represented by the reduced control parameter Dm~m ? {m, hence, to identify a nontrivial fixed point we set Dm~0 or equivalently n??, where the structure of the Feigenbaum graph turns to be completely self-similar under R. Therefore we conclude that G(1,?):G ? is the nontrivial fixed point of the RG flow, RfG ? g~G ? . In connection with this, let P t (k) be the degree distribution of a generic Feigenbaum graph G t in the period-doubling cascade after t iterations of R, and point out that the RG operation, RfG t g~G tz1 , implies a recurrence relation (1{P t (2))P tz1 (k)~P t (kz2), whose fixed point coincides with the degree distribution found in equation 2. This confirms that the nontrivial fixed point of the flow is indeed G ? .
Next, under the same RG transformation, the self-affine structure of the family of attractors yields RfG m (1,n)g~G m (1,n{1), generating a RG flow that converges to the Feigenbaum graph Figure 3. Aperiodic Feigenbaum graphs for mwm ? . A sequence of graphs associated with chaotic series after n chaotic-band reverse bifurcations, starting at m~4 for n~0, when the attractor extends along a single band and the degree distribution does not present any regularity (red links). For nw0 the phase space is partitioned in 2 n disconnected chaotic bands and the n-th self-affine image of m~4 is the n-th Misiurewicz point m 2 n{1 {2 n . In all cases, the orbit visits each chaotic band in the same order as in the periodic region mvm ? . This order of visits induces an ordered structure in the graphs (black links) analogous to that found for the period-doubling cascade. doi:10.1371/journal.pone.0022411.g003 is coarse-grained with one of its neighbors (indistinctively) into a block node that inherits the links of both nodes. This process coarse-grains every associated to the 1st chaotic band, R (n) fG m (1,n)g~G m (1,0). Repeated application of R breaks temporal correlations in the series, and the RG flow leads to a 2nd trivial fixed point R (?) fG m (1,0)g~G rand~R fG rand g, where G rand is the HV graph generated by a purely uncorrelated random process. This graph has a universal degree distribution P(k)~(1=3)(2=3) k{2 , independent of the random process underlying probability density (see [1,20]).
Finally, let us consider the RG flow inside a given periodic window of initial period m. As the renormalization process addresses nodes with degree k~2, the initial applications of R only change the core structure of the graph associated with the specific value m (see figure 4.b for an illustrative example). The RG flow will therefore converge to the 1st trivial fixed point via the initial path R (p) fG(m,n)g~G(1,n), with pƒm, whereas it converges to the 2nd trivial fixed point for G m (m,n) via R (p) fG m (m,n)g~G m (1,n). In the limit of n?? the RG flow proceeds towards the nontrivial fixed point via the path R (p) fG(m,?)g~G (1,?). Incidentally, extending the definition of the reduced control parameter to Dm(m)~m ? (m){m, the family of accumulation points is found at Dm(m)~0. A complete schematic representation of the RG flows can be seen in figure 4.c.
Interestingly, and at odds with standard RG applications to (asymptotically) scale-invariant systems, we find that invariance at Dm~0 is associated in this instance to an exponential (rather than power-law) function of the observables, concretely, that for the degree distribution. The reason is straightforward: R is not a conformal transformation (i:e: a scale operation) as in the typical RG, but rather, a translation procedure. The associated invariant functions are therefore non homogeneous (with the property g(ax)~bg(x)), but exponential (with the property g(xza)~cg(x)).
Finally, we derive, via optimization of an entropic functional for the Feigenbaum graphs, all the RG flow directions and fixed points directly from the information contained in the degree distribution. Amongst the graph theoretical entropies that have been proposed we employ here the Shannon entropy of the degree distribution P(k), that is h~{ P ? k~2 P(k) log P(k). By making use of the Maximum Entropy formalism, it is easy to prove that the degree distribution P(k) that maximizes h is exactly P(k)~(1=3)(2=3) k{2 , which corresponds to the distribution for the 2nd trivial fixed point of the RG flow G rand . Alternatively, with the incorporation of the additional constraint that allows only even values for the degree (the topological restriction for Feigenbaum graphs G(1,n)), entropy maximization yields a degree distribution that coincides with equation 2, which corresponds to the nontrivial fixed point of the RG flow G ? . Lastly, the degree distribution that minimizes h trivially corresponds to G 0 , i.e. the 1st trivial fixed point of the RG flow. Remarkably, these results indicate that the fixed-point structure of the RG flow are obtained via optimization of the entropy for the entire family of networks, supporting a suggested connection between RG theory and the principle of Maximum Entropy [31].
The network entropy h can be calculated exactly for G(1,n) (mvm ? or T~2 n ), yielding h(n)~log 4 : (1{2 {n ). Because increments of entropy are only due to the occurrence of bifurcations h increases with m in a step-wise way, and reaches asymptotically the value h(?)~log 4 at the accumulation point m ? . For Feigenbaum graphs G m (1,n) (in the chaotic region), in general h cannot be derived exactly since the precise shape of P(k) is unknown (albeit the asymptotic shape is also exponential [20]). Yet, the main feature of h can be determined along the chaoticband splitting process, as each reverse bifurcation generates two self-affine copies of each chaotic band. Accordingly, the decrease of entropy associated with this reverse bifurcation process can be described as h m (n)~log 4zh m (0)=2 n , where the entropy h m (n) after n reverse bifurcations can be described in terms of the entropy associated with the first chaotic band h m (0). In figure 1 we observe how the chaotic-band reverse bifurcation process takes place in the chaotic region from right to left, and therefore leads in this case to a decrease of entropy with an asymptotic value of log 4 for n?? at the accumulation point. These results suggest that the graph entropy behaves qualitatively as the map's Lyapunov exponent l, with the peculiarity of having a shift of log 4, as confirmed in figure 5. This unexpected qualitative agreement is reasonable in the chaotic region in view of the Pesin theorem [5], that relates the positive Lyapunov exponents of a map with its Kolmogorov-Sinai entropy (akin to a topological entropy) that for unimodal maps reads h KS~l , Vlw0, since h can be understood as a proxy for h KS . Unexpectedly, this qualitative agreement seems also valid in the periodic windows (lv0), since the graph entropy is positive and varies with the value of the associated (negative) Lyapunov exponent even though h KS~0 , hinting at a Pesin-like relation valid also out of chaos which deserves further investigation. The agreement between both quantities lead us to conclude that the Feigenbaum graphs capture not only the period-doubling route to chaos in a universal way, but also inherits the main feature of chaos, i.e. sensitivity to initial conditions.

Discussion
In summary, we have shown that the horizontal visibility theory combines power with straightforwardness as a tool for the analytical study of nonlinear dynamics. As an illustration we have established how the families of periodic and chaotic attractor bifurcation cascades of unimodal maps transform into families of networks with scale-invariant limiting forms, whose characterization can be deduced from two basic and universal properties of unimodal maps: ordering of consecutive positions in the attractors and self-affinity. Further, we have demonstrated that these networks and their associated degree distributions comply with renormalization group and maximum entropy principles, filtering out irrelevant variables and finding fixed-point networks which are independent of the map's nonlinearity. The entire Feigenbaum scenario is therefore fully described. The potential of the theory for revealing new information is indicated by the ability of the network entropy to emulate the Lyapunov exponent for both periodic and chaotic attractors. Extensions of this approach to other complex behavior, such as dynamical complexity associated to vanishing Lyapunov exponents, intermittency, quasiperiodic node with degree k~2 present at each renormalization step. (b) Example of an iterated renormalization process in a sample Feigenbaum graph at a periodic window with initial period m~9 after n~2 period-doubling bifurcations (an orbit of period T~m : 2 n~3 6). (c) RG flow diagram, where m identifies the periodic window that is initiated with period m and ñ designates the order of the bifurcation, ñ~nz1 for period-doubling bifurcations and ñ~{(nz1) for reverse bifurcations. Dm(m)~m ? (m){m denotes the reduced control parameter of the map, and m ? (m) is the location of the accumulation point of the bifurcation cascades within that window. Feigenbaum graphs associated with periodic series (Dm(m)w0, ñ w0) converge to G(1,0):G 0 under the RG, whereas those associated with aperiodic ones (Dm(m)v0, ñ v0) converge to G rand . The accumulation point m ? :m ? (1) corresponds to the unstable (nontrivial) fixed point G(1,?):G ? of the RG flow, which is nonetheless approached through the critical manifold of graphs G(m,?) at the accumulation points m ? (m). In summary, the nontrivial fixed point of the RG flow is only reached via the family of the accumulation points, otherwise the flow converges to trivial fixed points for periodic or chaotic regions. doi:10.1371/journal.pone.0022411.g004 routes to chaos, etc., are still open questions. Finally, observe that in symbolic dynamics [32] one usually defines a phase-space partition (Markov partition) in order to create a symbolic representation of the dynamics. This partition tiles phase space in a non-overlapping manner: every value of the series has a univocally associated symbol. While a Feigenbaum graph also symbolizes the series data (incidentally, without the need of defining an ad hoc partition), each series datum is not associated univocally to a symbol (the degree of the node): this symbol is a function, in principle, of the complete series, and incorporates global information. Furthermore, note that besides the symbolization that converts a time series into a series of node degrees, the Feigenbaum graphs also store the connectivity pattern amongst nodes -i.e. the topological structure of the graph. On this respect, the possible connections of HV theory with kneading theory [33] and symbolic dynamics [32] are of special interest. Figure 5. Horizontal visibility network entropy h and Lyapunov exponent l for the Logistic map. We plot the numerical values of h and l for 3:5vmv4 (the numerical step is dm~5 : 10 {4 and in each case the processed time series have a size of 2 12 data). The inset reproduces the same data but with a rescaled entropy h{ log (4). The surprisingly good match between both quantities is reminiscent of the Pesin identity (see text). Unexpectedly, the Lyapunov exponent within the periodic windows (lv0 inside the chaotic region) is also well captured by h. doi:10.1371/journal.pone.0022411.g005