# Analysis and Visualization of High-Dimensional Dynamical Systems’ Phase Space Using a Network-Based Approach

Authors: Shane St Luce; Hiroki Sayama
Journal: Complexity
 (2022)

---

## Abstract

The concept of attractors is considered critical in the study of dynamical systems as they represent the set of states that a system gravitates toward. However, it is generally difficult to analyze attractors in complex systems due to multiple reasons including chaos, high-dimensionality, and stochasticity. This paper explores a novel approach to analyzing attractors in complex systems by utilizing networks to represent phase spaces. We accomplish this by discretizing phase space and defining node associations with attractors by finding sink strongly connected components (SSCCs) within these networks. Moreover, the network representation of phase space facilitates the use of well-established techniques of network analysis to study the phase space of a complex system. We show the latter by introducing a new node-based metric calledattractivity which can be used in conjunction with the SSCC as they are highly correlated. We demonstrate the proposed method by applying it to several chaotic dynamical systems and a large-scale agent-based social simulation model.


---

## Body

## 1. Introduction

Practical applications of dynamical systems analysis can be found in various subjects of study including transportation networks, food systems, cancer, and game theory [1–5]. A commonality stemming from analyses in these subjects is an attempt to understand the evolution and long-term states of these systems. One way to do this is by studying the trajectories of these systems towards attractors in order to understand their long-term behavior. However, several definitions exist for what constitutes an attractor, leading to several different ways to formulate an analysis [6]. Many of these ways involve finding solutions through analytical computation. This makes it fundamentally difficult to analyze systems with nonlinear dynamics [7]. The problems become amplified in higher dimensions as well with many variables. Additionally, little research has been done on analyzing the dynamics within an attractor consisting of many states.In this paper, we offer a different approach to analyzing complex dynamical systems based on two significant attributes associated with attractors, states that are reached after an extended period of time and the inability of a system to escape an attractor. Networks are then used as the basis of algorithmical analyzing of the trajectories in complex dynamical systems instead of analytical methods. This allows networks to be used both as a dimension reduction technique and for analyzing the attractors. This leads to two computational approaches to analyze long-term behavior of a system with a network-based phase space. One approach focuses on finding the set of states associated with the attractor of the system by finding sink strongly connected components in the network-based phase space. The other approach focuses on calculating a node-specific value that we refer to as “node attractivity.” This metric can be viewed as an indicator of how strongly a node in the network-based phase space is associated with an attractor in the system. We will then present the results from applying this method to ODE-based systems as well as show its applicability to high-dimensional systems. The results will also be visualized with details on how modifications can be incorporated without changing the essence of the method.


## 2. Literature Review

System dynamics are important for understanding the long-term behavior of a system. Long-term behavior that lead to a single steady state, or fixed point attractor, have an extensive research history that includes their properties [8] and use cases such as in ecological and neurophysiological models [9, 10]. This form of analysis bodes well for determining stability and the dynamics of the system in the region surrounding the fixed point. However, there are other types of attractors, including periodic, quasi-periodic, and chaotic cycles [11]. With the rise in popularity of complexity science, other forms of attractors have also been used to identify more complex structures and situations such as health sciences [12, 13], racial identity [14], thermodynamics [15], and human interaction. As these attractors are more complex, existing approaches for analyzing these attractors have also become more complex in comparison to fixed points. Along with the increased awareness of chaotic behavior in the real world, this has led to research into how we can extract information out of chaotic (strange) attractors [17, 18].The embrace of complexity has fueled adoption of networks for understanding complex structures and computer simulations for modeling complex phenomena, which in turn has been fueled by the rise of computing power [7]. Networks provide a nearly universal conceptual structure to be applied for any complex system. Some information often calculated for complex systems such as entropy and fractal dimension can be calculated for networks as well [19, 20]. However, there are also network-specific information such as closeness and harmonic centralities, and clusters and strong/weak components that can be gathered algorithmically.Attractors can be useful in a network representation of a system as well. In neural networks, attractors are identified through recurrent neural connections [21]. As an example, this definition has been used to study the dynamics of perceptual bistability without assumptions used in other models [10]. This combination of attractors and networks was mostly found in neuroscience for fixed point attractors.


## 3. Methods

Recognizing that attractors are a critical component of complex systems, we will provide two complementary methods for finding a attractor-like structure within a network representation of a system. These methods function regardless of attractor type and rely on algorithmic processes to avoid potential complexity from analytical computation. Both approaches involve first computing a network representation of the phase space. We will show that minimal adjustment is made if the phase space is continuous rather than discrete. Afterward, the first approach will lead to the discovery of a network substructure that we consider as the attractor represented in a network-embedding. The second approach will provide a node-specific metric that can be calculated with the expectation of distinguishing between attractor nodes and nonattractor nodes as time passes.


### 3.1. Network Representation of Phase Space

A continuous phase space of a system can be represented as a network by first dividing the phase space into nonoverlapping (disjoint) regions. Each region can then be represented by a unique node in a network. For discrete states in a phase space, each state can be directly represented by a node. Next, we can add directed edges between nodes where the weight represents the frequency that a transition occurred from a state represented by one node to a state represented by another. The resulting network can then be analyzed algorithmically for attractors. Note the resulting network is the state transition graph, or state diagram, of the system when the states are discrete and each node uniquely represents a single state.


### 3.2. Identifying Attractors in Networks

#### 3.2.1. Sink Strongly Connected Component

Consider two major properties of an attractor: (1) no exit exists once a system enters an attractor and (2) any position in an attractor should be approachable from any other position in the same attractor. The first property coincides with asink in a network. The second property can be interpreted as a connectivity requirement for an attractor. In a directed network, this interpretation coincides with the definition of a strongly connected component (SCC). A set of nodes that meets both of these properties in a network would then be a sink strongly connected component (SSCC) and can represent an attractor for the associated system. All of the nodes with a path leading to a SSCC would then be a sink weakly connected component of the network to correlate with the basin of attraction. Stochastic state transitions are also covered by the possibility of a node existing in multiple sink weakly connected components. The resulting analogous terms are shown in Table 1.Table 1
Analogous terminology for attractors in network-based phase space.
Traditional conceptsNetwork-based concepts(Observed) phase spaceFull networkAttractorSink strongly connected componentBasin of attractionSink weakly connected componentA simulation approach involves executing many simulations with different initial states. The observed network can be reconstructed after every time step across all simulations to incorporate all states and state transitions from the beginning of the simulation to the most recent time step. The result is a series of networks where each network has at least the same number of nodes as the previous and the sum of state transitions increases by the number of simulations. The SSCC may then reappear with another new transition soon after. To reduce the possibility of SSCCs incorrectly disappearing once identified, we introduce the concept of “SSCC efficiency” for a component and a corresponding “SSCC efficiency threshold.” The SSCC efficiency of a component is defined as the percentage of nodes in the component that are also part of an SCC. With this metric, we can define the algorithm to finding the SSCC of a network as the following:(1)
Find all strongly connected components (SCCs) in the network. Each SCC becomes a candidate component.(2)
Check if the candidate component is a sink.(a)
If the candidate component is a sink, then it is the SSCC.(b)
Otherwise,(i)
create a new candidate component from the set of candidate component nodes and their outgoing neighbors;(ii)
if the candidate component has an SSCC efficiency greater than the thresholdt, then go to step 2a with the new candidate component;(iii)
otherwise, an SSCC is not found.(3)
Go back to step 2 for next SCC.This approach loosens the definition of an SSCC. In our case, we found it sufficient to use a threshold of 95% when one attractor was expected and a lower threshold of 65% was used when expecting multiple attractors. A visual example of the overall approach is represented in Figure1.Figure 1
A visual example of the proposed SSCC method. SSCC efficiency threshold is defined ast and set to 0.5. SSCC efficiency is defined as s. It is updated throughout the process by dividing the number of nodes in the original SCC (shown in red) by the number of nodes in the candidate component (shown in red and orange). Since t=0.5, an SSCC was found after a couple iterations. If t was 0.7, then the process would have ended at 5 instead of 7. If t was 0.9, then it would have stopped by 3.


#### 3.2.2. Node Attractivity

The SSCC approach provides a discrete, Boolean answer for whether or not a node should be associated with an attractor. To provide a continuous answer between 0 and 1, we offer an alternative approach utilizing harmonic centrality.Harmonic centrality utilizes harmonic mean in place of the arithmetic mean as used for calculating closeness centrality, givendv,u represents the distance from v to u [22](1)Hrawu=∑v≠u1dv,uHnormalizedu=1N−1∗∑v≠u1dv,u.This has led to harmonic centrality being used in a few different ways, including calculating network topology scores in large social networks [23] and transportation networks [1] as well as a hybridization of node-weighted centralities [24]. It contains two major properties for us to consider: (1) if a path does not exist from node u to node v, then there will be no contribution to the harmonic centrality and (2) as the shortest path between nodes u and v shrinks, the contribution of this path to the harmonic centrality increases.Now, consider our network-based phase space. If we use the inverse of each edge weight for our harmonic centrality calculation, then this would result in an inverse weighted harmonic centrality (IWHC). These values for nodes in a network can also be easily normalized into the range of0,1. As transition frequency increases, corresponding IWHC measurements also increase before normalization for affected nodes. Nonattractor nodes would maintain lower IWHC, creating a clear distinction between attractor nodes and nonattractor nodes. We would expect these results to eventually match the SSCC approach. For this reason, we can also call this metric as “node attractivity.” The computational complexity calculations of both methods can be found in Appendix A.Both the SSCC and IWHC approaches were used in Ref. [25] for analyzing the El Farol bar problem. Attractor (SSCC) nodes and nonattractor (nonSSCC) nodes were identified based on the SSCC approach. We then used logistic regression to compare a series of node-based metrics with the SSCC results. Results confirmed that the SSCC nodes contained high node attractivity while nonSSCC nodes contained significantly low node attractivity. Figure 2 provides a side-by-side comparison of the results from the chosen node-based metrics. We can see that node attractivity was the best performing one. F1 score is used over accuracy, precision, or recall due to the wide range of attractor sizes. More details related to this analysis can be found in Appendix B.Figure 2
Box plot of F1 scores from the results of a trained logistic regression model on a given metric versus the SSCC.


## 3.1. Network Representation of Phase Space

A continuous phase space of a system can be represented as a network by first dividing the phase space into nonoverlapping (disjoint) regions. Each region can then be represented by a unique node in a network. For discrete states in a phase space, each state can be directly represented by a node. Next, we can add directed edges between nodes where the weight represents the frequency that a transition occurred from a state represented by one node to a state represented by another. The resulting network can then be analyzed algorithmically for attractors. Note the resulting network is the state transition graph, or state diagram, of the system when the states are discrete and each node uniquely represents a single state.


## 3.2. Identifying Attractors in Networks

### 3.2.1. Sink Strongly Connected Component

Consider two major properties of an attractor: (1) no exit exists once a system enters an attractor and (2) any position in an attractor should be approachable from any other position in the same attractor. The first property coincides with asink in a network. The second property can be interpreted as a connectivity requirement for an attractor. In a directed network, this interpretation coincides with the definition of a strongly connected component (SCC). A set of nodes that meets both of these properties in a network would then be a sink strongly connected component (SSCC) and can represent an attractor for the associated system. All of the nodes with a path leading to a SSCC would then be a sink weakly connected component of the network to correlate with the basin of attraction. Stochastic state transitions are also covered by the possibility of a node existing in multiple sink weakly connected components. The resulting analogous terms are shown in Table 1.Table 1
Analogous terminology for attractors in network-based phase space.
Traditional conceptsNetwork-based concepts(Observed) phase spaceFull networkAttractorSink strongly connected componentBasin of attractionSink weakly connected componentA simulation approach involves executing many simulations with different initial states. The observed network can be reconstructed after every time step across all simulations to incorporate all states and state transitions from the beginning of the simulation to the most recent time step. The result is a series of networks where each network has at least the same number of nodes as the previous and the sum of state transitions increases by the number of simulations. The SSCC may then reappear with another new transition soon after. To reduce the possibility of SSCCs incorrectly disappearing once identified, we introduce the concept of “SSCC efficiency” for a component and a corresponding “SSCC efficiency threshold.” The SSCC efficiency of a component is defined as the percentage of nodes in the component that are also part of an SCC. With this metric, we can define the algorithm to finding the SSCC of a network as the following:(1)
Find all strongly connected components (SCCs) in the network. Each SCC becomes a candidate component.(2)
Check if the candidate component is a sink.(a)
If the candidate component is a sink, then it is the SSCC.(b)
Otherwise,(i)
create a new candidate component from the set of candidate component nodes and their outgoing neighbors;(ii)
if the candidate component has an SSCC efficiency greater than the thresholdt, then go to step 2a with the new candidate component;(iii)
otherwise, an SSCC is not found.(3)
Go back to step 2 for next SCC.This approach loosens the definition of an SSCC. In our case, we found it sufficient to use a threshold of 95% when one attractor was expected and a lower threshold of 65% was used when expecting multiple attractors. A visual example of the overall approach is represented in Figure1.Figure 1
A visual example of the proposed SSCC method. SSCC efficiency threshold is defined ast and set to 0.5. SSCC efficiency is defined as s. It is updated throughout the process by dividing the number of nodes in the original SCC (shown in red) by the number of nodes in the candidate component (shown in red and orange). Since t=0.5, an SSCC was found after a couple iterations. If t was 0.7, then the process would have ended at 5 instead of 7. If t was 0.9, then it would have stopped by 3.


### 3.2.2. Node Attractivity

The SSCC approach provides a discrete, Boolean answer for whether or not a node should be associated with an attractor. To provide a continuous answer between 0 and 1, we offer an alternative approach utilizing harmonic centrality.Harmonic centrality utilizes harmonic mean in place of the arithmetic mean as used for calculating closeness centrality, givendv,u represents the distance from v to u [22](1)Hrawu=∑v≠u1dv,uHnormalizedu=1N−1∗∑v≠u1dv,u.This has led to harmonic centrality being used in a few different ways, including calculating network topology scores in large social networks [23] and transportation networks [1] as well as a hybridization of node-weighted centralities [24]. It contains two major properties for us to consider: (1) if a path does not exist from node u to node v, then there will be no contribution to the harmonic centrality and (2) as the shortest path between nodes u and v shrinks, the contribution of this path to the harmonic centrality increases.Now, consider our network-based phase space. If we use the inverse of each edge weight for our harmonic centrality calculation, then this would result in an inverse weighted harmonic centrality (IWHC). These values for nodes in a network can also be easily normalized into the range of0,1. As transition frequency increases, corresponding IWHC measurements also increase before normalization for affected nodes. Nonattractor nodes would maintain lower IWHC, creating a clear distinction between attractor nodes and nonattractor nodes. We would expect these results to eventually match the SSCC approach. For this reason, we can also call this metric as “node attractivity.” The computational complexity calculations of both methods can be found in Appendix A.Both the SSCC and IWHC approaches were used in Ref. [25] for analyzing the El Farol bar problem. Attractor (SSCC) nodes and nonattractor (nonSSCC) nodes were identified based on the SSCC approach. We then used logistic regression to compare a series of node-based metrics with the SSCC results. Results confirmed that the SSCC nodes contained high node attractivity while nonSSCC nodes contained significantly low node attractivity. Figure 2 provides a side-by-side comparison of the results from the chosen node-based metrics. We can see that node attractivity was the best performing one. F1 score is used over accuracy, precision, or recall due to the wide range of attractor sizes. More details related to this analysis can be found in Appendix B.Figure 2
Box plot of F1 scores from the results of a trained logistic regression model on a given metric versus the SSCC.


## 3.2.1. Sink Strongly Connected Component

Consider two major properties of an attractor: (1) no exit exists once a system enters an attractor and (2) any position in an attractor should be approachable from any other position in the same attractor. The first property coincides with asink in a network. The second property can be interpreted as a connectivity requirement for an attractor. In a directed network, this interpretation coincides with the definition of a strongly connected component (SCC). A set of nodes that meets both of these properties in a network would then be a sink strongly connected component (SSCC) and can represent an attractor for the associated system. All of the nodes with a path leading to a SSCC would then be a sink weakly connected component of the network to correlate with the basin of attraction. Stochastic state transitions are also covered by the possibility of a node existing in multiple sink weakly connected components. The resulting analogous terms are shown in Table 1.Table 1
Analogous terminology for attractors in network-based phase space.
Traditional conceptsNetwork-based concepts(Observed) phase spaceFull networkAttractorSink strongly connected componentBasin of attractionSink weakly connected componentA simulation approach involves executing many simulations with different initial states. The observed network can be reconstructed after every time step across all simulations to incorporate all states and state transitions from the beginning of the simulation to the most recent time step. The result is a series of networks where each network has at least the same number of nodes as the previous and the sum of state transitions increases by the number of simulations. The SSCC may then reappear with another new transition soon after. To reduce the possibility of SSCCs incorrectly disappearing once identified, we introduce the concept of “SSCC efficiency” for a component and a corresponding “SSCC efficiency threshold.” The SSCC efficiency of a component is defined as the percentage of nodes in the component that are also part of an SCC. With this metric, we can define the algorithm to finding the SSCC of a network as the following:(1)
Find all strongly connected components (SCCs) in the network. Each SCC becomes a candidate component.(2)
Check if the candidate component is a sink.(a)
If the candidate component is a sink, then it is the SSCC.(b)
Otherwise,(i)
create a new candidate component from the set of candidate component nodes and their outgoing neighbors;(ii)
if the candidate component has an SSCC efficiency greater than the thresholdt, then go to step 2a with the new candidate component;(iii)
otherwise, an SSCC is not found.(3)
Go back to step 2 for next SCC.This approach loosens the definition of an SSCC. In our case, we found it sufficient to use a threshold of 95% when one attractor was expected and a lower threshold of 65% was used when expecting multiple attractors. A visual example of the overall approach is represented in Figure1.Figure 1
A visual example of the proposed SSCC method. SSCC efficiency threshold is defined ast and set to 0.5. SSCC efficiency is defined as s. It is updated throughout the process by dividing the number of nodes in the original SCC (shown in red) by the number of nodes in the candidate component (shown in red and orange). Since t=0.5, an SSCC was found after a couple iterations. If t was 0.7, then the process would have ended at 5 instead of 7. If t was 0.9, then it would have stopped by 3.


## 3.2.2. Node Attractivity

The SSCC approach provides a discrete, Boolean answer for whether or not a node should be associated with an attractor. To provide a continuous answer between 0 and 1, we offer an alternative approach utilizing harmonic centrality.Harmonic centrality utilizes harmonic mean in place of the arithmetic mean as used for calculating closeness centrality, givendv,u represents the distance from v to u [22](1)Hrawu=∑v≠u1dv,uHnormalizedu=1N−1∗∑v≠u1dv,u.This has led to harmonic centrality being used in a few different ways, including calculating network topology scores in large social networks [23] and transportation networks [1] as well as a hybridization of node-weighted centralities [24]. It contains two major properties for us to consider: (1) if a path does not exist from node u to node v, then there will be no contribution to the harmonic centrality and (2) as the shortest path between nodes u and v shrinks, the contribution of this path to the harmonic centrality increases.Now, consider our network-based phase space. If we use the inverse of each edge weight for our harmonic centrality calculation, then this would result in an inverse weighted harmonic centrality (IWHC). These values for nodes in a network can also be easily normalized into the range of0,1. As transition frequency increases, corresponding IWHC measurements also increase before normalization for affected nodes. Nonattractor nodes would maintain lower IWHC, creating a clear distinction between attractor nodes and nonattractor nodes. We would expect these results to eventually match the SSCC approach. For this reason, we can also call this metric as “node attractivity.” The computational complexity calculations of both methods can be found in Appendix A.Both the SSCC and IWHC approaches were used in Ref. [25] for analyzing the El Farol bar problem. Attractor (SSCC) nodes and nonattractor (nonSSCC) nodes were identified based on the SSCC approach. We then used logistic regression to compare a series of node-based metrics with the SSCC results. Results confirmed that the SSCC nodes contained high node attractivity while nonSSCC nodes contained significantly low node attractivity. Figure 2 provides a side-by-side comparison of the results from the chosen node-based metrics. We can see that node attractivity was the best performing one. F1 score is used over accuracy, precision, or recall due to the wide range of attractor sizes. More details related to this analysis can be found in Appendix B.Figure 2
Box plot of F1 scores from the results of a trained logistic regression model on a given metric versus the SSCC.


## 4. ODE System Results

For systems of differential equations, we can use a similar set-up and simulate system movement across Euclidean space, with a series of simulations starting from different initial states. Since these phase spaces are boundless, we will set a boundary for the initial positions of the simulations based on known positioning of the systems’ attractors.The phase space is also continuous, so we will discretize it by establishing a unit size and rounding down every element in a coordinate to a multiple of the unit size. This will determine the identity of the corresponding node in the network and the corresponding Euclidean space associated with it. This unit size will be called the “granularity” of our nodes. For example, if a new state is found at10.4,11.3,2.9 and the granularity is set to 1, the corresponding node in the network will be at 10,11,2. This also means every node coincides with an integer coordinate representing an n-cube with a length of 1 unit along each of n dimensions. If the granularity was set to 2, then the corresponding node in the network will be at 10,10,2, and every node would coincide with an even integer coordinate and represent an n-cube with a length of 2 units along each dimension. After determining the region for initial positions and the granularity to analyze, we can then proceed as we did for the El Farol bar problem simulation, and then visualize the resulting network in both 3D Euclidean space and/or in a unit-less 2D network-embedding space. All ODE-based systems use a step size of 0.01 during calculations. Metrics of the resulting networks can be found in Appendix C.


### 4.1. Lorenz System

The Lorenz system is a system of equations presented by meteorologist Lorenz as a simplified model of atmospheric convection flow [26]. A recognizable butterfly-like shaped attractor arises from this system in Euclidean space. The system of equations is given in Equation (1). For our experiment, we use the same values for the parameters used in the original paper, ρ=28, σ=10 and β=8/3(2)x˙=−σy−x,y˙=ρx−y−xz,z˙=xy−βz.For the Lorenz system, we used the region defined below in Equation (2) and generated initial coordinates at every even integer coordinate with a granularity of 2 for node positions. We then ran the simulation for 20 steps with no transient period(3)−20≤x≤20,−30≤y≤30,0≤z≤50.Figure3(a) shows the results of our proposed method in Euclidean space. Figure 3(b) shows the same results as Figure 3(a) except visualized as a network using the ForceAtlas 2 layout algorithm. We can call this visualization of the attractor as the “Lorenz Attractor Network.”Figure 3
Lorenz attractor visualized with similar node coloring based on X coordinate from low (blue) to high (yellow): (a) in Euclidean space and (b) as a network resulting from proposed approach.
(a)(b)Nodes in both figures are colored similarly. Those with highX values are yellow while nodes with low X values are blue/purple. This is highlighted the most by the tips of each wing of the attractor as they are distinctly in different colors. The wing on the positive X (X+) side is yellow while the wing on the negative X (X-) side is more or less purple.TheX+ wing and X-wing are clearly distinguishable in the Lorenz attractor network based on color where the X+ wing contains yellow nodes and the X-wing contains blue/purple nodes. This visualization successfully projects the attractor into two dimensions by replacing the Euclidean space-based spatial relationship between nodes with a layout based solely on system dynamics. Nodes representing states that are visited in close succession are in relatively close proximity to each other as well.


### 4.2. Yu’s Chaotic System with Multiple Attractors

Next, we will test with a system of equations known to produce multiple attractors. Equation (3) is the system of differential equations for Yu’s chaotic system [27] using a value of 2.5 for a. The region for initial states is defined below in Equation (4) with initial positions at every integer coordinate and a granularity of 0.5 is used. The simulation ran for 40 steps with no transient period(4)x˙=−x+0.5xz+yz,y˙=ay−1.2xz,z˙=xy−6z,a=2.5,(5)−18≤x≤18,−15≤y≤15,−10≤z≤10.Figure4 shows the full-generated SSCC network in 3D Euclidean space and in the network-embedding space. Both visuals show two attractors. In Figure 4(a), we can see that the attractors are circling about the z-axis. Two components are seen in the network visualization in Figure 4(b). The components have incidental overlap but are disconnected. We can also tell from Figure 4(b) that the node attractivity is higher inside the attractor, which is not as clearly seen in Figure 4(a). It was noticed that a large granularity with a long-running simulation can lead to merged attractors and basins due to nodes straddling between multiple basins.Figure 4
Attracting portions of Yu’s chaotic system visualized with similar node coloring based on node attractivity from low (blue/purple) to high (green/yellow): (a) in Euclidean space and (b) as a network resulting from proposed approach.
(a)(b)


### 4.3. Four-Dimensional Chaotic System

Lastly, we will test our method against a system with more than three dimensions. Equation (5) is the system of differential equations found in [28] containing different kinds of attractors. However with our results, we argue that the different kinds of attractors are all part of one multidimensional attractor.We set the parameter values to the same ones used in the aforementioned paper,a=6, b=11, and c=5, and define the region as shown in Equation (6). Since the region is rather large due to the additional dimension, the initial coordinates are created at coordinates with X, Y, Z, and W values at multiples of 4. The resulting network uses a granularity of 1, just as we did for the Lorenz network. The simulation ran for 70 steps with no transient period(6)x˙=ay−x,y˙=xz+w,z˙=b−xy,w˙=yz−cw,(7)−16≤x≤16,−20≤y≤20,−20≤z≤20,−20≤w≤20.Figure5 shows the attractor nodes of the identified SSCC with nodes colored based on attractivity from low (blue) to high (yellow). We found the dynamics which suggested a single two-winged attractor with reflection symmetry within the system. Visualizing the entire network including nonattractor nodes with the same ForceAtlas 2 layout algorithm led to a different apparent SSCC subgraph. We discuss this finding and more related to this system in Appendix E.Figure 5
SSCC of the four-dimensional chaotic system as a network. Nodes colors are based on node attractivity from blue (low) to yellow (high).


## 4.1. Lorenz System

The Lorenz system is a system of equations presented by meteorologist Lorenz as a simplified model of atmospheric convection flow [26]. A recognizable butterfly-like shaped attractor arises from this system in Euclidean space. The system of equations is given in Equation (1). For our experiment, we use the same values for the parameters used in the original paper, ρ=28, σ=10 and β=8/3(2)x˙=−σy−x,y˙=ρx−y−xz,z˙=xy−βz.For the Lorenz system, we used the region defined below in Equation (2) and generated initial coordinates at every even integer coordinate with a granularity of 2 for node positions. We then ran the simulation for 20 steps with no transient period(3)−20≤x≤20,−30≤y≤30,0≤z≤50.Figure3(a) shows the results of our proposed method in Euclidean space. Figure 3(b) shows the same results as Figure 3(a) except visualized as a network using the ForceAtlas 2 layout algorithm. We can call this visualization of the attractor as the “Lorenz Attractor Network.”Figure 3
Lorenz attractor visualized with similar node coloring based on X coordinate from low (blue) to high (yellow): (a) in Euclidean space and (b) as a network resulting from proposed approach.
(a)(b)Nodes in both figures are colored similarly. Those with highX values are yellow while nodes with low X values are blue/purple. This is highlighted the most by the tips of each wing of the attractor as they are distinctly in different colors. The wing on the positive X (X+) side is yellow while the wing on the negative X (X-) side is more or less purple.TheX+ wing and X-wing are clearly distinguishable in the Lorenz attractor network based on color where the X+ wing contains yellow nodes and the X-wing contains blue/purple nodes. This visualization successfully projects the attractor into two dimensions by replacing the Euclidean space-based spatial relationship between nodes with a layout based solely on system dynamics. Nodes representing states that are visited in close succession are in relatively close proximity to each other as well.


## 4.2. Yu’s Chaotic System with Multiple Attractors

Next, we will test with a system of equations known to produce multiple attractors. Equation (3) is the system of differential equations for Yu’s chaotic system [27] using a value of 2.5 for a. The region for initial states is defined below in Equation (4) with initial positions at every integer coordinate and a granularity of 0.5 is used. The simulation ran for 40 steps with no transient period(4)x˙=−x+0.5xz+yz,y˙=ay−1.2xz,z˙=xy−6z,a=2.5,(5)−18≤x≤18,−15≤y≤15,−10≤z≤10.Figure4 shows the full-generated SSCC network in 3D Euclidean space and in the network-embedding space. Both visuals show two attractors. In Figure 4(a), we can see that the attractors are circling about the z-axis. Two components are seen in the network visualization in Figure 4(b). The components have incidental overlap but are disconnected. We can also tell from Figure 4(b) that the node attractivity is higher inside the attractor, which is not as clearly seen in Figure 4(a). It was noticed that a large granularity with a long-running simulation can lead to merged attractors and basins due to nodes straddling between multiple basins.Figure 4
Attracting portions of Yu’s chaotic system visualized with similar node coloring based on node attractivity from low (blue/purple) to high (green/yellow): (a) in Euclidean space and (b) as a network resulting from proposed approach.
(a)(b)


## 4.3. Four-Dimensional Chaotic System

Lastly, we will test our method against a system with more than three dimensions. Equation (5) is the system of differential equations found in [28] containing different kinds of attractors. However with our results, we argue that the different kinds of attractors are all part of one multidimensional attractor.We set the parameter values to the same ones used in the aforementioned paper,a=6, b=11, and c=5, and define the region as shown in Equation (6). Since the region is rather large due to the additional dimension, the initial coordinates are created at coordinates with X, Y, Z, and W values at multiples of 4. The resulting network uses a granularity of 1, just as we did for the Lorenz network. The simulation ran for 70 steps with no transient period(6)x˙=ay−x,y˙=xz+w,z˙=b−xy,w˙=yz−cw,(7)−16≤x≤16,−20≤y≤20,−20≤z≤20,−20≤w≤20.Figure5 shows the attractor nodes of the identified SSCC with nodes colored based on attractivity from low (blue) to high (yellow). We found the dynamics which suggested a single two-winged attractor with reflection symmetry within the system. Visualizing the entire network including nonattractor nodes with the same ForceAtlas 2 layout algorithm led to a different apparent SSCC subgraph. We discuss this finding and more related to this system in Appendix E.Figure 5
SSCC of the four-dimensional chaotic system as a network. Nodes colors are based on node attractivity from blue (low) to yellow (high).


## 5. Conclusion

In this paper, we provided a novel approach for discovering attractors in complex systems by utilizing a network representation of the phase space. Our approach focuses on two key properties of an attractor. First, an attractor should be a component containing no (or exceptionally rare) outgoing edges. Second, each node should be approachable from every other node within the identified component. This definition leads us to the concept of a sink strongly connected component (SSCC) as the network-based entity associated with an attractor in a complex system. With this baseline, we proposed a novel metric named “node attractivity.” Although the computation is more time-consuming, node attractivity provides more details about the system behavior than the SSCC approach. However, these two approaches complement each other well as “node attractivity” aligns with SSCC results over time as shown by our logistic regression results.The generality of both approaches were seen by analyzing simulations of the El Farol bar problem as well as three systems of differential equations with varying dimensionality and number of attractors. The results from the El Farol bar problem showed how node attractivity performed against in degree, weighted in degree, and weighted out degree in predicting the SSCC. Applying these approaches to the Lorenz system revealed the behavioral tendencies of the Lorenz attractor in a unit-agnostic space, exposing a cyclic nature in a tornado-like structure in the “Lorenz attractor network.” We also showed that these approaches can successfully identify the number of attractors in a system given a set of initial positions and can clarify structural components of the attractor for systems of multiple dimensions by applying this approach to Yu’s chaotic system and a four-dimensional chaotic system.Future possibilities for research include testing this approach for systems with more than four dimensions. It is possible that the network becomes visually more difficult to understand, but computational analysis may remain consistent as the nature of this approach is to capture all of the behavior of system trajectories within a network structure, regardless of dimensions. We can also find possible research avenues by focusing just on node attractivity. One key benefit that can arise from node attractivity is to develop a proxy of the probability of long-term system state. It may be sufficient to sum-normalize the node attractivity since node attractivity is already nonnegative.Further research can be performed to see how well this works, or if an intermediate function can be included to model the long-term system state. We previously mentioned potential issues that can appear when using our approach on systems with multiple attractors due to straddling nodes between basins. It would be worthwhile to explore ways to mitigate this issue, such as condensing the space of initial coordinates while adding a transitory period. Using different condensed initial regions would likely improve identifying attractor nodes and nodes straddling between basins since they would presumably appear in multiple SSCCs separately. Lastly, we can take these approaches and apply them to other domains such as transportation systems, voting behavior, or other decision-making systems in search of any commonalities in the SSCC or node attractivity of different systems.



---

*Source: 3937475-2022-06-29.xml*