- Research article
- Open Access

# Constructing network topologies for multiple signal-encoding functions

- Lili Wu
^{1}, - Hongli Wang
^{1, 2}Email authorView ORCID ID profile and - Qi Ouyang
^{1, 2, 3}Email author

**Received:**29 May 2018**Accepted:**28 December 2018**Published:**11 January 2019

## Abstract

### Background

Cells use signaling protein networks to sense their environment and mediate specific responses. Information about environmental stress is usually encoded in the dynamics of the signaling molecules, and qualitatively distinct dynamics of the same signaling molecule can lead to dramatically different cell fates. Exploring the design principles of networks with multiple signal-encoding functions is important for understanding how distinct dynamic patterns are shaped and integrated by real cellular networks, and for building cells with targeted sensing–response functions via synthetic biology.

### Results

In this paper, we investigate multi-node enzymatic regulatory networks with three signal-encoding functions, i.e., dynamic responses of oscillation, transient activation, and sustained activation upon step stimulation by three different inducers, respectively. Taking into account competition effects of the substrates for the same enzyme in the enzymatic reactions, we searched for robust subnetworks for each signal-encoding function by three-node-network enumeration and then integrated the three subnetworks together via node-merging. The obtained tri-functional networks consisted of four to six nodes, and the core structures of these networks were hybrids of the motifs for the subfunctions.

### Conclusions

The simplest but relatively robust tri-functional networks demonstrated that the three functions were compatible within a simple negative feedback loop. Depending on the network structure, the competition effects of the substrates for the same enzyme within the networks could promote or hamper the target functions, and can create implicit functional motifs. Overall, the networks we obtained could in principle be synthesized to construct dynamic control circuits with multiple target functions.

## Keywords

- Enzymatic networks
- Signal-encoding
- Network motifs
- Design principle

## Background

The relationship between structure and function in biological circuits is a focus of systems biology [1, 2]. Structure reflects the topology of interactions within a circuit, and function is usually considered a static quantity (a steady state) or a temporal behavior (dynamics) of the circuit’s output [3]. In various signaling systems, cells transmit information via dynamic control of signaling molecules [3]. The information about the identity and quality of a stimulus can be encoded in temporal patterns of a signaling molecule by modulating its amplitude, duration, frequency and so on. Cells mediate gene expression programs and induce cellular responses according to the sequential dynamics of signaling molecules [3–13]. Furthermore, signaling dynamics usually have greater information transmission capacities compared to non-dynamics responses [14–16]. For example, the tumor suppressor p53 shows a series of pulses with number modulation under γ-radiation, resulting in cell-cycle arrest. By contrast, UV radiation triggers a single p53 pulse with a dose-dependent amplitude and duration, leading to apoptosis [3, 17, 18]. In yeast, the transcription factor Msn2 plays a major role in the general stress response. Glucose limitation stress induces nucleocytoplasmic shuttling of Msn2. The duration of the initial peak and frequency of the sequential bursts in Msn2 increase with the strength of the stress [6, 19, 20]. Under osmotic stress, Msn2 exhibits a similarly initial nuclear peak and returns to normal. The duration of the peak increases with the osmotic pressure [3, 21–25]. In addition, oxidative stress leads to prolonged nuclear Msn2 accumulation whose amplitude is dependent on the concentration of H_{2}O_{2} [6, 26, 27]. Other systems, such as ERK, NF-κB and Crz all drive distinct responses via dynamic control [3, 28–30].

These observations in real signaling systems raise fascinating questions about how distinct input information is encoded into different dynamic patterns of the signaling molecule and the design principles for such networks with multiple signal-encoding functions. Instead of considering the complex molecular interaction details in real signaling networks, we here intend to perform a theoretical investigation of multi-node enzymatic regulatory networks with three signal-encoding functions, i.e., dynamic responses of oscillation, transient activation, and sustained activation when a respective input node is stimulated. These dynamic patterns are similar to the observed dynamic response of Msn2 in yeast under glucose limitation, osmotic stress, and oxidative stress, respectively.

General computational search strategies for detecting network topologies for a target function include evolutionary algorithms and the network enumeration approach [2]. These approaches have been used to explore simple network solutions for dynamic functions such as oscillations [31–33], switch-like responses [34, 35], adaptation [36–38] and dose-response alignment [39]. Core motifs and design principles for these individual dynamic functions have been investigated previously. By contrast, multi-functional networks that can execute several distinct functions have rarely been considered [1], despite the accumulation of experimental observations of multiple signal-encoding processes. The stumbling blocks for searching multi-functional networks are mainly due to the increased network scale and the functional complexity, which require enormous computing power. As the multi-functional networks we consider have at least four nodes, with three nodes for inputting three external signals and one node for output, it is practically infeasible to enumerate all possibilities to detect multi-functional network topologies. Evolutionary approaches are also inapplicable for deducing the architectural landscape of tri-functional networks. Accordingly, here we adopt the scenario of module combination for searching multi-functional networks capable of the desired multiple signal-encoding function. The module combination approach has recently been used to generate bi-functional networks by hybridizing simpler mono-functional modules [1]. It was also applied successfully to detect network structures with Pavlovian-like functions capable of learning and recalling [40].

In our calculations, the external stimulus, which is a temporal step function, acts on three different input nodes, and generates on the output node sustained oscillations (stable limit cycle), transient activation (adaptation), and sustained activation, respectively. This mimics the signal-encoding of Msn2 in yeast in response to glucose limitation, osmotic stress, and oxidative stress, respectively. To obtain the target multi-functional networks, we search first for robust subnetworks (sub-functional modules) for each target input-output response by enumerating all three-node network topologies and then meld different modules through node merging. All computations are performed with enzymatic regulatory networks. In addition, we take competition effects among different substrates into account in situations when two or more types of substrates bind competitively to the same enzyme. This ubiquitous hidden interaction has been considered in gene regulatory networks [41], but is commonly omitted in most theoretical or experimental models in enzymatic circuits [42]. We find that competition effects can increase the diversity of functional modules by creating implicit cross couplings within the circuits, whereas in combined networks, competition effects can promote or hamper the overall functionality. These findings may be helpful for understanding multi-functional networks in real biological systems and also provide a guide to construct information-encoding circuits in synthetic biology.

## Results

### Module selection for each signal-processing function

We use the method of three-node-network enumeration to explore the functional networks of oscillation (F_{1}) and adaptation (F_{2}), respectively. In a two- or three-node network, one node serves as the input node (Node R) and one as the output node (Node O), and a third node will be the intermediate node (Node M). For each topology, we sample 10,000 parameter sets using the Latin hypercube sampling method [43]. For each parameter set, the ordinary differential equations that describe the network dynamics in response to one of the two stimuli are then solved numerically. The topologies with at least one or more parameter sets that produce the targeted dynamical function are considered functional circuits. The robustness of sub-functions F_{1} and F_{2} is measured by the number of parameter sets, i.e., Q_{1}- and Q_{2}-values, that generate the target dynamic functions. For functions F_{1} and F_{2}, we select the network topologies with relatively high Q-values and fewer edges as functional F_{1} and F_{2} module pools.

The oscillation (F_{1}) pool contains 73 networks with 3 to 6 edges with Q_{1} > 50 and 8 simplest networks with 3 edges with 10 < Q_{1} < 50 (as marked in Fig. 2a). Clustering analysis of these networks (Fig. 2b) shows that the basic oscillation motifs for achieving function F_{1} are negative feedback loops between the three nodes. Motifs a_{1} and a_{3} are simply two forms of the repressilator, and motifs a_{2} and a_{4} are the delayed feedback loop. All 73 networks contain explicitly one of the two negative feedback loops which are well-known motifs for oscillations. This is true even for networks with edge numbers greater than 6. Additionally, the F_{1} functional pool also contains the 8 simplest networks. They do not explicitly have the structure of a repressilator or delayed feedback loop (Fig. 2c). These topologies achieve self-sustained oscillations because one of the two basic motifs of the repressilator and delayed feedback loop is implicitly contained in their topologies due to the hidden interactions. Two such examples are illustrated in Fig. 2g, where the oscillations are achieved through the implicit a_{3} and a_{4} motifs. Thus, the modular pool of F_{1} consists of 81 networks as marked in Fig. 2a.

Similarly, as marked in Fig. 2d, the F_{2} functional pool also consists of 81 selected networks, i.e., 76 networks with 3 to 6 edges and Q_{2} > 55, and 5 networks with 3 edges and 10 < Q_{2} < 55. Clustering analysis of the 76 networks reveals the backbones of the functional circuits (Fig. 2e), two negative feedback loops with a buffer node (motifs b_{1}, b_{3}), a delayed feedback loop (motif b_{2}), and an incoherent feedforward loop with a proportional node (motif b_{4}). The motifs for function F_{2} are similar to Ma’s results [31], but different in that motifs b_{3} and b_{4} (for which Q_{2} < 55) are less robust than motifs b_{1} and b_{2} (for which Q_{2} > 55). This difference is attributable to the inclusion of competitive regulations in both motifs b_{3} and b_{4}: node M_{3} and node O_{3} compete for enzyme R_{3} to convert from the inactive to active state. Consequently, there exit implicit positive interactions between nodes M_{3} and O_{3}, and the regulatory interactions are more complex than those discussed in [36], in which the competitive effect was ommitted. Figure 2f illustrates the 5 networks with three nodes and lower Q_{2} value as marked in Fig. 2d. These simple functional circuits all have implicit interactions, of which two are motifs b_{3} and b_{4}. The left three motifs have very different topologies from the ordinary adaptation motifs. Checking the implicit interactions reveals that these seemingly very different adaptation motifs all have implicitly the core structures discussed in detail in [36]. Two such examples are illustrated in Fig. 2h where adaptation is achieved through the core structures of b_{1} and b_{3}.

The functional circuits with sustained activation (F_{3}) are not searched by enumeration. We prefer to adopt the simplest structure, i.e., a link of activation from the input node R_{3} to the output node O_{3}. Thus, we have constructed the module pools for the sub-functions of oscillation (F_{1}, with 81 topologies), transient activation (F_{2}, with 81 topologies), and sustained activation (F_{3}, with 1 topology), respectively.

### Module combination for bi-functional networks of F_{1} and F_{2}

_{1}and F

_{2}module pools, we next construct F

_{1}-F

_{2}bi-functional networks that encode the external signal received by nodes R

_{1}and R

_{2}respectively into oscillations and adaptation in the output node O. The network topologies from the two sub-functional pools are united by merging their nodes. As demonstrated in Fig. 3, we adopt the following node-merging rules:

- (i).
Output O nodes have to be merged with each other as we assume a unique output node.

- (ii).
M nodes may be merged with each other (unless there are contradictory linkages).

- (iii).
R nodes may be merged with M nodes from the other pathway (unless there are contradictory linkages).

- (iv).
R nodes for receiving signals cannot be merged with each other.

- (v).
R nodes are not permitted to merge with O nodes.

Here we do not consider the merging between R and O nodes and node merging within pathway. As there are two input nodes for receiving the signals and a unique output node in a bi-functional network, there are three possible ways to combine the F_{1} and F_{2} circuits as illustrated in the upper left of Fig. 3a-c. One approach is to simply merge the output nodes O_{1} and O_{2}. Another is two-node merging of M_{1}-M_{2}, R_{1}-M_{2}, and R_{2}-M_{2} in addition to the O_{1}-O_{2} merging. The third approach is the three-node merging of O_{1}-O_{2}, R_{1}-M_{2}, and R_{2}-M_{1}.

We enumerate all possible networks combined from the two module pools as candidate topologies for bi-functional circuits. For each candidate network, we randomly sample 10,000 parameter sets and check whether it is bifunctional. The number of parameter sets that can achieve both oscillation and adaptation is recorded as the Q_{1,2}-value to estimate the robustness of the bi-function. In the bi-functional networks with Q_{1,2} > 0, we obtain 3215 five-node, 1337 four-node, and 41 three-node networks by one-, two-, and three-node merging, respectively. The Q_{1,2} distributions of these three classes of bi-functional combined networks are illustrated in the lower left of Fig. 3a-c. The networks with relatively high Q_{1,2} are highlighted with red dashed boxes in the histograms of Fig. 3. These robust networks obtained by one-, two- and three-node merging are then clustered (middle panels in Fig. 3a-c). Structure skeletons from the clustering analyses are listed alongside (right panels in Fig. 3a-c). The robust bi-functional networks obtained via one-node and two-node merging are composed of hybrid motifs from the F_{1} and F_{2} module pools, such as hybrids of a_{3} + b_{1}, a_{3} + b_{2}, a_{4} + b_{2}, and a_{1} + b_{3} in five-node combined networks and a_{3} + b_{1} and a_{3} + b_{4} in four-node combined networks. We find that the networks containing the hybrid motif a_{3} + b_{1} have the largest proportion; the mechanism for achieving the bi-function is given in Additional file 1. As three-node merging is full merging of the network, it is difficult to find very robust bi-functional networks with large Q_{1,2} values (no more than 5). Due to regulation conflicts, the F_{1} motifs (a_{1}, a_{2}, a_{3}, and a_{4}) and F_{2} motifs (b_{1}, b_{2}, b_{3}, and b_{4}) cannot be fully merged to form a bifunctional motif. Instead, the bi-functional networks of combination are derived from non-core motifs that can achieve both oscillation and adaptation, such as the structure skeletons c_{1} and c_{2} depicted in Fig. 3c.

By combining the nodes of functional F_{1} and F_{2} networks, the merged nodes function as a joint uniting the sub-functional networks. This results in competing effects between substrates and additional implicit interactions within the networks. This side effect of node merging could promote, hamper, or have little influence on the circuit’s function, depending on the strength of the competing strength. In Fig. 3d-f, we list several examples of bi-functional networks that contain implicit interactions as denoted by the dashed linkages. It is noted that the three-node bi-functional networks belong to the intersection of three-node networks in F1 and F2 pools, and that only a portion of networks in the intersection can achieve the bi-function. The node merging procedure does not discover new three-node networks that are bifunctional. For instance, all compatibly combined networks of oscillation motifs a_{1}-a_{4} and adaptation motifs b_{1}-b_{4} by full merging of the three nodes have a zero Q_{1,2} value.

### Extending the bi-functional networks into tri-functional networks

_{3}) for receiving the third signal. One approach is to adopt an existing regulatory node (i.e., a node other than R

_{1}, R

_{2}and O) as the receptor for input of the third signal. The other strategy is to introduce a new node as R

_{3}into a bifunctional network for receiving the third signal. As the third function requires sustained activation of the output node O, an activation from node R

_{3}to the output node O is simply added. For the circuit in Fig. 4a, the node M

_{1}or M

_{2}can be used as the third receptor node R

_{3}; otherwise, an additional node R

_{3}is added that activates the output node. Both strategies add no new competitive regulations to the networks.

For each candidate tri-functional network obtained from the above extension strategies, we sample 1000 random parameter sets for the newly introduced parameters while keeping the parameters in the original circuits unchanged. The number of parameter sets that can achieve the target tri-functional behavior is recorded as the Q(F_{3}|F_{1,2})-value, which manifests the robustness for the triple signal-encoding function (see Additional file 1). The distribution of *Q*(*F*_{3}| *F*_{1, 2}) for all tri-functional networks obtained is shown in Fig. 4b. There are 3190 tri-functional networks with six nodes, 2092 tri-functional networks with five nodes, and 255 tri-functional networks with four nodes. The networks with *Q*(*F*_{3}| *F*_{1, 2}) > 500 are highlighted in the red frame in Fig. 4b, with 4 six-node networks, 134 five-node networks, and 3 four-node networks.

*Q*(

*F*

_{3}|

*F*

_{1, 2}) > 100 [44]. A dot denotes a tri-functional network, and networks with

*Q*(

*F*

_{3}|

*F*

_{1, 2}) > 500 are represented by large dots. All networks are distributed with respect to the number of links and nodes in the networks. We have linked any two networks that could be transformed by the gain or removal of one regulatory interaction. Most of the tri-functional networks we obtained are connected, and complex tri-functional networks can evolve from simpler core topologies. Several core networks are illustrated at the bottom in Fig. 5c, i.e., the networks labeled from N

_{1}to N

_{7}. Network N

_{4}is constructed by uniting motifs a

_{3}and b

_{1}and adding the R

_{3}node. The network N

_{3}is devised on the bi-functional networks with hybrid motifs of a

_{1}+ b

_{3}via the first strategy. It has two possible nodes (R

_{3}and M) used as a buffer node for F

_{2}and can be reduced to the network N

_{5}by removing node M. The robustness of network N

_{5}is comparable with that of network N

_{3}. Thus, network N

_{5}is the simplest robust tri-functional network in our calculation. In the following section, we perform an analysis of the dynamics of network N

_{5}.

### Dynamics analysis for the simplest tri-functional network

_{5}in Fig. 4c, we demonstrate its multiple signal-encoding dynamics in detail for illustrative purposes. The network dynamics in concentrations of active forms of enzymes R

_{1}, R

_{2}, R

_{3}and O is described by the following ordinary differential equations,

In this description, the total concentrations of enzymes are normalized, and the variables vary from 0 to 1. In the normal condition, the receptor nodes R_{1}, R_{2}, R_{3} sense the surroundings with a base signal *I*_{1} = *I*_{2} = *I*_{3} = 0.1. With the parameter values in Fig. 5, the system rests at its stable steady state (with *R*_{1} = 0.1179, *R*_{2} = 0.3431, *R*_{3} = 0.0024, *O* = 0.1690). When node R_{1} is stimulated with *I*_{1} ≥ 0.2 and *I*_{2} = *I*_{3} = 0.1, the system enters into the regime of self-sustained oscillation (Fig. 5a). Oscillations are generated due to the R_{1}-R_{3}-O negative feedback loop within the circuit. The step signal is thus interpreted by the system into stable limit cycle oscillation in the output node. As depicted in Fig. 5d, the network dynamics undergoes a subcritical Hopf bifurcation at *I*_{1C} = 0.25 when the strength of stimulus I_{1} is tuned as the control parameter.

*I*

_{2}> 0.1,

*I*

_{1}=

*I*

_{3}= 0.1), the circuit’s output responds to the stimulus I

_{2}transiently and returns to its pre-stimulated level. This transient activation results from the saturated regulations of R

_{3}from both signal I

_{3}and output node O, with

*K*

_{3}≪ 1 −

*R*

_{3},

*K*

_{7}≪

*R*

_{3}. In this case, Eq. 1c for R

_{3}is simplified as,

*I*

_{3},

The independence of *O*_{SS} on *I*_{2} results in a perfect adaptation. Thus, the necessary condition for adaptation under the second stimulus requires that the node R_{3} functions as a buffer node within the negative feedback loop in the circuit. Figure 5b shows the responses of the output variable when node R_{2} is triggered by *I*_{2} of different strengths. The steady state of the network dynamics as a function of the signal amplitude *I*_{2} is depicted in Fig. 5e.

When node R_{3} of the circuit receives the third stimulus (*I*_{3} > 0.1, *I*_{1} = *I*_{2} = 0.1) of different amplitudes, the output node is always activated (see Fig. 5c), and the equilibrium state of the output increases with I_{3} strength (Fig. 5f). Therefore, the three signal-encoding processes are mutually compatible and completed in the simplest and robust topology.

## Discussion

Biochemical reactions in cellular signaling processes form complex molecular networks analogous to electrical circuits. There has been a long tradition of theoretical studies that looked into modules of small biochemical networks and signaling motifs that are viewed as simple building blocks of complex networks. Proteins functioning as the transfer and processing of information can be linked into motifs of reoccurring biochemical circuits that perform specific tasks such as to amplify signal, integrate and store information et al. [51, 52]. Motifs of common features that may help formalize the mapping from biochemical reactions to pathway block diagrams in signaling processes have been previously examined [53, 54]. In signaling and transcription regulation networks, the main classes of motifs that can carry out specific dynamic functions such as bistability, adaptation, oscillations, excitable pulse and et al. have been reviewed in references such as [55] and [2]. Compared to previous studies that considered mainly mono-functional network motifs, we here report a systematic approach for extracting multi-functional motifs, i.e., depending on different external signals, the simple motifs can output dynamic responses of oscillation, transient activation, and sustained activation respectively.

The second novel point of this study lies in that competition between substrates for the same enzyme has been taken into account in our mathematical description of signaling circuits. Substrate competition that was usually omitted in previous theoretical studies can have pronounced effects. It has been reported that ultra-sensitivity can be generated through substrate competition between two sets of phosphorylation sites [56] or through negative cooperativity in which the binding of ligand to a multimeric receptor makes it more difficult for subsequent ligands to bind [57]. Competitive binding can lead to high-dose hook effect (also known as prozone effect) in reactions where one protein acts as a linker between parts of a complex [58]. In the interaction of Ca^{2+} and calcium binding proteins calmodulin (CaM), the competitive binding among CaM’s seven partners can tune the Ca^{2+}/CaM binding frequency [59]. In this study, we demonstrate that competitive effects could add implicit interactions within the networks that can promote or hamper the anticipated behaviors, and even create new dynamics unexpected from the explicit network topologies. For simple network structures without explicit functional motifs, the additional hidden regulations within the competitive nodes, combined with the rest of the linkages in the networks, could form implicit functional motifs that foster the networks to achieve target functions. For complex network topologies with dense linkages within the network, such as in the bi-functional networks obtained by two-node merging, the competitive effects could produce more complicated connections within the circuits that can drastically influence the network dynamics.

The multiple signal-encoding networks we consider involve at least four nodes. It is time-consuming to exhaust computationally all possible network topologies. We have adopted a compromise approach to efficiently construct multi-functional networks by enumerating modules of three nodes that can achieve each sub-function and then combing them by node merging. As the exploration is not exhaustive, the tri-functional networks we obtained may have occupied only a part of all functional topologies. As demonstrated, we have constructed tri-functional networks with four to six nodes via module combinations. The functions of oscillation and adaptation in the tri-functional networks could share a negative feedback loop, as long as the circuit could be driven to oscillate by the first signal and a node within as a buffer node for adaptation in the presence of the second signal. Alternatively, the function of adaptation could also be achieved via approaches such as the negative feedback loop with a buffer node (M_{2}) in networks N_{1}, N_{2}, or N_{3}. The third function of sustained activation, it requires at least one activation regulation directly or indirectly from node R_{3} to the output node O.

Compared to the functional networks we designed, the real signal-encoding network of Msn2 in *Saccharomyces cerevisiae* is much more complex [45–50]. A recent report that integrated experiments and modelling revealed that, the stimulus-dependent Msn2 dynamics is controlled by coupled feedback loops [60]. At the core of this Msn2 pathway is the Ras–cAMP–PKA negative feedback loop which is crucial in shaping Msn2 dynamics. It is coupled to the signal-dependent mutual inhibition between protein kinase PKA and the yeast AMP-activated protein kinase Snf1. The negative feedback regulating PKA activity at the core of Msn2 pathway is abundant as sub-module in the functional networks that we identified here. As can be seen in Fig. 3 and Fig. 4, the negative feedback appears in forms of activator-inhibitor, delayed negative feedback, and repressilator. In the networks we found, coupled negative and positive feedback loops are not apparent. But when implicit interactions due to substrate competing are taken into account, such as those demonstrated in Fig. 2g, Fig. 2h, and in Fig. 3d, Fig. 3e, Fig. 3f, coupled negative and positive feedbacks are abundantly found as sub-modules in the networks identified here.

In p53 signaling pathway, the information of hundreds of cellular stress stimuli is encoded through the single node of p53 protein, and stimulus-specific p53 dynamics then triggers distinct cellular outcomes. While the real p53 signaling pathway is very complex, modelling studies [61] revealed that the flexibility of p53 activity is rooted in the core p53-Mdm2 negative feedback loop in the pathway: the transcription factor p53 activates the expression of Mdm2, and Mdm2 in turn inhibits p53 by either turning down its transcription or triggering its degradation. In the multifunctional networks that we find here, the negative feedback loop featuring the p53 pathway is also a fundamental sub-module. Negative feedback loops involving the output node are ubiquitous in bi- and tri- functional circuits we find here. In single cell study of p53 dynamics [17], γ-radiation and UV light cause DNA double strand breaks (DSBs) and single-stranded DNA (ssDNA) damage, respectively. It was uncovered that the mechanism of stimulus-specific response is relevant to two core coupled negative feedback loops. The upstream kinases ATM responding to DSBs and ATR responding to ssDNA both relay the damage signal to p53. This activates p53-Mdm2 and p53-Wip1 negative feedback loops. The core topology of two coupled negative feedback loops in p53 pathway in response to DNA damages is closely related to some of the networks that we identify here. As depicted in Fig. 3a, two coupled negative feedback loops are characteristic of the combined a_{3} + b_{2} and a_{4} + b_{2} networks, and apparently feature the combined networks in Fig. 3d and Fig. 3e. Although the tri-functional networks we constructed are simple and omit the real aspects of dynamic control, they could still be helpful for understanding the mechanism of real signal-encoding processes. Furthermore, the networks we obtained and the method we adopted could also be used in synthetic biology to construct multiple functional dynamic circuits for practical purposes.

## Conclusions

The three signal-processing functions are compatible in a well-designed multi-functional network, in which the output dynamics is switched between oscillation, adaptation, and sustained activation by changing the input node for the signals, with no need to change the internal links in the networks. Competition between substrates for limited enzymes plays as a framework that renders implicit interactions between the nodes that are not explicitly linked. The substrate competition brings virtual links between nodes in the network and can create new implicit motifs that are seemingly nonfunctional.

## Methods

We aim to search for networks capable of encoding the identity of distinct stimuli into the dynamics of the circuit’s output as a whole. The dynamics of oscillation, adaptation (transient activation) and sustained activation are the three target input-output responses considered. We assume that the functional network has three independent receptors each receiving a respective stimulus signal. Under normal conditions, all input signals are set as a low baseline, and the functional network remains in the resting state. When facing one of the three stimuli, the input signal is represented by a step function (Fig. 1a). The information about the stimulus is processed through the unknown components of the networks (Fig. 1b), and then the identity of the stimulus is encoded in the dynamics of the circuit’s output. The target dynamic outputs are illustrated in Fig. 1c as oscillations (F_{1}), transient activation (F_{2}), and sustained activation (F_{3}). We intend to find the unknown networks in Fig. 1b that generate the target output when each receptor node receives its stimulus from an external signal.

In our calculations, the networks we consider are limited to enzymatic regulatory interactions. Each node in the network represents an enzymatic protein with a fixed total concentration (normalized to 1) and can be interconverted between its active and inactive forms by other enzymes that are active. Two types of links are considered. The arrow in the network denotes activation of the target enzyme, that is, the enzyme catalytically transforms its substrate from the inactive form to active form. The blunt-headed arrow indicates catalytic inactivation of the active substrate. In our model, we assume that nodes in the network can possibly be auto-activated or auto-inhibited. If a node in the network has no activation or deactivation from other nodes, it is assumed to be regulated by basally available enzymes in the background.

_{1}in the network converts both the enzymes M

_{1}and O

_{1}to active states M

^{*}

_{1}and O

^{*}

_{1}(Fig. 1d), the activation rate (RA) of M

_{1}and O

_{1}take the following forms (see additional file 1 for the full dynamical equations):

where *R*_{1}, *M*_{1} and *O*_{1} are the concentrations of active enzymes R^{*}_{1}, M^{*}_{1} and O^{*}_{1}, *K*_{1} and *K*_{2} denote the Michaelis-Menten constants, and *V*_{1} and *V*_{2} denote the regulation rate constants. Due to competition between enzymes M_{1} and O_{1} for binding enzyme R^{*}_{1}, an increase in the concentration of O^{*}_{1} and therefore a decrease in inactive O_{1} (due to the constant total concentration) would increase the concentration of enzymes available for the inactive M_{1}. This leads to an increased rate of M^{*}_{1} production and vice versa. The situation is similar when enzyme R^{*}_{1} deactivates both active enzymes M^{*}_{1} and O^{*}_{1}. In these two cases, the competition effect results in implicit positive interactions between the two nodes as illustrated in Fig. 1d and Fig. 1e (note the dashed arrows in the circuits). Figure 1f illustrates the situation when enzyme R^{*}_{1} activates M_{1} but deactivates O^{*}_{1}, resulting in implicit negative interactions of suppression. These competitive effects inducing implicit positive or negative interactions in enzymatic circuits have commonly been omitted in most theoretical or experimental models.

For a concrete demonstration of the competition effect between substrates, Fig. 1g and Fig. 1h demonstrate with a simple three-node network the distinct temporal behaviors generated with modified Michaelis-Menten equations and the normal Michaelis-Menten equations, respectively. With the same set of parameter values and same initial conditions, the description that takes into account competitive bindings depicts self-sustained oscillations (Fig. 1g) while a stationary steady state is generated when the normal Michaelis-Menten mechanism is adopted (Fig. 1h). The ordinary differential equations for the circuit in Fig. 1g and Fig.1 h are listed in additional file 1. One sees that the competition effect generates implicit mutual activation interactions as demonstrated with dashed links in Fig. 1g and sustains a completely different dynamical behavior. From the circuit that is added with two more edges, the oscillations not expected in the original topology are supported due to that a negative feedback loop forms with the aid of the implicit interactions coming from the competition effects. The competition effect could thus play significant roles in the topology of functional circuits and need to be considered in the design principles of multiple signal-encoding functional networks.

We adopt the method of module combination to construct the multi-functional networks. The scenario consists of three steps. First, we search separately for functional network topologies capable of oscillation and adaptation by enumerating all possible networks with three nodes. Typical and simple functional networks, i.e., those topologies that have less links and are relatively more robust, are then selected to form module pools for these two signal-encoding sub-functions. The combined network topologies are then constructed from the two sub-functional pools by merging nodes. Those that are capable of achieving both oscillation and adaptation in response to external stimuli are selected as bi-functional circuits. The final topologies for our target multiple signal-encoding functions are completed by adding to the bi-functional circuits an extra node that directly activates the output node upon receiving the signal, or alternatively by adopting an existing regulatory node for input of the third signal that directly or indirectly activates the output node. Ultimately, we obtain tri-functional networks having four to six nodes.

The modified Michaelis-Menten equations for network dynamics are integrated numerically. The code is written in computer C language with Visual Studio 2008. The ODE solver *rkf45* of Runge-Kutta algorithm is chosen. The scripts are available upon request. Rate constant values ranging from 0.1 to 10 are considered. The Michaelis-Menten constant changes in the range 10^{−3}~10. In order to determine computationally whether oscillations emerge out or not, the variance of network output variable in a time window is monitored after the transient behavior dies out. An apparently not zero variance indicates oscillations arise. The oscillations with both amplitude and oscillation period larger than 0.1 are considered in this study. To numerically determine adaptation, the transient peak value O_{peak} and asymptotic end value O_{end} of the output variable are monitored and calculate their differences from the initial value O_{ini}. The numerical criteria for adaptation are that *O*_{peak} − *O*_{ini} > 0.2 and ∣*O*_{end} − *O*_{ini} ∣ < 0.1 are both satisfied. Sustained activation is determined by the criterion *O*_{end} − *O*_{ini} > 0.2.

In order to analyze topological properties of three-node networks with oscillations and adaptation, we perform clustering analyses for the functional networks as in [39]. In a three-node network, there are totally nine possible links because there are at most three links from each node to other nodes and to itself. Values of 1, − 1, and 0 are assigned for links of activation, inhibition, and no regulation, respectively. A network topology can be therefore described with a digital number of length nine. The topological difference between two networks having three nodes can then be measured by calculating the hamming pair-wise distance from the digital numbers. The clustering analyses of the functional networks is performed by using *clustergram* in Matlab. The clustering for bi-functional networks is similarly analyzed.

The original c code files used in our computation for enumeration of sub-functions (Additional file 2), and merging of functional modules (Additional files 3, 4, 5) as well as a model in SBML format (Additional file 6) for the typical network (N5 in Fig. 4c) can be found in Additional files.

## Declarations

### Acknowledgements

The authors thank the financial support of the Ministry of Science and Technology of China (2015CB910300 to HW) and National Natural Science Foundation of China (11434001 to QO).

### Funding

HW was supported by the Ministry of Science and Technology of China (2015CB910300). QO was supported by the National Natural Science Foundation of China (11434001).

### Availability of data and materials

The codes and SBML file for generating the data are available in Additional files.

### Consent to publish

Not applicable.

### Authors’ contributions

QO initiated the work; HW, QO and LW designed the research; LW performed simulations; LW and HW wrote the paper.

### Ethics approval and consent to participate

Not applicable.

### Competing interests

The authors declare that they have no competing interests.

### Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Open Access**This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.

## Authors’ Affiliations

## References

- Jimenez A, Cotterell J, Munteanu A, Sharpe J. A spectrum of modularity in multi-functional gene circuits. Mol Syst Biol. 2017;13(4):925.PubMedPubMed CentralView ArticleGoogle Scholar
- Lim WA, Lee CM, Tang C. Design principles of regulatory networks: searching for the molecular algorithms of the cell. Mol Cell. 2013;49(2):202–12.PubMedPubMed CentralView ArticleGoogle Scholar
- Purvis JE, Lahav G. Encoding and decoding cellular information through signaling dynamics. Cell. 2013;152:954–6.View ArticleGoogle Scholar
- Hao N, Budnik BA, Gunawardena J, O'Shea EK. Tunable signal processing through modular control of transcription factor translocation. Science. 2013;339(6118):460–4.PubMedPubMed CentralView ArticleGoogle Scholar
- Levine JH, Lin Y, Elowitz MB. Functional roles of pulsing in genetic circuits. Science. 2013;342(6163):1193–200.PubMedPubMed CentralView ArticleGoogle Scholar
- Hao N, O'Shea EK. Signal-dependent dynamics of transcription factortranslocation controls gene expression. Nat Struct Mol Biol. 2012;19:31–40.View ArticleGoogle Scholar
- Locke JCW, Young JW, Fontes M, Jiménez MJH, Elowitz MB. Stochastic pulse regulation in bacterial stress response. Science. 2011;334(6054):366.PubMedPubMed CentralView ArticleGoogle Scholar
- Sadeh A, Movshovich N, Volokh M, Gheber L, Aharoni A. Fine-tuning of the Msn2/4–mediated yeast stress responses as revealed by systematic deletion of Msn2/4 partners. Mol Biol Cell. 2011;22(17):3127–38.PubMedPubMed CentralView ArticleGoogle Scholar
- Mitchell A, Wei P, Wendell AL. Oscillatory stress stimulation uncovers an Achilles' heel of the yeast MAPK signaling network. Science. 2015;350:1379–83.PubMedPubMed CentralView ArticleGoogle Scholar
- AkhavanAghdam Z, Sinha J, Tabbaa OP, Hao N. Dynamic control of gene regulatory logic by seemingly redundant transcription factors. eLife. 2016;5:e18458.PubMedPubMed CentralView ArticleGoogle Scholar
- Li Y, Roberts J, AkhavanAghdam Z, Hao N. Mitogen-activated protein kinase (MAPK) dynamics determine cell fate in the yeast mating response. J Biol Chem. 2017;292:20354.PubMedPubMed CentralView ArticleGoogle Scholar
- Benzinger D, Khammash M. Pulsatile inputs achieve tunable attenuation of gene expression variability and graded multi-gene regulation. Nat Commun. 2018;9:3521.PubMedPubMed CentralView ArticleGoogle Scholar
- Nandagopal N, Santat LA, LeBon L, Sprinzak D, Bronner ME, Elowitz MB. Dynamic ligand discrimination in the notch signaling pathway. Cell. 2018;172:869–80.PubMedView ArticleGoogle Scholar
- Selimkhanov J, Taylor B, Yao J, Pilko A, Albeck J, Hoffmann A, et al. Accurate information transmission through dynamic biochemical signaling networks. Science. 2014;346(6215):1370.PubMedPubMed CentralView ArticleGoogle Scholar
- O'Shaughnessy EC, Palani S, Collins JJ, Sarkar CA. Tunable signal processing in synthetic MAP kinase cascades. Cell. 2011;144(1):119–31.PubMedPubMed CentralView ArticleGoogle Scholar
- Hansen AS, O'Shea EK. Promoter decoding of transcription factor dynamics involves a trade-off between noise and control of gene expression. Mol Syst Biol. 2013;9:704.PubMedPubMed CentralView ArticleGoogle Scholar
- Batchelor E, Loewer A, Mock C, Lahav G. Stimulus-dependent dynamics of p53 in single cells. Mol Syst Biol. 2011;7(1):488.PubMedPubMed CentralView ArticleGoogle Scholar
- Purvis JE, Karhohs KW, Mock C, Batchelor E, Loewer A, Lahav G. p53 dynamics control cell fate. Science. 2012;336(6087):1440–4.PubMedPubMed CentralView ArticleGoogle Scholar
- Garmendia-Torres C, Goldbeter A, Jacquet M. Nucleocytoplasmic oscillations of the yeast transcription factor Msn2: evidence for periodic PKA activation. Curr Biol. 2007;17(12):1044–9.PubMedView ArticleGoogle Scholar
- De Wever V, Reiter W, Ballarini A, Ammerer G, Brocard C. A dual role for PP1 in shaping the Msn2-dependent transcriptional response to glucose starvation. EMBO J. 2005;24(23):4115–23.PubMedPubMed CentralView ArticleGoogle Scholar
- Hohmann S. Osmotic stress signaling and Osmoadaptation in yeasts. Microbiol Mol Biol Rev. 2002;66(2):300–72.PubMedPubMed CentralView ArticleGoogle Scholar
- Garcia R, Rodriguez-Pena JM, Bermejo C, Nombela C, Arroyo J. The high osmotic response and cell wall integrity pathways cooperate to regulate transcriptional responses to zymolyase-induced cell wall stress in Saccharomyces cerevisiae. J Biol Chem. 2009;284(16):10901–11.PubMedPubMed CentralView ArticleGoogle Scholar
- Mettetal JT, Muzzey D, Gomez-Uribe C, van Oudenaarden A. The frequency dependence of osmo-adaptation in Saccharomyces cerevisiae. Science. 2008;319(5862):482–4.PubMedPubMed CentralView ArticleGoogle Scholar
- Gat-Viks I, Shamir R. Refinement and expansion of signaling pathways: the osmotic response network in yeast. Genome Res. 2007;17(3):358–67.PubMedPubMed CentralView ArticleGoogle Scholar
- Gennemark P, Nordlander B, Hohmann S, Wedelin D. A simple mathematical model of adaptation to high osmolarity in yeast. In Silico Biol. 2006;6(3):193–214.PubMedGoogle Scholar
- Charizanis C, Juhnke H, Krems B, Entian KD. The oxidative stress response mediated via Pos9/Skn7 is negatively regulated by the Ras/PKA pathway in Saccharomyces cerevisiae. Mol Gen Genet. 1999;261(4–5):740–52.PubMedView ArticleGoogle Scholar
- Hasan R, Leroy C, Isnard AD, Labarre J, Boy-Marcotte E, Toledano MB. The control of the yeast H2O2 response by the Msn2/4 transcription factors. Mol Microbiol. 2002;45(1):233–41.PubMedView ArticleGoogle Scholar
- Tay S, Hughey JJ, Lee TK, Lipniacki T, Quake SR, Covert MW. Single-cell NF-κB dynamics reveal digital activation and analog information processing in cells. Nature. 2010;466(7303):267.PubMedPubMed CentralView ArticleGoogle Scholar
- Sasagawa S, Y-i O, Fujita K, Kuroda S. Prediction and validation of the distinct dynamics of transient and sustained ERK activation. Nat Cell Biol. 2005;7(4):365–73.PubMedView ArticleGoogle Scholar
- Huang S, Ingber DE. Shape-dependent control of cell growth, differentiation, and apoptosis: switching between attractors in cell regulatory networks. Exp Cell Res. 2000;261(1):91–103.PubMedView ArticleGoogle Scholar
- Wagner A. Circuit topology and the evolution of robustness in two-gene circadian oscillators. Proc Natl Acad Sci U S A. 2005;102(33):11775–80.PubMedPubMed CentralView ArticleGoogle Scholar
- Castillo-Hair SM, Villota ER, Coronado AM. Design principles for robust oscillatory behavior. Syst Synth Biol. 2015;9(3):125–33.PubMedPubMed CentralView ArticleGoogle Scholar
- Kobayashi Y, Shibata T, Kuramoto Y, Mikhailov AS. Evolutionary design of oscillatory genetic networks. Euro Phys J B. 2010;76(1):167–78.View ArticleGoogle Scholar
- Shah NA, Sarkar CA. Robust network topologies for generating switch-like cellular responses. PLoS Comput Biol. 2011;7(6):e1002085.PubMedPubMed CentralView ArticleGoogle Scholar
- Zeng W, Du P, Lou Q, Wu L, Zhang HM, Lou C, et al. Rational Design of an Ultrasensitive Quorum-Sensing Switch. ACS Synth Biol. 2017;6(8):1445–52.PubMedView ArticleGoogle Scholar
- Ma W, Trusina A, El-Samad H, Lim WA, Tang C. Defining network topologies that can achieve biochemical adaptation. Cell. 2009;138(4):760–73.PubMedPubMed CentralView ArticleGoogle Scholar
- Shi W, Ma W, Xiong L, Zhang M, Tang C. Adaptation with transcriptional regulation. Sci Rep. 2017;7:42648.PubMedPubMed CentralView ArticleGoogle Scholar
- Paul F, Eric DS. A case study of evolutionary computation of biochemical adaptation. Phys Biol. 2008;5(2):026009.View ArticleGoogle Scholar
- Yan L, Ouyang Q, Wang H. Dose-response aligned circuits in signaling systems. PLoS One. 2012;7(4):e34727.PubMedPubMed CentralView ArticleGoogle Scholar
- Xi JY, Ouyang Q. Using sub-network combinations to scale up an enumeration method for determining the network structures of biological functions. PLoS One. 2016;11(12):e0168214.PubMedPubMed CentralView ArticleGoogle Scholar
- Bintu L, Buchler NE: Garcia HG, Gerland U, Hwa T, Kondev J, et al.: Transcriptional regulation by the numbers: models. Curr Opin Genet Dev 2005;15(2):116–124.Google Scholar
- Rondelez Y. Competition for catalytic resources alters biological network dynamics. Phys Rev Lett. 2012;108(1):018102.PubMedView ArticleGoogle Scholar
- Iman RL, Davenport JM, Zeigler DK: Latin hypercube sampling (program user's guide).[LHC, in FORTRAN]. Sandia Labs., Albuquerque, NM (USA); 1980.Google Scholar
- Cotterell J, Sharpe J. An atlas of gene regulatory networks reveals multiple three-gene mechanisms for interpreting morphogen gradients. Mol Syst Biol. 2010;6:425.PubMedPubMed CentralView ArticleGoogle Scholar
- Gutin J, Sadeh A, Rahat A, Aharoni A, Friedman N. Condition-specific genetic interaction maps reveal crosstalk between the cAMP/PKA and the HOG MAPK pathways in the activation of the general stress response. Mol Syst Biol. 2015;11(10):829.PubMedPubMed CentralView ArticleGoogle Scholar
- Chasman D, Ho YH, Berry DB, Nemec CM, MacGilvray ME, Hose J, et al. Pathway connectivity and signaling coordination in the yeast stress-activated signaling network. Mol Syst Biol. 2014;10(11):759.PubMedPubMed CentralView ArticleGoogle Scholar
- Jacquet M, Renault G, Lallet S, De Mey J, Goldbeter A. Oscillatory nucleocytoplasmic shuttling of the general stress response transcriptional activators Msn2 and Msn4 in Saccharomyces cerevisiae. J Cell Biol. 2003;161(3):497–505.PubMedPubMed CentralView ArticleGoogle Scholar
- Natalia Petrenkoa RzVC, Megan N. McCleanc, Alexandre V. Morozovb,d, James R. Broacha: Noise and interlocking signaling pathways promote distinct transcription factor dynamics in response to different stresses. Mol Biol Cell 2013;24:2045–2057.Google Scholar
- Morano KA, Grant CM, Moye-Rowley WS. The response to heat shock and oxidative stress in Saccharomyces cerevisiae. Genetics. 2012;190(4):1157–95.PubMedPubMed CentralView ArticleGoogle Scholar
- Boisnard S, Lagniel G, Garmendia-Torres C, Molin M, Boy-Marcotte E, Jacquet M, et al. H2O2 activates the nuclear localization of Msn2 and Maf1 through thioredoxins in Saccharomyces cerevisiae. Eukaryot Cell. 2009;8(9):1429–38.PubMedPubMed CentralView ArticleGoogle Scholar
- Bray D. Protein molecules as computational elements in living cells. Nature. 1995;376:307–12.PubMedView ArticleGoogle Scholar
- Bhalla US, Iyengar R. Emergent properties of networks of biological signaling pathways. Science. 1999;283:381–7.PubMedView ArticleGoogle Scholar
- Bhalla US. The chemical organization of signaling interactions. Bioinformatics. 2002;18:855–63.PubMedView ArticleGoogle Scholar
- Tyson JJ, Chen KC, Novak B. Sniffers, buzzers, toggles and blinkers: dynamics of regulatory and signaling pathways in the cell. Curr Opin Cell Biol. 2003;15(2):221–31.PubMedPubMed CentralView ArticleGoogle Scholar
- Shoval O, Alon U. SnapShot: network motifs. Cell. 2010;143:326–e1.PubMedView ArticleGoogle Scholar
- Kim SY, Ferrell JJE. Substrate competition as a source of ultrasensitivity in the inactivation of Wee1. Cell. 2007;128:1133–45.PubMedView ArticleGoogle Scholar
- Ha SH, Ferrell JJE. Thresholds and ultrasensitivity from negative cooperativity. Science. 2016;352:990–3.PubMedPubMed CentralView ArticleGoogle Scholar
- Roy RD, Rosenmund C, Stefan MI. Cooperative binding mitigates the high-dose hook effect. BMC Syst Biol. 2017;11:74.PubMedPubMed CentralView ArticleGoogle Scholar
- Romano DR, Pharris MC, Patel NM, Kinzer-Ursem TL. Competitive tuning: Competition’s role in setting the frequency- dependence of Ca
^{2+}-dependent proteins. PLoS Comput Biol. 2017;13(11):e1005820.PubMedPubMed CentralView ArticleGoogle Scholar - Jiang Y, AkhavanAghdam Z, Tsimring ALS, Hao N. Coupled feedback loops control the stimulus-dependent dynamics of the yeast transcription factor Msn2. J Biol Chem. 2017;292:12366–72.PubMedPubMed CentralView ArticleGoogle Scholar
- Hunziker A, Jensen MH, Krishna S. Stress-specific response of the p53-Mdm2 feedback loop. BMC Syst Biol. 2010;4:94.PubMedPubMed CentralView ArticleGoogle Scholar