## Abstract

A key goal in synthetic biology is the construction of molecular circuits that robustly adapt to perturbations. Although many natural systems display perfect adaptation, whereby stationary molecular concentrations are insensitive to perturbations, its *de novo* engineering has proven elusive. The discovery of the antithetic control motif was a significant step toward a universal mechanism for engineering perfect adaptation. Antithetic control provides perfect adaptation in a wide range of systems, but it can lead to oscillatory dynamics due to loss of stability, and moreover, it can lose perfect adaptation in fast growing cultures. Here, we introduce an extended antithetic control motif that resolves these limitations. We show that molecular buffering, a widely conserved mechanism for homeostatic control in nature, stabilises oscillations and allows for near-perfect adaptation during rapid growth. We study multiple buffering topologies and compare their performance in terms of their stability and adaptation properties. We illustrate the benefits of our proposed strategy in exemplar models for biofuel production and growth rate control in bacterial cultures. Our results provide an improved circuit for robust control of biomolecular systems.

## I. INTRODUCTION

Synthetic biology promises to revolutionise many sectors such as healthcare, chemical manufacture and materials engineering^{7}. A number of such applications require precise control of biomolecular processes in face of environmental perturbations and process variability^{22}. An important requirement in such control systems is *perfect adaptation*, a property whereby chemical concentrations remain insensitive to perturbations^{20,23}. The molecular mechanisms that can produce perfect adaptation has been extensively studied in natural systems^{3,14,20,23}. In these systems, perfect adaptation can be produced by a range of feedforward and feedback mechanisms^{3,20}. Such natural systems have been shaped by evolutionary processes, but it remains unclear if they are sufficiently robust and tune-able for *de novo* engineering of perfect adaptation in synthetic circuits.

One approach to engineer perfect adaptation relies on the use of feedback control. As illustrated in Figure 1A, this strategy requires circuits that sense the output and act upon the inputs of a biomolecular process. The groundbreaking work by Briat and colleagues^{1,4} identified *antithetic feedback* as a promising candidate for engineering perfect adaptation in living systems. Antithetic control involves a feedback mechanism with two molecular components that sequester and annihilate each other (see Figure 1B). It enables a system output to robustly follow an input signal and remain insensitive to various types of perturbations, akin to what integral feedback achieves in classic control engineering strategies^{2}.

The original antithetic control motif, however, has two weaknesses that can limit its applicability: it is often not effective when cells are growing rapidly, and the feedback mechanism can cause unwanted oscillations under a range of conditions^{27}. Specifically, dilution effects caused by cell growth cause “leaky integration” - so called because integration is a form of memory and dilution causes that memory to leak over time^{29}. This prevents perfect adaptation from occurring during rapid growth. Although in some motif configurations, the loss of perfect adaptation can be partly mitigated with a stronger feedback^{29}, in general the use of strong feedback results in the loss of stability and undesirable oscillations^{27}. Such oscillations can be stabilised in specific motifs^{15}, and in more general cases the combination of antithetic control with classic Proportional-Integral-Derivative (PID) control has been shown to improve temporal regulation^{5,9}. Yet to date, there is no general strategy to avoid oscillations and prevent the loss of adaptation during rapid growth.

Here we propose an extended antithetic control system that resolves the above limitations. We show that the addition of molecular buffers improves stability and suppresses undesirable oscillations, and moreover it can allow for near-perfect adaptation in fast growth regimes. Molecular buffering is a widespread regulatory mechanism in nature (e.g. ATP, calcium & pH buffers^{17,21,32}) that has received modest attention in the literature as compared to other regulatory mechanisms. Recent work found that the combination of buffering and feedback is often critical for robust regulation^{16,17}. Buffering has the ability to attenuate fast disturbances and stabilise feedback control^{17} (see Figure 1C), and can also be essential for the control of multiple coupled outputs^{18}. Here, we first show that a number of buffering topologies can stabilise the original antithetic control system and preserve perfect adaptation. We then show that buffering can allow increased feedback strength or ‘gain’ without producing oscillations, which in turn reduces the steady state error even in fast growth conditions. To illustrate the utility of this new antithetic control strategy, we examine two case studies that involve the control of biofuel production and growth rate in microbes.

## II. BACKGROUND

### A. Perfect adaptation and antithetic control

Antithetic control employs a feedback mechanism with two molecular components that sequester and annihilate each other (see Figure 1B). In its most basic formulation, an antithetic system contains a two-species molecular process to be controlled, and a two-species antithetic controller. The two species of the controlled process (*x*_{1} and *x*_{2}) can represent a variety of molecular systems, including e.g. mRNA and protein as in Figure 1A. The goal of the antithetic control system is to desensitise the steady state concentration of *x*_{2} with respect to external perturbations. Such perturbations include, for example, insults of molecular species coming from upstream or downstream processes, changes in cellular growth conditions, or alterations to binding affinities between species.

In the absence of stochastic effects, the feedback system can be modelled by the ODEs:
where *z*_{1} and *z*_{2} are the concentrations of species in the antithetic controller, and *θ*_{1}, *k*, and *θ*_{2} are positive parameters representing first-order kinetic rate constants. The parameter *µ* describes a zero-order influx of controller species *z*_{1}, while *η* is a second-order kinetic rate constant. We further assume that molecular species are diluted by cellular growth, degraded by other molecular components, or consumed by downstream cellular processes, all of which we model as a first-order clearance with rate constant *γ*_{p}. The controller species *z*_{1} and *z*_{2}, on the other hand, are assumed to be diluted by cellular growth with a rate constant *γ*_{c}.

In the absence of dilution effects (*γ*_{c} = 0), from the model in (1) we can write:
which after integration becomes
The above equation means that, if the system has a stable equilibrium, the steady state concentration of *x*_{2} is *x*_{2} = *µ*/*θ*_{2}, and hence independent of all model parameters except *µ* and *θ*_{2}. Therefore the antithetic control system displays perfect adaptation because the steady state of *x*_{2} is robust to perturbations in parameters *k, θ*_{1}, *θ*_{2} and *η*.

A caveat of antithetic feedback is that it can have a destabilising effect. When the controller species are not diluted (*γ*_{c} = 0), it can be shown that a parametric condition for stability is^{27}:
As shown by the stability diagram in Figure 1B, strong antithetic feedback can cause the system to lose perfect adaptation and display oscillatory dynamics. Moreover, in the presence of dilution of the controller species (*γ*_{c} > 0), the antithetic controller is unable to produce perfect adaptation^{27}. The adaptation error can be reduced with stronger feedback, for example by increasing the rate constants *θ*_{1}, *k*, or *θ*_{2}. Yet as mentioned above, stronger feedback can cause unwanted oscillations^{27}. These caveats are particularly relevant in bioproduction applications that require fast culture growth^{10,35}.

### B. Molecular buffering

Buffering is the use of molecular reservoirs to maintain the concentration of chemical species^{17}. It is a widespread regulatory mechanism found across all do-mains of life, with common examples including pH, ATP and calcium buffering^{21,32}. Molecular buffering can have a number of regulatory roles^{17,18}, including acting as a stabilising mechanism for other molecular feedback systems^{17}.

To provide a background on buffering models, we consider the simple case of a chemical species (*x*) that is subject to feedback regulation, as shown in Figure 1C. A general model for such process is:
where *w* is a molecular buffer for the regulated species *x*, and *p*(*x*) is a feedback-regulated production rate of *x*. The parameters *γ*_{x} and *γ*_{w} are first-order clearance rate constants of *x* and *w*, respectively. The terms *g*_{w} and *g*_{x} describe the reversible binding of species *x* and the buffer *w*. The steady state occurs when production matches degradation (i.e. ) and conversion from *x* to *w* matches the reverse conversion plus removal (i.e. ).

It can be shown that after linearisation and assuming that the buffering reactions rapidly reach quasi-equilibrium, the model (3) can be simplified to (see SI1):
where is the deviation of x from the steady state is the linearised feedback gain and *B* is the buffer equilibrium ratio:
where is the deviation of *w* from the steady state . The parameter *B* is buffer-specific and compares the change in the concentration of a regulated species, *x*, to the change in the concentration of a buffering species, *w*, when the buffering reactions are at quasi-equilibrium^{17,18}.

From (4) we observe that buffering slows down the rate of change of the output *x* by a factor of (1 + *B*). This slowed rate allows for attenuation of fast disturbances and stabilisation of unwanted oscillations^{17} (see SI1). In the next section, we study the ability of buffering to stabilise oscillations in antithetic control systems.

## III. BUFFERING CAN STABILISE ANTITHETIC INTEGRAL FEEDBACK

In this section we study a modified version of the antithetic feedback controller that includes buffering of its molecular components. We consider a number of architectures (Figure 2A) and identify those that suppress oscillations caused by the instability illustrated in Figure 1B.

To study the impact of buffering on the antithetic controller, we consider a mathematical model in which species *z*_{1}, *z*_{2}, and *x*_{2} are buffered by molecules *w*_{1}, *w*_{2}, and *w*_{x}, respectively. For simplicity we use linear buffering reaction rates on the basis that the most important nonlinearity is the mutual annihilation of controller species *z*_{1} and *z*_{2}^{9}. As in Eq. (4), we assume that the buffers rapidly reach quasi-equilibrium to obtain the following extended model (see SI2):
where *B*_{x}, *B*_{1} and *B*_{2} are the equilibrium ratios for each buffer, and *γ*_{x} represents the degradation rate of buffer *w*_{x}. The extended model in (6) reduces to the original antithetic system in (1) if *B*_{x} = 0, *B*_{1} = 0 and *B*_{2} = 0.

In the extended antithetic controller, the parameters (*B*_{x}, *B*_{1}, *B*_{2}, *γ*_{x}) are additional tuning knobs that can be used to shape the closed-loop dynamics. In Figure 2A we show the three considered buffering architectures. As a result, for rapid buffering we can compute the buffer concentrations as *w*_{2} = *B*_{2}*z*_{2}, *w*_{1} = *B*_{1}*z*_{1} and *w*_{x} = *B*_{x}*x*_{2}.

We first show that buffered antithetic feedback pre-serves perfect adaptation. From (6) we write
which after integrating becomes:
The above integral ensures that if the system is stable then the steady state of *x*_{2} is *x*_{2} = *µ*/*θ*_{2}, hence independent of all parameters except *µ* and *θ*_{2}. The steady state of *x*_{2} is thus robust to perturbations in the original parameters *k, θ*_{1}, *θ*_{2}, and *η*, as well as the additional parameters introduced by the buffering mechanism *B*_{x}, *B*_{1}, *B*_{2}, and *γ*_{x}. This means that the buffered antithetic feedback displays perfect adaptation as in the original formulation in Eq. (1).

As shown by the stability condition in (2), the original antithetic control system becomes oscillatory when the feedback gain is too strong or the degradation of *x*_{1} and *x*_{2} is too slow^{27}. To analyse the stabilising role of buffering, we first consider the system in the absence of degradation of the *w*_{x} buffer (i. e. *γ*_{x} = 0). Assuming rapid buffering and strong integral binding (large *η*), we found that the system is stable when (see SI2):
From the condition (7) we observe that increasing *B*_{1} reduces the lower bound for *γ*_{p} and improves stability; this suggests that buffering of *z*_{1} provides a route to suppress oscillations. The stability condition also shows that buffering of *z*_{2} has no impact on stability, whereas buffering of *x*_{2} can destabilise the system and produce oscillations. As shown in the stability diagram in Figure 2B, this finding means that molecular buffering in topologies 1 and 3 can counteract each other.

We next sought to identify variations of topology 3 that can ameliorate its destabilising effect. We found that degradation of the *x*_{2} buffer (*w*_{x}) can indeed stabilise the closed-loop and eliminate the unwanted oscillations; this new topology is shown in Figure 3A. Under the same assumptions as condition (7) (rapid buffering and strong binding rate constant *η*), the stability conditions are (see SI3):

Condition (8) reduces to the one in (7) if *γ*_{x} = 0 and there is no buffering of *z*_{1} and *z*_{2}, i.e. *B*_{1} = *B*_{2} = 0.

As shown by the stability diagram in Figure 3B, topology 3 with degradation provides an effective solution to stabilise the closed-loop. The condition in (8) suggests that the ratio *γ*_{x}/*γ*_{p} has a key role in stability. Buffers with shorter half-lives (larger *γ*_{x}) tend to almost completely remove the instability, even for low buffer equilibrium ratios *B*_{x}. As we show in SI3, under the condition *γ*_{x}/*γ*_{p} > 1/3, increases in *B*_{x} tend to stabilise the closed-loop. This includes the important case where both *x*_{2} and its buffer are degraded at the same rate, i.e. *γ*_{x} = *γ*_{p}. For large values of *B*_{x}, the stabilisation effect is even stronger and becomes independent of the half-life of *w*_{x}. We found a similar stabilisation effect in systems with *x*_{1} buffering that include degradation of the buffer (see SI4).

## IV. ACHIEVING NEAR PERFECT ADAPTATION IN FAST GROWTH

In this section we investigate the effect of buffering on antithetic control with dilution. It is well-known that dilution by cell growth can disrupt perfect adaptation in the original antithetic control system^{1}, and thus here we explore the impact of dilution in the proposed topologies with molecular buffering. To study the effect of dilution, we modify (6) to include dilution terms for the control species *z*_{1} and *z*_{2}, as well as the buffer for *z*_{1}:
where *γ*_{c} represents the dilution rate constant of the control species *z*_{1}, *z*_{2} and buffer species at *z*_{1}. As in the previous model in (6), the model (9) can be obtained under the assumption of rapid equilibrium of the buffering reaction. Moreover, we further assume that dilution of *x*_{1} and *x*_{2} and the buffer at *x*_{2} can be lumped into their first-order degradation rates.

### A. Topology 3

For the case of dilution with a single buffer at *x*_{2}, we set *B*_{1} = 0 in (9). The resulting steady state is (see SI5)
where
and . The second term in (10) is always smaller than unity. Therefore the steady state of the out-put is *x*_{2} *< µ*/*θ*_{2} and the system loses perfect adaptation. Moreover, the deviation of the steady state of *x*_{2} from the reference point *µ*/*θ*_{2} is (see SI5):
where *x*_{2n} = *µ*/*θ*_{2} is the reference input. Increases to *B*_{x}, *γ*_{c} or *γ*_{x} in (12) thus amplify the steady state error, while increases to the feedback strength *k θ*_{1} *θ*_{2} brings the system closer to perfect adaptation. We further obtained parametric conditions for stability (see SI5):
Taken together, the relations in (12)–(13) define an upper bound for the best possible steady state error. Specifically, in (12) we see that stronger feedback gain can increase Ω _{x} and so reduce the steady state error. Buffering of *x*_{2} tends to stabilise the oscillations and, at the same time, allows the steady state error to be reduced by stronger feedback gain, without the risk of instability observed in the original formulation. This phenomenon is illustrated in Figure 4A, which shows the stability condition (13). Notably, we observe that increasing *B*_{x} improves stability only in regions for low and high values of *γ*_{x}, and not intermediate values. Figure 4B shows simulations of the stabilising effect of molecular buffering for the case of *γ*_{x} = 10, which enables a decrease of steady state error by means of stronger feedback gain. We also found that topology 3 improves stability even without buffer degradation (*γ*_{x} = 0), unlike the case when there is no dilution in (7) and Figure 2.

### B. Topology 1

We found that buffering at *z*_{1} can similarly reduce steady state error via increases to the feedback gain. However, unlike Topology 3, this compensatory effect only occurs for systems that utilise slow buffers. To examine this result in detail, we set *B*_{x} = 0 in (9) and compute the resulting steady state (see SI7):
where
and . As in the previous case, the steady state satisfies *x*_{2} *< µ*/*θ*_{2} and thus the system loses perfect adaptation (see Figure 5). Moreover, in this case the steady state error is (see SI7):

where *x*_{2n} = *µ*/*θ*_{2} is the reference input. Increasing *B*_{1} or *γ*_{c} in (16) increases the steady state error, while increasing the feedback strength *kθ*_{1}*θ*_{2} brings the system closer to perfect adaptation. Assuming rapid buffering, the conditions for stability is (see SI7):
From condition (17) we observe that rapid *z*_{1} buffering results in no net change to the steady state error; this effect can be observed in simulations in Figure 5A. Increases to *B*_{1} enable higher gain feedback without destabilising the system, but they also worsens the steady state error in (16) and the two effects cancel each other.

As we show in SI8, slow buffering of *z*_{1} also enables an increase in the feedback gain without producing oscillations, which is similarly at the expense of an increased adaptation error. But in this case the net effect of buffering and increased feedback gain is positive and there is an overall reduction of the adaptation error, as shown in the simulations in Figure 5B.

## V. CASE STUDIES

### A. Model for biofuel production

To illustrate the potential of the proposed control topologies, here we employ an existing model for biofuel production that incorporates antithetic control^{6} (see also^{12}), shown in Figure 6A. The synthetic system produces biofuel from sugars through a metabolic pathway. The biofuel product can be toxic to the cell and so efflux pump proteins are expressed to remove the toxic metabolic product. However, at large concentrations the efflux protein pump can also be toxic. A feedback mechanism can help robustly regulate these two competing toxic products. The antithetic feedback mechanism senses the biofuel concentration to control the expression of efflux pump protein. An increase in the pump protein then reduces the biofuel concentration, completing the loop. Stability is known to be a major issue for the system, as it susceptible to oscillation for large *η*, which is the typical design case^{6}.

The extended model of the biofuel circuit with antithetic feedback and the addition of a protein buffer is:
where *n* is the normalized cell density, which is assumed to follow logistic growth with additional death rates due to toxicity of intracellular biofuel concentration *b*_{i} and efflux protein pump *p*. The variables *z*_{1} and *z*_{2} are the controller species, while the production of the protein pump *p* is assumed to be proportional to controller species *z*_{2}. The variable *w* is the buffering species which buffers *z*_{1} through a reversible reaction via chemical species *I* that inhibits *z*_{1} sequestering when bound to *z*_{1}. The variable *b*_{e} is the extracellular concentration of biofuel.

Buffering of *z*_{1} can be seen to stabilise the process in the simulations of the model (see Figure 6A). These simulations show that the oscillations, which occur when *η* is large, quickly settle to the steady state when buffering is introduced. This stabilising effect resembles impact of buffering in Topology 1 in Section III. It also shows the stabilising effect for models that are destabilised by strengthening the antithetic binding mechanism via in-increased h. This example thus illustrates that the stabilising effect of buffering also occurs in more complex systems than those considered in the previous section.

### B. Model for growth control

For the synthetic growth control case study, we use an existing model of the synthetic growth control circuit, which includes the new addition of buffering27 (see Figure 6). The variable *N* represents the population size and is assumed to follow logistic growth, with an additional death rate due to toxicity that is proportional to the concentration of *CcdB* per cell. *Ccdb* is a protein that is toxic to the cell. *mRNA* is messenger RNA while *asRNA* is a short antisense RNA that has a complementary sequence to the mRNA, which enables sequestration between the two. *mRNA* and *asRNA* form the antithetic integral controller. The transcription of mRNA is induced by a quorum-sensing ligand. The term *G*_{a} represents the gain between *N* and mRNA induction resulting from the quorum-sensing molecule AHL. *W* represents a buffer of Ccdb, which consists of an inactivated form of Ccdb that can reversibly bind to an inhibitor molecule *I*. The adapted model of the genetic circuit with the new addition of a protein buffer is
where [·] represents intracellular concentrations for each species and the last line indicating the rate of change of *W* is new to the model.

The buffering of *CcdB* as shown above is equivalent to *x*_{1} buffering in the model (6) and Figure 2, as *N* is the output and equivalent to *x*_{2}. Buffering at *x*_{1} provides a similar benefit as buffering at *x*_{2} and so can also enable near-perfect adaptation, where we omitted these results in the previous section for brevity (see SI6).

Buffering of *CcdB* in conjuction with increased feedback gain can be shown to reduce steady state error in the simulations in Figure 6. Increased feedback gain is implemented in these simulations by increasing the translation rate of *CcdB*.

## VI. DISCUSSION

Perfect adaptation has been subject of intense study in the synthetic biology community. Although perfectly adapting systems are ubiquituous in nature, their implementation has proven particularly elusive. The antithetic control motif, first discovered by Briat et al4 and implemented by Aoki et al1, provides a new molecular mechanism to build perfect adaptation into a wide range of synthetic gene circuits. A number of works have sought to find alternative circuits that provide adaptation properties similar to antithetic control. For example, several authors have shown that ultrasensitive feedback can display some of the features of perfect adaptation25,28, and the idea was recently extended in great detail for synthetic gene circuits30. Other works have sought to devise molecular implementations of Proportional-Integral-Derivative control9, as this is a widely adopted strategy for perfect adaptation in engineered control systems.

Here we have addressed caveats of the original antithetic control system with an extended architecture that has improved stability properties. The proposed circuit combines an antithetic motif with a molecular buffering mechanism. Molecular buffering is widely conserved in natural systems, and common examples include the ATP buffering by creatine phosphate, pH buffering and calcium buffering. In all these examples, a molecular buffer sequesters a target molecule into an inactive form, resulting in a system with improved ability to mitigate fast perturbations. In the case of antithetic control, the addition of buffering results in the stabilisation of unwanted oscillations and, moreover, provides near-perfect adaptation even in rapid growth conditions where the performance of antithetic control is known to be particularly poor.

After detailed examination of mathematical models for various circuit architectures (Fig. 2), we found two candidate systems with improved stability properties, either by buffering species of the antithetic motif it-self, or by buffering and degrading a target species to be controlled. The first circuit, called topology 1 in Fig. 2, provides stability over a large range of parameters values than classic antithetic control and can generally stabilise the oscillations produced by antithetic control. Moreover, topology 1 requires buffering of a molecular species of the antithetic motif itself, and therefore it provides a promising strategy to stabilise variables that are not easily buffered directly, such as population size or metabolite species as illustrated by the example in Fig 6A.

The second circuit, termed topology 3 in Fig. 2, requires buffering the molecular output of the process to be controlled. We found that, although in principle this topology can have a destabilising effect, when coupled with buffer degradation it provides an effective way to mitigate oscillations in fast growth regimes (Fig. 3). Interestingly, there is a similar effect when applied to inter-mediate species instead of the output species in the controlled process, such as *x*_{1} in the original circuit shown in Fig. 2 or CcdB in the growth control case study in Fig. 6. Buffering an intermediate species also provides an alternative location when the output is not easily buffered. A drawback of topology 3 is that degradation of the buffer may require the implementation of additional molecular mechanisms. If the degradation mechanism requires expression of heterologous proteins, this can increase the genetic burden on the host cell and impair its physiology^{8,26}.

The effect of buffering on perfect adaptation is strikingly similar to a strategy employed in industrial process control, where buffer tanks are employed to regulate and smooth out the impact of disturbances^{13}. In our case, the specific implementation of the molecular buffers is a subject of future study, as this will largely depend on the type of biomolecular process to be controlled. For example, buffers for gene expression may require gene products to be sequestered, which can be achieved through several mechanisms such as reversible protein-protein binding^{24}, phosphorylation^{19}, small molecule inhibitors^{11}, or DNA decoy sites^{34}. In metabolism and signalling systems, ubiquitous examples are the interconversion between a target species and a buffer (e. g. reversible catalysis between ATP and creatine phosphate^{17,32}) or sequestering by dedicated proteins (e. g. Ca^{2+}or H^{+} ions^{17,21,32}).

Our main goal in this paper was to show that molecular sequestration can improve perfect adaptation in the antithetic control motif. Since buffering is known to stabilise a much wider range of molecular networks^{16}, it also has the potential to improve other circuits implementing perfect adaptation, e. g. those that rely on ultrasensitive behaviour^{30}. Another promising line of inquiry is investigating production feedback mechanisms with similar kinetic effects to degradation^{16}, which may enable topology 3 type buffers to stabilise the systems without an increase in burden. Further, the effect of nonlinear buffering reactions also requires investigation as these can produce an effective increase of the buffer equilibrium ratio without increasing the concentration of the buffer itself^{18}. For simplicity, we have focused exclusively on deterministic dynamics, but the analysis of stochastic effects emerging from the interplay between molecular buffering and antithetic control are particularly attractive, as it is known that buffering does not amplify stochastic fluctuations^{17} yet some phenomena are known to emerge only in the presence of molecular noise^{4,31,33}.

As synthetic gene circuits grow in size and complexity, there is a growing need for mechanisms that can enhance their robustness in a range of operational conditions. In the longer term, this will require the availability of a catalogue of gene circuits that can produce perfect adaptation in response to perturbations. In this work we have presented one such architecture, and thus laid theoretical groundwork for the discovery of new biomolecular systems with improved functionality.

## METHODS

All mathematical models are based on systems of ordinary differential equations. The stability conditions in Eqs. (7), (8), (13) and (17) were obtained using frequency domain transformations (Laplace and Fourier) of the linearised models, along with detailed examination of the magnitude and phase equations of the resulting characteristic polynomials for the closed-loop systems^{27}. Simulations were carried out using standard ODE solvers in MATLAB. All calculations and model descriptions can be found in the Supplementary Material.

## AUTHOR CONTRIBUTIONS

Conceptualization, E.J.H.; Methodology, E.J.H.; Formal Analysis: E.J.H., and D.A.O.; Investigation, E.J.H., and D.A.O.; Writing - Original Draft, E.J.H.; Writing - Review & Editing, E.J.H., and D.A.O.; Funding Acquisition, E.J.H.; Supervision, E.J.H., and D.A.O.

## SUPPLEMENTARY MATERIAL

### SI1. BUFFERING

In this section, we provide a background on buffering, including methods for analysing models with buffering. We start with the simple model in (3) of a single regulated species that is being buffered, such that
where *x* is the output species being regulated, *w* is the buffering species, *p* is the production rate of *x, γ*_{x} is the removal kinetic rate of *x, g*_{w} is the forward buffering reaction rate and *g*_{x} is the reverse buffering reaction rate. Incorporation of feedback is represented by the *x* dependence of production. The nominal steady state (*d* = 0) occurs when production matches degradation and the buffer is at steady state .

To analyse (S1), we reduce the two state model to one state by assuming that the buffering reactions rapidly reach equilibrium. To carry this out, we first linearise (S1), which results in Where is the linearised feedback gain, and and are the linearised kinetic rates for the forward and reverse buffering reaction.

If we assume that the buffer rapidly reaches equilibrium then *w* is at quasi-steady state and so . This steady state results in
.

We set the slow variable as Δ*x*_{T} = Δ*w* + Δ*x* where . Using the definition of *B* above, we have Δ*x*_{T} = (1 + *B*)Δ*x*. Thus and so
which is a reduced one state model, where the second state can be determined from Δ*w* = *B*Δ*x*. In technology, integral feedback is often paired with proportional and derivative feedback (PID control)^{2}. In biology, rapid buffering without degradation is equivalent to negative derivative feedback and rapid buffering with degradation is equivalent to PD feedback with degradation. These equivalences can be observed in^{17}
where the buffer equilibrium ratio *B* corresponds to the derivative feedback ‘gain’.

To study the effect of buffering on stability, we can also modify the model in (S1) to include a delay of the production feedback term
where *p*(*x*(*t − τ*)) represents the production feedback with a delay of time *τ*. The reduced model is
where Δ*x*_{τ} = Δ*x*(*t − τ*).

It can be observed in Figure S1 that buffering can stabilise the oscillations that result from feedback delay.

### SI2. BUFFERING CAN STABILISE ANTITHETIC INTEGRAL FEEDBACK: ALL SPECIES

In this section, we analyse the stabilising effect of buffering (without degradation) on antithetic integral feedback. We base our studies on a simple model involving the antithetic integral feedback (without buffering)^{27}. Consider
where *x*_{2} is the output concentration being controlled, *x*_{1} is another concentration in the process being controlled, and *z*_{1} and *z*_{2} represent the molecular species involved in the perfect adaptation mechanism. We next introduce buffering to (S3). We show how the model reduction method described in SI1 can be used to simplify the model for one case. We then use the same method for all buffers. As a first case, we introduce buffering at the controlled variable *x*_{2} and simplify the model by assuming that the buffering reactions rapidly reach equilibrium. With buffering at *x*_{2}, we have
where *w*_{x} is the buffering species at *x*_{2} and *b*_{x}, *b*_{w} are the kinetic rate constants for the buffering reactions. Although the buffering equilibrium ratio is defined in terms of the linearised model and deviations from steady state, we can use the same notation here for the nonlinear model as the buffering reactions are linear. If the buffering reaction is at equilibrium then
If we assume that the buffer rapidly reaches equilibrium then *w* is at quasi-steady state and so
We note that although the buffer equilibrium ratio is defined in (5) in terms deviations from steady state, we can use the buffer equilibrium ratios in the nonlinear model as the reaction rates are linear. We set *x*_{T} = *w*_{x} + *x* as the slow variable and so *x*_{T} = (1 + *B*_{x})*x*_{2}. Thus and so
If we include buffering on *z*_{1}, *z*_{2}, *x*_{2} and apply a similar model reduction by assuming rapid buffering, then we have
where *B*_{1}, *B*_{2} are the buffer equilibrium ratios of the buffers at *z*_{1} and *z*_{2} respectively.

#### A. Steady State

We first determine the steady states of the system, which is useful both for determining perfect adaptation and as a prerequisite for stability analysis. Using
we can see that the steady state for the output is
The correspondence of this steady state with perfect adaptation is discussed further in Section II A and III. The species *x*_{1}, *z*_{1} and *z*_{2} have the corresponding steady states

### B. Stability

We next determine the effect of buffering on the stability of the model. We follow the methodology used by Olsman and colleagues27 to study the stability of antithetic integral feedback with the addition of buffering. Linearising the system (S5) about the steady states, we have
where are the deviations from steady state and
Taking the Laplace transform of (S6), we have
where *Z*_{1} =ℒ {Δ*z*_{1}}, *Z*_{2} = ℒ{Δ*z*_{2}}, *X*_{1} = ℒ{Δ*x*_{1}}, *X*_{2} =ℒ {Δ*x*_{2}} are the Laplace transforms of the time-domain concentration deviations. Substituting, we have
Simplifying, we have
Substituting from (S8), we have
The characteristic equation used to analyse stability is thus

#### 1. Characterisation of roots

Following the methodology used by Olsman and colleagues27, we first characterise the roots of (S9). If we substitute *s* = *γ*_{p}*σ* then
Taking the limit of strong binding for the sequestration process in antithetic integral feedback then
It can be observed that there is a large, real root at . We next examine the region , where we have
The magnitude and phase constraints are
If we assume that *σ* is real and positive, then the LHS of the phase constraint is
which contradicts. Thus unstable roots are not purely real.

If we assume that *σ* is real and *−*1/(1 + *B*_{x}) *< σ <* 0 then the LHS of the phase constraint is
and so it is possible to have stable real roots. If we set *f* = |*σ*|| 1 + *σ* ||1 + (1 + *B*_{x})*σ*| from the magnitude constraint, then there is a local maxim of *f* at
Thus there are two stable real roots between *−*1/(1 + *B*_{x}) *< σ <* 0 if
as *f* (*σ*_{max}) is larger than the RHS of the magnitude constraint. There is a bifurcation at the boundary resulting in a pair of complex conjugate roots if
For *−*1 *< σ < −*1/(1 + *B*_{x}) then
which violates the phase constraints. For *σ < −*1 then
which meets the phase constraints. Thus for the stability boundary with strong binding there is a negative real root and a complex pair of roots in the region , as well as one large negative root.

#### 2. Stability Conditions

To determine the stability boundary, we next determine the conditions for which the complex roots are purely imaginary. Substituting we have
From these equations, the magnitude and phase constraint are
for some integer *k*. Using the strong binding assumption then from above and so
Rewriting the phase constraint, we have
Applying tan(·) and trigonometric identities to both sides, we have
Solving, we have
Thus the stability boundary occurs at
From above, we know that all roots are real and stable if *α*/(*γ*_{p}(1 + *B*_{1})) is sufficiently small, and so the stability condition is
or
Thus increasing *B*_{1} improves stability and increasing *B*_{2} has no effect on stability. Further, increasing *B*_{x} worsens stability, although this effect saturates as *B*_{x} increases.

## SI3. BUFFERING CAN STABILISE ANTITHETIC INTEGRAL FEEDBACK: *x*_{2} BUFFERING WITH DEGRADATION

In this section, we analyse the stabilising effect of buffering at the output species *x*_{2} on antithetic integral feedback, where the buffering can be degraded. Consider the model (S4) with buffering at *x*_{2} that is degraded
where *w*_{x} is the buffering species of *x*_{2}, and *b*_{x}, *b*_{w} are the kinetic rates for the buffering reactions. If we assume rapid buffering such that *w*_{x} + *x*_{2} is the slow variable, *w*_{x} = *B*_{x}*x*_{2} and , and use the methodology from SI1 and SI2 then we have the reduced model

### Steady State Analysis

We next analyse the steady state. We have and so the steady state of the output is We also have corresponding steady states where .

### Stability Analysis

We next analyse the stability of the system. If we linearise about the steady states, we have
This system can be rewritten as
where and *β* = *ηµ*. Taking the Laplace transforms, we have
where *Z*_{1} = ℒ{Δ*z*_{1}}, *Z*_{2} = ℒ{Δ*z*_{2}}, *X*_{1} = ℒ{Δ*x*_{1}}, *X*_{2} = ℒ{Δ*x*_{2}}are the Laplace transforms of the time-domain concentration deviations. Substituting, we have
Rewriting and substituting, we have
Simplifying and taking the limit of strong binding for the sequestration process in the antithetic integral feedback then
Thus we have the characteristic equation
Substituting *s* = *γ*_{p}*s*, we have
Using the same argument as that in SI2, for the stability boundary with strong binding there is a negative real and complex pair of roots in the region , as well as one large negative root.

To determine the boundary of stability, we next determine the conditions for which the roots are purely imaginary. Substituting *s* = *iωγ*_{p}, we have
Taking the strong binding limit where , we have
The phase and magnitude constraints are
for some integer *k*. Solving the phase constraint, we have
For this constraint we require
which reduces to
Substituting into the magnitude equation, we have
Rearranging, we have the stability condition
We can differentiate the right hand side with respect to *B*_{x} to determine whether increasing *B*_{x} has a stabilising effect. If we set
then
Thus for small *B*_{x} then buffering stabilises antithetic integral feedback if . For large *B*_{x} then increasing *B*_{x} improves stability if *γ*_{x} *>* 0.

### SI4. RAPID *x*_{1} BUFFERING WITH DEGRADATION CAN STABILISE ANTITHETIC INTEGRAL FEEDBACK

In this section, we analyse the stabilising effect of buffering at the intermediate species *x*_{1} on antithetic integral feedback, where the buffering can be degraded. This section uses identical methodology and obtains equivalent results to SI3.

Consider the model (S3) with buffering at *x*_{1} that is degraded
where *w*_{i} is the buffering species of *x*_{1}, and *b*_{i}, *b*_{w} are the kinetic rates for the buffering reactions. If we assume rapid buffering such that *w*_{i} + *x*_{1} is the slow variable, *w*_{i} = *B*_{i}*x*_{1} and , and use the methodology from SI3 then

### Steady State Analysis

We next analyse the steady state. We have and so the steady state of the output is The corresponding steady states of the other species are where .

### Stability Analysis

We next study the stability of the system. If we linearise about the steady states, we have
This system can be rewritten as
where and *β* = *ηµ*. Taking the Laplace transforms, we have
where *Z*_{1} = ℒ{Δ*z*_{1}}, *Z*_{2} = ℒ{Δ*z*_{2}}, *X*_{1} = ℒ{Δ*x*_{1}}, *X*_{2} = ℒ{Δ*x*_{2}} are the Laplace transforms of the time-domain concentration deviations. Substituting, we have
Rewriting and substituting, we have
The above equation is equivalent to *x*_{2} buffering (see SI4), and so for strong integral binding we have the stability condition
which is equivalent to the case of *x*_{2} buffering.

### SI5. RAPID *x*_{2} BUFFERING WITH DEGRADATION CAN ENABLE NEAR-PERFECT ADAPTATION DESPITE LEAKY INTEGRATION

In this section, we analyse the ability of buffering at *x*_{2} to enable near perfect adaptation by stabilising antithetic integral feedback. Consider the model (S4) with buffering at *x*_{2} that is degraded
where *w*_{x} is the buffering species of *x*_{2}, and *b*_{x}, *b*_{w} are the kinetic rates for the buffering reactions. Using the methodology from SI1 and SI2, we have the reduced model

### Steady State Analysis

We next analyse the steady state. We have
and so the steady state of the output is
We also have the steady state
and so
The other species have the steady states
where . Substituting, we have
We assume strong binding of the sequestration mechanism in antithetic integral feedback (*η* large), which for steady state is the condition
With this assumption, we have
We ignore the zero solution, which corresponds to a negative solution in (S12), and so the steady state concentrations are
The steady state error of *x*_{2} is

### Stability Analysis

We next analyse the stability of the system. If we linearise (S11) about the steady states, we have
This system can be rewritten as
where and *β* = *ηµ*. Taking the Laplace transforms, we have
where *X*_{1}, *X*_{2,} *Z*_{1} and *Z*_{2} are the Laplace transforms of Δ*x*_{1,} Δ*x*_{2,} Δ*z*_{1} and Δ*z*_{2}. Substituting, we have
Rewriting and substituting, we have
Simplifying and taking the limit of strong binding for the sequestration mechanism in antithetic integral feedback
then
Thus we have the characteristic equation
Substituting *s* = *γ*_{p}*σ*, we have
Using the same argument as that in SI2, for the stability boundary with strong binding there is a negative real and complex pair of roots in the region , as well as one large negative root.

To determine the boundary of stability, we next determine the conditions for which the roots are purely imaginary. Substituting *s* = *iωγ*_{p}, we have
Taking the strong binding limit where , we have
The phase and magnitude constraints are
for some integer *k*. Solving the phase constraint, we have
For this constraint we require
which reduces to
Substituting into the magnitude equation, we have
Simplifying, we have
Multiplying by , we have
and so the stability condition is
where
From earlier, we know that the steady state error of *x*_{2} is

### SI6. RAPID *x*_{1} BUFFERING WITH DEGRADATION CAN ENABLE NEAR-PERFECT ADAPTATION DESPITE LEAKY INTEGRATION

In this section, we analyse the ability of buffering at *x*_{1} to enable near perfect adaptation by stabilising antithetic integral feedback. This section uses identical methodology and obtains equivalent results to SI5.

Consider the model with *x*_{1} buffering and dilution
where *w*_{x} is the buffering species of *x*_{1}, and *b*_{i}, *b*_{w} are the kinetic rates for the buffering reactions. Assuming rapid buffering and using the same methodology as previous sections, the reduced model is
where is the buffer equilibrium ratio.

### Steady State Analysis

We next analyse the steady state of the system. We have
and so the steady state of the output is
We also have the steady state
and so
Now at steady state we have
where . Substituting, we have
Assuming strong binding
we have
and so, ignoring the solution *x*_{1} = 0, the steady state is
The steady state error of *x*_{2} is

### Stability Analysis

We next study the stability of the system. If we linearise about the steady state, we have
This system can be rewritten as
where and *β* = *ηµ*. Taking the Laplace transforms, we have
where *X*_{1}, *X*_{2}, *Z*_{1} and *Z*_{2} are the Laplace transforms for Δ*x*_{1}, Δ*x*_{2}, Δ*z*_{1} and Δ*z*_{2}. Substituting, we have
Rewriting and substituting, we have
The above equation is equivalent to that for *x*_{2} buffering (see SI5), and so for strong integral binding we have the equivalent stability constraint and steady state error

### SI7. RAPID *z*_{1} BUFFERING WITH DEGRADATION HAS A TRADE-OFF DUE TO LEAKY INTEGRATION

In this section, we analyse the trade-offs for rapid buffering at *z*_{1} on stability and the steady state error from perfect adaptation. For buffering at *z*_{1} with dilution, we use the model
where *x*_{2} is the output concentration being controlled, *x*_{1} is another concentration in the process being controlled, and *z*_{1} and *z*_{2} represent the molecular species involved in the perfect adaptation mechanism. Assuming that the buffer is rapid then *w* is at quasi-steady state then
If *x*_{T} = *w* + *z*_{1} is the slow variable then *x*_{T} = (1 + *B*_{1})*z*_{1}. Thus and so
Thus we have

### Steady State Analysis

We next determine the steady state and and any error from perfect adaptation. For the case of dilution, we have
resulting in
We also have
and so
Now at steady state we have
where . Substituting, we have
Assuming strong binding of the sequestration mechanism
we have
and so, ignoring the zero solution, the steady state is
The steady state error of *x*_{2} is
We can see that increasing *B*_{1} increases the steady state error of *x*_{2} when there is degradation/dilution of *z*_{1} and *z*_{2}.

### Stability Analysis

We next study the stability of the system with degradation. If we linearise about the steady states, we have
This system can be rewritten as
where
Taking the Laplace transforms, we have
where *X*_{1}, *X*_{2}, *Z*_{1} and *Z*_{2} are the laplace transforms for Δ*x*_{1}, Δ*x*_{2}, Δ*z*_{1} and Δ*z*_{2}. Substituting, we have
Rewriting and substituting, we have
Taking the limit of strong binding then
Thus we have the characteristic equation
Substituting *s* = *γ*_{p}*σ*, we have
Using the same argument as above, for the stability boundary with strong binding there is a negative real and complex pair of roots in the region , as well as one large negative root.

To determine the boundary of stability, we next determine the conditions for which the roots are purely imaginary. Substituting *s* = *iωγ*_{p}, we have
Taking the strong binding limit where , we have
The phase and magnitude constraints are
for some integer *k*. Solving the phase constraint, we have
For this, we require
which reduces to
Substituting into the magnitude equation, we have
As a consequence, the stability constraint is
We can observe that increasing *B*_{1} improves the stability constraint. However, the steady state error of *x*_{2} is
Thus there is a steady state error constraint that is independent of *B*_{1}, and so increasing *B*_{1} does not enable the removal of leaky integration. This result is a consequence of the added degradation of the buffer, which cancels the stabilising effect of the buffer.

### SI8. NON-RAPID BUFFERING CAN ALLOW NEAR PERFECT ADAPTATION WITH LEAKY INTEGRATION

In this section, we analyse the ability of non-rapid buffering at *z*_{1} to enable near perfect adaptation by stabilising antithetic integral feedback. We use the model
where the buffer *w* is not assumed to rapidly reach equilbrium. As a result, the model cannot be reduced in a similar manner to previous sections. The steady state for (S15) is identical to the rapid case in SI7. The linearisation is
This system can be rewritten
where
Taking the Laplace transform of , we have
where *W* and *Z*_{1} are the Laplace transforms of *w* and *z*_{1}. We have
where
Thus
Combining, we have
Simplifying, we have
Taking the strong binding limit of the sequestration mechanism (*β* ≫ max {(*α* + *γ*_{c})(*α* + *γ*_{c}(1 + *B*_{1})), *γ*_{p}(*α* + *γ*_{c}(1 + *B*_{1}))}), we have
and so
Rewriting *C*_{b}, we have
or
Substituting *s* = *iωγ*_{p}, we have
Taking the strong binding limit where , we have
The magnitude constraint is
and the phase constraint is
for some integer *k*. Using trigonometric identities, we have
This can be simplified to
Ignoring the trivial solution *ω* = 0, we have
If for simplicity we assume that *B*_{1} is large, then we have
and so
Thus *ω*^{2} is monotonically increasing wrt *τ* for .

From the magnitude constraint S16, we have the stability constraint
We can rewrite Ω _{1} as
We can observe that Ω _{1} is monotonically increasing with *ω* as *τ* = 1/(*γ*_{c} + *b*_{w}) *<* 1/*γ*_{c}. Combining, Ω _{1} is monotonically increasing with respect to *τ*. Thus increasing *τ* improves the stability constraint for large *B*_{1}. As the steady state error is
then increasing *τ* decreases the steady state error constraint.

### SI9. BODE INTEGRAL OF BUFFERING AND ANTITHETIC INTEGRAL FEEDBACK

Feedback is a highly effective method of robust regulation, but this mechanism is subject to fundamental limits. The Bode integral describes one of these fundamental limit, where improving the regulation at one frequency of a disturbance will worsen regulation of disturbances at other frequencies.

To observe this mathematically, we let *r*(*t*) be the reference signal (corresponding to *µ*). We can mathematically decompose *r*(*t*) and *e* = *r*(*t*) − *x*_{2}(*t*) into their ‘fast’ and ‘slow’ components (via a Fourier transform), which we write as *R*(*iω*) and *E*(*iω*) respectively. A useful measure of regulation is the sensitivity function . The model of the system is
We can write the open-loop model of the two state
where *u*_{a} is the process input for antithetic integral feedback and *u*_{b} is the input for buffering at *x*_{2}, and *z*_{1}, *z*_{2} and *w*_{x} are controller variables.

Buffering at *z*_{1} and antithetic integral feedback act through the input *u*_{a} while buffering at *x*_{2} acts through *u*_{b}. The bode integral is a fundamental constraint on the effectiveness of feedback in any system. It provides a constraint on the overall regulatory effectiveness in terms of the sensitivity function. With output buffering, Bode’s integral is^{16}
where the integral of *S*(*iω*) represents an overall measure of regulation, *b*_{x} is the kinetic rate of the forward buffering reaction and it is assumed that the system without feedback is stable. The integral of *S*(*iω*) sums the regulation of disturbances at different ‘speeds’. Without buffering, if regulation is improved at one ‘speed’ of regulation, it worsens at other ‘speeds’. However, increasing *b*_{x} reduces the integral. Thus the trade-off does not occur with buffering, which can uniformly improve regulation.

In contrast, control species buffering is part of the feedback regulation mechanism and so Bode’s integral is^{16}
Thus control buffering does not remove fundamental constraints, despite stabilising buffering. The tradeoff remains such that improving regulation at one frequency will worsen regulation of disturbances at other frequencies.

## ACKNOWLEDGEMENTS

The authors would like to thank Morgan Kelly for her assistance. EJH was would like to thank the contribution of Judith and David Coffey.

The authors declare that there are no competing interests.