Adaptive Responses Limited by Intrinsic Noise

Sensory systems have mechanisms to respond to the external environment and adapt to them. Such adaptive responses are effective for a wide dynamic range of sensing and perception of temporal change in stimulus. However, noise generated by the adaptation system itself as well as extrinsic noise in sensory inputs may impose a limit on the ability of adaptation systems. The relation between response and noise is well understood for equilibrium systems in the form of fluctuation response relation. However, the relation for nonequilibrium systems, including adaptive systems, are poorly understood. Here, we systematically explore such a relation between response and fluctuation in adaptation systems. We study the two network motifs, incoherent feedforward loops (iFFL) and negative feedback loops (nFBL), that can achieve perfect adaptation. We find that the response magnitude in adaption systems is limited by its intrinsic noise, implying that higher response would have higher noise component as well. Comparing the relation of response and noise in iFFL and nFBL, we show that whereas iFFL exhibits adaptation over a wider parameter range, nFBL offers higher response to noise ratio than iFFL. We also identify the condition that yields the upper limit of response for both network motifs. These results may explain the reason of why nFBL seems to be more abundant in nature for the implementation of adaption systems.


Introduction
Sensory systems have a number of mechanisms to keep track of their environment and to respond appropriately to any changes. Adaptation is one such mechanism where the system shows a transient response to a temporal perturbation followed by a return back to its original state [1][2][3][4][5]. Such adaptive response contributes to detecting stimulus over a wide dynamic range of environmental conditions [6]. By responding to perturbations, adaptation systems can also perform sensing of temporal changes in stimulus. Such a temporal sensing capability of adaptation system is used in bacterial chemotaxis [7][8][9][10][11]. An adaptive response has also been shown to occur in many eukaryotic systems, such as yeast osmoregulation [12], olfactory receptor cells [13] and calcium homeostasis [14]. In chemotactic cell Dictyostelium, an adaptive response is used to sense the direction of wave propagation in the chemoattractant concentration field [15][16][17]. While these diverse systems show adaptive response, common principles underlie their mechanism. It has been reported that the perfect adaptation can be achieved by two types of network motifs [18]. One is incoherent feedforward loop (iFFL), and the other is negative feedback loop (nFBL). For instance, the chemotaxis receptor in bacteria uses feedback mechanism to achieve perfect adaptation [19], whereas the adaptive response in chemotaxis signaling pathway of Dictyostelium cells is operated by a feed-forward mechanism [20].
Adaptive responses are known to be extremely sensitive to small changes in the stimulus [21]. Such a sensitive response can be achieved by amplification mechanisms, such as molecular cooperativities among subunits and a push-pull type reaction of antagonistic enzymes [22,23]. Perfect adaptation systems have a larger tendency towards more sensitive response as compared with imperfect adaptation [24]. Generally, sensitive amplification mechanisms are also sensitive to stochastic variation in the stimulus, and they may further generate large stochastic variations in their activities [25]. Because the cellular reactions are often operated by a small number of molecules, such stochastic fluctuations are inevitable. For a relatively simple reaction that works under the thermal equilibrium conditions and reactions which work essentially with one degree of freedom, a relation between fluctuations and responses has been known. For non-equilibrium condition, relations between fluctuation and responses at steady state values has been also studied (see [26] and references cited therein). Adaptive mechanism works far from thermodynamic equilibrium and thus adaptive response requires free energy input [27,28]. Therefore, a simple relation between response and fluctuation may not be expected as in the situations close to thermodynamic equilibrium.
Previously, a linear relationship between behavioral variability and response time in bacterial chemotaxis has been found experimentally [29]. A theoretical study indicates that such a linear relationship can hold in a limited parameter range, otherwise a combination of fluctuations and nonlinearity generally produces a complex response with many time scales, which deviates from the linear relationship [30]. A relation between noise suppression and responsiveness in feedback loops for non adaptive systems has been reported previously [31,32]. For adaptive networks, their noise filtering [33,34], noise generation and noise propagation characteristic [35,36] have been studied. However, how the relation between response in adaptation system, and its steady state fluctuations before perturbations has not been explored yet. Such a relation can be important to understand the limitations of adaptive response in a noisy enviornment. While with a high response magnitude, smaller perturbation can be detected by the adaptive system, it becomes difficult to distinguish between small random fluctuations and deterministic changes in the presence of noise. Therefore, we study the quantitative relation between the peak response and fluctuations before perturbation in adaptation reactions. Here, we focus on the intrinsic noise of the adaptive response, that is the noise generated in the adaptation reaction, and ignore the fluctuations in the input perturbation. Thus, in this paper, we consider the adaptive network motifs of iFFL and nFBL that consist of two signaling molecules as shown in Fig 1(B) and 1(C). Using numerical simulation and theoretical analysis, we first show that the peak amplitude in the adaptive response is limited by noise. We then study the conditions when the peak amplitude achieves its upper limit, and also study the key differences between iFFL and nFBL in this context.
can take active and inactive states. We denote the active form by A a and B a and the inactive form by A i and B i , respectively. These states can inter convert into each other by enzymatic reactions. Thus, the total numbers of molecules of A and B, denoted by N A and N B , respectively, are constant in time. For the sake of convenience, the numbers of A a , B a and S are denoted by A, B and S, respectively.
In iFFL, the input node S activates both A and B, and the active form B a inactivates A a (incoherent feedforward loop, Fig 1B). Thus, S regulates A through two regulatory pathway with opposing regulatory effects. The inactivation of B is carried out by a basal enzyme whose number remains fixed, and whose concentration we assume to be unity. This reaction scheme can be represented by where the rates of reaction per molecule for activation and deactivation, for the species A and B, are given by G A a Â S, G A d Â B, G B a Â S and G B d , respectively. In nFBL, the input node S regulates only the output node A (Fig 1B). While A activates B, node B inactivates A, completing the negative feedback loop. Node B decays to the inactive form by the basal enzyme whose number is fixed with unity concentration. The reaction scheme for the nFBL can thus be written as: These reaction rates are described by the Michaelis-Menten kinetics. For instance, the total activation rate for species A for the iFFL can be written as where V SA is the maximum reaction rate, and K SA is the Michaelis constant giving the half maximum concentration. The reaction rate for the deactivation can be written similarly. The total reaction rates are given by For the activation of B, the rate for iFFL is given by whereas for nFBL We shall denote the steady state numbers of A and B by A Ã and B Ã respectively.
This makes the number of node A independent of the input, as it is activated by one pathway and simultaneously deactivated by the other. In the nFBL given by Eqs (2), (4), (5), (6), and (8), the perfect adaptation of A a is achieved when the reaction rates for B, Eqs (6) and (8), do not depend on the number of B at the reaction steady state. Such a situation is obtained if K B ( B and K AB ( N B − B in Eqs (6) and (8), respectively. Then, the reaction equilibrium in B gives V B = AV AB , leading towards A independent of S, as The correlation between gain and noise in adaptation reactions When a step stimulus is applied, the adaptive network shows a response by changing its activity A, which is followed by adaptation through returning towards the prestimulus level. The magnitude of this response varies for different parameters. In Using the time series of the responses above, the response magnitude can be quantified by gain g, which is defined as a ratio between the fractional changes in the input signal S and in output signal A, that is, where ΔA is the change in A from the steady state value A Ã before external stimulus to the maximum response A max (see Fig 2A), when the signal intensity is changed from S to S + ΔS. For this quantification, we used A Ã and A max obtained by the differential equation without stochastic effects. Noise can be quantified by the variance s 2 a of time series of A at the steady state using the Gillespie algorithm. The variance is usually proportional to the average, we consider the ratio between the variance and the average, which is called Fano factor f defined as We note here that the average and the variance are measured for the number of molecules instead of concentrations. The exhibited response for a given parameter can be classified into three groups using (A Ã , A max , A ÃÃ ), where A ÃÃ is the steady state after the stimulus as shown in Fig 1A. The first group shows perfect adaptation with jA max − A Ã j > 0 and A Ã % A ÃÃ . The second group shows adaptive response which is not perfect, that is jA max − A Ã j > 0 and A ÃÃ 6 ¼ A Ã . The last group does not show adaptive response with A max % A ÃÃ . The tendency of network to return to its prestimulus state after a response is quantified by the adaptation tendency α, given by Here, we classify the response into the first two groups with adaptive tendency when α ! 0.2, whereas the response is non adaptive when α < 0.2. The precision of perfect adaptation is quantified by the adaptation error , given by Here, we consider that the system shows the perfect adaptation when < 0.01.
To study the relation between response and noise, we quantified the gain and Fano factor over 10 6 parameters for both iFFL and nFBL. Here, we focus on the parameters that show only adaptive responses (the first and second groups), and exclude the non adaptive parameters for which (α < 0.2). These networks have eight different kinetic parameters, as introduced in the last section. Each of them were sampled from a uniform distribution on a logarithmic scale from the range [10 −3 ,10 3 ] using random sampling(see Methods for detail). We also excluded those parameters from our analysis for which either the steady state was too small or too large, ,000 is the system volume. We also removed those parameters for which the gain was too small, i.e. g <.001. For iFFL, 8.3% parameters among 10 6 randomly chosen parameters showed adaptive response, among which 41% showed perfect adaptation and 59% were imperfect. For nFBL, 12.8% showed adaptive response, among which 5.3% showed perfect adaptation while 94.7% were imperfect. Among the parameters with adaptive responses, 56.64% showed non oscillatory response, and 43.35% showed damped oscillation. Among the parameters with perfect adaptation, 43.56% was non-oscillatory while 56.44% showed damped oscillation.
The joint histograms of the gain g and Fano factor f is shown for the parameters with adaptive responses and with perfect adaptation for iFFL and for nFBL in  nFBL, we observe that the gain seem to be not larger than the variance, given by indicating that the gain is limited by noise. We note that the average and the variance on the right hand side are considered with respect to the number of molecules rather than concentration. When these values are considered in the unit of concentration, an appropriate pre-factor is necessary to convert the unit of right hand side to a dimensionless number. For iFFL shown in Fig 3(A), the parameters with adaptive responses are distributed in a wide range of both gain and noise strengths from 10 −3 to 10 1 satisfying inequality Eq (15). The same tendency is seen even for the parameters with perfect adaptation (Fig 3(B)). For nFBL shown in Fig 3(C), the parameters with adaptive responses are distributed in a wide range with inequality Eq (15). In contrast to iFFL, however, they are rarely distributed in the range with smaller values of both gain and noise. Furthermore, when the adaptation is perfect, as shown in Fig 3(D), the parameters lie relatively closer to the line which satisfies equality in Eq (15). Such a tendency is even more pronounced when we focus on the parameters with only non-oscillatory responses as shown in Fig 3(

E) and 3(F).
To see that with the perfect adaptation condition the points in nFBL ( Fig 3F) are distributed closer to the line of equality than in iFFL (Fig 3B), we consider the ratio between the gain and Fano factor, g/f. In Fig 4, using the samples in Fig 3A and 3E, we plot the average of the ratio g/f calculated over the samples for which the adaptation error is less than a given threshold. This threshold decreases from 10 − 1 to 10 − 3, implying the selected points show increasing more perfect adaptation. In the plot, the average of ratio g/f of nFBL is larger than that of iFFL for any value of , indicating that nFBL is closer to the equality line than iFFL. Furthermore, as the perfect adaptation becomes more strict as the adaptation error decreases, the average value of ratio g/f for nFBL increases approaching to the equality line, whereas that for iFFL decreases. This indicates that in nFBL the perfect adaptation condition takes the gain closer to the limit of Fano factor.

The adaptive response is limited by noise
Here, we derive the gain-noise inequality given by Eq 15 analytically within a linear range. We first introduce the linear chemical Langevin equation. The perfect adaptation condition is given in terms of the linear coefficients. Then, the gain and noise are derived and compared.
We consider small deviations in A and B at the steady state values, denoted by a and b, respectively. By performing the linearization for the chemical Langevin equation Eq (35) shown in Appendix A and applying the linear noise approximation [38] at the steady state, we obtain the linearized Langevin equation given by with the linear coefficient or regression matrix, where γ a and γ b are the coupling with the input signal S, and ξ a and ξ b are the Gaussian white noise with hξ a i = hξ b i = 0 and hξ i (t)ξ j (t 0 )i = D ij δ(t − t 0 ), (i, j = a, b). Here, K, γ a , γ b and D = {D ij } are can be written with the reaction constants, which are explicitly given in Appendix B. Using the eigenvalues of K, λ 1 > 0 and λ 2 > 0, and the corresponding eigenvectors, t (1, p 1 ) and t (p 2 , 1), K can be rewritten as K = PΛP −1 with where χ = λ 2 /λ 1 . The regression matrix K for iFFL is given with p 1 = 0, while the nFBL is given with γ b = 0. The perfect adaptation condition that a vanishes (a = 0) at steady state is given as In iFFL, with this perfect adaptation condition, the regression matrix reduces to Because the eigenvalues λ 1 and λ 2 are always positive real number, the ratio χ for iFFL can take any positive value. In the case of nFBL, with the perfect adaptation condition, the regression matrix is given by Here, the eigenvalues can be real or complex. For complex eigenvalues, damping oscillations takes place in the response. For real eigenvalues and non-oscillatory response, the ratio χ is constrained to be always not larger than unity, i.e., χ 1. Next, we derive the gain g using the linearized equation Eq (16) and condition Eq (19), which is given by The variance of a, s 2 a , is obtained by using the generalized fluctuation-response relation, Kσ 2 + where σ 2 is the variance-covariance matrix, K is the regression matrix, and D is the noise strength. When D ab = 0, the variance of a is calculated as Now, we compare the gain and noise by using Eqs (22) and (23). By comparing the gain and noise, we find that the sufficient condition for the inequality Eq (15) is given by As the noise strengths D aa for both iFFL and nFBL can be written as the above condition Eq (24) is always satisfied in both network motifs. Therefore, the gain is always less than the intrinsic noise, and it is a direct consequence of the relation between noise strength and input coupling constant for the adapting variable. This relation can be also be extended to the case with D ab 6 ¼ 0.

Equality conditions restricts gain around unity in iFFL
We next provide insight into how the equality between gain and noise in Eq (15) are achieved for iFFL. With the condition for iFFL, p 1 = 0, and Eq (25), the difference between gain and noise is given by The second term on the right hand side indicates the contribution of noise generated in reactions of A. Because it does not pass through reactions in B, the second term is independent of the ratio χ that includes the property of reaction in B. We also notice that since the second term is independent of B, it does not depend on whether the adaptation is perfect or not. In Eq (26), as χ decreases to zero, the second and third terms are canceled with each other, while the first term vanishes. As a result, the deviation between gain and noise vanishes, and the equality in Eq (15) holds. In contrast, as χ increases, the first term increases, while the deviation between the second and third terms increases.
We verified this dependence on χ numerically, as shown in Fig 5(A) and 5(B). As χ decreases, both gain and noise approach to unity, reaching the equality. In contrast, as χ increases, the gain systematically decreases as w w 1Àw according to Eq (22), while the Fano factor f is maintained around unity or larger than unity for many parameter sets (Fig 5(A)) independent of χ, because the noise contribution from the reactions of A is independent of χ. As a result, the deviation between noise and gain increases with a larger value in Eq (26).
For iFFL, the number A Ã can take any values up to N A without any constraint by the perfect adaptation condition. When A Ã is large and close to N A , the noise (the first two terms in Eq (26)) and the gain (the third term) decrease as A approaches to N A , because the response coefficient g a ¼ G A a ðN A À A Ã Þ decreases as A Ã approaches to N A . This dependence of gain and noise on A Ã is also found numerically as shown in the Fig 6(A) and 6(B). Thus, as A Ã approaches to N A , the gain approaches to the noise, though the gain is small.
Another way to achieve the equality in Eq (15) is to decrease χ. When A Ã /N A ( 1, the ratio χ can be approximated as For this condition, the coefficient of gain in Eq (22) and the second term in Eq (26) are approximately given by which can increase up to unity as the first term in the denominator decreases. By increasing G A d in Eqs (27) and (28), χ can decrease systematically approaching the equality condition in Eq

Perfect adaptation decreases the noise in nFBL
For nFBL with non-oscillatory responses, we consider the conditions under which the gain approaches to the equality in Eq (15). With the condition γ b = 0 for nFBL, and Eq (25), the difference between gain and noise is given by In contrast to iFFL, the second term on the right hand side is a function of χ, because the noise contribution from reactions of A is affected by reactions of B. Because χ 1 is always satisfied as mentioned above, the second and third terms are easier to be close to each other than iFFL. In fact, as shown in Fig 5(C) and 5(D), the gain and noise exhibit similar tendency for any value of χ. Thus, the equality in Eq (15) is achieved easier than iFFL. When A is sufficiently smaller than the total number N A , the coefficient in the second and third terms is approximately given by which can be larger than unity. This indicates that even for sufficiently large χ approaching to zero on the right hand side in Eq (29), a higher gain is possible. For non oscillatory parameters and for adaptive response, the dependence of gain and noise on the steady state is shown in Fig  6(C) and 6(D). In contrast to iFFL, steady state values close to N A don't occur in nFBL. In Fig 3(E) without perfect adaptation conditions, the points are present far below the equality line. Among all points shown in Fig 3(E), only points that exhibit perfect adaptation are shown in Fig 3(F). As mentioned above, those points with perfect adaptation are distributed closely to the equality line. To study this behavior, we introduce the deviation from the perfect adaptation, ω, given by From Eq (19), the perfect adaptation condition is given by ω = 0. Then, the gain and the noise in A are respectively written as and In Fig 7(A), the gain in Eq (32) (blue line) and the second term in the noise in Eq (33) (green line) are plotted as functions of ω. When ω = 0 and the perfect adaptation condition is satisfied, these two values are close to each other. As ω increases, the gain decreases while the noise increases. Therefore, in nFBL the perfect adaptation condition can contribute to increase the gain and decrease the noise at the same time. In Fig 7(B) and 7(C), we show the effect of adaptation error on the noise and gain for nFBL. As the adaptation error decreases, the gain and noise both approach same value, thus leading to high signal to noise ratio in the limit of perfect adaptation.

Discussion
In this paper, we have studied the relation between response and noise in the network motifs that show perfect adaptation. By numerical simulation, we found that the gain measured at the peak response is not larger than the intrinsic noise, measured as the Fano factor. For the case with perfect adaptation, we showed analytically that the inequality in Eq (15) holds in the linear regime. In iFFL, modulating the gain towards its limiting value of Fano factor constraints both of them around a range of unity. Thus, for iFFL, parameters showing perfect adaptation, higher gain and unity gain to noise ratio is difficult to achieve. In contrast, for nFBL, while imperfect adaptation itself shows higher gain to noise ratio, perfect adaptation narrows down the range of gain around its limit of Fano factor. We also note that in nFBL the parameters that show perfect adaptation property also exhibit higher gains. In contrast, in iFFL the parameters with perfect adaptation are distributed in a range from low values to around unity, though the frequency of perfect adaptation among randomly chosen parameters are higher than that in nFBL.
We can have some insight about the differences between iFFL and nFBL for perfect and imperfect adaptation by looking at the analytical terms of the response function and the correlation function. Between the correlation function C(t) and the response function h(t) = exp(−K t), the equation dCðtÞ dt ¼ ÀhðtÞ Dþa 2 holds, where α characterizes the degree of non-equilibrium defined as Kσ 2 − σ 2 t K − α = 0 [39]. By integrating this equation with respect to time t, and introducing a step response function R, which is given by RðtÞ ¼ R t 0 hðsÞds, we have at time t Ã = log χ/λ 1 (χ − 1), which is the time when the response in A reaches the maximum value. The response in A to a step increase in S is given by the first row of R(t Ã ) Á γ, i.e., 1Þ, the gain g can be calculated from the first term on the left hand side in Eq (34) by applying t (1,1) to Eq (34) from the right. When the ratio of eigenvalues, χ, is sufficiently small, the (1,2) component in R vanishes as evident from the regression matrix K in Eq (20). Thus, only (1,1) component in R(t Ã ) contributes to the gain g. We also note that the diagonal components in α is always zero by the definition, indicating that the second term does not contribute to the gain g. On the right hand side in Eq (34), C(t Ã ) vanishes as χ approaches to zero. Consequently, we see that the equality holds in Eq (15), when χ is sufficiently small for iFFL. This does not depend on whether the adaptation is perfect. In a similar way, for nFBL, considering γ ¼ 1 2S D Á t ð1; 0Þ, the (1,1) components in Eq (34) are enough to calculate the gain g. As seen in Fig 5(C) and 5(D), even when χ is not very small, the gain and Fano factor are close with each other. This indicates that (1,1) component in the first terms in Eq (34) are close to each other while the (1,1) components in second terms are almost canceled at t = t Ã . The condition that this occurs in Eq (34) and the gain approaches to the limit of Fano factor remains to be a future problem. We have studied the noise generated by the adaption reaction (intrinsic noise). For the noise in an input stimulus (extrinsic noise), filtering properties of adaptive response has been studied [33,34]. Fast fluctuations in the input are averaged out by the low pass filtering effect, while slow fluctuations are attenuated by adaptation. Thus, adaptive reactions behaves as a bandpass filter. Because of this filtering property, the magnitude of the extrinsic noise is dominated by the intrinsic noise.
In the bacterial chemotaxis, the adaptive responses is well studied. The mechanism of adaptation is considered to be negative feedback loop, in which the receptor activity shows response to ligand-binding, which is then attenuated through a modulation of methylation level by methylesterase CheB, the activity of which depends on the receptor activity. Examples of iFFL has been reported for eukaryotic cells. In a cultured mammalian cell, epidermal growth factor (EGF) induces transient activation of extracellular-signal-regulated kinase (ERK). To explain this adaptive response, an iFFL mechanism upstream of ERK has been proposed, in which phosphorylated EGF receptors activate both Ras GEF and GAP with different time constants [40]. Similarly, for the adaptive response of chemotaxis signaling pathway in Dictyostelium cell, iFFL of Ras and its GAP and GEF has been proposed [20].
Previously, it has been reported that iFFL is more robust than nFBL with respect to the perfect adaptation [18]. Our analysis also confirms that the frequency of perfect adaptation is higher in iFFL than nFBL. However, our result indicates that the frequency to have higher gain is larger for nFBL. Moreover, in nFBL the gain is close to the upper limit of the Fano factor, while in iFFL the parameters are broadly distributed away from the upper limit. Also, our analysis shows that the perfect adaptation condition itself can contribute to increase the gain up to the limit of Fano factor for nFBL, which could be of great advantage. In particular, in small cells with smaller number of molecules, such as in bacteria, nFBL should work more robustly against molecular noise. We hope our results would help to understand the mechanism and occurrence of adaptive networks in biological systems.

Sampled parameter space
To quantify the gain and noise for adaptive networks, the Michaelis-Menten parameters for the iFFL and nFBL were sampled randomly, where each of the network needs 8 parameters. The Michaelis-Menten parameters V m and K m were sampled from a uniform distribution on a logarithmic scale, from the range (.001 − 1000), where V m is the dimensionless maximal velocity, and K m is the Michaelis constant, in the units of moles.
For calculating the deterministic gain using the differential equations, the simulations were performed with the concentration of each species. For computing the Fano factor, and for observing the stochastic response, the simulations were performed with the actual number of molecules with a uniform system size or volume V = 10,000 for each parameter. The Michaelis constant in terms of number of molecules, K M is computed by K M = VK m . The maximum velocity in terms of number of molecules is the same as in the concentration space, that is V M = V m .

Quantifying the response
For a given parameter, first the steady state of the system was obtained. For iFFL, the steady state number can be obtained analytically in a closed form, while for nFBL, the steady state number was obtained using numerically solving the nonlinear differential equation using the Runge-Kutta method. If the steady state number was too small or large, that is (A/V < 0.001, B/V < 0.001), or ((N A − A)/V < 0.001, (N B − B)/V < 0.001), where V = 10,000 then we call it an ill-behaved parameter and discard it.
For parameters which were not ill-behaved, a small perturbation ΔS was applied to it in the steady state, and the deterministic part of the response was obtained using the differential equation. For solving the differential equation, we used the fourth order Runge-kutta method with a default step size of h = .001 and the global error tolerance, that is the difference between successive values (y n+1 − y n ), to be less than h/50, where y i is the numerical solution of the differential equation at the i th step. If the solver failed to reach the steady state with this step size, we successively used a smaller step size of h = 10 −4 or 10 −5 . A response was considered to be perfectly adaptive if the adaptation error < 0.01. The response was calculated in the concentration space and then converted into numbers using the system volume information.
From the obtained time series of the response, the quantities A Ã , the initial steady state, A max , the maximum number, and A ÃÃ , the final number of the node A was obtained. The maximum change in the output is ΔA = (A max − A Ã ). Using these quantities, we compute the following: A max ÀA Ã , which quantifies the tendency of the network to get back to its pre stimulus state from its maximum state.

DS=S
, which quantifies the deviation from perfect adaptation.

Gain g ¼ jDAj=A Ã
DS=S , which gives a scale free estimate of the response for a perturbation of ΔS in the input signal.
For nFBL, the oscillatory and non oscillatory parameters were differentiated by observing the response time series, where the total number of times the slope of the response changed its sign was calculated. If the total sign change happened at least three times, then we considered the response as oscillatory response.

Quantifying the Noise
To quantify the Fano factor, the system was simulated by using the stochastic simulation algorithm described by Gillespie [37]. All the simulations were performed with a system volume of V = 10,000, and reaction was ran for a minimum of 100 time units and 10 5 reaction events. The steady state of the system itself was chosen as the initial condition. The resultant time series was used to calculate the variance σ 2 of the output. From this output, the Fano factor can be estimated by f ¼ s 2 A Ã as defined in Eq (12) where A Ã is the steady state value.  The total number of molecules was N A = N B = 200,000, S = 10,000, and a large perturbation of ΔS = 10000 was provided at the steady state, where the volume is V = 10,000. Fig 3 and following figures. For obtaining the joint histogram of Fano factor and gain in Fig 3, a total of 10 6 random parameters were sampled for each of iFFL and nFBL, and both Fano factor and gain was calculated for them. In each of the above case, the total number of molecules was N A = N B = 200,000, S = 10,000, and a 1% perturbation of ΔS = 100 was provided at the steady state, where the volume is V = 10,000. The same parameters were also used for the Figs (4-7) as well.

Making the log histogram plot
To represent the distribution of Fano factor and gain or for other pair of variables, we used log of the histogram plot, as there was large variation in the parameter density in low density and high density regions. To obtain the log histogram plot, first the relevant pair of variables were converted on a log base(base 10). Subsequently, the histogram was obtained by counting the number of parameters falling in each square bin of chosen size. This count was again taken on a log 10 base, and the resultant matrix was plotted on a pseudo color plot using the function pcolor of matlab, with maximum density being assigned the red color, and the minimum density as blue color. The binwidths were typically chosen to be of size.1 for both x and y coordinates.

Appendix A
Here, we show the chemical Langevin equations for iFFL and nFBL, which are given by and respectively, where G A a , G A d , G B a and G B d are the reaction rates as defined in Eq (3), and ξ a and ξ b are the Gaussian white noise with hξ a i = hξ b i = 0 and hξ i (t)ξ j (t 0 )i = D ij δ(t − t 0 ), (i, j = a, b). The noise strength D = {D ij } can be derived from the property of Poisson process in individual reactions [38,41]. Thus, D ij are given by D aa ¼ 2G A d A Á B, D bb ¼ 2G B d B, and D ab = D ba = 0.

Appendix B
Here, we show the explicit forms of the linear coefficients K, γ and D in the linearized Langevin equations Eq (16). For iFFL, the liner coefficients K in Eq 35 are given by The coupling coefficients γ with the input S also are obtained in the above step, which are given by g a For the nFBL, using the Eq (35), K is given as The coupling coefficient γ with the input S is obtained and given by: The noise strength D can be obtained by the linear noise approximation [38].