Attractors in Neural Network Circuits: Beauty and Chaos

is one thing in common between memories, oscillating chemical reactions and double pendulums? All these systems have a basin of attraction for possible states, like a magnet that draws the system towards certain trajectories. Complex systems with multiple inputs usually evolve over time, generating intricate and sometimes chaotic behaviors. Attractors represent the long-term behavioral pattern of dynamical systems — a pattern to which a system converges over time regardless of its initial conditions. 

Neural networks have become ubiquitous in our current Artificial Intelligence era, typically serving as powerful tools for representation extraction and pattern recognition. However, these systems can also be viewed through another fascinating lens: as dynamical systems that evolve and converge to a manifold of states over time. When implemented with feedback loops, even simple neural networks can produce strikingly beautiful attractors, ranging from limit cycles to chaotic structures.

Neural Networks as Dynamical Systems

While neural networks in general sense are most commonly known for embedding extraction tasks, they can also be viewed as dynamical systems. A dynamical system describes how points in a state space evolve over time according to a fixed set of rules or forces. In the context of neural networks, the state space consists of the activation patterns of neurons, and the evolution rule is determined by the network’s weights, biases, activation functions, and other tricks.

Traditional NNs are optimized via gradient descent to find its endstate of convergence. However, when we introduce feedback — connecting the output back to the input — the network becomes a recurrent system with a different kind of temporal dynamic. These dynamics can exhibit a wide range of behaviors, from simple convergence to a fixed point to complex chaotic patterns.

Understanding Attractors

An attractor is a set of states toward which a system tends to evolve from a wide variety of starting conditions. Once a system reaches an attractor, it remains within that set of states unless perturbed by an external force. Attractors are indeed deeply involved in forming memories [1], oscillating chemical reactions [2], and other nonlinear dynamical systems. 

Types of Attractors

Dynamical Systems can exhibit several types of attractors, each with distinct characteristics:

• Point Attractors: the simplest form, where the system converges to a single fixed point regardless of starting conditions. This represents a stable equilibrium state.
• Limit Cycles: the system settles into a repeating periodic orbit, forming a closed loop in phase space. This represents oscillatory behavior with a fixed period.
• Toroidal (Quasiperiodic) Attractors: the system follows trajectories that wind around a donut-like structure in the phase space. Unlike limit cycles, these trajectories never really repeat but they remain bound to a specific region.
• Strange (Chaotic) Attractors: characterized by aperiodic behavior that never repeats exactly yet remains bounded within a finite region of phase space. These attractors exhibit sensitive dependence on initial conditions, where a tiny difference will introduce significant consequences over time — a hallmark of chaos. Think butterfly effect.

Setup

In the following section, we will dive deeper into an example of a very simple NN architecture capable of said behavior, and demonstrate some pretty examples. We will touch on Lyapunov exponents, and provide implementation for those who wish to experiment with generating their own Neural Network attractor art (and not in the generative AI sense).

We will use a grossly simplified one-layer NN with a feedback loop. The architecture consists of:

1. Input Layer:
   • Array of size D (here 16-32) inputs
   • We will unconventionally label them as y₁, y₂, y₃, …, y_D to highlight that these are mapped from the outputs
   • Acts as a shift register that stores previous outputs
2. Hidden Layer:
   • Contains N neurons (here fewer than D, ~4-8)
   • We will label them x₁, x₂, …, x_N
   • tanh() activation is applied for squashing
3. Output Layer
   • Single output neuron (y₀)
   • Combines the hidden layer outputs with biases — typically, we use biases to offset outputs by adding them; here, we used them for scaling, so they are factually an array of weights
4. Connections:
   • Input to Hidden: Weight matrix w[i,j] (randomly initialized between -1 and 1)
   • Hidden to Output: Bias weights b[i] (randomly initialized between 0 and s)
5. Feedback Loop:
   • The output y₀ is fed back to the input layer, creating a dynamic map
   • Acts as a shift register (y₁ = previous y₀, y₂ = previous y₁, etc.)
   • This feedback is what creates the dynamical system behavior
6. Key Formulas:
   • Hidden layer: u[i] = Σ(w[i,j] * y[j]); x[i] = tanh(u[i])
   • Output: y₀ = Σ(b[i] * x[i])

The critical aspects that make this network generate attractors:

• The feedback loop turns a simple feedforward network into a dynamical system
• The nonlinear activation function (tanh) enables complex behaviors
• The random weight initialization (controlled by the random seed) creates different attractor patterns
• The scaling factor s affects the dynamics of the system and can push it into chaotic regimes

In order to investigate how prone the system is to chaos, we will calculate the Lyapunov exponents for different sets of parameters. Lyapunov exponent is a measure of the instability of a dynamical system…

\[delta Z(t)| approx e^lambda t |delta (Z(0))|\]

\[lambda = n_t sum_k=0^n_t-1 ln frac\]

…where n_t​ is a number of time steps, Δy_k ​is a distance between the states y(x_i) and y(x_i+ϵ) at a point in time; ΔZ(0) represents an initial infinitesimal (very small) separation between two nearby starting points, and ΔZ(t) is the separation after time t. For stable systems converging to a fixed point or a stable attractor this parameter is less than 0, for unstable (diverging, and, therefore, chaotic systems) it is greater than 0.

Let’s code it up! We will only use NumPy and default Python libraries for the implementation.

import numpy as np
from typing import Tuple, List, Optional


class NeuralAttractor:
    """
    
    N : int
        Number of neurons in the hidden layer
    D : int
        Dimension of the input vector
    s : float
        Scaling factor for the output

    """
    
    def __init__(self, N: int = 4, D: int = 16, s: float = 0.75, seed: Optional[int] = 
None):
        self.N = N
        self.D = D
        self.s = s
        
        if seed is not None:
            np.random.seed(seed)
        
        # Initialize weights and biases
        self.w = 2.0 * np.random.random((N, D)) - 1.0  # Uniform in [-1, 1]
        self.b = s * np.random.random(N)  # Uniform in [0, s]
        
        # Initialize state vector structures
        self.x = np.zeros(N)  # Neuron states
        self.y = np.zeros(D)  # Input vector

We initialize the NeuralAttractor class with some basic parameters — number of neurons in the hidden layer, number of elements in the input array, scaling factor for the output, and random seed. We proceed to initialize the weights and biases randomly, and x and y states. These weights and biases will not be optimized — they will stay put, no gradient descent this time.

    def reset(self, init_value: float = 0.001):
        """Reset the network state to initial conditions."""
        self.x = np.ones(self.N) * init_value
        self.y = np.zeros(self.D)
        
    def iterate(self) -> np.ndarray:
        """
        Perform one iteration of the network and return the neuron outputs.
        
        """
        # Calculate the output y0
        y0 = np.sum(self.b * self.x)
        
        # Shift the input vector
        self.y[1:] = self.y[:-1]
        self.y[0] = y0
        
        # Calculate the neuron inputs and apply activation fn
        for i in range(self.N):
            u = np.sum(self.w[i] * self.y)
            self.x[i] = np.tanh(u)
            
        return self.x.copy()

Next, we will define the iteration logic. We start every iteration with the feedback loop — we implement the shift register circuit by shifting all y elements to the right, and compute the most recent y₀ output to place it into the first element of the input.

    def generate_trajectory(self, tmax: int, discard: int = 0) -> Tuple[np.ndarray, 
np.ndarray]:
        """
        Generate a trajectory of the states for tmax iterations.
        
        -----------
        tmax : int
            Total number of iterations
        discard : int
            Number of initial iterations to discard

        """
        self.reset()
        
        # Discard initial transient
        for _ in range(discard):
            self.iterate()
        
        x1_traj = np.zeros(tmax)
        x2_traj = np.zeros(tmax)
        
        for t in range(tmax):
            x = self.iterate()
            x1_traj[t] = x[0]
            x2_traj[t] = x[1]
            
        return x1_traj, x2_traj

Now, we define the function that will iterate our network map over the tmax number of time steps and output the states of the first two hidden neurons for visualization. We can use any hidden neurons, and we could even visualize 3D state space, but we will limit our imagination to two dimensions.

This is the gist of the system. Now, we will just define some line and segment magic for pretty visualizations.

import numpy as np
import matplotlib.pyplot as plt
import matplotlib.collections as mcoll
import matplotlib.path as mpath
from typing import Tuple, Optional, Callable


def make_segments(x: np.ndarray, y: np.ndarray) -> np.ndarray:
    """
    Create list of line segments from x and y coordinates.
    
    -----------
    x : np.ndarray
        X coordinates
    y : np.ndarray
        Y coordinates

    """
    points = np.array([x, y]).T.reshape(-1, 1, 2)
    segments = np.concatenate([points[:-1], points[1:]], axis=1)
    return segments


def colorline(
    x: np.ndarray,
    y: np.ndarray,
    z: Optional[np.ndarray] = None,
    cmap = plt.get_cmap("jet"),
    norm = plt.Normalize(0.0, 1.0),
    linewidth: float = 1.0,
    alpha: float = 0.05,
    ax = None
):
    """
    Plot a colored line with coordinates x and y.
    
    -----------
    x : np.ndarray
        X coordinates
    y : np.ndarray
        Y coordinates

    """
    if ax is None:
        ax = plt.gca()
        
    if z is None:
        z = np.linspace(0.0, 1.0, len(x))
    
    segments = make_segments(x, y)
    lc = mcoll.LineCollection(
        segments, array=z, cmap=cmap, norm=norm, linewidth=linewidth, alpha=alpha
    )
    ax.add_collection(lc)
    
    return lc


def plot_attractor_trajectory(
    x: np.ndarray,
    y: np.ndarray,
    skip_value: int = 16,
    color_function: Optional[Callable] = None,
    cmap = plt.get_cmap("Spectral"),
    linewidth: float = 0.1,
    alpha: float = 0.1,
    figsize: Tuple[float, float] = (10, 10),
    interpolate_steps: int = 3,
    output_path: Optional[str] = None,
    dpi: int = 300,
    show: bool = True
):
    """
    Plot an attractor trajectory.
    
    Parameters:
    -----------
    x : np.ndarray
        X coordinates
    y : np.ndarray
        Y coordinates
    skip_value : int
        Number of points to skip for sparser plotting

    """
    fig, ax = plt.subplots(figsize=figsize)
    
    if interpolate_steps > 1:
        path = mpath.Path(np.column_stack([x, y]))
        verts = path.interpolated(steps=interpolate_steps).vertices
        x, y = verts[:, 0], verts[:, 1]
    
    x_plot = x[::skip_value]
    y_plot = y[::skip_value]
    
    if color_function is None:
        z = abs(np.sin(1.6 * y_plot + 0.4 * x_plot))
    else:
        z = color_function(x_plot, y_plot)
    
    colorline(x_plot, y_plot, z, cmap=cmap, linewidth=linewidth, alpha=alpha, ax=ax)
    
    ax.set_xlim(x.min(), x.max())
    ax.set_ylim(y.min(), y.max())
    
    ax.set_axis_off()
    ax.set_aspect('equal')
    
    plt.tight_layout()
    
    if output_path:
        fig.savefig(output_path, dpi=dpi, bbox_inches='tight')

    return fig

The functions written above will take the generated state space trajectories and visualize them. Because the state space may be densely filled, we will skip every 8th, 16th or 32th time point to sparsify our vectors. We also don’t want to plot these in one solid color, therefore we are coding the color as a periodic function (np.sin(1.6 * y_plot + 0.4 * x_plot)) based on the x and y coordinates of the figure axis. The multipliers for the coordinates are arbitrary and happen to generate nice smooth color maps, to your liking.

N = 4
D = 32
s = 0.22
seed=174658140

tmax = 100000
discard = 1000

nn = NeuralAttractor(N, D, s, seed=seed)

# Generate trajectory
x1, x2 = nn.generate_trajectory(tmax, discard)

plot_attractor_trajectory(
    x1, x2,
    output_path='trajectory.png',
)

After defining the NN and iteration parameters, we can generate the state space trajectories. If we spend enough time poking around with parameters, we will find something cool (I promise!). If manual parameter grid search labor is not exactly our thing, we could add a function that checks what proportion of the state space is covered over time. If after t = 100,000 iterations (except the initial 1,000 “warm up” time steps) we only touched a narrow range of values of the state space, we are likely stuck in a point. Once we found an attractor that is not so shy to take up more state space, we can plot it using default plotting params:

One of the stable types of attractors is the limit cycle attractor (parameters: N = 4, D = 32, s = 0.22, seed = 174658140). It looks like a single, closed loop trajectory in phase space. The orbit follows a regular, periodic path over time series. I will not include the code for Lyapunov exponent calculation here to focus on the visual aspect of the generated attractors more, but one can find it under this link, if interested. The Lyapunov exponent for this attractor (λ=−3.65) is negative, indicating stability: mathematically, this exponent will lead to the state of the system decaying, or converging, to this basin of attraction over time.

If we keep increasing the scaling factor, we are more likely to tune up the values in the circuit, and perhaps more likely to find something interesting.

Here is the toroidal (quasiperiodic) attractor (parameters: N = 4, D = 32, s = 0.55, seed = 3160697950). It still has an ordered structure of sheets that wrap around in organized, quasiperiodic patterns. The Lyapunov exponent for this attractor has a higher value, but is still negative (λ=−0.20).

As we further increase the scaling factor s, the system becomes more prone to chaos. The strange (chaotic) attractor emerges with the following parameters: N = 4, D = 16, s = 1.4, seed = 174658140). It is characterized by an erratic, unpredictable pattern of trajectories that never repeat. The Lyapunov exponent for this attractor is positive (λ=0.32), indicating instability (divergence from an initially very close state over time) and chaotic behavior. This is the “butterfly effect” attractor.

As we further increase the scaling factor s, the system becomes more prone to chaos. The strange (chaotic) attractor emerges with the following parameters: N = 4, D = 16, s = 1.4, seed = 174658140. It is characterized by an erratic, unpredictable pattern of trajectories that never repeat. The Lyapunov exponent for this attractor is positive (λ=0.32), indicating instability (divergence from an initially very close state over time) and chaotic behavior. This is the “butterfly effect” attractor.

Just another confirmation that aesthetics can be very mathematical, and vice versa. The most visually compelling attractors often exist at the edge of chaos — think about it for a second! These structures are complex enough to exhibit intricate behavior, yet ordered enough to maintain coherence. This resonates with observations from various art forms, where balance between order and unpredictability often creates the most engaging experiences.

An interactive widget to generate and visualize these attractors is available here. The source code is available, too, and invites further exploration. The ideas behind this project were largely inspired by the work of J.C. Sprott [3]. 

References

[1] B. Poucet and E. Save, Attractors in Memory (2005), Science DOI:10.1126/science.1112555.

[2] Y.J.F. Kpomahou et al., Chaotic Behaviors and Coexisting Attractors in a New Nonlinear Dissipative Parametric Chemical Oscillator (2022), Complexity DOI:10.1155/2022/9350516.

[3] J.C. Sprott, Artificial Neural Net Attractors (1998), Computers & Graphics DOI:10.1016/S0097-8493(97)00089-7.