Predicting oscillations in complex networks with delayed feedback

TL;DR

It is revealed that oscillations emerge from the interplay of structural complexity and delay, with reduced models uncovering their critical thresholds and showing that greater connectivity lowers the delay required for their onset.

Abstract

Oscillatory dynamics are common features of complex networks, often playing essential roles in regulating function. Across scales from gene regulatory networks to ecosystems, delayed feedback mechanisms are key drivers of system-scale oscillations. The analysis and prediction of such dynamics are highly challenging, however, due to the combination of high-dimensionality, non-linearity and delay. Here, we systematically investigate how structural complexity and delayed feedback jointly induce oscillatory dynamics in complex systems, and introduce an analytic framework comprising theoretical dimension reduction and data-driven prediction. We reveal that oscillations emerge from the interplay of structural complexity and delay, with reduced models uncovering their critical thresholds and showing that greater connectivity lowers the delay required for their onset. Our theory is empirically tested in an experiment on a programmable electronic circuit, where oscillations are observed once structural complexity and feedback delay exceeded the critical thresholds predicted by our theory. Finally, we deploy a reservoir computing pipeline to accurately predict the onset of oscillations directly from timeseries data. Our findings deepen understanding of oscillatory regulation and offer new avenues for predicting dynamics in complex networks.

Figures (4)

• Figure 1: Complex oscillations in dynamical systems.a. Schematic illustration of real-world systems exhibiting oscillatory behavior: forest ecosystems (trees), gene regulatory networks, and power grids. Blue arrows denote interspecific interactions, while orange lines represent intraspecific time-delay effects. b. General formulation of a complex time-delayed system. Left: Temporal evolution of species abundances in system (\ref{['e1']}), exhibiting sustained oscillatory dynamics. Right: Linear stability analysis of system (\ref{['e1']}). The system is stable when all eigenvalues of the Jacobian matrix have negative real parts, i.e., when they lie entirely in the left half of the complex plane. Eigenvalues associated with species are marked by black dots. The time delay destabilizes the system, as indicated by the emergence of purely imaginary eigenvalues. The underlying structure follows an ER network with $N=10$, connectivity $C=0.5$, and parameters $r_i=0.4$, $s_i=0.08$, and $\tau=0.3$.
• Figure 2: Critical delay for oscillations induced by network complexity and time delay.a. The time series of system (\ref{['e1']}) under various connectivity and delays. For $C = 0.2$, the species abundance remains stable at $\tau = 0.21$ and $\tau = 0.24$, but exhibits oscillatory behavior at $\tau = 0.26$. As connectivity increases to $C = 0.22$, the system also transitions from a stable to an oscillatory regime, with oscillations already appearing at $\tau = 0.24$. When $C = 0.25$, oscillations emerge as early as $\tau = 0.21$, and their amplitude further increases with larger delays. b. The colored regions correspond to the dynamical states of the original high-dimensional system (\ref{['e1']}) on ER random networks, whereas the solid lines depict the theoretical estimates (\ref{['e4']}) derived from the reduced one-dimensional system (\ref{['e2']}). As $\tau$ and $C$ increase, the system states shift from stability to oscillatory behavior, as shown by the green (stable) and yellow (oscillatory, OS) regions. c. The reduced one-dimensional model (\ref{['e2']}) predicts the critical delay of the original high-dimensional system (\ref{['e1']}). Increasing $\beta_{\mathrm{eff}}$ decreases the critical delay of the system. Both theoretical (solid lines) and numerical (dotted lines) results are shown for ER networks with varying connectivity $C$ and network size $N = 100$, while different markers represent results from ten empirical networks displayed around the panel. d. The reduced 1D model accurately reproduces the abundance dynamics of the original high-dimensional model for both ER networks and empirical mutualistic networks. The network structures are the same as those used in the left panel. In panel (a), $N = 100$; for panels (b) and (c), $r_i=0.4$, $s_i=0.08$, $D_i=1$, $E_i=1$, and the networks are ER random networks.
• Figure 3: Prediction and validation of oscillatory behavior.a. The design block diagram of the oscillation control circuit experimental system consists of the MCU chip, the H-bridge driver circuit, the subtractor circuit, and the integrator circuit connected in series. The network nodes and their topological information are stored in the registers of the MCU. b. The physical PCB circuit, realized based on the design block diagram, is capable of simultaneously detecting the state waveforms mapped from up to four nodes. c. Experimental voltage u for various delays and ER random networks with different connectivity: $C = 0.08$, $C = 0.12$, and $C = 0.16$. The dashed line represents the predictions obtained from reservoir computing. All networks have $N = 100$ nodes, and model parameters are set to $r_i = 0.4$, $s_i = 0.08$, and $D_i = E_i = 1$. d. Schematic diagram of the reservoir computing architecture, which includes three layers: input, reservoir (hidden), and output. The reservoir is trained using time series corresponding to three delays ($\tau = 0.23$, $0.24$, and $0.25$) and demonstrates generalization beyond the training range by adjusting the input delay parameter. e. Prediction of oscillatory dynamics in system (\ref{['e1']}). Time series of $x_1(t)$ from system (\ref{['e1']}) and its prediction at $\tau = 0.28$ are shown to validate the accuracy and robustness of the training. f. Effective system abundance $x_{\mathrm{eff}}(t)$ as a function of $\tau$ for three cases. Theoretical results are obtained using the GBB method (system (\ref{['e4']})), real data are generated from system (\ref{['e1']}), and data-driven results are produced by the reservoir computing model.
• Figure 4: Comparison of three dimension reduction methods.a-c. The critical delay $\tau^*$ of the original high-dimensional system (\ref{['e1']}) is compared with that obtained from three reduced representations: the GBB system (\ref{['e2']}), the SDR system (\ref{['e13']}), and the MFA system (\ref{['e14']}). Panel (a) shows results for ER random networks with $N=100$ and varying connectivity $C$, panel (b) for empirical networks, and panel (c) for SW networks with $N=100$, varying nearest neighbors $K$ and $C=0.15$. Numerical results are obtained from direct simulations, while theoretical results represent the predicted critical delays from Eq. (\ref{['e4']}) for the reduced systems. The solid line with slope 1 serves as a reference for perfect agreement. d-f. Time series of the effective system abundance $x_{\mathrm{eff}}(t)$ are shown for the same four systems (the original system and the three reduced forms) under representative ER, empirical, and SW networks, with varying delay $\tau$, respectively. These panels highlight how the dynamics of each network type are reproduced across the different dimension reduction methods. Model parameters are fixed at $r_i = r= 0.4$, $s_i = s = 0.08$, $D_i = D = 1$, and $E_i = E = 1$.