## Abstract

Recent advances in synthetic posttranslational protein circuits are substantially impacting the landscape of cellular engineering and offer several advantages compared to traditional gene circuits. However, engineering dynamic phenomena such as oscillations in protein-level circuits remains an outstanding challenge. Few examples of biological posttranslational oscillators are known, necessitating theoretical progress to determine realizable oscillators. We construct mathematical models for two posttranslational oscillators, using few components that interact only through reversible binding and phosphorylation/dephosphorylation reactions. Our designed oscillators rely on the self-assembly of two protein species into multimeric functional enzymes that respectively inhibit and enhance this self-assembly. We limit our analysis to within experimental constraints, finding (i) significant portions of the restricted parameter space yielding oscillations and (ii) that oscillation periods can be tuned by several orders of magnitude using recent advances in computational protein design. Our work paves the way for the rational design and realization of protein-based dynamic systems.

## INTRODUCTION

Protein oscillators play a major regulatory role in organisms ranging from prokaryotes to humans. In most known biological cases, the oscillation is realized through transcription/translation cycles. Few examples of purely posttranslational oscillators have been found in biology (*1*, *2*). At the same time, posttranslational protein circuits are increasingly sought after for synthetic applications, since they have the potential to exhibit faster response to environment changes, allow for more direct control over the circuit behavior, be directly coupled to a functional output, and can be used in contexts that do not include the vast genetic apparatus (*3*–*5*). While significant recent work has enabled the design of synthetic posttranslational protein–based logic gates (*4*, *5*), engineering tunable dynamic phenomena such as oscillations in a synthetic posttranslational context remains an outstanding challenge (*6*, *7*).

The best-studied example of biological posttranslational protein oscillators is the KaiABC system in cyanobacteria (*8*). By placing only the proteins KaiA, KaiB, and KaiC in a test tube, along with abundant adenosine triphosphate, the KaiC proteins collectively get sequentially phosphorylated and dephosphorylated, forming an oscillatory cycle (*9*, *10*). While the KaiC proteins generally exist in a hexameric state, monomers are shuffled among the hexamers during only a certain phase of the oscillatory cycle (*11*). The KaiABC system demonstrates that protein oscillators need not use transcription/translation cycles or large numbers of components to achieve oscillatory behavior.

Motivated by the KaiABC system, we set out to design a protein-based oscillator that could be reconstituted in vitro using only a small number of components at relatively high copy numbers, so that any resulting oscillations are not stochastic. To facilitate the future translation of this theoretical study to an experimental system, we base the architecture of our system on biochemical constraints and on a design space navigable through computational protein design. We constrain the kinetic reaction network to only include three protein species and to only allow reversible binding and phosphorylation/dephosphorylation enzyme reactions.

As in KaiABC, the oscillating system will cycle through periods of a protein species being phosphorylated or not. When the target protein is phosphorylated, the phosphorylation must induce a change to propel the global state along the oscillatory cycle. This change can be achieved by affecting the enzyme kinetics of the kinase and phosphatase in one of two ways: either by altering the conformation of the target protein or by directly modifying the kinase and phosphatase enzymes. We first consider systems that do not modify the kinase or phosphatase.

The simplest such system, a protein with one phosphorylation site being modified by a kinase and a phosphatase, cannot yield oscillations regardless of parameter choices (*7*). When two phosphorylation sites are included, oscillations are possible only under the assumption that each of the four possible phosphorylation states has significantly different rates of subsequent phosphorylations and dephosphorylations (*7*). While biology seems to have designed a system in KaiABC capable of undergoing the many conformational changes necessary to implement this form of oscillations (*9*), the design of even two (let alone several) protein structures from the same sequence remains a significant challenge for the field of computational protein design (*12*).

These challenges are not unique to molecules with two phosphorylation sites. For example, since oscillations for molecules with two phosphorylation sites are effected by enzyme sequestration (*7*), we consider a molecule containing a single phosphorylation site alongside a kinase- or phosphatase-sequestering domain (or a binding domain for an external compound that itself contains an enzyme-binding domain). These systems are capable of producing oscillations—but only if phosphorylation and binding accompany a significant conformational change in the molecule that modifies the rate constants of subsequent reactions. Even assuming that such a conformational change were designed, we have found no evidence of sustained oscillations in these systems within the parameter regimes of typical binding/unbinding rate constants and typical kinase and phosphatase activity (i.e., the catalytic rate and Michaelis constants *k*_{cat} and *K*_{M}, discussed further below). See section S1 for further discussion.

Systems that focus on modifications to the enzymes themselves are therefore more likely candidates for the production of experimentally realizable oscillations. Biology has found several ways to tie phosphorylation to enzymatic activity. The most straightforward conceptually, having the activity of an enzyme dependent on its own phosphorylation state (*13*) remains a challenge to implement in the context of computational protein design (*14*). However, the field has achieved remarkable success in the design of protein-protein interactions (*15*) which can be modified by phosphorylation (*16*, *17*).

Bootstrapping off of this success, we consider proteins that self-assemble into multimeric functional kinases and phosphatases (*18*, *19*). We are motivated, in part, by the success of using split proteases to implement posttranslational protein–based logic gates (*4*, *5*). In our design, when the proteins’ binding interfaces are phosphorylated, their self-assembly is impeded, reducing the concentration of functional enzymes available in the system.

Oscillations are thus achieved by the push and pull of two opposing factors: Self-assembled kinases inhibit the self-assembly of new proteins by phosphorylating both kinase and phosphatase monomers; meanwhile, self-assembled phosphatases counteract this inhibition. Incoherent inputs such as these are known to enhance the robustness of oscillations (*20*). Our overall oscillator design is motivated by analogy to known successful oscillators, particularly the dual-feedback genetic oscillator (*1*, *21*–*23*). Just as that oscillator relies on the interplay between the inhibiting effects of LacI and the activating effects of AraC, our oscillations rely on the interplay between kinase and phosphatase multimers, which respectively inhibit and activate their own self-assembly (Fig. 1A).

The rest of our manuscript is organized as follows. First, we describe the oscillatory circuits and their experimental constraints. Next, we develop simplified mathematical models for two distinct protein-based oscillators: In one, multimers are designed to form closed, bounded structures; in the other, they form unbounded fibers. We show that simple analytical formulae describing the first oscillator can predict both the regions of parameter space admitting oscillations and the oscillations’ resulting frequencies. We then demonstrate that the second oscillator design is complementary to the first in that it can admit robust and experimentally realizable oscillations in a regime of parameter space where the first cannot. Finally, we discuss the significance of our findings.

## RESULTS

### Self-assembly–based protein oscillators are designed within experimental constraints

*Designed synthetic oscillators rely on few protein species with specified interactions*. The main components of our oscillators are two proteins, which we call κ and ρ. Each individual protein of type κ (ρ) has two complementary parts of a split kinase (phosphatase) and a phosphorylation site. When the respective sites are dephosphorylated, copies of protein κ (ρ) can self-assemble into a functional kinase (phosphatase), which we call *K* (*P*). Thus, self-assembled kinases inhibit the self-assembly of new proteins, while self-assembled phosphatases counteract the inhibition. In addition to the proteins κ and ρ, we include a constitutive phosphatase

The resulting circuit topology (Fig. 1A) is analogous to that used in the dual-feedback genetic oscillator (*1*, *23*). The multimeric kinase plays an analogous role to the LacI protein in the genetic oscillator, repressing the production of new multimers; the multimeric phosphatase plays an analogous role to that of the AraC protein, activating the production (or more precisely, counteracting the kinase inhibition).

Two related networks based on these proteins can be designed. In the first, self-assembly is into closed symmetric homomultimers of specified size; in the second, the monomers are designed such that they self-assemble into one-dimensional unbounded fibers.

*Experimental realizability constrains parameter sets*. Because we are motivated by experimental feasibility, we consider only physically realizable parameters for our models. Binding rates *k*_{b} are typically in the range 10^{−2} to 10^{0} μM^{−1}s^{−1} (*24*) with dissociation constants *k*_{d} typically in the 10^{−3} to 10^{3} μM range (*25*). Both of these quantities can be tuned on the basis of the geometry, energy, and symmetry of the binding interface between the proteins, which we assume here to be designed de novo. Less straightforward to design are the Michaelis constants and catalytic rates of the kinase and phosphatase, especially since these depend not only on the enzyme but on the substrate. Mutational screens can be used to adjust these parameters, but predicting the effect of a mutation on *k*_{cat} or *K*_{M} is highly nontrivial (*26*). We were unable to find studies measuring kinase and phosphatase rates on the same substrate. Instead, as a standard to demonstrate physical realizability, we consider the parameters for sample Ser/Thr enzymes: wild-type λ-phosphatase acting on para-nitrophenylphosphate (*k*_{cat} = 2.0 × 10^{3} s^{−1}; *K*_{M} = 1.0 × 10^{4} μM) and wild-type MST4 (kinase) acting on the short peptide chain NKGYNTLRRKK (*k*_{cat} = 3.1 s^{−1}; *K*_{M} = 14 μM) (*26*, *27*). We assume throughout that the constitutive phosphatase *P*. We also treat the self-assembled enzyme as only one functional protein because the copies of the enzyme are all colocalized.

### Bounded self-assembly can yield oscillations whose behavior is well-predicted by analytical formulae

*The onset of oscillations for the bounded self-assembly system is well-predicted by two dimensionless parameters*. We first consider a system where unphosphorylated kinase and phosphatase monomers self-assemble in an all-or-nothing manner into closed symmetric homomultimers. We denote by *n* the number of kinase monomers κ in the functional kinase multimer *K* and by *m* the analogous number of phosphatase monomers ρ in the multimer *P*. The reaction network is shown in Fig. 1B.

Since the full equations describing this bounded self-assembly system (eq. S1) are too complex to directly tackle analytically, we numerically integrate them within the parameter ranges outlined above (section S5). Our results, shown in Fig. 2, demonstrate a significant portion of parameter space within experimental constraints capable of admitting sustained oscillations.

To simplify these equations to an analytically tractable form, we make the Michaelis-Menten approximation that enzymatic intermediates are in quasi-steady state. We then make the approximation that the enzyme and substrate concentrations are low compared to the Michaelis constants, such that the concentration of enzymatic intermediates can be entirely neglected within our analytical approximations (section S2). This approximation, like others that we will consider, is not obeyed by all oscillating solutions found numerically (Fig. 2) but is nonetheless useful in clarifying the fundamentals of a large swath of the oscillations. We find that, in contrast to well-known examples from other systems that rely on enzyme sequestration to achieve oscillations (*6*, *7*), neglecting enzyme sequestration does not preclude oscillations for our systems.

To reduce our systems further to only two differential equations, we assume a separation of timescales between the self-assembly and the enzymatic activity. In particular, we assume that phosphorylation/dephosphorylation reactions equilibrate much faster than self-assembly. The opposite separation–of–timescales limit yields oscillations only for extremely large values of *m*, which are infeasible to realize experimentally (see further discussion in section S4). We thereby arrive at the following two-dimensional system of equations

We briefly define the parameters: *k*_{bκ} is the binding rate for κ into its multimeric state, *k*_{uκ} is the respective unbinding rate, and *k*_{dκ} is the inverse ratio of the two; η_{Kκ} is the specificity constant *k*_{cat}/*K*_{M} for the kinase *K* acting on κ, η_{Pκ} is the same for the phosphatase *P*, and η_{κ} = η_{Kκ}/η_{Pκ}; κ_{tot} is the total concentration of monomeric κ added to the system, a conserved quantity. Similar quantities are defined for ρ. We assume the concentration of the constitutive phosphatase

To describe the oscillatory behavior of the system, we seek the eigenvalues of the Jacobian in the vicinity of a fixed point (*K*^{⋆}, *P*^{⋆}). Oscillations require coupling between the equations, motivating the approximations that in the oscillatory regime, _{ρ}), κ_{tot} ≫ *nK*^{⋆}, and ρ_{tot} ≫ *mP*^{⋆}. Defining the dimensionless concentrations

We constrain ourselves to *m* ≤ *n* + 1 so that, within our approximations, there is no more than one physical fixed point in the system as long as

Sustained oscillations in the system typically correspond to complex eigenvalues of the Jacobian with positive real parts. However, following the Poincaré-Bendixson theorem, as long as our system has a single fixed point, instability of the fixed point must imply oscillations even if they are beyond the linear regime. Translated into constraints on *P*^{⋆}, instability of the fixed point corresponds to

We now proceed to express Eq. 4 only in terms of the input parameters. Because Eq. 2 cannot be solved for *P*^{⋆} for general *n*, *m*, we consider the approximation that *m* < *n* + 1 to γ ≫ 1. This approximation is most accurate for small values of *m*, since fewer terms are neglected. The approximation is motivated by the intuition that oscillations require *P* to be non-negligible compared to

By simplifying Eq. 4 within this limit where

The implication of Eq. 5, that oscillations require a larger dissociation rate of the phosphatase monomers compared to the kinase monomers (at least, for *n* + 2 ≥ *m*), agrees with intuition found by visualizing individual oscillation cycles (Fig. 3A). In *K-P* space, oscillations proceed in a counterclockwise fashion: Starting from the unphosphorylated state, the monomers self-assemble into multimers (top right). The larger dissociation rate of phosphatase multimers leads those to dissociate first and get phosphorylated by the abundant kinase multimers (bottom right). Next, the kinase multimers slowly dissociate, enabling the gradual dephosphorylation and self-assembly of the phosphatase monomers by the constitutive and self-assembled phosphatases (bottom left). Once kinase multimer levels have decreased and enough phosphatase multimers have formed, the latter quickly dephosphorylate the remaining monomers (top left), and the system returns to its initial state (top right).

We verify that the approximate formula given by Eq. 5 is valid in describing Eq. 1 by comparing it to oscillations found by random parameter searches in Fig. 3B. We numerically integrate Eq. 1 with random parameters chosen to satisfy the experimental constraints described previously (including setting η_{κ} = η_{ρ}) and with *n* = *m* = 2. We constrain concentrations κ_{tot} and ρ_{tot} to be within 10^{−3} to 10^{2} μM, while we set the bounds of ^{−8} and 10^{−2} μM. For each parameter set, we numerically estimate the fixed point using Python’s scipy.optimize.root function. We only show parameter sets estimated to agree with the approximations described before Eq. 2 (with >5× substituted for ≫). We found no oscillations in ∼2.5 × 10^{4} parameter sets for which η_{K} or η_{P} is less than *K*, *P*) = (0,0). Oscillations are almost exclusively found in the quadrant γ > 1, ν > 1. Values of γ slightly less than unity are also found to produce oscillations, as shown in the figure.

These results show that oscillations are found when the dissociation rate of the phosphatase multimer is in a middle range between two extremes. To see oscillations, the phosphatase multimer must dissociate significantly faster than the kinase multimer (ν > 1); however, the dissociation rate must concurrently be small enough such that the fixed point concentration of phosphatase multimer is larger than the constitutive phosphatase concentration (γ > 1). In contrast to intuition from other systems, which signifies that higher-order nonlinearities increase the parameter range producing oscillations (*6*), here, we found that more nonlinear self-assembly (i.e., higher values of *n* and *m*) makes oscillations less frequent. We find that oscillations are robust to even order-of-magnitude variations in most other parameters (section S3 and fig. S2).

Our results suggest that oscillations in the concentrations of unphosphorylated monomers or of phosphatase multimers may be most straightforward to visualize experimentally, as these concentrations typically vary by several orders of magnitude across an oscillatory cycle (fig. S8). Parameter sets yielding largest oscillation amplitudes are also typically farthest from the bifurcation point (fig. S10).

*The frequency of resulting oscillations can be well-predicted by assuming that kinase multimer dissociation is rate limiting*. We next seek to predict how system parameters tune the frequency of resulting oscillations when they appear. Within the linear regime around the fixed point, in the limit of Eq. 5, the frequency of oscillations ω is predicted to be

While Eq. 6 agrees well with the true squared frequency for those parameter sets where it is positive, oscillations are frequently found in the nonlinear regime in which it is not applicable (fig. S7A). However, the intuition given by Eq. 6, that oscillation frequency is determined by the unbinding rates of the kinase and phosphatase multimers, may still be valid outside the linear regime of the fixed point. This intuition is reasonable given the separation–of–timescales limit in which we are operating, of enzymatic reactions equilibrating much faster than self-assembly. Since oscillations for the *n* = *m* = 2 system that we considered numerically require that *k*_{uρ} > *k*_{uκ} (i.e., ν > 1), Eq. 6 implies that the limiting reaction in the oscillations is the unbinding of the kinase multimers. To test that implication, we compare *k*_{uκ} to the frequency of oscillations found through numerical integration in Fig. 3C. (We make no constraints on the fixed points of the parameter sets considered here.) We find a strong correlation (*R*^{2} = 0.66; *R*^{2} = 0.93 in log space) and a root mean square relative error of ∼3.9, demonstrating that *k*_{uκ} can accurately predict the frequency of oscillations in this system.

Oscillations found within experimental constraints for Eq. 1 have periods ranging from fractions of a second to >1 day (fig. S7B). For oscillations found for the full system of equations plotted in Fig. 2, we find periods within a slightly more constrained range than for the simplified system but still spanning orders of magnitude, between ~1 min and >1 day.

### Unbounded self-assembly can yield oscillations within experimental constraints in the limit of fast self-assembly compared to enzymatic activity

We now consider a second system in which individual species κ and ρ can self-assemble incrementally into one-dimensional unbounded fibers (Fig. 4A). Unlike the bounded case in which no phosphorylation sites are accessible in the multimeric state, in this system, one is (corresponding to the final protein in the fiber). An *n*-mer of species *X* (where *X* is either κ or ρ), *X _{n}*, can be created either from binding two smaller molecules

*X*and

_{m}*X*

_{n−m}or from the spontaneous breaking of a bond of a larger molecule. The concentration of

*X*decreases when an

_{n}*X*molecule either binds to any other molecule or breaks any of its

_{n}*n −*1 bonds. The equations for self-assembly of species

*X*are therefore given by

As in the first system, each protein of type κ (ρ) includes a split kinase (phosphatase). While in the bounded system, enzymes required exactly *n* or *m* monomers to self-assemble, here, an enzyme is created by any group of more than one monomer (i.e., *X _{n}* is a functional enzyme as long as

*n*≥ 2). We assume that when a multimer is phosphorylated, its final monomer dissociates from the fiber and cannot reassociate in its phosphorylated state. A less stringent assumption, that the phosphorylated monomer does not dissociate automatically but merely prevents new monomers from binding to that end of the molecule, appears to be incompatible with oscillations, at least in both separation–of–timescales limits.

The full equations describing the system are given in the Supplementary Materials (eq. S8). As previously, we search for a two-dimensional set of simplified equations by considering a separation of timescales between enzymatic reactions and self-assembly, along with the same simplifying assumptions as considered for the bounded self-assembly system (eq. S9). The limit considered for the first system, of fast phosphorylation/dephosphorylation compared to self-assembly, would disallow multimers from forming in this system, since, here, phosphorylation is accompanied by dissociation of the final multimer in the chain. Therefore, the separation–of–timescales limit considered in the bounded self-assembly system is no longer applicable for this system. Instead, we consider the opposite limit, of fast self-assembly compared to enzymatic activity. At steady state, *X _{n}* is given by

*x*= 2

*k*

_{bX}

*X*

_{tot}/

*k*

_{uX}= 2

*X*

_{tot}/

*k*. The same steady state is reached even if self-assembly involves binding and unbinding only a single monomer at a time.

_{dX}The concentrations of phosphorylated monomers as a function of time are given by κ* ^{p}* and ρ

*. The total amount of kinase present is given by*

^{p}*k*= 2(κ

_{tot}− κ

*)/*

^{p}*k*

_{dκ}and

*p*= 2(ρ

_{tot}− ρ

*)/*

^{p}*k*

_{dρ}

Unlike in the previous system for which the frequency and onset of oscillations can be determined by simple formulae by taking a limit of the two-dimensional system, no such limits give similarly straightforward results for Eq. 9. Instead, we analyze Eq. 9 through random parameter searches (Fig. 4B). As shown in the figure, we find that although these equations assume an opposite separation–of–timescales limit to that yielding experimentally realizable oscillations for the system of bounded self-assembly, they can nevertheless yield oscillations within a significant region of parameter space consistent with experimental constraints.

We find that increased values of κ_{tot}/ρ_{tot} and decreased values of

## DISCUSSION

In summary, we have presented two posttranslational protein–based oscillators motivated by the biological KaiABC system and by the synthetic dual-feedback genetic oscillator. Both systems that we present rely on split kinase and phosphatase self-assembling to form functional enzymes and on that self-assembly being inhibited by phosphorylation of the split monomers. The two systems differ mainly in the nature of the self-assembly as all-or-nothing into bounded structures of specified size or incremental into unbounded one-dimensional fibers.

Both systems are capable of producing oscillations within experimental constraints, using experimentally determined wild-type values for kinase and phosphatase activity and for a range of designed self-assembly rates. We have shown that neither complex reactions nor large number of species are necessary to achieve oscillations: Both networks that we present use only three protein species interacting only through reversible binding and phosphorylation/dephosphorylation reactions, and the resulting oscillations can be understood as arising from a minimal system of two differential equations in both cases.

Although the systems that we described shared much in common, they produced robust oscillations in opposite separation–of–timescales limits from one another: The first primarily oscillates when self-assembly is much slower than enzymatic reactions; the second when it is much faster. These two networks are thus complementary: Depending on the parameter regime most easily accessible to an experimentalist, one or the other network might be preferable to implement. Nevertheless, the conditions giving rise to oscillations shared similarities in the two systems. Smaller values of constitutive phosphatase (equivalent to larger values of γ in the case of bounded self-assembly) lead to more robust oscillations in both systems, as do smaller values of *k*_{dκ} (i.e., larger values of ν).

Our work paves the way toward the rational design and experimental realization of protein-based far-from-equilibrium dynamic systems. The models described here were designed to be feasible to synthesize experimentally and are guiding an implementation in the test tube that is currently under way.

## SUPPLEMENTARY MATERIALS

Supplementary material for this article is available at http://advances.sciencemag.org/cgi/content/full/6/51/eabc1939/DC1

This is an open-access article distributed under the terms of the Creative Commons Attribution license, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

## REFERENCES AND NOTES

**Acknowledgments:**We thank A. Murugan for helpful discussions.

**Funding:**This work was supported by the Harvard Materials Research Science and Engineering Center grant no. DMR-1420570 and ONR grant no. N00014-17-1-3029 (to M.P.B.). O.K. acknowledges funding from the DoD through NDSEG Fellowship 32 CFR 168a, the NSF-Simons Center for Mathematical and Statistical Analysis of Biology at Harvard University award no. 1764269, and the Harvard Quantitative Biology Initiative. A.C. is a recipient of the Human Frontier Science Program Long-Term Fellowship and a Washington Research Foundation Senior Fellow. M.P.B. is an investigator of the Simons Foundation. M.P.B. is a Research Scientist at Google Research.

**Author contributions:**O.K., C.P.G., A.I.C., and M.P.B. designed and implemented the numerical analyses. O.K. and A.I.C. designed and performed the analytical analyses. A.C., N.B.W., and D.B. developed the experimental feasibility. O.K. and A.C. implemented the figures. All authors developed the oscillator designs, designed the figures, and wrote the paper.

**Competing interests:**The authors declare that they have no competing interests.

**Data and materials availability:**All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials. The code used in this study is available at https://github.com/ofer-kimchi/protein-oscillator. Additional data related to this paper may be requested from the authors.

- Copyright © 2020 The Authors, some rights reserved; exclusive licensee American Association for the Advancement of Science. No claim to original U.S. Government Works. Distributed under a Creative Commons Attribution License 4.0 (CC BY).