# A thermodynamic limit on molecular computation

Easun Arunachalam1,2,3, Milo M. Lin1,2,4,5,∗

1Green Ctr. for Systems Biology, University of Texas Southwestern Medical
Ctr., Dallas, TX, USA
2Lyda Hill Dept. of Bioinformatics, University of Texas Southwestern Medical
Ctr.
3Dept. of Molecular and Cellular Biology, Harvard University, Cambridge, MA,
USA
4Dept. of Biophysics, University of Texas Southwestern Medical Ctr.
5Ctr. for Alzheimer’s and Neurodegenerative Diseases, University of Texas
Southwestern Medical Ctr.
∗To whom correspondence should be addressed; E-mail:
<EMAIL_ADDRESS>
> Molecular systems can perform discrete computations by converting between
> different steady-state probability distributions over system states. Rather
> than using equilibrium steady-states for computation, biology uses
> thermodynamic forces to access non-equilibrium steady-states, which require
> constant energy expenditure to maintain even when no computations are being
> performed. The information-theoretic advantage of this costly paradigm is
> poorly understood. Here, we derive a concise expression relating the
> thermodynamic force to the information gain achievable in any molecular
> system performing an arbitrary computational step, thereby quantifying the
> maximum advantage derived from non-equilibrium driving. Numerical sampling
> demonstrates the tightness of this universal bound, and as examples, we show
> that Ras signaling and microtubule assembly can closely approach it.

#### One-sentence summary:

Information gain in any biomolecular system is bounded by the thermodynamic
force driving the system away from equilibrium.

Molecular systems can stochastically interconvert between different states,
and the probability distribution over these states represents the information
content of the system. In cases such as protein folding (?), modulation of
enzyme activity (?), or signaling (?, ?), the probability of a single output
state is biologically most relevant. In other cases, such as self-assembly (?,
?), biologically relevant information is encoded by the entire probability
distribution over all states. A discrete computational step is a modification
of system parameters, which alters the steady-state probability distribution
(Figure 1A). The molecular transformations underpinning life can be modeled as
a series of such computational steps.

If a computational step converts between equilibrium steady states, energy is
not required to maintain the initial or final probability distributions. Such
equilibrium computation has been synthetically realized, e.g. by programming
DNA nano-machines (?, ?). In contrast, living systems perform non-equilibrium
computation by harnessing thermodynamic forces, often in the form of an excess
concentration of an energy currency such as adenosine or guanosine
triphosphate (ATP or GTP). In the cell, these forces are typically kept
constant while the system changes its other parameters (which we refer to as
passive parameters) from $\mathbf{c}\to\mathbf{c}^{\prime}$ in order to
convert between non-equilibrium steady states. Examples of such computational
steps include mutations, changes in the concentrations of ligands or
cofactors, or chemical alterations such as phosphorylation.

Thermodynamic forces maintain non-equilibrium steady states in between
computational steps, resulting in continuous energy consumption even when no
computations are being performed. It is therefore natural to ask whether this
energetically costly form of computation confers quantifiable advantages over
equilibrium computation. Here we provide a universal relation showing how
information gain in any molecular computation is limited by the thermodynamic
force. Information gain and suppression of noise arising from parameter
fluctuations are shown to be equally constrained by the thermodynamic force;
they are fundamentally symmetric sides of the same coin.

Figure 1: Equilibrium versus non-equilibrium computation. (A) Computation is
defined as a change in the steady-state probability distribution of system
states. (B) Protein-receptor binding is an example of equilibrium computation.
(C) Protein phosphorylation-gated receptor binding is an example of non-
equilibrium computation. (D) The rate constant and circuit representations of
a Markov model are equivalent to one another. The transition, rate constant,
and voltage shown in red are equivalent representations of the single non-
equilibrium element in this system. The vector $\mathbf{c}$ contains all the
passive parameters for this three-state system. (E) Information gain between
systems differing by 1-4 passive parameters as a function of the applied non-
equilibrium driving force.

#### Statistical mechanics of non-equilibrium computation.

We consider discrete-state Markovian systems coupled to a bath at temperature
$T$. In the rate constant representation, systems at equilibrium are
parameterized by passive rate constants $\kappa_{ij}$ associated with
transitions between each pair of connected states $i$ and $j$. The probability
of state $i$, denoted by $\tilde{p}_{i}$, evolves according to
$\frac{d\tilde{p}_{i}}{dt}=\kappa_{ij}\tilde{p}_{i}-\kappa_{ji}\tilde{p}_{j}$;
the steady-state probability of $i$ is denoted by $p_{i}$.

At equilibrium, the net probability current between $i$ and $j$, denoted by
$I_{ij}=\kappa_{ij}p_{i}-\kappa_{ji}p_{j}$, is zero for all $ij$ (detailed
balance). Non-equilibrium systems can be constructed by starting from a
reference equilibrium system and adding a nonzero active rate constant
$\alpha_{ij}$ to the passive rate constant $\kappa_{ij}$ between any two
connected states $i$ and $j$, which can break detailed balance. Non-
equilibrium systems are thus parameterized by passive and active rate
constants $\kappa_{ij}$ and $\alpha_{ij}$ for all $ij$. We shall frequently
compare a non-equilibrium system to its reference equilibrium system, which
shares the same $\kappa_{ij}$s, but for which $\alpha_{ij}=0$ for all $ij$. We
define the applied driving force $\theta_{ij}$ in terms of the active and
passive rate constants between states $i$ and $j$ (Box 1):

$\theta_{ij}=kT\ln\left(1+\dfrac{\alpha_{ij}}{\kappa_{ij}}\right)$ (1)

where $kT$ is the product of the Boltzmann constant and the temperature, which
sets the energy scale for the system. We denote by $\bm{\theta}$ the set of
all active parameters $\theta_{ij}$, and by $\mathbf{c}$ the set of all
passive parameters $\kappa_{ij}$ for all $ij$; the probability of state $i$ is
thus represented by $p_{i}(\bm{\theta},\mathbf{c})$. To compute steady-state
probabilities for non-equilibrium systems, we employ an alternative
parameterization (a circuit representation (?)) that makes analytic solutions
tractable and enables efficient numerical sampling of different reference
equilibrium systems (Fig. 1D and Supplementary Note 1).

In many contexts, biomolecular computation corresponds to changing the passive
parameters $\mathbf{c}$, while keeping the non-equilibrium driving forces
$\bm{\theta}$ fixed. We formally define a computational step as a transition
in passive parameters $\mathbf{c}\to\mathbf{c}^{\prime}$. The Kullback-Leibler
divergence $D$ between the steady-state probability distribution of the post-
computation system and the probability distribution of the pre-computation
system, also known as the information gain, is:

$D(\bm{\theta};\mathbf{c}^{\prime}\|\mathbf{c})=\sum_{i=1}^{N}p_{i}(\bm{\theta},\mathbf{c}^{\prime})\ln{\Big{[}\dfrac{p_{i}(\bm{\theta},\mathbf{c}^{\prime})}{p_{i}(\bm{\theta},\mathbf{c})}\Big{]}}.$
(2)

To illustrate the difference between equilibrium and non-equilibrium
computation, consider simple examples of each case. Protein complex formation
is a prototypical example of equilibrium computation (Figure 1B). In a three-
state model of this process, protein A can reversibly bind protein B, the
A$\cdot$B complex can reversibly bind a receptor R, and the A$\cdot$B$\cdot$R
complex can reversibly dissociate. The relative populations of the different
states of $A$ can be changed by increasing the number of B proteins or
mutating them to increase their binding affinity to A, which will increase the
probability of finding A in the complexed states. Operation of this system
requires energy expenditure only during the computational steps, when copies
of B are produced or degraded (Figure 1B).

Phosphorylation-gated binding is a prototypical example of non-equilibrium
computation (Figure 1C). An active kinase can consume energy to phosphorylate
protein A, and a phosphatase can dephosphorylate A; phosphorylated A can
reversibly bind a receptor R; and reversible dephosphorylation not mediated by
an enzyme can lead to receptor unbinding. In contrast to equilibrium
computation, the power consumption of the system has two distinct components:
energy is consumed during the computational steps to activate kinases, and is
also continuously consumed to maintain the steady-state net clockwise
probability flux in between computational steps (Figure 1C).

Both examples can be represented by specific parameterizations of a general
three-state model (Figure 1D). To illustrate the difference between
equilibrium versus non-equilibrium computation in this family of three-state
models, we randomly generated a large number of sets of passive parameters
$\mathbf{c}$, performed a computational step by changing 1-4 individual
parameters to yield a new set $\mathbf{c}^{\prime}$, applied a fixed non-
equilibrium driving force $\theta_{12}$ between states 1 and 2 in each
network, and calculated the information gain as a function of the number of
parameters changed and $\theta_{12}$. We then examined the maximum information
gain, over the different sets of $\mathbf{c},\mathbf{c}^{\prime}$ pairs, as a
function of $\theta_{12}$ (Figure 1E). At equilibrium ($\theta_{12}=0$),
information gain increases modestly as more passive parameters are allowed to
change. Regardless of how many parameters are changed, the presence of an
applied driving force ($\theta_{12}>0$) enhances the maximum information gain
compared to the same computational step at equilibrium ($\theta_{12}=0$). The
information gain is sensitive to the applied driving force within a limited
range, beyond which it saturates. These data suggest that non-equilibrium
computation may be needed to achieve a target level of information gain if the
number or dynamic range of passive parameters that can be modulated is
limited. Below, we derive a mathematical relation that precisely establishes
this principle and applies to arbitrarily complex systems.

Box 1: What is the applied driving force $\theta$ ? Applied driving forces can
induce a net current $I$ in the system (as in Fig. 1C). The associated power
consumption $\mathcal{P}$ of the system, which is also the rate of heat
dissipation at steady state, is (?): $\mathcal{P}=\sum_{mn}\theta_{mn}I_{mn}.$
(B1) where the sum is performed over each pair of connected states $m$ and
$n$, $I_{mn}$ is the current between $m$ and $n$, and $\theta_{mn}$ is related
to the rate constants between $m$ and $n$ as per Eq. 1. Because power is
generically the product of a force and a velocity, and the current $I_{mn}$ is
a “probability velocity”, we can interpret $\theta_{mn}$ as the applied
driving force. The applied driving force in a system is also related to
chemical potential differences: the sum of the net applied driving force along
any loop is equal to the net chemical potential difference along that loop,
which is also called the cycle affinity (?, ?) (Supplementary Note 1):
$\sum_{mn\>\in\>\text{loop}}\theta_{mn}-\theta_{nm}=\sum_{mn\>\in\>\text{loop}}\Delta\mu_{mn}.$
(B2) In contrast to the rate of dissipation $\mathcal{P}$ or the chemical
potential difference $\Delta\mu$ between two states, which are generally
emergent properties of the system, the applied driving force $\theta$ between
two states is an explicit function of the rate constants between those states.
Therefore, $\theta$ serves as a directly controllable design parameter. For
example, in biological systems, $e^{\theta/kT}$ is often the fold-increase of
ATP concentration over its concentration at equilibrium. This driving leads to
breaking of detailed balance and net entropy production, requiring the
constant external supply of energy to regenerate ATP. In Eq. 3 we introduce
the notion of the total applied driving force $\Theta$, which is the maximum
possible sum of $\theta_{ij}$s along any non-intersecting path in state space.
This work demonstrates how $\Theta$ dictates the information processing limits
of a system.

#### Applied driving forces bound amplification and suppression of individual
states.

The fundamental components of information gain are the probabilities of
individual states. Thus, to obtain our main result, we first proved bounds on
how the probabilities of individual states can change under non-equilibrium
driving.

The steady-state probability of any state can be expressed in terms of rate
constants present in distinct directed spanning trees of a network. The
relevant object then becomes non-self-intersecting paths between two states;
the presence of many driven transitions in a single path enables extreme
probability amplification and suppression. We define the total applied driving
force $\Theta$ to be the maximum value of the applied driving force along any
non-self-intersecting path in the network (Box Statistical mechanics of non-
equilibrium computation.):

$\Theta\equiv\max_{\text{n.s.i. path }u}\left(\sum_{\text{edge
}mn\>\in\>u}\theta_{mn}\right)$ (3)

In Supplementary Note 3, we use Kirchhoff’s matrix-tree theorem (?) to show
that the steady-state probability of any state $i$ in the presence and absence
of the applied driving force ($p_{i}(\Theta,\mathbf{c})$ and
$p_{i}(0,\mathbf{c})$, respectively) is constrained by:

$\displaystyle p_{i}(0,\mathbf{c})e^{-\Theta/kT}\leq
p_{i}(\Theta,\mathbf{c})\leq p_{i}(0,\mathbf{c})e^{\Theta/kT}.$ (4)

Hence, $\Theta$ limits the maximal amplification or suppression of any state.
This bound is general, and applies to all steady-state systems regardless of
topological complexity. This result differs from previous work (?, ?, ?) which
bounded the relative probabilities between two different states of the same
system. In contrast, Eq. 4 bounds the change in the probability of the same
state if the system is altered by thermodynamic driving, and therefore
requires treatment of the denominator in the matrix-tree theorem
(Supplementary Note 3, Eq. S8).

To illustrate the tightness of this bound, we generated a range of state-space
networks with diverse topologies and passive parameters $\mathbf{c}$ (Fig.
2A). In each system, we randomly selected four transitions along which to
apply a non-equilibrium driving force ($\theta=2\>kT$) and compared the
equilibrium probability $p_{i}^{*}$ and non-equilibrium probability $p_{i}$ of
each state $i$. We observe that $p_{i}$ cannot exceed
$p_{i}^{*}e^{4\theta/kT}$, and cannot be lower than
$p_{i}^{*}e^{-4\theta/kT}$. Moreover, there are many states for which $p_{i}$
lies along sub-bounds given by $p_{i}^{*}e^{\pm\theta/kT}$, $p_{i}^{*}e^{\pm
2\theta/kT}$, and $p_{i}^{*}e^{\pm 3\theta/kT}$. Our proof highlights the
common origin of the overall bound and sub-bounds as follows. In the example
set of systems containing four driven transitions analyzed in Fig. 2, there
may exist arrangements in which all four driven transitions point in the same
direction along a single non-intersecting path. In this case,
$\Theta=4\theta$, and some states may be maximally suppressed or amplified.
However, for other arrangements, there may not exist any non-intersecting
paths containing all of the driven transitions. If at most one, two or three
driven transitions are part of a non-intersecting path, we have
$\Theta=\theta$, $2\theta$ or $3\theta$, respectively. Thus, even though there
may exist four driven transitions in the overall network, no state probability
can be suppressed or amplified to the extent possible when $\Theta=4\theta$
(Fig. 2C).

This result also explains why the maximum amplification or suppression depends
on whether driven transitions are arranged in a series or parallel
configuration. When arranged in series, a single non-intersecting path can
contain all transitions, and hence amplification depends on the number of
driven transitions (Fig. 2D). In contrast, when arranged in parallel (e.g.
transitions from several different states are directed towards the same
state), no path can contain more than one driven transition. Hence, the
maximum amplification is determined by the largest driving force among the
different transitions, regardless of the total number of driven transitions
(Fig. 2E). Overall, the bound for each value of $\Theta$ (Eq. 4) is tight.

Figure 2: Total applied driving force tightly bounds probability amplification
and suppression in arbitrary networks. (A) Examples of network topologies
whose probability amplification and suppression properties were studied. Four
driven transitions were added in different locations in each network, and the
equilibrium and non-equilibrium probabilities of each state were calculated.
(B) Equilibrium and non-equilibrium probabilities of individual states in each
network. Each dot represents a single state. (C) Different arrangements of the
four driven transitions lead to differences in the maximum number of driven
transitions along any non-intersecting path (i.e. different total applied
driving forces), and hence different bounds on probability amplification and
suppression. (D) Equilibrium and non-equilibrium probabilities of individual
states in networks with different numbers of transitions driven in series. (E)
Equilibrium and non-equilibrium probabilities of individual states in networks
with different numbers of transitions driven in parallel.

#### Maximum signal amplification and suppression limit.

Eq. 4 bounds the change in probability of a system state if the driving force
$\Theta$ is applied. We next used this relation to derive bounds relevant for
biomolecular computation, in which the passive parameters $\mathbf{c}$ are
changed while $\Theta$ is held constant.

The maximum amplification or suppression of the probability of any state $i$
after a computational step $\mathbf{c}\to\mathbf{c}^{\prime}$, relative to the
amplification of the same computational step at equilibrium, is given by:

$\Bigg{|}\ln{\left(\dfrac{p_{i}(\Theta,\mathbf{c})}{p_{i}(\Theta,\mathbf{c}^{\prime})}\right)}-\ln{\left(\dfrac{p_{i}(0,\mathbf{c})}{p_{i}(0,\mathbf{c}^{\prime})}\right)}\Bigg{|}\leq\dfrac{2\Theta}{kT}.$
(5)

This amplification/suppression relation applies to the change in probability
of any single biologically relevant state $i$ of a system (e.g. a signaling
state). An input event, such as phosphorylation, changes the passive
parameters of the signaling protein from $\mathbf{c}\to\mathbf{c}^{\prime}$,
thereby altering the probability of the signaling state, which conveys the
signal to downstream components of the signaling pathway. We consider the
specific case of the Ras signaling protein below.

#### Maximum and minimum information gain limit.

In other cases such as self-assembly, the overall probability distribution
over all system states is biologically meaningful. Here, the relevant measure
of computation is the information gain rather than the fold-change of an
individual state. We used Eq. 4 to derive bounds on non-equilibrium
information gain ($D(\Theta;\mathbf{c}\|\mathbf{c}^{\prime})$) as a function
of the equilibrium information gain ($D(0;\mathbf{c}\|\mathbf{c}^{\prime})$)
and the total applied driving force $\Theta$ (Supplementary Note 4):

$\displaystyle\begin{split}\text{upper bound:
}&D(\Theta;\mathbf{c}^{\prime}\|\mathbf{c})\leq\dfrac{2\Theta}{kT}+D(0;\mathbf{c}^{\prime}\|\mathbf{c})\mathrm{e}^{\Theta/kT}\\\
\text{lower bound:
}&D(0;\mathbf{c}^{\prime}\|\mathbf{c})\leq\dfrac{2\Theta}{kT}+D(\Theta;\mathbf{c}^{\prime}\|\mathbf{c})\mathrm{e}^{\Theta/kT}\end{split}$
(6)

The bound on non-equilibrium information gain
$D(\Theta;\mathbf{c}^{\prime}\|\mathbf{c})$ has two components. The first
term, which is independent of the equilibrium information gain and scales
linearly with the thermodynamic force, is relevant to computational steps such
as the addition of a catalyst. Such steps lead to zero information gain at
equilibrium but non-zero information gain in the presence of the thermodynamic
force. This type of computation can be implemented by introducing auxiliary
molecules or modifications that are sub-stoichiometric relative to the system,
and can potentially be achieved rapidly and without overcrowding. For example,
the concentrations of some kinases are more than three orders of magnitude
lower than their target substrates (?). The second term is proportional to the
equilibrium information gain and scales exponentially with the thermodynamic
force. However, this exponential amplification can only act on information
gain already present in the equilibrium computation, and must therefore be
implemented via stoichimetric modifications.

Some systems can consume energy to generate larger changes in their state
probability distributions than would be possible through modification of
passive parameters $\mathbf{c}$ alone, resulting in
$D(\Theta;\mathbf{c}^{\prime}\|\mathbf{c})>D(0;\mathbf{c}^{\prime}\|\mathbf{c})$.
Such systems are capable of additional information gain that is enabled by
non-equilibrium driving. Other systems consume energy to avoid changes in
their state probability distributions that would occur when the underlying
passive parameters were modified from $\mathbf{c}\to\mathbf{c}^{\prime}$, i.e.
$D(\Theta;\mathbf{c}^{\prime}\|\mathbf{c})<D(0;\mathbf{c}^{\prime}\|\mathbf{c})$.
In such cases, non-equilibrium driving enables noise suppression, making
systems robust against fluctuations in $\mathbf{c}$. Eq. 6 provides upper and
lower bounds on additional information gain and noise suppression,
respectively.

We tested the tightness of the bound in Eq. 6 for a total applied driving
force of $\Theta=2kT$ by constructing a series of state-space networks with
different topologies. For each network, we considered two pairs of passive
parameters $\mathbf{c}$ and $\mathbf{c}^{\prime}$. First, we generated random
combinations of passive parameters $\mathbf{c}$ and $\mathbf{c}^{\prime}$, and
calculated their information gain (shown in purple in Fig. 3A). While these
randomly generated sets of passive parameters approached the bound in some
regions, they did not reach it for low or high values of the equilibrium
information gain. To systematically test the tightness of the bound in these
regimes, we optimized $\mathbf{c}$ and $\mathbf{c}^{\prime}$ to maximize or
minimize $D(\Theta;\mathbf{c}^{\prime}\|\mathbf{c})$ given a particular value
of the equilibrium information gain $D(0;\mathbf{c}^{\prime}\|\mathbf{c})$
(shown in gray in Fig. 3A) (Supplementary Note 5). While these systems cannot
approach the bound arbitrarily tightly (Supplementary Note 4), optimization
brings them quite close to the upper and lower bounds.

Figure 3:  Numerical optimization illustrates the tightness of the information
gain bound both near and far from equilibrium. (A) Equilibrium and non-
equilibrium information gain between pairs of systems (denoted by parameter
sets $\mathbf{c}$ and $\mathbf{c}^{\prime}$). Randomly chosen $\mathbf{c}$ and
$\mathbf{c}^{\prime}$ are shown in purple; optimized $\mathbf{c}$ and
$\mathbf{c}^{\prime}$ are shown in gray. The
$D(\Theta;\mathbf{c}^{\prime}\|\mathbf{c}),D(0;\mathbf{c}^{\prime}\|\mathbf{c})$
plane was divided into narrow bins of equal width (on a logarithmic scale).
Points from optimization of network pairs and simulation pairs were placed
into these bins. For visibility, a small subset of the bins with at least one
point are shown here with equal probability. Bounds are indicated by dashed
lines. (B) Achievable values of equilibrium and non-equilibrium information
gain are shown for different values of the total applied driving force
$\Theta$. Each filled point corresponds to a pair of passive parameters, and
dashed lines denote the corresponding bound for each value of $\Theta$. The
locations of representative systems shown in (C) are indicated by L and E .
(C) Representative pairs of four-state completely connected systems with an
applied driving force of $\Theta=5\>k_{B}T$ in the linear L and exponential E
regimes.

#### The information gain bound is relevant far from equilibrium.

Having verified the tightness of the bound for a modest value of the applied
driving force, we sought to test whether the bound was also relevant far from
equilibrium. Focusing on a completely connected four-state network, we
optimized $\mathbf{c}$ and $\mathbf{c}^{\prime}$ for a range of applied
driving forces up to $15\>kT$, comparable to the driving force supplied by the
hydrolysis of ATP or GTP (?). We found that the bound was similarly tight even
for systems driven far from equilibrium (Fig. 3B).

We examined representative sets of passive parameters $\mathbf{c}$ for the
four-state complete network that generated information gain values in the
regions where either the linear or exponential terms in the bound dominate for
$\Theta=5\>kT$ (Fig. 3C). In the linear region, the steady-state probability
distributions for the $\mathbf{c}$ and $\mathbf{c}^{\prime}$ parameter sets
are nearly identical for the equilibrium case, and are different in the non-
equilibrium case. The probability currents also change qualitatively from
$\mathbf{c}$ to $\mathbf{c}^{\prime}$. In the exponential region, the
probability distribution differs slightly between $\mathbf{c}$ and
$\mathbf{c}^{\prime}$ in the equilibrium case, and changes drastically in the
non-equilibrium case; again, the pattern of currents changes. Taken together,
these data indicate that the changes in passive parameters needed to closely
approach the bound are non-trivial, leading to qualitative variation in
probability and current distributions in different information gain regimes.

#### Biological systems can closely approach signal amplification and
information gain bounds.

Having shown that synthetic systems can closely approach the bounds on signal
amplification and information gain, we compared the performance of real
biological systems to these bounds.

We first considered the Ras signaling switch, which transduces signals
involved in cell growth and differentiation. Ras can be modeled as a four-
state system with one driven transition, consisting of a cycle of GTP binding
and hydrolysis, and GDP and phosphate unbinding (Fig. 4A). The presence of a
guanine nucleotide exchange factor (GEF) affects the rate of GTP and GDP
(un)binding, and the presence of a GTPase activating protein (GAP) affects the
rate of GTP hydrolysis. Thus, the presence or absence of GEF and GAP change
the passive parameters of the system, and hence the probabilities of
individual states (Fig. S1).

The probability of the GTP-bound“on” state, in which Ras can activate
downstream machinery, is the biologically meaningful quantity. Hence, the
relevant bound is that on the amplification of the probability of that single
state, or the signal ratio. While Eq. 5 represents the most general form of
the bound on the probability of a single state, which is valid for any
arbitrarily driven network, we can derive a tighter bound for the specific
case of Ras:
$\big{|}\ln{\left(\frac{p_{i}(\Theta,\mathbf{c})}{p_{i}(\Theta,\mathbf{c}^{\prime})}\right)}-\ln{\left(\frac{p_{i}(0,\mathbf{c})}{p_{i}(0,\mathbf{c}^{\prime})}\right)}\big{|}\leq\frac{\Theta}{kT}$
(Supplementary Note 6). Holding each rate constant fixed at previously
measured values (Table S1), we ask to what extent GAP or GEF can alter the
probability of the “on” state. In Fig. 4B, we show that the signal ratio in
the presence of GEF versus the presence of GAP is tightly bounded by
$e^{\Theta/kT}$ for modest values of the applied driving force, and reaches a
plateau at higher values.

We next considered the polymerization of microtubules, which are central to
cytoskeletal dynamics and cell division. The appropriate measure of
differences between microtubule length distributions is their KL divergence,
or information gain. The elementary steps in our simplified model of
microtubule self-assembly are shown in Fig. 4C. The reversible addition of
GTP-bound tubulin dimers occurs with forward and reverse rate constants
$\kappa_{f}+\alpha$ and $\kappa_{r}$, respectively (Supplementary Note 7). At
any length, the growing microtubule can undergo rapid and complete
disassembly, known as “catastrophe without rescue”, with frequency
$f_{\mathrm{cat}}$ (?). We mapped this process to the circuit framework to
obtain a closed-form expression for the steady-state length distribution of
microtubules as a function of all elementary parameters (Supplementary Note 7
and Fig. S2-S3).

The rate of catastrophe varies between different phases of the cell cycle (?),
and is known to be modulated by applied load and binding of regulatory
proteins (?). Given the fixed driving force associated with GTP hydrolysis in
vivo and the measured rate constants for the (de)polymerization of tubulin
(Table S2), we asked to what extent changes in $f_{\mathrm{cat}}$ alter the
length distribution of microtubules. At equilibrium, the microtubule length
distribution is independent of $f_{\mathrm{cat}}$; hence, the equilibrium
information gain is zero, and the upper bound for non-equilibrium information
gain is $2\Theta/kT$.

We computed the information gain between length distributions for typical
values of $f_{\mathrm{cat}}$ during interphase and during mitosis (Fig. 4D).
The information gain increases with increasing driving force for modest values
of $\Theta$, then reaches a plateau for physiological values, before finally
decreasing at high values, where all systems are dominated by microtubules
approaching the maximum length. Given that numerous factors can increase
catastrophe frequency (?), we also considered the limit of
$f_{\mathrm{cat}}\to\infty$ (Fig. 4D); the information gain between the value
of $f_{\mathrm{cat}}^{\mathrm{interphase}}$ and $f_{\mathrm{cat}}\to\infty$
closely approached the bound for large $\Theta$. A notable feature of this
behavior is that the relative distance to the bound decreases as the system is
driven more strongly.

Taken together, the analyses of Ras and microtubule polymerization show that
biological systems can closely approach signal amplification and information
gain bounds within particular ranges of the applied driving force, even with a
limited set of control parameters.

Figure 4: Comparison of signaling and self-assembly performance to
amplification and information gain bounds. (A) Structure of the Ras signaling
network. Cofactors are omitted from reverse reactions for clarity. GAP and GEF
each alter rate constants of specific processes; their presence or absence are
the control parameters of the system. (B) The ratio of probability of the “on”
state in the presence of GEF versus the presence of GAP, or the “signal
ratio”, as a function of applied driving force. (C) Structure of the state
space for microtubule polymerization. The catastrophe rate $k_{\text{cat}}$ is
the single control parameter. (D) Nonequilibrium information gain between
systems with different catastrophe rates.

#### Summary and implications

Previous work has highlighted how non-equilibrium driving affects information
processing in biological systems. For example, it has been shown that energy
expenditure generates the sharp morphogen profiles necessary for eukaryotic
development which are inaccessible at equilibrium (?) and enables error
reduction in copying processes (?, ?). Considerable theoretical work has
focused on understanding how energy expenditure constrains the fluctuations or
dynamic response of a system (?, ?, ?, ?, ?, ?).

However, it remains unclear how applied driving forces constrain changes in
probability distributions of arbitrary systems in response to changes in
passive system parameters. This fundamental question about the thermodynamic
limits on biomolecular computation is not addressed by the Landauer energy
requirement for irreversible computation (?), which does not apply to the
class of logically reversible transitions between different steady states (?,
?, ?) for which the probability distribution is completely determined by the
system parameters $\mathbf{c}$ and $\bm{\theta}$.

Here we have derived bounds on how much systems can change the probability of
individual states and alter their overall information content depending on the
strength of the total applied driving force. These bounds apply to all
Markovian systems at steady state, regardless of size or complexity. We
demonstrate that the probability amplification and suppression bound is tight,
and that diverse molecular systems can closely approach the information gain
bound, which is relevant for systems both near and far from equilibrium.

Our work places a lower bound on the thermodynamic force necessary to increase
information gain beyond the extent possible by equilibrium computation. In
practice, information processing in the cell is also limited by noisy
integration of signals and sub-sampling of probability distributions (?).
Pioneering work by Berg and Purcell (?) showed that for measurement of an
input signal, the minimum noise for a fixed mean occupancy of the receptor is
limited by the integration time. In a similar spirit, we show that the change
in the output signal corresponding to any fixed molecular computation is
limited by the thermodynamic force.

Finally, we show that canonical examples of biological signaling and self-
assembly can function near the limits of signal amplification and information
gain derived here. The plateau of the signal ratio and information gain near
the biologically relevant values of $\Theta\approx 15-20\>kT$ (?, ?) suggests
a robustness to natural fluctuations of intracellular GTP and GDP levels,
which provide a driving force that can vary by about $1\>kT$ in vivo (?).

Given the energetic efficiency of equilibrium computation, it is perhaps
surprising that biological systems rely so heavily on non-equilibrium
processes. Our analysis offers an highly general explanation for this
strategy, and examples of systems that can take advantage of it. Future work
will reveal the extent to which living systems consistently exploit the
additional information gain allowed by the bound.

## References

  * 1. A. Sali, E. Shakhnovich, M. Karplus, Nature 369, 248 (1994).
  * 2. R. Otten, et al., Science 370, 1442 (2020).
  * 3. N. R. Latorraca, A. Venkatakrishnan, R. O. Dror, Chemical Reviews 117, 139 (2017).
  * 4. D. K. Simanshu, D. V. Nissley, F. McCormick, Cell 170, 17 (2017).
  * 5. J. Brugués, V. Nuzzo, E. Mazur, D. J. Needleman, Cell 149, 554 (2012).
  * 6. B. H. Morrow, G. F. Payne, J. Shen, Journal of the American Chemical Society 137, 13024 (2015).
  * 7. B. Yurke, A. J. Turberfield, A. P. Mills Jr, F. C. Simmel, J. L. Neumann, Nature 406, 605 (2000).
  * 8. L. Qian, E. Winfree, Science 332, 1196 (2011).
  * 9. M. M. Lin, Physical Review Letters 125, 218101 (2020).
  * 10. J. Schnakenberg, Reviews of Modern Physics 48, 571 (1976).
  * 11. U. Seifert, Reports on Progress in Physics 75, 126001 (2012).
  * 12. C. Maes, K. Netočnỳ, Annales Henri Poincaré (Springer, 2013), vol. 14, pp. 1193–1202.
  * 13. U. Çetiner, J. Gunawardena, Physical Review E 106, 064128 (2022).
  * 14. S. Liang, P. D. L. Rios, D. M. Busiello, arXiv preprint arXiv:2212.12074 (2022).
  * 15. B. M. Martins, P. S. Swain, PLoS Computational Biology 9, e1003175 (2013).
  * 16. D. L. Nelson, A. L. Lehninger, M. M. Cox, Lehninger Principles of Biochemistry (Macmillan, 2008).
  * 17. F.-A. Piedra, et al., Molecular Biology of the Cell 27, 3515 (2016).
  * 18. L. Belmont, A. Hyman, K. Sawin, T. Mitchison, Cell 62, 579 (1990).
  * 19. N. B. Gudimchuk, J. R. McIntosh, Nature Reviews Molecular Cell Biology 22, 777 (2021).
  * 20. J. Estrada, F. Wong, A. DePace, J. Gunawardena, Cell 166, 234 (2016).
  * 21. J. J. Hopfield, Proceedings of the National Academy of Sciences 71, 4135 (1974).
  * 22. Y. Song, C. Hyeon, The Journal of Physical Chemistry Letters 11, 3136 (2020).
  * 23. A. C. Barato, U. Seifert, Physical Review Letters 114, 158101 (2015).
  * 24. T. R. Gingrich, J. M. Horowitz, N. Perunov, J. L. England, Physical Review Letters 116, 120601 (2016).
  * 25. I. Lestas, G. Vinnicombe, J. Paulsson, Nature 467, 174 (2010).
  * 26. J. M. Horowitz, T. R. Gingrich, Nature Physics 16, 15 (2020).
  * 27. N. Ohga, S. Ito, A. Kolchinsky, Physical Review Letters 131, 077101 (2023).
  * 28. J. M. Horowitz, K. Zhou, J. L. England, Physical Review E 95, 042102 (2017).
  * 29. R. Landauer, IBM Journal of Research and Development 5, 183 (1961).
  * 30. E. Fredkin, T. Toffoli, International Journal of Theoretical Physics 21, 219 (1982).
  * 31. C. H. Bennett, International Journal of Theoretical Physics 21, 905 (1982).
  * 32. C. H. Bennett, Studies In History and Philosophy of Science Part B: Studies In History and Philosophy of Modern Physics 34, 501 (2003).
  * 33. M. Razo-Mejia, et al., Physical Review E 102, 022404 (2020).
  * 34. H. C. Berg, E. M. Purcell, Biophysical Journal 20, 193 (1977).
  * 35. D. Zala, et al., F1000Research 6 (2017).
  * 36. E. Ahn, P. Kumar, D. Mukha, A. Tzur, T. Shlomi, Molecular Systems Biology 13, 953 (2017).
  * 37. S. E. Neal, J. F. Eccleston, A. Hall, M. R. Webb, Journal of Biological Chemistry 263, 19718 (1988).
  * 38. C. Allin, M. R. Ahmadian, A. Wittinghofer, K. Gerwert, Proceedings of the National Academy of Sciences 98, 7754 (2001).
  * 39. T. Rudack, F. Xia, J. Schlitter, C. Kötting, K. Gerwert, Proceedings of the National Academy of Sciences 109, 15295 (2012).
  * 40. J. O. Park, et al., Nature Chemical Biology 12, 482 (2016).
  * 41. S. Kunzelmann, M. R. Webb, Biochemical Journal 440, 43 (2011).
  * 42. A. Desai, T. J. Mitchison, Annual Review of Cell and Developmental Biology 13, 83 (1997).
  * 43. L. Brun, B. Rupp, J. J. Ward, F. Nédélec, Proceedings of the National Academy of Sciences 106, 21173 (2009).
  * 44. T. Mitchison, M. Kirschner, Nature 312, 237 (1984).
  * 45. H. Bowne-Anderson, M. Zanic, M. Kauer, J. Howard, BioEssays 35, 452 (2013).
  * 46. D. L. Gard, M. W. Kirschner, The Journal of Cell Biology 105, 2191 (1987).

Acknowledgements: We thank members of the Lin lab for useful discussions;
Peter Foster, Jordan Horowitz, and Daniel Needleman for feedback; and Yong
Hyun Song for a critical reading of the manuscript. We acknowledge the Texas
Advanced Computing Center (TACC) at the University of Texas at Austin for
providing HPC resources that have contributed to the research results reported
within this paper. Some graphical elements were adapted from BioRender.

Funding: This work was supported by the Chan Zuckerberg Initiative Theory in
Biology grant (M.M.L.) and a National Science Foundation Graduate Research
Fellowship (E.A.).

Author contributions: E.A. and M.M.L. derived results, performed calculations,
and wrote the paper.

Competing interests: Authors declare that they have no competing interests.

Data and materials availability: Code used for calculations has been deposited
on GitHub (link to be added).

Supplementary Materials

Supplementary Text

Figs. S1 to S3

Tables S1 to S2

References 37-46

## Supplementary Information

##  1 Circuit representation

For systems of nontrivial complexity, obtaining general insights into steady-
state behavior using the rate constant framework is typically analytically
intractable (?). It is convenient to transform variables from rate constants
to an equivalent circuit representation described by the parameters of the
equilibrium system and the thermodynamic forces that drive the system away
from equilibrium (?).

In this representation, each state has an associated free energy; the free
energies of states $i$ and $j$ are related to equilibrium rate constants by:

$\displaystyle\kappa_{ij}e^{-G_{i}/kT}=\kappa_{ji}e^{-G_{j}/kT},$ (S1)

where $kT$ is the product of the Boltzmann constant and the temperature, which
sets the energy scale for the system. The free energies provide a complete
description of systems at equilibrium because the probability of such a system
visiting any state $i$ relative to any other state $j$ is given by the
Boltzmann distribution: $p_{i}=p_{j}e^{(G_{j}-G_{i})/kT}$. The currents that
emerge in response to non-equilibrium driving are tuned by the resistance $R$
associated with each transition, defined by:

$\displaystyle R_{ij}=e^{G_{i}/kT}/\kappa_{ij}=R_{ji}$ (S2)

Intuitively, the resistance between two states is inversely proportional to
the forward flux (equal to the reverse flux) in the equilibrium system. Note
that the forward flux from state $i$ to $j$ at equilibrium (equal to the
reverse flux from $j$ to $i$) is equal to $Z/R_{ij}$ where the normalization
constant $Z$ is the partition function.

The applied driving force between states $i$ and $j$, denoted by
$\theta_{ij}$, depends on the ratio of active and passive rate constants
between states $i$ and $j$. $\theta_{ij}$ takes the form of a battery in the
circuit representation. The voltage $\mathcal{E}$ between $i$ and $j$, which
is a function of the applied driving force between those states, is given by:

$\displaystyle\begin{split}\mathcal{E}_{ij}=e^{G_{i}/kT}\left(e^{\theta_{ij}/kT}-1\right)p_{i}.\end{split}$
(S3)

The circuit equation, which is an Ohm’s law generalization of the Boltzmann
distribution to non-equilibrium systems (?), relates the applied driving force
to the induced currents:

$\displaystyle\begin{split}p_{i}e^{G_{i}/kT}-p_{j}e^{G_{j}/kT}&=\sum_{mn}\left(\mathcal{E}_{mn}-I_{mn}R_{mn}\right),\end{split}$
(S4)

where the sum is taken along any path connecting $i$ to $j$ parameterized by a
series of adjacent states $m$ and $n$ along the path.

The circuit representation naturally decomposes system parameters into
independently tunable categories: passive parameters which control equilibrium
behavior (the free energies), passive parameters which tune probability flows
when the system is not at equilibrium (the resistances), and active parameters
(applied driving forces). In the main text, we denote by $\mathbf{c}$ the set
of all passive parameters in the system (the free energy $G_{i}$ of each state
$i$ and the resistance $R_{ij}$ between each pair of adjacent states $i$ and
$j$). This parameterization ensures that all free energies and resistances can
be independently tuned while preserving detailed balance in the absence of
applied driving forces.

##  2 Relationship between applied driving forces and chemical potential
differences

Consider a loop in a state-space network where edges are indexed by $ij$
(indicating an edge between vertex $i$ and vertex $j$). The sum of the
chemical potential differences along each edge in the loop is given by:

$\displaystyle\sum_{ij\in\text{ loop}}\Delta\mu_{ij}$
$\displaystyle=\sum_{ij}kT\ln\left(\dfrac{J_{ij}}{J_{ji}}\right)=kT\sum_{ij}\ln\left(\dfrac{p_{i}\left(\kappa_{ij}+\alpha_{ij}\right)}{p_{j}\left(\kappa_{ji}+\alpha_{ji}\right)}\right)$
(S5)
$\displaystyle=kT\left[\sum_{ij}\ln\left(\dfrac{1+\alpha_{ij}/\kappa_{ij}}{1+\alpha_{ji}/\kappa_{ji}}\right)+\ln\left(\prod_{ij}\dfrac{p_{i}}{p_{j}}\right)+\ln\left(\prod_{ij}\dfrac{\kappa_{ij}}{\kappa_{ji}}\right)\right]$

We note that for a sum over the edges in a loop, the probability of each state
enters once in the numerator and once in the denominator, and that for every
pathway connecting states $m$ and $n$, the ratios of the forward and reverse
product of equilibrium rate constants is unity:

$\displaystyle\dfrac{\kappa_{m,m+1}\kappa_{m+1,m+2}\cdots\kappa_{n-2,n-1}\kappa_{n-1,n}}{\kappa_{m+1,m}\kappa_{m+2,m+1}\cdots\kappa_{n-1,n-2}\kappa_{n,n-1}}$
$\displaystyle=1$

Therefore,

$\displaystyle\sum_{ij\in\text{ loop}}\Delta\mu_{ij}$
$\displaystyle=kT\sum_{ij\in\text{ loop}}\left(\theta_{ij}-\theta_{ji}\right)$
(S6)

##  3 Proof of the probability amplification limit

Let us consider an arbitrary network (or graph) $G$ with $N$ states. Following
(?), we define a maximal tree $T(G)$ of our network $G$ to have the following
properties:

  1. 1.

$T(G)$ is a covering subgraph of $G$, that is, all edges of $T(G)$ are edges
of $G$, and $T(G)$ contains all vertices of $G$;

  2. 2.

$T(G)$ is connected;

  3. 3.

$T(G)$ is acyclic.

The $M$ maximal trees are indexed by $\mu$: $T^{\mu}(G)$,
$\mu=\left\\{1..M\right\\}$. Using each $T^{\mu}(G)$ we now construct a series
of directed spanning trees $T^{\mu}_{i}(G)$, $\left\\{1..N\right\\}$ by
assigning a direction to each of its edges so that they all point towards
vertex $i$. Therefore, each edge is associated with a directed rate constant
pointing in the direction of the edge.

We now define a function $\Pi(T^{\mu}_{i}(G))$, which gives the product of all
the associated directed rate constants in $T^{\mu}_{i}(G)$, which we call the
tree product.

From Kirchoff’s matrix-tree theorem we have that the steady-state probability
of state $i$ is given by

$\displaystyle p_{i}$ $\displaystyle=\dfrac{S_{i}}{S},\qquad\text{ where
}S_{i}\equiv\sum_{\mu=1}^{M}\Pi(T_{i}^{\mu}(G))\text{ and
}S\equiv\sum_{i=1}^{N}S_{i}$ (S7)

Define $A_{i}^{\mu}$ to be the reference equilibrium tree product of tree
$\mu$ rooted at $i$, and $A_{i}^{\mu}\Gamma_{i}^{\mu}$ to be the corresponding
non-equilibrium tree product. (i.e. $\Gamma_{i}^{\mu}$ is the ratio of the
non-equilibrium tree product to the equilibrium tree product.) Then, the non-
equilibrium probability of state $i$ can be expressed as:

$p_{i}=\dfrac{\sum_{\mu}A_{i}^{\mu}\Gamma_{i}^{\mu}}{\sum_{\mu}\sum_{j}A_{j}^{\mu}\Gamma_{j}^{\mu}}$
(S8)

For a fixed spanning tree $\mu$, changing the vertex from $i$ to $j$ can only
change the tree product in the terms along the path within the spanning tree
that connects $i$ and $j$. Therefore,

$\dfrac{\Gamma_{j}^{\mu}}{\Gamma_{i}^{\mu}}=\dfrac{\prod_{mn\in(i\to
j\mid\mu)}\left(1+\frac{\alpha_{mn}}{\kappa_{mn}}\right)}{\prod_{nm\in(j\to
i\mid\mu)}\left(1+\frac{\alpha_{nm}}{\kappa_{nm}}\right)}$

where the product over $mn\in(i\to j\mid\mu)$ indicates the edges (indexed by
$mn$) present in the path within spanning tree $\mu$ that connects vertex $i$
to vertex $j$. Thus,

$e^{-\frac{\Theta}{kT}}\leq\dfrac{\Gamma_{j}^{\mu}}{\Gamma_{i}^{\mu}}\leq
e^{\frac{\Theta}{kT}}$ (S9)

$\Theta$ is the maximum sum of $\ln\left(1+\frac{\alpha_{mn}}{k_{mn}}\right)$
along any path across all spanning trees. Therefore,

$p_{i}\leq
e^{\frac{\Theta}{kT}}\dfrac{\sum_{\mu}A_{i}^{\mu}\Gamma_{i}^{\mu}}{\sum_{\mu}\sum_{j}A_{j}^{\mu}\Gamma_{i}^{\mu}}$

Note that

$\dfrac{A_{i}^{\mu}}{A_{j}^{\mu}}=e^{\Delta G_{ij}/kT},$

where $\Delta G_{ij}$ is the equilibrium free energy difference between states
$i$ and $j$, and is thus independent of $\mu$. Therefore,

$p_{i}\leq
e^{\frac{\Theta}{kT}}\dfrac{\sum_{\mu}A_{i}^{\mu}\Gamma_{i}^{\mu}}{\sum_{\mu}\sum_{j}e^{-\Delta
G_{ij}/kT}A_{i}^{\mu}\Gamma_{i}^{\mu}}$

Furthermore, noting that $\sum_{j}e^{-\Delta G_{ij}/kT}=1/p_{i}^{*}$, where
$p_{i}^{*}$ is the equilibrium probability of state $i$, we finally have:

$p_{i}\leq
p_{i}^{*}e^{\frac{\Theta}{kT}}\dfrac{\sum_{\mu}A_{i}^{\mu}\Gamma_{i}^{\mu}}{\sum_{\mu}A_{i}^{\mu}\Gamma_{i}^{\mu}}=p_{i}^{*}e^{\frac{\Theta}{kT}}.$

Following the same logic in reverse, we have that

$\displaystyle p_{i}\geq
e^{-\frac{\Theta}{kT}}\dfrac{\sum_{\mu}A_{i}^{\mu}\Gamma_{i}^{\mu}}{\sum_{\mu}\sum_{j}A_{j}^{\mu}\Gamma_{i}^{\mu}}=e^{-\frac{\Theta}{kT}}\dfrac{\sum_{\mu}A_{i}^{\mu}\Gamma_{i}^{\mu}}{\sum_{\mu}\sum_{j}e^{-\Delta
G_{ij}/kT}A_{i}^{\mu}\Gamma_{i}^{\mu}}=p_{i}^{*}e^{-\frac{\Theta}{kT}}$

Thus, we have upper and lower bounds for the change in the probability of a
single state $i$ upon driving:

$\displaystyle e^{-\frac{\Theta}{kT}}\leq\dfrac{p_{i}}{p_{i}^{*}}\leq
e^{\frac{\Theta}{kT}}$ (S10)

This bound is fully general and applies to arbitrary systems. It considers the
probability of a given state in isolation, and therefore is purely
multiplicative. However, there is an additional constraint on probabilities
both in the presence and absence of driving: the probabilities of all states
must sum to unity. This can, for specific choices of rate constants, further
restrict probability amplification and suppression.

##  4 Proof of the information gain limit

Starting from the definition of the non-equilibrium information gain, we use
the probability amplification limit and add and subtract identical terms:

$\displaystyle\small\begin{split}D(\Theta;\mathbf{c}\|\mathbf{c}^{\prime})&=\sum_{i=1}^{N}p_{i}(\Theta,\mathbf{c})\ln{\big{[}{p_{i}(\Theta,\mathbf{c})}\big{]}}-\sum_{i=1}^{N}p_{i}(\Theta,\mathbf{c})\ln{\big{[}{p_{i}(\Theta,\mathbf{c}^{\prime})}\big{]}}\\\
&\qquad-\sum_{i=1}^{N}p_{i}(\Theta,\mathbf{c})\ln{\big{[}{p_{i}(0,\mathbf{c})}\big{]}}+\sum_{i=1}^{N}p_{i}(\Theta,\mathbf{c})\ln{\big{[}{p_{i}(0,\mathbf{c})}\big{]}}\\\
&\qquad+\sum_{i=1}^{N}p_{i}(\Theta,\mathbf{c})\ln{\big{[}{p_{i}(0,\mathbf{c}^{\prime})}\big{]}}-\sum_{i=1}^{N}p_{i}(\Theta,\mathbf{c})\ln{\big{[}{p_{i}(0,\mathbf{c}^{\prime})}\big{]}}\\\
&=\underbrace{\sum_{i=1}^{N}p_{i}(\Theta,\mathbf{c})\ln{\Big{[}\dfrac{p_{i}(\Theta,\mathbf{c})}{p_{i}(0,\mathbf{c})}\Big{]}}}_{<\Theta/kT}-\underbrace{\sum_{i=1}^{N}p_{i}(\Theta,\mathbf{c})\ln{\Big{[}\dfrac{p_{i}(\Theta,\mathbf{c}^{\prime})}{p_{i}(0,\mathbf{c}^{\prime})}\Big{]}}}_{>-\Theta/kT}+\sum_{i=1}^{N}p_{i}(\Theta,\mathbf{c})\ln{\Big{[}\dfrac{p_{i}(0,\mathbf{c})}{p_{i}(0,\mathbf{c}^{\prime})}\Big{]}}\\\
&<\dfrac{2\Theta}{kT}+\sum_{i=1}^{N}p_{i}(\Theta,\mathbf{c})\ln{\Big{[}\dfrac{p_{i}(0,\mathbf{c})}{p_{i}(0,\mathbf{c}^{\prime})}\Big{]}}\\\
&<\dfrac{2\Theta}{kT}+\sum_{i=1}^{N}e^{\frac{\Theta}{kT}}p_{i}(0,\mathbf{c})\ln{\Big{[}\dfrac{p_{i}(0,\mathbf{c})}{p_{i}(0,\mathbf{c}^{\prime})}\Big{]}}\end{split}$

The maximum information gain per computational step of the machine therefore
scales exponentially with the total applied driving force $\Theta$:

$D(\Theta;\mathbf{c}\|\mathbf{c}^{\prime})\leq\dfrac{2\Theta}{kT}+D(0;\mathbf{c}\|\mathbf{c}^{\prime})e^{\frac{\Theta}{kT}}$
(S11)

Similarly, we find a lower bound for the non-equilibrium information gain
(equivalent to an upper bound for the equilibrium information gain in terms of
the non-equilibrium information gain):

$\displaystyle\small\begin{split}D(0;\mathbf{c}\|\mathbf{c}^{\prime})&=\underbrace{\sum_{i=1}^{N}p_{i}(0,\mathbf{c})\ln{\Big{[}\dfrac{p_{i}(0,\mathbf{c})}{p_{i}(\Theta,\mathbf{c})}}\Big{]}}_{\leq\Theta/kT}-\underbrace{\sum_{i=1}^{N}p_{i}(\Theta,\mathbf{c})\ln{\Big{[}\dfrac{p_{i}(0,\mathbf{c}^{\prime})}{p_{i}(\Theta,\mathbf{c}^{\prime})}\Big{]}}}_{\geq-\Theta/kT}+\sum_{i=1}^{N}p_{i}(0,\mathbf{c})\ln{\Big{[}\dfrac{p_{i}(\Theta,\mathbf{c})}{p_{i}(\Theta,\mathbf{c}^{\prime})}\Big{]}}\\\
&\leq\dfrac{2\Theta}{kT}+\sum_{i=1}^{N}p_{i}(0,\mathbf{c})\ln{\Big{[}\dfrac{p_{i}(\Theta,\mathbf{c})}{p_{i}(\Theta,\mathbf{c}^{\prime})}\Big{]}}\\\
&\leq\dfrac{2\Theta}{kT}+\sum_{i=1}^{N}e^{-\frac{\Theta}{kT}}p_{i}(\Theta,\mathbf{c})\ln{\Big{[}\dfrac{p_{i}(\Theta,\mathbf{c})}{p_{i}(\Theta,\mathbf{c}^{\prime})}\Big{]}}\end{split}$

Hence, the minimum information gain per computational step of the machine
similarly scales exponentially with the total applied driving force $\Theta$:

$D(0;\mathbf{c}\|\mathbf{c}^{\prime})\leq\dfrac{2\Theta}{kT}+D(\Theta;\mathbf{c}\|\mathbf{c}^{\prime})e^{\frac{\Theta}{kT}}$
(S12)

Note that in this derivation, we have assumed that the only bound on the
amplification and suppression of individual probabilities is the bound in Eq.
S10. However, as noted in Sec. 3, there is a further constraint that the
probabilities of all states must sum to unity. Thus, not all states can be
amplified or suppressed to the maximal extent ($e^{\pm\Theta/kT}$). Despite
the fact that we do not consider this complication, our optimizations approach
the information gain bound reasonably closely – but the bound is not tight.

##  5 Monte Carlo optimization of state-space networks

The set of states (vertices) in each state-space network (graph) is denoted by
$V$, and the set of edges is denoted by $E$. The steady-state probabilities of
each state were expressed analytically in terms of the free energy of each
state and the resistances associated with transitions between states. Here, a
“system” denotes a particular set of passive parameters (free energies and
resistances) for a state-space network. We sought to identify pairs of systems
$\mathbf{c}$ and $\mathbf{c}^{\prime}$ that yielded equilibrium and non-
equilibrium steady-state probability distributions that differed from one
another by a specific amount.

In particular, we minimized a harmonic penalty associated with the distance of
a particular system’s probability distributions in the
$D(\Theta;\mathbf{c}^{\prime}\|\mathbf{c})$ vs.
$D(0;\mathbf{c}^{\prime}\|\mathbf{c})$ plane to a target point
$\left(D^{*}(0;\mathbf{c}^{\prime}\|\mathbf{c}),D^{*}(\Theta;\mathbf{c}^{\prime}\|\mathbf{c})\right)$:

$\displaystyle L$
$\displaystyle=\left[D(0;\mathbf{c}^{\prime}\|\mathbf{c})-D^{*}(0;\mathbf{c}^{\prime}\|\mathbf{c})\right]^{2}+\left[D(\Theta;\mathbf{c}^{\prime}\|\mathbf{c})-D^{*}(\Theta;\mathbf{c}^{\prime}\|\mathbf{c})\right]^{2}$
(S13)

Monte Carlo optimization over the vectors of parameters
$\mathbf{c}=\left\\{G_{i}\right\\}_{i\in V}\cup\left\\{R_{ij}\right\\}_{ij\in
E}$ and $\mathbf{c}^{\prime}=\left\\{G_{i}^{\prime}\right\\}_{i\in
V}\cup\left\\{R_{ij}^{\prime}\right\\}_{ij\in E}$ was performed with an
inverse temperature of $\beta=4\times 10^{4}-1\times 10^{5}/kT$ for $4\times
10^{3}-4\times 10^{4}$ time steps in Wolfram Mathematica 12.3. Optimization to
target points in inaccessible regions beyond the information gain and noise
suppression bounds only yielded points that lay within the accessible region.

##  6 Analysis of Ras signaling

Parameters for the four-state model of Ras are shown in Table S1. The rate
constants for the Ras cycle were obtained from Refs. (?, ?, ?) and nucleotide
concentrations were obtained from (?). The equilibrium GTP concentration was
calculated from the non-equilibrium reaction quotient and the physiological
free energy of GTP hydrolysis, or phosphorylation potential, $\Delta
G_{\mathrm{P}}$. The phosphate association rate constant $\kappa_{43}$ was set
such that the system reached detailed balance at the equilibrium ratio of
[GTP]/[GDP]. GTPase-activating proteins (GAP) increase the rate of GTP
hydrolysis $\kappa_{32}$ by a factor of $10^{5}$ (?). To implement detailed
balance, we assume that this lowering of the free energy barrier causes a
concomitant increase in $\kappa_{23}$. Guanine nucleotide exchange factors
(GEF), such as SOS, stabilize the apo state relative to GTP and GDP-bound
states, leading to a $10^{4}-10^{5}$-fold increase in the rates of GTP and GDP
dissociation (?).

Table S1: Kinetic and thermodynamic parameters for the Ras model Parameter | Value | Reference
---|---|---
$\kappa_{12}$ (M-1 s-1) | $1.40\times 10^{8}$ |  (?)
$\kappa_{21}$ (-GEF) (s-1) | $1.00\times 10^{-4}$ |  (?)
$\kappa_{21}$ (+GEF) (s-1) | $7.78\times 10^{-1}$ |  (?, ?)
$\kappa_{23}$ (-GAP) (s-1) | $3.40\times 10^{-4}$ |  (?)
$\kappa_{23}$ (+GAP) (s-1) | $3.40\times 10^{1}$ |  (?, ?)
$\kappa_{32}$ (-GAP) (s-1) | $0.16\times 10^{-4}$ |  (?)
$\kappa_{32}$ (+GAP) (s-1) | $1.6$ |  (?, ?)
$\kappa_{34}$ (s-1) | $0.102$ |  (?)
$\kappa_{43}$ (s-1) | $1.98\times 10^{-6}$ | n/a
$\kappa_{41}$ (-GEF) (s-1) | $4.20\times 10^{-4}$ |  (?)
$\kappa_{41}$ (+GEF) (s-1) | $3.27$ |  (?, ?)
$\kappa_{14}$ (M-1 s-1) | $0.51\times 10^{8}$ |  (?)
$\Delta G_{\mathrm{P}}$ (kcal/mol) | $\approx-12$ |  (?)
$\left[\mathrm{GTP}\right]$ (M) | $6.77\times 10^{-4}$ |  (?)
$\left[\mathrm{GTP}\right]_{\mathrm{eq}}$ (M) | $2.39\times 10^{-12}$ |  (?, ?)
$\left[\mathrm{GDP}\right]$ (M) | $3.02\times 10^{-5}$ |  (?)
$\left[\mathrm{P}_{\mathrm{i}}\right]$ (M) | $5\times 10^{-3}$ |  (?)

Using the paramters above, we calculated the probabilities of individual
states the presence and absence of GAP and GEF (Fig. S1.

Figure S1: (A) Ras can be modeled as a four-state system: (1) apo, with no
nucleotide bound; (2) GTP-bound, which is the ”on” state, or signaling state;
(3) GDP + phosphate bound; and (4) GDP-bound. (B) Probabilities of the
individual states are shown as a function of GTP concentration.

We next considered the limits on amplification or suppression of the
probability of a single state possible when changing system parameters. Eq. 5
is the most general form of the bound, and it is valid for any arbitrarily
driven network. However, we can derive a tighter bound for the specific case
of the “on” state of Ras:

$\displaystyle\bigg{|}\ln{\left(\frac{p_{i}(\Theta,\mathbf{c})}{p_{i}(\Theta,\mathbf{c}^{\prime})}\right)}-\ln{\left(\frac{p_{i}(0,\mathbf{c})}{p_{i}(0,\mathbf{c}^{\prime})}\right)}\bigg{|}\leq\frac{\Theta}{kT}.$
(S14)

Note that this is because the driven transition leads directly to the “on”
state, so there does not exist a non-intersecting path in which the driven
transition does not point towards the “on” state. Thus, its probability cannot
be suppressed by driving, i.e. $1\leq
p_{\mathrm{on}}(\Theta,\mathbf{c})/p_{\mathrm{on}}(0,\mathbf{c})\leq
e^{\Theta/kT}$.

##  7 Analysis of microtubule self-assembly

We consider a simplified model of microtubule self-assembly, which includes
addition and subtraction of tubulin dimers (note that a dimer of
$\alpha$-tubulin and $\beta$-tubulin is the building block for the
microtubule), and catastrophe. We neglect rescue from catastrophe.

We define the elementary parameters as follows. If GTP were allowed to
equilibrate with GDP, the reversible assembly of GTP-bound tubulin dimers
would occur with forward and reverse rate constants $\kappa_{f}$ and
$\kappa_{r}$, respectively. Thus, $\kappa_{f}/\kappa_{r}=e^{-\Delta
G_{\mathrm{P}}^{\text{\lx@text@degree}}/kT}$, where $\Delta
G_{\mathrm{P}}^{\text{\lx@text@degree}}$ is the equilibrium dimer binding free
energy. In cells, GTP levels are maintained at a concentration far higher than
the equilibrium level, effectively giving rise to an additional forward rate
constant $\alpha$, which is proportional to excess [GTP].

Catastrophe is triggered by the stochastic disruption of the growing
microtubule cap (?), which allows cap-modifying substrates to act as sub-
stoichimetric catalysts of microtubule shrinkage. The dynamic instability
steady state is reached when the rate of catastrophe balances that of net
dimer addition, in contrast to the detailed balance steady state of an
equilibrium system (?). We next solve for the steady-state microtubule length
distribution in terms of these parameters.

###  7.1 Calculation of length distribution

Figure S2: Circuit diagram for microtubule assembly

Here we employ the circuit representation (?) of the microtubule
polymerization state space shown in Fig. S2. In terms of the $n$th mesh
current ($I_{n}$, indicated by gray arrows), the voltage equation taken along
the path of the $n$th battery (with voltage
$\mathcal{E}\propto\alpha/\kappa_{f}$) is:

$P_{n+1}e^{\beta G_{n+1}}-P_{n}e^{\beta
G_{n}}=\dfrac{\alpha}{\kappa_{f}}P_{n}e^{\beta G_{n}}-R_{n}I_{n}.$ (S15)

where $R_{n}=e^{\beta G_{n}}/\kappa_{f}$ and $G_{n}=nG$. Note that $e^{\beta
G}=k_{b}/\kappa_{f}$, where $\kappa_{f}$ and $\kappa_{r}$ are the equilibrium
forward and reverse rates, respectively. Using these definitions, we can solve
for the probability of the $(n+1)$th state in terms of the previous state
probability and current:

$P_{n+1}=\left(1+\dfrac{\alpha}{\kappa_{f}}\right)e^{-\beta
G}P_{n}-I_{n}\dfrac{e^{-\beta G}}{\kappa_{f}}$

Taking the potential difference from state $n+1$ and state 1 along the
catastrophe path:

$P_{1}e^{\beta G_{1}}-P_{n+1}e^{\beta
G_{n}+1}=-R_{\mathrm{cat},n}(I_{n}-I_{n+1}),$

Where $R_{\mathrm{cat},n}=e^{\beta
G_{n+1}}/f_{\mathrm{cat}}=R_{n+1}(\dfrac{\kappa_{f}}{f_{\mathrm{cat}}})$.
Therefore, the $(n+1)$th current is:

$I_{n+1}=I_{n}-f_{\mathrm{cat}}P_{n+1}+f_{\mathrm{cat}}P_{1}e^{-\beta G_{n}}$

In vector notation, the recursive probability and current equations become:

$\begin{bmatrix}1&0\\\
f_{\mathrm{cat}}&1\end{bmatrix}\begin{bmatrix}P_{n+1}\\\
I_{n+1}\end{bmatrix}=\begin{bmatrix}\left(1+\frac{\alpha}{\kappa_{f}}\right)e^{-\beta
G}&-\frac{e^{-\beta G}}{\kappa_{f}}\\\ 0&1\end{bmatrix}\begin{bmatrix}P_{n}\\\
I_{n}\end{bmatrix}+\begin{bmatrix}0\\\ \frac{f_{\mathrm{cat}}P_{1}}{e^{\beta
G_{n}}}\end{bmatrix}$

Left-multiplying both sides by the inverse of the left-hand-side matrix, the
recursion relation is:

$\begin{bmatrix}P_{n+1}\\\ I_{n+1}\end{bmatrix}=M\begin{bmatrix}P_{n}\\\
I_{n}\end{bmatrix}+\begin{bmatrix}0\\\ {\frac{f_{\mathrm{cat}}P_{1}}{e^{\beta
G_{n}}}}\end{bmatrix}$

where the transition matrix $M$ is given by:

$M=\begin{bmatrix}\left(1+\frac{\alpha}{\kappa_{f}}\right)e^{-\beta
G}&-\frac{e^{-\beta G}}{\kappa_{f}}\\\
-f_{\mathrm{cat}}\left(1+\frac{\alpha}{\kappa_{f}}\right)e^{-\beta
G}&f_{\mathrm{cat}}\frac{e^{-\beta G}}{\kappa_{f}}+1\end{bmatrix}$

The probability and current of state $n$ in terms of those of state 1 is thus:

$\begin{bmatrix}P_{n+1}\\\ I_{n+1}\end{bmatrix}=M^{n}\begin{bmatrix}P_{1}\\\
I_{1}\end{bmatrix}+\sum_{k=0}^{n-1}(e^{\beta G}M)^{k}\begin{bmatrix}0\\\
{{f_{\mathrm{cat}}P_{1}}{e^{-\beta Gn}}}\end{bmatrix}$

Diagonalizing $M$:

$M=V\begin{bmatrix}\lambda_{-}&0\\\ 0&\lambda_{+}\end{bmatrix}V^{-1}$

Where the columns of $V$ are the eigenvectors of $M$ and $\lambda_{-}$ and
$\lambda_{+}$ are the eigenvalues of $M$:

$\lambda_{\pm}=\dfrac{e^{-\beta
G}}{2\kappa_{f}}\Bigg{(}\alpha+\kappa_{f}(1+e^{\beta
G})+f_{\mathrm{cat}}\pm\sqrt{(\alpha+\kappa_{f}-\kappa_{f}e^{\beta
G})^{2}+f_{\mathrm{cat}}(2\alpha+2(1+e^{\beta
G})\kappa_{f}+f_{\mathrm{cat}})}\Bigg{)}$

Which simplifies to:

$D=-\ln{\lambda_{-}}=-\ln{\left[1-\dfrac{\sqrt{(\alpha+f_{\mathrm{cat}}+\kappa_{f}-\kappa_{r})^{2}+4f_{\mathrm{cat}}\kappa_{r}}-(\alpha+f_{\mathrm{cat}}+\kappa_{f}-\kappa_{r})}{2\kappa_{r}}\right]}$

Note that $\lambda_{-}\leq 1$ whereas $\lambda_{+}\geq 1$.

The transfer matrix equation is then

$\begin{bmatrix}P_{n+1}\\\
I_{n+1}\end{bmatrix}=V\begin{bmatrix}\lambda_{-}^{n}&0\\\
0&\lambda_{+}^{n}\end{bmatrix}V^{-1}\begin{bmatrix}P_{1}\\\
I_{1}\end{bmatrix}+\sum_{k=0}^{n-1}e^{\beta
kG}V\begin{bmatrix}\lambda_{-}^{k}&0\\\
0&\lambda_{+}^{k}\end{bmatrix}V^{-1}\begin{bmatrix}0\\\
{{f_{\mathrm{cat}}P_{1}}{e^{-\beta Gn}}}\end{bmatrix}$

Expanding this expression and taking the geometric sum yields $P_{n}$:

$P_{n}=P_{1}e^{-\beta
G(n-1)}{\dfrac{f_{\mathrm{cat}}}{f_{\mathrm{cat}}-\alpha(e^{\beta
G}-1)}}+A_{1}\lambda_{-}^{n-1}+A_{2}\lambda_{+}^{n-1},$

where the $A_{1}$ and $A_{2}$ are explicit functions of the elementary
parameters. For nonzero $f_{\mathrm{cat}}$ the probability monotonically
decreases for larger $n$, thus the coefficient $A_{2}$ must be zero. Solving
this boundary condition for $I_{1}$ and substituting into the expression for
$A_{1}$, we obtain the length distribution:

$P_{n}={\dfrac{f_{\mathrm{cat}}}{f_{\mathrm{cat}}-\alpha(e^{\beta
G}-1)}}P_{1}e^{-\beta G(n-1)}+{\dfrac{\alpha(e^{\beta G}-1)}{\alpha(e^{\beta
G}-1)-f_{\mathrm{cat}}}}P_{1}e^{-D(n-1)},$ (S16)

where $D=-\ln{\lambda_{-}}$. For microtubule assembly, the physiologically
relevant parameters are shown in Table S2.

The first-order equilibrium rate constant for microtubule polymerization was
calculated from the measured second-order non-equilibrium rate constant
($\kappa_{f}^{\prime}$), the tubulin concentration, and the applied driving
force $\Delta G_{\mathrm{P}^{\text{\lx@text@degree}}}$:

$\displaystyle\kappa_{f}$
$\displaystyle=\kappa_{f}^{\prime}\left[\mathrm{tubulin}\right]\exp{\left(-\Delta
G_{\mathrm{P}}/RT\right)}$

Table S2: Kinetic and thermodynamic parameters for microtubule assembly Parameter | Value | Reference
---|---|---
$\kappa_{f}^{\prime}$ (M-1 s-1) | $1.40\times 10^{8}$ |  (?)
$\left[\mathrm{tubulin}\right]$ (M) | $1.00\times 10^{-4}$ |  (?)
$\Delta G_{\mathrm{P}}$ (kcal/mol) | $\approx-12$ |  (?)
$\kappa_{f}$ (s-1) | $1.00\times 10^{-4}$ |  (?, ?, ?)
$\kappa_{r}$ (s-1) | $7.78\times 10^{-1}$ |  (?)
$f_{\mathrm{cat}}^{\mathrm{(mitosis)}}$ (s-1) | $0.12\pm 0.015$ |  (?)
$f_{\mathrm{cat}}^{\mathrm{(interphase)}}$ (s-1) | $0.018\pm 0.011$ |  (?)
$\langle L\rangle^{\mathrm{(mitosis)}}$ ($\text{\,}\mathrm{\SIUnitSymbolMicro m}$) | $6.2\pm 3.0$ |  (?)
$\langle L\rangle^{\mathrm{(interphase)}}$ ($\text{\,}\mathrm{\SIUnitSymbolMicro m}$) | $23\pm 11$ |  (?)

Using the result in Eq. S16 and the parameters in Table S2, we calculate the
mean length of microtubules as a function of catastrophe rate and GTP
concentration (Fig. S3). This is a coarse measure of how the length
distribution is controlled by a single passive parameter and the applied
driving force. In the absence of the driving force, manipulating
$f_{\text{cat}}$ does not change the length distribution; only out of
equilibrium does this have an effect. Thus, for changes in this passive
parameter, the equilibrium information gain
$D(0;f_{\text{cat}}^{\prime}\|f_{\text{cat}})=0$. The bound on the non-
equilibrium information gain is therefore
$D(\Theta;f_{\text{cat}}^{\prime}\|f_{\text{cat}})\leq 2\Theta/kT$.

Figure S3: Mean length of microtubules as a function of catastrophe rate
$f_{\mathrm{cat}}$ and GTP concentration (applied driving force).

###  7.2 Accounting for finite length

Physical constraints such as cell size impose an upper bound on the length of
microtubules, which we take to be 100 $\text{\,}\mathrm{\SIUnitSymbolMicro
m}$. We denote this maximum length by $n_{\text{max}}$. While Eq. S16 was
derived for unbounded $n$, we account for the effect of finite
$n_{\text{max}}$ using the approximate probability distribution
$P_{n}^{\prime}$:

$\displaystyle P^{\prime}_{n}$ $\displaystyle=\begin{cases}P_{n}&1\leq
n<n_{\text{max}}\\\
1-\sum_{n=1}^{n_{\text{max}}-1}P_{n}&n=n_{\text{max}}\end{cases}$ (S17)

Note that this is formally equivalent to leaving the length unrestricted and
coarse-graining all microstates with length greater than or equal to
$n_{\text{max}}$ into a single state. The total applied driving force $\Theta$
is the product of the applied driving force per driven transition $\theta$ and
the number of driven transitions $L_{\mathrm{max}}$.

###  7.3 Simplification in the case of high catastrophe rate

For $f_{\mathrm{cat}}\to\infty$, we obtain a simplified expression for the
length distribution:

$\displaystyle P^{\prime}_{n}$
$\displaystyle=\begin{cases}\left(\frac{\kappa_{r}}{\kappa_{f}}-1\right)\left(\frac{\kappa_{f}}{\kappa_{r}}\right)^{n}&1\leq
n<n_{\text{max}}\\\
\left(\frac{\kappa_{f}}{\kappa_{r}}\right)^{n_{\text{max}}}&n=n_{\text{max}}\end{cases}$
(S18)