Tractable nonlinear memory functions as a tool to capture and explain dynamical behaviours (2005.04751v2)

Abstract: Mathematical approaches from dynamical systems theory are used in a range of fields. This includes biology where they are used to describe processes such as protein-protein interaction and gene regulatory networks. As such networks increase in size and complexity, detailed dynamical models become cumbersome, making them difficult to explore and decipher. This necessitates the application of simplifying and coarse graining techniques in order to derive explanatory insight. Here we demonstrate that Zwanzig-Mori projection methods can be used to arbitrarily reduce the dimensionality of dynamical networks while retaining their dynamical properties. We show that a systematic expansion around the quasi-steady state approximation allows an explicit solution for memory functions without prior knowledge of the dynamics. The approach not only preserves the same steady states but also replicates the transients of the original system. The method also correctly predicts the dynamics of multistable systems as well as networks producing sustained and damped oscillations. Applying the approach to a gene regulatory network from the vertebrate neural tube, a well characterised developmental transcriptional network, identifies features of the regulatory network responsible dfor its characteristic transient behaviour. Taken together, our analysis shows that this method is broadly applicable to multistable dynamical systems and offers a powerful and efficient approach for understanding their behaviour.

• The paper’s main contribution is deriving explicit, tractable nonlinear memory functions that capture the dynamic influence of bulk components on key subnetwork behaviors.
• It introduces a self-consistent Zwanzig–Mori projection method around the QSS approximation to reduce complexity in simulating high-dimensional nonlinear systems.
• The approach accurately reproduces multistability, oscillations, and transient dynamics in models like biological regulatory networks, enhancing mechanistic insights.

This paper introduces a method for reducing the complexity of nonlinear dynamical systems, particularly relevant for systems like biological regulatory networks, where large size and complex interactions make traditional analysis and simulation difficult. The core idea is to use a Zwanzig-Mori (ZM) projection approach to explicitly track a chosen subset of the system's components (the "subnetwork") while replacing the remaining components (the "bulk") with memory functions that capture their dynamic influence on the subnetwork.

The key contribution is the derivation of an explicit, tractable form for these nonlinear memory functions. Unlike previous ZM applications which often rely on linear approximations or specific system structures, this method systematically expands around the Quasi-Steady State (QSS) approximation for the bulk. The method assumes that, for any given state of the subnetwork, the bulk has a unique steady state. This allows defining a projection operator based on setting the bulk variables to these QSS values.

The time evolution of the subnetwork species $x_s$ is then described by:

$$\frac{d}{d t}x_s=v_s(x(t))+\int_0^tdt'\ M_s(x(t'),t-t')$$

Here, $v_s(x)$ is the QSS drift, which is simply the rate of change of $x_s$ from the full system's equations, but with the bulk species fixed at their QSS values $x^*(x)$. The integral term contains the memory function $M_s(x, \tau)$, which corrects the QSS dynamics by accounting for the fact that the bulk is not instantaneously at QSS.

The paper derives $M_s(x, \tau)$ as:

$M_s(x,\tau) = \sum_{b'}c_{b'}(x)f_{b's}(x,\tau)$

where $c_{b'}(x)$ depends on derivatives of the rate functions evaluated at the QSS, and $f_{b's}(x,\tau)$ describes how initial deviations from QSS propagate in the linearized bulk dynamics, influenced by the flow of the subnetwork state under the QSS drift $v(x)$. This results in a memory function that is nonlinear in the past subnetwork state $x(t')$ and depends on the time difference $\tau = t-t'$.

Implementation Considerations and Self-Consistent Approximation (ZMs)

Directly simulating the integro-differential equation above (referred to as ZMn) can be computationally expensive due to the memory integral requiring storage of past states and evaluation of $M_s(x(t'), \tau)$. To address this, the paper proposes a self-consistent approximation (ZMs). This approximation replaces the QSS flow $\bm{\phi}_v(x(t'), t-t')$ in the expression for $f_{b's}$ with the actual state $x(t)$. This seemingly small change allows the memory integral term $$\tilde{\mathcal{M}_s(t) = \sum_b m_b(t) f_{bs}^0(x(t))}$$ to be calculated by introducing auxiliary variables $m_b(t)$. These $m_b(t)$ satisfy a system of linear ordinary differential equations (ODEs) driven by the current subnetwork state $x(t)$:

$$\frac{d}{dt} m_b(t) = c_b(x(t)) + \sum_{b'} m_{b'}(t)l_{b'b}(x(t))$$

with initial conditions $m_b(0)=0$. The matrix $l_{b'b}(x)$ depends on derivatives of the rate functions evaluated at the QSS.

Practical Implementation of ZMs:

Implementing the ZMs approach requires solving a coupled system of ODEs: the original subnetwork ODEs augmented by one auxiliary ODE for each bulk species:

$$\frac{d}{dt}x_s(t)=v_s{(x(t))}+\sum_{b'}m_{b'}(t)f^{0}_{sb'}(x(t))$$

$$\frac{d}{dt} m_b(t) = c_b(x(t)) + \sum_{b'} m_{b'}(t)l_{b'b}(x(t))$$

1. Define System: Express the original system $\frac{d\bm{x}}{dt}=R(\bm{x})$ and choose the subnetwork ($x$) and bulk ($x$) species.
2. Calculate QSS: Solve $R^b(x, x^*(x))=0$ for $x^*(x)$. This is a crucial step and requires that $x^*(x)$ is unique for the chosen bulk. This might involve numerical root finding.
3. Calculate QSS Drift: Compute $v_s(x) = R_s(x, x^*(x))$.
4. Calculate $f^0_{bs}(x)$, $c_b(x)$, and $l_{bb'}(x)$: Compute the derivatives $\frac{\partial R_s}{\partial x_b}$, $\frac{\partial R_b}{\partial x_{s'}}$, $\frac{\partial R_{b'}}{\partial x_b}$, and the Jacobian $J_{bb'} = \frac{\partial R_{b'}}{\partial x_b}$. Evaluate these and $R_s$ at $(x, x^*(x))$. Compute the inverse Jacobian $(\bm{J}^{-1})_{b'b''}$. Use these to calculate $$f^0_{bs}(x) = \frac{\partial R_s}{\partial x_b}(x, x^*(x))$$, $$l_{bb'}(x) = J_{b'b}(x) + \sum_{s''b''}(\bm{J}^{-1})_{b'b''}(x)\frac{\partial R_{b''}}{\partial x_{s''}}(x)\frac{\partial R_{s''}}{\partial x_b}(x)$$, and $$c_{b'}(x) = \sum_{s'} \sum_{b''}(\bm{J}^{-1})_{b'b''}(x)\frac{\partial R_{b''}}{\partial x_{s'}}(x) R_{s'}(x)$$. Note that these functions $v_s$, $f^0_{bs}$, $c_b$, and $l_{bb'}$ depend nonlinearly on the subnetwork state $x$.
5. Set up Augmented ODEs: Define the system of ODEs for $[x_s(t), m_b(t)]$.
6. Integrate: Numerically integrate this augmented system of ODEs starting from desired initial conditions $[x_s(0), m_b(0)]$. Note that $m_b(0)=0$ corresponds to the bulk initially being at its QSS relative to the initial subnetwork state.

The computational cost is primarily in calculating the functions $v_s$, $f^0_{bs}$, $c_b$, and $l_{bb'}$ at each step of the ODE integration (which involves solving for $x^*$ and matrix inversion) and in integrating the augmented system of ODEs (number of variables = $|x| + |x|$). This is generally more efficient than the ZMn method for time course simulations, especially for long times or when high accuracy is needed.

Practical Applications and Insights:

The paper demonstrates the method's utility on several example systems:

1. Multistable Systems (Hill functions): Applying the method to bistable and tetrastable systems shows it accurately captures complex transient dynamics and basin boundaries, which are crucial for understanding state transitions or "cell fate decisions" in biological contexts. The QSS approximation alone fails to replicate the timing of transients and distorts the basin shapes.
2. Oscillations (Repressilator): The ZMs method correctly reproduces both damped and sustained oscillations and the associated Hopf bifurcation point. The QSS approximation predicts only damped oscillations. Analysis of the memory amplitude reveals how its sign changes drive the oscillatory behavior.
3. Neural Tube Gene Regulatory Network [Cohen2014]: This real-world biological network exhibits multistability and complex transients dependent on a spatial parameter (Shh signaling level). The ZMs method successfully captures the tristability basin boundaries and non-monotonic transients observed in the full system, where the QSS approximation fails qualitatively.
   • Key Insight (Memory Decomposition): By decomposing the memory functions into "channels" (signal flow from subnetwork species to bulk species and back), the authors identify specific molecular interactions (e.g., involving Pax6 and Irx3) that are responsible for generating crucial transient behaviors like delayed gene expression (Nkx2.2 delay relative to Olig2). This provides a mechanistic explanation for observed dynamics from the reduced model.
   • Robustness: The memory effects are shown to increase the robustness of the transient dynamics to variations in initial conditions, ensuring consistent behavior despite perturbations, a vital property for biological development.

Advantages and Limitations:

• Advantages:
  • Provides tractable, explicit nonlinear memory functions.
  • Accurately captures complex dynamics (multistability, oscillations, transients) where QSS fails.
  • The ZMs version is computationally efficient for time course simulation via augmenting ODEs.
  • Allows decomposition of memory effects into channels for mechanistic insight.
  • Applicable to generic nonlinear ODE systems (assuming unique bulk QSS).
• Limitations:
  • Requires the bulk to have a unique steady state for any given subnetwork state. This might restrict the choice of subnetwork/bulk split or applicability to certain network structures.
  • Relies on linearizing deviations around QSS, which might introduce inaccuracies if bulk dynamics deviate significantly from the QSS manifold.
  • The random force term is neglected, meaning it captures deterministic effects of the bulk not being at QSS, but not stochastic effects from initial conditions far from QSS.

In summary, this paper offers a powerful, generalizable, and practical method for coarse-graining complex dynamical systems while retaining their key behaviors. The ability to derive explicit nonlinear memory functions and simulate them efficiently via augmented ODEs, coupled with the interpretation offered by memory decomposition, makes it a valuable tool for analyzing and understanding the dynamics of high-dimensional nonlinear systems, particularly in biology.