Control of Stochastic and Induced Switching in Biophysical Networks

Noise caused by fluctuations at the molecular level is a fundamental part of intracellular processes. While the response of biological systems to noise has been studied extensively, there has been limited understanding of how to exploit it to induce a desired cell state. Here we present a scalable, quantitative method based on the Freidlin-Wentzell action to predict and control noise-induced switching between different states in genetic networks that, conveniently, can also control transitions between stable states in the absence of noise. We apply this methodology to models of cell differentiation and show how predicted manipulations of tunable factors can induce lineage changes, and further utilize it to identify new candidate strategies for cancer therapy in a cell death pathway model. This framework offers a systems approach to identifying the key factors for rationally manipulating biophysical dynamics, and should also find use in controlling other classes of noisy complex networks.

Popular Summary

Noise is a fundamental aspect of many physical and biophysical systems and is capable of inducing large-scale shifts in system behavior. For example, noise can induce a change in the global gene expression state of a cell, potentially altering its behavior drastically. Such rare but profound state transitions are important for organismal development and the cellular transition to cancer. Controlling these transitions requires controlling how different states respond to noise, whereas current control methods focus primarily on manipulating system states directly. To address this issue, we present here a highly scalable algorithm to rationally identify parameter interventions to control the response to noise in complex nonlinear networks.

This method, termed “optimal least action control,” works by modulating the heights of the barriers separating different stable states, and can also be used to manipulate the stability of states even in the absence of noise. A central insight of our approach is that state transition dynamics on a network can be mapped to a greatly distilled network of state transitions; in this way, a network can be used to solve a network problem. We apply this method to numerous biophysical models, including a model of the cell death pathway in which we identify top candidate gene targets for potentially new cancer treatment strategies.

Our algorithm represents a fundamentally new direction in the control of large dynamical systems: Rather than controlling the system state directly, one can leverage noise and rationally manipulate the dynamics of that system indirectly. We expect that our results will be applicable to a range of network systems in biological, infrastructural, and ecological realms.

Figure 1

Control of the response to noise illustrated for a multiwell quasipotential. (a) Without noise the state of the system is fixed and will not change over time. (b) In the presence of noise the state of the system wanders within the attraction basin and will be eventually ejected into the attraction basin of another stable state. For small noise such transitions are exponentially more likely to occur through the lowest barrier, even though other transitions are also possible in principle. (c) Using OLAC on a set of tunable parameters, we can alter the quasipotential to produce a desired response to noise, in this case tailored to increase transitions from state 2 to state 1 and to reduce transitions to state 3. (d) If desired, OLAC can also find combinations of tunable parameters that can alter the quasipotential in order to eliminate a stable state through a bifurcation; in this illustration state 2 is eliminated, in favor of state 1. (e)–(h) NEST for each of the quasipotentials in (a)–(d), respectively, where the nodes represent stable states and the continuous edges represent transition rates. A wider continuous edge indicates a higher transition rate. The dotted edges in (e) indicate transitions that could occur in the presence of noise.

Figure 2

OLAC applied to a two-dimensional model of VPC differentiation. (a), (b) State space representation of stable states (nodes) and optimal transition paths (edges) before intervention (a) and after OLAC is applied to maximize the limiting occupancy of lineage 1° (b). Node size indicates occupancy and edge width indicates the (negative of the) log transition rate along the corresponding minimum action path; the color code on the transition paths indicates the derivative of the quasipotential. The background shows the velocity field for the given parameters. For equal EGF and Notch signaling (a), lineage 2° is the one with highest occupancy. The OLAC solution (b) indicates that high EGF signaling and low Notch signaling will lead to the maximum occupancy of lineage 1°. (c) Trajectory in the parameter space for one realization of the optimization procedure, where the contour plot indicates the limiting occupancy of lineage 1°. Note that each step of the optimization routine increases the occupancy of this state. (d), (e) NESTs for the initial (d) and optimized (e) systems. In all panels, the noise intensity is assumed to be $ϵ=0.007$ (as used previously [24]).

Figure 3

Controlled stochastic lineage switching in a multidimensional model of HPC differentiation. (a), (b) NESTs of the model for the unmodified parameters (a) and for OLAC applied to optimize the limiting occupancy of the $\beta$ cell state (red) (b), for $ϵ=0.01$. In both panels, node size indicates the limiting occupancy of the state and edge width indicates the (negative of the) log transition rate. OLAC identifies that only three (out of ten) transcription factors need to be tuned in the model to optimize $\beta$-cell occupancy: MafA (increased), Brn4 (decreased), and $\delta$-gene (decreased). (c) Transition hierarchy into the $\beta$ cell state for the optimized system. Two states ($\delta$ and $\alpha /\beta$) transition directly; two others ($\alpha$ and $\alpha /\delta$) require passage through an intermediate state—in this case $\alpha /\beta$.

Figure 4

Optimal therapeutic interventions for the cell death regulatory network model. The network has 22 genes (circles), 3 input parameters (top), and 42 target tunable parameters (edges), where the edge heads distinguish between excitatory (arrow) and inhibitory (bar) regulatory interactions. OLAC is applied to induce a bifurcation transition from survival to apoptotic states. The parameters involved in one or more optimal intervention are consistently upregulated (green) or downregulated (red) by the interventions; those not recruited by any optimal intervention are shown in black. The edge width indicates the prevalence of that parameter in optimal interventions for the various cell types and intervention strengths considered (up to the elimination of the survival state). The size of each circle represents the sensitivity of the gene to proapoptotic interventions, defined as the change in gene expression between the uncontrolled scenario and the smallest intervention scenario at which the survival state is eliminated; the color distinguishes up and down expression. Genes in white are those that change expression substantially between the survival and apoptotic states. These results are based on 6 simulated cell types and 9 different intervention strengths (see Supplemental Material [20], Fig. S3).

Figure 5

Transition time for the cell death model as a function of barrier height and noise strength. Transition time is calculated using Eq. (6). The black line separates those barrier height–noise strength combinations that occur over a therapeutically relevant time scale (bottom right) and those that do not (top left).