Feedback Regulation and Its Efficiency in Biochemical Networks

Intracellular biochemical networks fluctuate dynamically due to various internal and external sources of fluctuation. Dissecting the fluctuation into biologically relevant components is important for understanding how a cell controls and harnesses noise and how information is transferred over apparently noisy intracellular networks. While substantial theoretical and experimental advancement on the decomposition of fluctuation was achieved for feedforward networks without any loop, we still lack a theoretical basis that can consistently extend such advancement to feedback networks. The main obstacle that hampers is the circulative propagation of fluctuation by feedback loops. In order to define the relevant quantity for the impact of feedback loops for fluctuation, disentanglement of the causally interlocked influences between the components is required. In addition, we also lack an approach that enables us to infer non-perturbatively the influence of the feedback to fluctuation in the same way as the dual reporter system does in the feedforward networks. In this work, we address these problems by extending the work on the fluctuation decomposition and the dual reporter system. For a single-loop feedback network with two components, we define feedback loop gain as the feedback efficiency that is consistent with the fluctuation decomposition for feedforward networks. Then, we clarify the relation of the feedback efficiency with the fluctuation propagation in an open-looped FF network. Finally, by extending the dual reporter system, we propose a conjugate feedback and feedforward system for estimating the feedback efficiency non-perturbatively only from the statistics of the system.


Dissecting fluctuation in biochemical networks
Biochemical molecules in a cell fluctuate dynamically because of the stochastic nature of intracellular reactions, fluctuation of the environment, and the spontaneous dynamics of intracellular networks [18,32,40].Some part of the fluctuation is noise that impairs or disturbs the robust operation of the intracellular networks.The other part, however, conveys information on complex dynamics of various factors inside and outside of the cell [20,31].Dissecting fluctuation into distinct components with different biological roles and meanings is crucial for understanding the mechanisms how a cell controls and harnesses the noise and how information is transferred over apparently noisy intracellular network [5,9].
The intracellular fluctuation is generated from the intracellular reactions or comes from the environment, and then propagates within the intracellular networks.The fluctuation of individual molecular species within the networks is therefore a consequence of the propagated fluctuation from the different sources.Decomposition of the fluctuation into the contributions from the sources is the indispensable step for understanding their biological roles and relevance.When fluctuation propagates from one component to another unidirectionally without circulation, the fluctuation of the downstream can be decomposed into two contributions.One is the intrinsic part that originates within the pathway between the components.The other is the extrinsic part that propagates from the upstream component.Such decomposition can easily be extended for the network with cascading or branching structures in which no feedback exists.This fact drove the intensive anatomical analysis of the intracellular fluctuation in the last decade.

Decomposition of sources of fluctuation
In order to dissect fluctuation into different components, two major strategies have been developed.One is to use the dependency of each component on different kinetic parameters in the network.By employing theoretical predictions on such dependency, we can estimate the relative contributions of different components from single-cell experiments with perturbations.Possible decompositions of the fluctuation were investigated theoretically for various networks such as single gene expression [3,26,27], signal transduction pathways [41,47], and cascading reactions [45].Some of them were experimentally tested [3,26,28].
The other is the dual reporter system in which we simultaneously measure a target molecule with its replica obtained by synthetically duplicating the target.From the statistics of the target and the replica, i.e, mean, variance, and covariance, we can discriminate the intrinsic and extrinsic contributions to the fluctuation because the former is independent between the target and the replica whereas the latter is common to them.The idea of this strategy was proposed and developed in [27,43], and verified experimentally for different species [11,23,33].Its applicability and generality were further extended [4,6,16,17,35].Now these strategies play the fundamental role to design single-cell experiments and to derive information on the anatomy of fluctuation from the experimental observations [8,15,29,34,37,44].

Feedback regulation and its efficiency
Even with the theoretical and the experimental advancement in decomposing fluctuation, most of works focused on the feedforward (FF) networks in which no feedback and circulation exist.As commonly known in the control theory [7,39], feedback (FB) loops substantially affect fluctuation of a network by either suppressing or amplifying it.Actually, the suppression of fluctuation in a single gene expression with a FB loop was experimentally tested in [2] earlier than the decomposition of fluctuation.While the qualitative and the quantitative impacts of the FB loops were investigated both theoretically and experimentally [1,21,24,25,35,42,46], we still lack a theoretical basis that can consistently integrate such knowledge with that on the fluctuation decomposition developed for the FF networks.
The main problem that hampers the integration is the circulation of fluctuation in the FB network.Because fluctuation generated at a molecular component propagates the network back to itself, we need to disentangle the causally interlocked influence between the components to define the relevant quantity for the impact of the FB loops.From the experimental point of view, in addition, quantification of the impact of FB loops by perturbative experiments is not perfectly reliable because artificial blocking of the FB loops inevitably accompanies the change not only in fluctuation but also in the average level of the molecular components involved in the loops.It is quite demanding and almost impossible for most cases to inhibit the loops by keeping the average level unchanged.We still lack an approach that enables us to infer the influence of the FB non-perturbatively as the dual reporter system does.

Outline of this work
In this work, we resolve these problems by extending the work on the fluctuation decomposition [27,41] and the dual reporter system [11,43].By using a single-loop FB network with two components and its linear noise approximation (LNA) [10,19,46], we first provide a definition of the FB loop gain as FB efficiency that is consistent with the fluctuation decomposition in [27,41].
Then, we clarify the relation of the FB efficiency with the fluctuation propagation in a corresponding open-looped FF network.Finally, by extending the dual reporter system, we propose a conjugate FB and FF system for estimating the feedback efficiency only from the statistics of the system nonperturbatively.We also give a fluctuation relation among the statistics that may be used to check the validity of the LNA for a given network.
The rest of this paper is organized as follows.In Sect.2, we review the decomposition of fluctuation for a simple FF system derived in [27,41] by using the LNA.In Sect.3, we extend the result shown in Sect. 2 to a FB network by deriving a decomposition of the fluctuation with feedback.Using this decomposition, we define the FB loop gain that is relevant for quantifying the impact of the FB to the fluctuation.In Sect.4, we give a quantitative relation of the loop gain in the FB network with the fluctuation propagation in a corresponding open-looped FF network.In Sect.5, we propose a conjugate FF and FB network as a natural extension of the dual reporter systems used mainly for the FF networks.We clarify that the loop gain can be estimated only from the statistics, i.e., mean, variances, and covariances, of the conjugate network.We also show that a fluctuation relation holds among the statistics, which generalizes the relation used in the dual reporter sys-tem.In Sect.6, we discuss a link of the conjugate network with the directed information, and give future directions of our work.

Biochemical Network
In this section, we summarize the result for the decomposition of fluctuation obtained in [27,41] by using the LNA, and also its relation with the dual reporter system that was employed in [11,23,33] to quantify the intrinsic and extrinsic contributions from the experimental measurements.

Stochastic chemical reaction and its linear noise approximation
Let us consider a chemical reaction network consisting of N different molecular species and M different reactions.We assume that the stochastic dynamics of the network is modeled by the following chemical master equation: where n = (n 1 , . . ., n N ) T ∈ N N ≥0 is the numbers of the molecular species, P(t, n) is the probability that the number of molecular species is n at t, and a k (n) ∈ R M ≥0 and s k are the propensity function and the stoichiometric vector of the kth reaction, respectively [12,14,19].The propensity function characterizes the probability of occurrence of the kth reaction when the number of the molecular species is n, and the stoichiometric vector defines the change in the number of the molecular species when the kth reaction occurs.
In general, it is almost impossible to directly solve Eq. (1) both analytically and numerically because it is a high-dimensional or infinite-dimensional differential equation.To obtain insights for the dynamics of the reaction network, several approximations have been introduced [12,19].Among others, the first-order approximation is the deterministic reaction equation that is described for the given propensity functions and the stoichiometric vector as where Ω is the system size, u ∈ R n ≥0 is the concentration of n as u = n/Ω, a(n) := (a 1 (n)), . . ., a M (n)) T , and S * ,k := s k is the stoichiometric matrix.Equation (2) was successfully applied for chemical reaction networks with large system size where the fluctuation of the concentration of the molecular species can be neglected.When the system size is not sufficiently large, however, Eq. ( 2) is not relevant for analyzing the fluctuation of the network.
The LNA is a kind of the second-order approximation of the Eq. ( 1) that characterizes the fluctuation around a fixed point, n, of Eq. ( 2) that satisfies Sa(n) = 0.The stationary fluctuation of the network around n is then obtained by solving the following Lyapunov equation [10,19,46] where and Σ is the covariance matrix of n.When the propensity function a(n) is affine with respect to n, the dynamics of Ωu(t) determined by Eq. ( 2) is identical to that of the first cumulant of n, i.e., the average of n, as where n := n nP(t, n).In addition, the second cumulant, i.e., the covariance matrix, follows the Lyapunov equation as Therefore, if the propensity function a(n) is affine, the stationary fluctuation of n is exactly described by Eq. ( 3).For a non-affine a(n), Eq. ( 5) and Eq. ( 6) can also be regarded as an approximation of the full cumulant equations by the cumulant closure [12] under which we ignore the influence of the second and the higher order cumulants to Eq. ( 5), and that of the third and the higher order ones to Eq. ( 6) 1 .Even though the propensity function a(n) is not affine, Eq. ( 3) (or Eq. ( 6)) can produce a good approximation of the fluctuation, provided that the fixed point, n, is a good approximation of the average, n , and that the local dynamics around the fixed point is approximated enough by its linearization.In addition, compared with other approximations, the LNA enables us to obtain an analytic representation of the fluctuation because Eq. ( 3) is a linear algebraic equation with respect to Σ. Owing to this property, the LNA and its variations played the crucial role to reveal the analytic representation of the fluctuation decomposition in biochemical reaction networks [27,46].As in these previous works, we employ the LNA to obtain an analytic representation for the feedback efficiency.

Decomposition of fluctuation
Starting from the LNA, Paulsson derived an analytic result on how noise is determined in a FF network with two components (Fig. 1 (A)) [27].Here, we briefly summarize the result derived in [27].Let n 1 = x and n 2 = y for notational simplicity, and consider the FF reaction network (Fig. 1 (A)) with the following propensity function and stoichiometric matrix, Because a ± x depends only on x, x regulates y unidirectionally.Then, for a fixed point (x, ȳ) of Eq. ( 2), that satisfies K and D in Eq. ( 3) becomes where āx := a + x (x) = a − x (x), āy = a + y (x, ȳ) = a − y (x, ȳ), and H i,j := for i, j ∈ {x, y}.d x and d y are the minus of the diagonal terms of K and represent the effective degradation rates of x and y. k yx is the off-diagonal term of K that represents the interaction from x to y.
H ij is the susceptibility of the ith component to the perturbation of the jth one.τ x := x/ā x and τ y := ȳ/ā y are the effective life-time of x and y, respectively.Except this section, we mainly use ds and ks as the representation of the parameters rather than Hs and τ s introduced in [27]2 .
By solving Eq. ( 3) analytically, the following fluctuation-dissipation relation was derived in [27] as x.The term (III) is further decomposed into (i), (ii), and (iii).(i) is the fluctuation of x, and therefore, identical to (I).(ii) and (iii) determine the efficiency of the propagation of the fluctuation from x to y, which correspond to the time-averaging and sensitivity of the pathway from x to y, respectively.This representation captures the important difference of the intrinsic and the extrinsic fluctation such that the intrinsic one, the term (II), can be reduced by increasing the average of ȳ whereas the extrinsic one, the term (III), cannot.
While Eq. ( 10) provides an useful interpretation on how the fluctuation propagates in the FF network, it is not appropriate for the extension to the FB network because the contribution of x to the fluctuation of y is described by the CV of x as the term (i).Because the fluctuation of x and y depend mutually if we have a FB between x and y, we need to characterize the fluctuation of y without directly using the fluctuation of x.To this end, we adopt the variances and covariances as the measure of the fluctuation and use the following decomposition of the fluctuation for the FF network: x x where we define G xx , G yy , and G yx as The terms (I), (II), and (III) in Eq. ( 11) correspond to those in Eq. ( 10).
The interpretation of the terms within (I), (II), and (III) is, however, different.In Eq. ( 11 x to y.If we use the notation in Eq. ( 10), G yx is described as The decomposition of the fluctuation of y into (II) and (III) is consistent between Eq. ( 10) and Eq. ( 11) while the further decompositions within (II) and (III) are different.It was recently revealed in [4,6,16,17] that the decomposition into (II) and (III) is linked to the variance decomposition formula in statistics as where X (t) := {x(τ ); τ ∈ [0, t]} is the history of x(t), and E and V are the expectation and the variance, respectively 4 .Because the variance decomposition formula holds generally, Eq. ( 14) is more fundamental than Eq.(10) and Eq. ( 11) as the decomposition of fluctuation.This decomposition was further analyzed in [4,6,16].

Dual reporter system
The decomposition, Eq. ( 10) or Eq. ( 11), guides us how to evaluate the intrinsic and the extrinsic fluctuation in y experimentally.When we can externally control the mean of y without affecting the term (III), we can estimate the relative contributions of (II) and (III) as the intrinsic and the extrinsic fluctuation by plotting σ 2 y as a function of the mean of y.When x and y correspond to a mRNA and a protein in the single gene expression, the translation rate works as such a control parameter [22].This approach was intensively employed to estimate the efficiency of the fluctuation propagation in various intracellular networks [3,26,28].
Another way to quantify the intrinsic and the extrinsic fluctuation is the dual reporter system adopted in [11,23,33] whose network structure is shown in Fig. 1 (B) 5 .In the dual reporter system, a replica of y is attached to the downstream of x as in Fig. 1 (B) where y denotes the molecular species of the replica.The replica, y , must have the same kinetics as y, and must be measured simultaneously with y.If y is a protein whose expression is regulated by another protein, x, as in Fig. 1 (C), then y can be synthetically constructed by duplicating the gene of y and attaching fluorescent probes with different colors to y and y as in Fig. 1 (D) [11].Under the LNA, the covariance between y and y can be described as Thus, by using only the statistics of the dual reporter system, the intrinsic and the extrinsic components in y can be estimated as where it is unnecessary to control any kinetic parameters externally.This approach is generalized as the conjugate reporter system together with the more general decomposition, Eq. ( 14), in [4,6,16].
3 Feedback Loop Gain in a Small Biochemical Network While we had substantial progress in the decomposition of the fluctuation and its experimental measurement for the FF networks in the last decade, its extension to FB networks is yet to be achieved.In this section, we extend the decomposition of the fluctuation in the simple FF network (Eq.( 11)) to the corresponding FB one depicted in Fig. 2 (A) .As in Eq. ( 9), K and D in the Lyapunov equation (Eq.( 3)) for the FB network can be described as By defining the path gain from y back to x as we can derive a decomposition of the fluctuation of x and y as where we define If the FB from y to x does not exist, i.e., k xy = 0, then L x = L y = G xy = 0 and Eq. ( 19) is reduced to Eq. ( 11).Thus, L x and L y account for the effect of the FB.L x and L y are denoted as the FB loop gains in this work.Eq. (19) clearly demonstrates that the representation in Eq. ( 10) does not work with the FB because σ 2 y cannot be described with σ 2 x any longer.In contrast, we can interpret the terms in Eq. ( 19) consistently with those in Eq. ( 11) because all the terms, (I), (II) and (III), are unchanged in Eq. (19).The new term, (IV ), in the expression for σ 2 x appears to account for the propagation of the fluctuation generated by the birth and death events of y back to x.

Equations (19) and (20) indicate how the FB affects the fluctuation of
x and y.First, the fluctuation is suppressed when L x and L y are negative whereas it is amplified when they are positive.When d x and d y are positive6 , the sign of L x and L y are determined by the sign of k xy k yx .Thus, the FB loop is negative when x regulates y positively and y does x negatively or vise versa.This is consistent with the normal definition of the sign of a FB loop.
Second, the efficiency of the FB depends on the source of the fluctuation.
For example, when L x L y , e.g., the time-scale of x is much faster than y, then Eq. ( 19) can be approximated as Thus, the FB does not work for the term (I) that is the part of fluctuation of x whose origin is the birth and death events of x itself.This result reflects the fact that the slow FB from x to itself via y cannot affect the fast component

Relation between Fluctuation Propagation and Feedback Gain
As shown in Sect.3, L x and L y are quantitatively related to the efficiency of the FB.Because L x L y = G xy G yx holds, the loop gains are also linked to the propagation of the fluctuation from x to y and from y back to x.However, the meaning of the individual gains, L x and L y , is still ambiguous.
To clarify the relation between the loop gains and the propagation of fluctuation, we consider a three-component FF network shown in Fig. 3 (A) and (B) that are obtained by opening the FB network in Fig. 2 (A).In the network shown in Fig. 3 (A), x is regulated not by y but by its replica denoted as y .We assume that y is not driven by x, and thereby, y → x → y forms a FF network.K and D in Eq. ( 3) for this network become Because y is the replica of y, we assume that all the kinetic parameters of y are equal to those of y as k y x = k yx , d y = d y , and āy = āy .
By solving Eq. ( 3), we can obtain the following decomposition of the fluctuation of y as where The term (IV ) for x is similar to the propagation of the fluctuation from y to x in the FB network, but it represents the propagation of the fluctuation from the replica y in this opened FF network.The new term (V ) accounts for the propagation of the fluctuation from y down to y.The gain of this propagation has two terms, G yxy and G yy .The first term, G yxy , is the total gain of the fluctuation propagation from y to x and from x to y because G yxy = G yx G xy holds.In order to see the meaning of the second term G yy , we need to rearrange Eq. ( 22) as This representation clarifies that G yy describes the propagation of the fluctuation from y to y that cannot be reflected to the fluctuation of the intermediate component, x.In addition, by solving Eq. ( 3), we can see that the gain G yy is directly related to the covariance between y and y as This implies that we have at least two types of propagation of the fluctuation.One described by G yxy is that the fluctuation of the upstream, i.e., y , is absorbed by the intermediate component, i.e., x, and then the absorbed fluctuation propagates to the downstream, i.e., y.This component does not convey the information of the upstream because that does not affect the covariance between the upstream and the downstream.The other described by G yy is that the fluctuation of the upstream propagates to the downstream without affecting the intermediate component.This fluctuation conveys the information on the upstream to the downstream because it is directly linked to their covariance.The fact that L y is related to the latter indicates that the FB efficiency is directly linked to the information transfer in the opened loop from y to y.By considering another opened loop network where the replica of x is introduced as in Fig. 3 (C), we can obtain the following result: where G xx := L 2 x /2, and we also have By combining Eq. ( 25) and Eq. ( 27), we can estimate the loop gains, L x and L y , by only measuring the averages and covariances of the opened networks as follows: While this strategy sounds eligible theoretically, it has an experimental difficulty in constructing the opened networks.In the opened network, the replica, e.g., y , must be designed so that it is free from the regulation of x by keeping all the other properties and kinetic parameters the same as those of y.For example, if x and y are regulatory proteins and if they regulate each other as transcription factors as Fig. 2 (B), the replica, y , must not be regulated by x but its expression rate must be equal to the average expression rate of y under the regulation of x as in Fig. 3 (B).This requires fine-tuning of the expression rate of y by modifying the DNA sequences relevant for the rate.In order to conduct this tuning, we have to measure several kinetic parameters of the original FB networks that undermines the advantage of the opened network that measurements of the kinetic parameters are unnecessary to estimate L x and L y via Eq.( 28).

Estimation of Feedback Loop Gain by a Conjugate FB & FF Network
As a more promising strategy for the measurement of the loop gains, we propose a conjugate FB and FF network that is an extension of the dual reporter system for the estimation of the intrinsic and the extrinsic components.In the conjugate network, we couple the original FB network with a replica that is opened as in Fig. 4 (A).x and y are the same as the original FB network.
The replica, y , is regulated by x as y is but does not regulate x back.Thus, x and the replica, y , form an FF network.If x and y are regulatory proteins as in Fig. 2 (B), the replica y can be engineered by duplicating the gene y with its promoter site, and by modifying the coding region of the replica so that y looses the affinity for binding to the regulatory region of x as in Fig. 4 (B).This engineering is much easier than that required for designing the opened network.

FF FB FB FF
Next, we show how to use this conjugate network to measure the loop gains.For this network, K and D in Eq. ( 3) become Because the replica y affects neither x nor y, the fluctuation of x and that of y are the same as those of the FB network in Eq. (19).The variance of the replica y can be decomposed as where Rearrangement of this equation leads to In addition, we have the following expression for the covariance between y and y as A similar result can be obtained by replicating x as in Fig. 4 (C) and (D).
By using these relations, we have where F x := σ 2 x / x , F y := σ 2 y / y , F x,x := σ x,x / x , and F y,y := σ y,y / y are the Fano factors of x and y and normalized covariances.This result

Verification of the relations by numerical simulation
We verify Eq. ( 34) and Eq. ( 35) by using numerical simulation.For the simulation, we use a conjugate network of gene regulation in which the replicas of both x and y are involved to measure L x and L y , simultaneously as in Fig. 5(A).For the variable n = (x, y, x , y ) T , the stoichiometric matrix and the propensity function are The simulation is conducted by the Gillespie's next reaction algorithm [13].
First, we test linear negative feedback regulation defined by f (y) := max f 0 [1 − y Ky ], 0 , and g(x) := max g 0 x Kx , 0 under which the LNA holds exactly as long as its trajectory has sufficiently small probability to reach the boundaries of x = 0 or y = 0.In Fig. 5(B), the distributions of p(x, y), p(x, y ), and p(x , y) are plotted for different parameter values of K x .The FB is strong for small K x whereas it is weak for large K x .(L x , L y ) estimated by Eq. ( 34 In order to test how nonlinearity affects the estimation, we also investigated non-linear negative feedback regulation defined by f (y) := f 0 1 1+(y/Ky) ny , and g(x) := g 0 (x/Kx) nx 1+(x/Kx) nx .We change the Hill coefficients, n x and n y , by keeping the fixed point unchanged as in Fig. 6(A).Compared with the linear case, the variability of the estimators is almost similar even though the feedback regulation is nonlinear (Fig. 6(B)).In addition, the estimators show a good agreement with the analytical value except for very large value of |L x | and |L y |.This suggests that Eq. ( 34) works as good estimators when the trajectories of the system are localized sufficiently near the fixed point as in Fig. 6(A All these results indicate that the estimators obtained in Eq. ( 34) can be used estimate the FB efficiency as long as the efficiency is moderate and the trajectories of the system are localized near the fixed point.

Discussion
In this work, we extended the fluctuation decomposition formula obtained for the FF network to the FB network.In this extension, the FB loop gains are naturally derived as the measure to quantify the efficiency of the FB.By considering the opened FF network obtained by opening the loop of the FB network, the relation between the loop gains and the fluctuation propagation in the FF network was clarified.In addition, we proposed the conjugate FB and FF network as a methodology to quantify the loop gains by showing that the loop gains are estimated only from the statistics of the conjugate network.By using numerical simulation, we demonstrated that the loop gains can actually be estimated by the conjugate network while we need more investigation on the bias and variance of the estimators.Furthermore, the fluctuation relation that holds in the conjugate network was also verified.
We think that our work gives a theoretical basis for the conjugate network as a scheme for experimental estimation of the FB loop gains.
As for the further problems on the efficiency of FB and its quantification, the generality of the FB loop gain should be clarified.In the case of the FF network, the fluctuation decomposition proposed initially by using the LNA was generalized as the variance decomposition formula with respect to the conditioning of the history of the upstream fluctuation [4,17].In the case of the FB network, similarly, the decomposition and the FB loop gains were obtained and defined via the LNA.However, its generalization is not straightforward because the contribution of the interlocked components cannot be dissected in the FB network because of the circulative flow of the fluctuation.The previous work on the relation between the fluctuation decomposition and the information transfer for the FF network [4] suggests that the FB loop gains can also be interpreted in terms of the information transfer.As a candidate for such information measure, we illustrate a connection of the conjugate FB and FF network with the directed information and Kramer's causal conditioning.Let us consider the joint probability of the histories of x(t) and y(t), P[X (t), Y(t)].From the definition of the conditional probability, we can decompose this joint probability as where P y||x [Y(t)||X (t)] and P x||y [X (t)||Y(t − 1)] are the Kramer's causal conditional probability [30].If no FB exists from y back to x, then this decom-

XFig. 1 (
Fig.1 (A)The structure of the two-component FF network.Interpretations of this network as single-gene expression[11,22,26,27], two-gene regulation[46], and signal transduction[47] are shown.(B) The structure of the dual reporter system [4, 11, 17, 27, 43].(C) A schematic diagram of the FF network for two-gene regulation.(D)A schematic diagram of the dual reporter network for the two-gene regulation.

. ( 10 )
This representation measures the intensity of the fluctuation by the coefficient of variation (CV) 3 , and describes how the fluctuation generates and propagates within the network.(I) is the intrinsic fluctuation of x that originates from the stochastic birth and death of x. 1/x reflects the Poissonian nature of the stochastic birth and death, and 1/H xx is the effect of auto-regulatory FB.Similarly, (II) is the intrinsic fluctuation of y that originates from the stochastic birth and death of y. (III), on the other hand, corresponds to the extrinsic contributions to the fluctuation of y due to the fluctuation of ), the terms (i) are interpreted as the fluctuation purely generated by the birth and death reactions of x and y by neglecting any contribution of the auto-FBs.Because the simple birth and death of a molecular species without any regulation follow the Poissonian statistics, the intensity of the fluctuation is equal to the means of the species.Thus, x ≈ x and ȳ ≈ y in the terms (i) represent the generation of the fluctuation by the birth and death of x and y, respectively.The fluctuation generated is then amplified or suppressed by the auto-regulatory FBs.G xx and G yy in the terms (ii) account for this influence, and are denoted as auto-FB gains in this work.Finally, G yx in the term (iii) quantifies the efficiency of the propagation of the fluctuation from x to y.We denote G yx as the path gain from

Fig. 2 (
Fig. 2 (A) The structure of the two-component FB network.(B) A schematic diagram of the FB network for two-gene regulation.

of the fluctuation of x 7 .
Finally, when L x and L y satisfy 1 − L x − L y = 1 − k xy k yx /d x d y = 0, the fluctuation of both x and y diverges due to the FB.Since this condition means that the determinant of K becomes 0 and K is the Jacobian matrix of Eq. (2), the fluctuation of x and y diverges due to the destabilization of the fixed point, n, by the FB.

Fig. 3 (
Fig. 3 (A) The structure of the the opened FF network obtained by replicating y in the FB network.(B) A schematic diagram of the network in (A) for gene regulation.(C) The structure of the the opened FF network obtained by replicating x in the FB network.(D) A schematic diagram of the network in (C) for gene regulation.

Fig. 4 (
Fig. 4 (A) The structure of the conjugate FB and FF network obtained by replicating y in the FB network.(B) A schematic diagram of the network in (A) for gene regulation.(C) The structure of the conjugate FB and FF network obtained by replicating x in the FB network.(D) A schematic diagram of the network in (C) for gene regulation.

Fig. 6 (
Fig. 6 (A) The distributions of p(x, y) , p(x, y ), and p(x , y) sampled from the simulation for different parameter values of nx and ny where nx = ny = |n|.The other parameters are f0 = g0 = 300, Kx = Ky = 150, and dx = dy = 1.The blue, red, and green dots are the distributions of p(x, y) , p(x, y ), and p(x , y), respectively.Solid and dashed curves are the nullclines of Eq. (2) defined as dx/dt = 0 and dy/dt = 0 .(B) (Lx, Ly) estimated with the relations, (a), (b), and (c), in Eq. (34) from the numerical simulation.The same parameter values were used as in (A).For each estimation, the means and the covariances required in Eq. (34) were calculated from 9 × 10 4 samplings.For each parameter value, we calculated the estimates 10 times to see their variation.Black markers show the analytically obtained values of (Lx, Ly) for each parameter value.(C) A plot of the estimates of Fx + F x − 2F x,x and Fy + F y − 2F y,y shown in Eq.(35).The parameter values used were the same as those in (B).
Fig. 6(A).For the very large value of |L x | and |L y | where |n| = 2 7/2 , all the estimators are slightly biased towards smaller values, and (b) in Eq. (34) has larger variance than the others.In addition, similarly to the linear case, the estimators require larger sampling when |L x | and |L y | are much less than 1.Even with the nonlinear regulation, the fluctuation relation holds robustly as shown in Fig. 6(C), which is consistent with the good agreement of the estimators of L x and L y with their analytical values as in Fig. 6(B).