Information Theoretical Study of Cross-Talk Mediated Signal Transduction in MAPK Pathways

Abstract

Biochemical networks having similar functional pathways are often correlated due to cross-talk among the homologous proteins in the different networks. Using a stochastic framework, we address the functional significance of the cross-talk between two pathways. A theoretical analysis on generic MAPK pathways reveals cross-talk is responsible for developing coordinated fluctuations between the pathways. The extent of correlation evaluated in terms of the information theoretic measure provides directionality to net information propagation. Stochastic time series suggest that the cross-talk generates synchronisation in a cell. In addition, the cross-interaction develops correlation between two different phosphorylated kinases expressed in each of the cells in a population of genetically identical cells. Depending on the number of inputs and outputs, we identify signal integration and signal bifurcation motif that arise due to inter-pathway connectivity in the composite network. Analysis using partial information decomposition, an extended formalism of multivariate information calculation, also quantifies the net synergy in the information propagation through the branched pathways. Under this formalism, signature of synergy or redundancy is observed due to the architectural difference in the branched pathways.

1. Introduction

The decision making processes at the cellular level are initiated by some specialised signalling networks [1,2]. These networks play a pivotal role in making robust and precise cellular response towards endogenic and exogenic perturbations. In addition, the process of decision making resolves cellular fate as well as survival strategies in diverse peripheral conditions. Although both prokaryotic and eukaryotic cells are comprised of several common signalling networks, few signalling networks are incorporated mostly in the eukaryotes [3].

As an evolutionary outcome, cells have developed an optimal protein-protein (cognate and noncognate) interaction within the signalling pathway to transduce extra-cellular information efficiently [4]. One such signalling network is the mitogen-activated protein kinase (MAPK) pathway that plays the central role to attune with extra-cellular signal in eukaryotic cells [5,6,7,8,9]. Although different MAPK pathways with diverse inputs and outputs belong to a higher living species, they are sometimes interconnected through overlapping sets of signalling components. Depending on the interconnections, MAPK pathways can be classified into different groups that use one or more than one common signalling components. Moreover, as a result of cross-interaction, a single regulon regulates multiple targets in addition to its cognate target. Such type of signal association is defined as cross-talk. Cross-coupling in the signalling network can modify the functionality of a network topology and can subsume errors compared to the uncoupled one. Cross-interactions have been identified not only in eukaryotes but also in prokaryotes, as observed in the bacterial two-component system [10,11,12,13,14]. In the eukaryotic system, cross-talk has been identified in numerous situations [15,16,17]. Furthermore, cross-talk and several of its variants have also been identified at different stages of gene regulation [18,19,20,21,22].

Since cross-talk is observed in a broad range of biological processes, one may interrogate the functional utility of such network coordination. The cross-coupling mechanism is conveyed through generations and observed in a significant number of evolutionary descent. This signature remains prominent in the course of evolution in spite of other modifications that are taking place at the cellular level. Such character indicates that cross-talk might have a definite functional role to build up synchronised cellular regulations by spending the stored energy. If this is true, how does a cell balance the trade-off between network association and potential cost? Few comprehensive experiments on the network connectivity suggest that networks of a well-delineated cluster are correlated with each other but are uncorrelated to the rest of the network [8]. Synchronisation is necessary to attain natural activity but needs to maintain a threshold value. Otherwise, too much synchronisation may lead to a physiological disorder like epilepsy [23]. Inter-pathway cross-talk becomes prominent due to the limitation of common resources, defined as the overloaded condition. However, cross-talk effect becomes faint in the underloaded condition, where the level of available resources are present in sufficient amount [19,20,22,24,25]. A key source of survival strategy under diverse environmental conditions is the generation of fluctuations which induces non-genetic variability in a cellular population. In such a situation, cells readjust to cope with the limited resources by introducing cross-correlation among a set of genes and thus implementing a successful bet-hedging program [22]. Cross-talk also facilitates synchronisation in different organs such as cardio-respiratory interaction, brain and tissues [23].

To address the functionality of cross-talk, we undertake a representative network comprised of well-characterized features of signal transduction. To be specific, we focus on MAPK pathway, a well-studied eukaryotic signalling machinery, conserved with three kinase cascades. In S. cerevisiae, five MAPK signalling pathways are present, out of which only three (pheromone response, filamentous growth response and osmostress adaptation) use a common kinase protein Ste11 [5,8]. In fact, pheromone response and filamentous growth pathways also use the same kinase Ste7. Pheromone MAPK cascade (Ste11 → Ste7 → Fus3) is activated by mating pheromone. Under low nutrient condition, filamentous growth MAPK cascade (Ste11 → Ste7 → Kss1) becomes active whereas high external osmolarity activates the osmoadaption cascade (Ste11 → Pbs2 → Hog1) [5,8]. Due to inter-cascade correlation among the three signalling pathways, one pathway can be activated by the signal of another pathway in the absence of its own signal. Several experimental results suggest that such cross-talk is filtered out by cross-pathway inhibition, kinetic insulation and formation of scaffold protein [5,8,26,27,28,29]. Although activation through inter-pathway cross-talk and cross-pathway inhibition compensates each other, information is exchanged among the pathways during these interactions. This leads to distinct queries (i) is it possible for an individual signalling pathway to convey its input signal reliably downstream without experiencing any influence from the other channels of signal? (ii) Since the inter-pathway connectivity is known not to allow the uniqueness of transduced signals - what are the physiological advantages of cross-talk? (iii) Is there any participation of pathway output in the cooperative regulation of a downstream target in a synchronised manner? (iv) How is it possible for correlated pathways to keep up static as well as dynamic synchronisation in a single cell environment that is prevalently stochastic in nature? (v) Do correlated fluctuations have any capability to control the variability in the correlation between two different proteins?

In the present manuscript, we study generic S. cerevisiae MAPK pathways to address the potential functionality of inter-pathway cross-talk within a stochastic framework. We consider two equivalent interacting MAPK pathways, each one consisting of a linear chain of three MAPK cascade proteins [5,8]. Both pathways get activated by their corresponding external signals propagating downstream through phosphorylation (activation) and dephosphorylation (deactivation) of the cascade proteins. In addition, due to cross-talk, phosphorylation of the intermediate components of the two pathways is influenced by the activated kinase of the other pathway along with the cognate one. As the population of each cascade protein is not sufficiently high within a single cell and experiences a fluctuating environment, we express all associated chemical reactions in terms of the stochastic differential equation. We solve the coupled set of nonlinear Langevin equations using linear noise approximation (LNA) [30,31] and calculate the variance of each kinase and covariance between two different kinases (see Section 2 and Appendix). Recent theoretical development [32] shows that LNA is not only limited to high copy number but also exact up to second moments of any chemical species involved in a second-order reaction. The fluctuations associated with at least one of the species participating in each of the second-order reaction are Poissonian and uncorrelated with the fluctuations of other species. Also, LNA remains valid for faster activation and deactivation (or synthesis and degradation) rates of the corresponding components compared to the coarse-grained (steady state) time scale [30,31,32,33,34,35,36,37,38,39,40].

To classify the signal transduction efficacy through two pathways in the presence of cross-association, we quantify two as well as multivariate mutual information. Distributions of all kinase proteins are approximately considered Gaussian, allowing us to adopt a reduced expression of mutual information [41,42]. The reduced equation mainly depends on the variance and the covariance of the corresponding kinase. We validate our analytical calculation by exact stochastic simulation [43]. In the first subsection, we quantify two variable mutual information under the influence of cross-talk parameter. We also investigate the mutual information between two non-cognate kinases and find causality of this coordination. Since causality leads to synchronisation [23,44], it is important to measure causality relation between the pathways, i.e., who regulates whom and to which extent. If both pathways interact with each other and transduce information of the corresponding input signal with different degrees, then it is challenging to characterize the magnitude and direction of signal propagation. To overcome such difficulty, we define a new measure, net information transduction, using the expressions of two cross mutual information, which satisfactorily quantifies the amount of net signal propagation. In the connection of measuring directionality, it is important to mention the concept of transfer entropy [45] which has been applied in several systems [46,47,48,49]. Transfer entropy quantifies the directed information transfer within a system. In the calculation of transfer entropy, one considers the time series of the random variable and makes use of the time lagged data. In our study, however, we make use of mutual information evaluated at steady state to define a relative dimensionless measure, net information transduction D (see Equation (1)). D is a normalized quantity and is bounded within a range $\pm 1$. It helps to diagnose the direction of net information flow between the two parallel cross interacting signalling pathways. We also verify inter-pathway synchronisation with the help of coordinated fluctuations of stochastic trajectories of two parallel kinases. This result implies how two kinases are synchronised within a cell. To understand this phenomenon further, we investigate how much correlation develops between the steady state levels of two different kinases in a hypothetical population of genetically identical cells.

In the second subsection, we quantify three variable mutual information (i.e., two sources to one target, or one source to two targets) when both the channels of information flow work separately. We make use of three variable mutual information along with two variable mutual information to define the measure of net synergy [50]. In the present work, we quantify net synergy using the theory of partial information decomposition (PID) [51] as the formalism deals with the calculation of multivariate mutual information. Furthermore, it defines unique, synergistic and redundant information among the variables of interest. Following the initial development by Williams and Beer [51], several modifications of the formalism have been proposed till date [52,53,54]. In the context of dynamical systems obeying Gaussian noise processes the formalism has been further extended by Barrett [55]. We note here that the measure of interaction information is the negative of the measure of net synergy [56]. Considering the interactions between two parallel MAPK pathway, we identify signal integration and signal bifurcation motif. In signal integration motif, two MAPKKK proteins transmit signal into a single MAPKK protein which amounts to a sub-motif with two inputs and one output. On the other hand, in signal bifurcation motif, a single MAPKKK protein transmits signal to two MAPKK protein thus leading into a sub-motif with one input and two outputs. Quantification of information transmission in these two sub-motifs requires calculation of both two variable and multivariate mutual information which in turn provides the metric of net synergy. We observe the sign of net synergy value changes depending on the signal integration as well as signal bifurcation and is mainly controlled by pathway architecture.

2. Results and Discussion

2.1. Two Variable Mutual Information

The parameters ${\epsilon }_2$ and ${\epsilon }_1$ control the signalling channel X and Y, respectively. In Figure 1, we show the mutual information profile as a function of ${\epsilon }_2$ for two different sets of parameters (see

Table 1. Reactions and corresponding parameter values for the MAPK network motif of S. cerevisiae [6,9,66], related to faster relaxation rate. Other Parameters are $s_x=s_y=10$ molecules/cell, $x_{T1}=x_1+x_{p1}=250$ molecules/cell, $x_{T2}=x_2+x_{p2}=1700$ molecules/cell, $x_{T3}=x_3+x_{p3}=5000$ molecules/cell, $y_{T1}=y_1+y_{p1}=250$ molecules/cell, $y_{T2}=y_2+y_{p2}=1700$ molecules/cell and $y_{T3}=y_3+y_{p3}=5000$ molecules/cell. The kinetic schemes adopted in the present work follows the model of Heinrich et al. [63].

Description	Reaction		Propensity Function	Rate Constant
Activation of $x_1$	$x_1+s_x\stackrel{k_x}{⟶}{x}_{p1}+s_x$		$k_x{s}_x{x}_1$	$k_x={10}^{-4}$ molecules${}^{-1}$ s${}^{-1}$
Deactivation of $x_{p1}$	$x_{p1}\stackrel{{\alpha }_1}{⟶}{x}_1$		${\alpha }_1{x}_{p1}$	${\alpha }_1=0.01$ s${}^{-1}$
Activation of $y_1$	$y_1+s_y\stackrel{k_y}{⟶}{y}_{p1}+s_y$		$k_y{s}_y{y}_1$	$k_y={10}^{-4}$ molecules${}^{-1}$ s${}^{-1}$
Deactivation of $y_{p1}$	$y_{p1}\stackrel{{\beta }_1}{⟶}{y}_1$		${\beta }_1{y}_{p1}$	${\beta }_1=0.01$ s${}^{-1}$
Activation of $x_2$	$x_2+x_{p1}\stackrel{k_{12x}}{⟶}{x}_{p2}+x_{p1}$		$k_{12x}{x}_{p1}{x}_2$	$k_{12x}={10}^{-4}$ molecules${}^{-1}$ s${}^{-1}$
Activation of $x_2$	$x_2+y_{p1}\stackrel{{\epsilon }_2}{⟶}{x}_{p2}+y_{p1}$		${\epsilon }_2{y}_{p1}{x}_2$	${\epsilon }_2=(0-1)\times {10}^{-4}$ molecules${}^{-1}$ s${}^{-1}$
Deactivation of $x_{p2}$	$x_{p2}\stackrel{{\alpha }_2}{⟶}{x}_2$		${\alpha }_2{x}_{p2}$	${\alpha }_2=0.05$ s${}^{-1}$
Activation of $y_2$	$y_2+y_{p1}\stackrel{k_{12y}}{⟶}{y}_{p2}+y_{p1}$		$k_{12y}{y}_{p1}{y}_2$	$k_{12y}={10}^{-4}$ molecules${}^{-1}$ s${}^{-1}$
Activation of $y_2$	$y_2+x_{p1}\stackrel{{\epsilon }_1}{⟶}{y}_{p2}+x_{p1}$		${\epsilon }_1{x}_{p1}{y}_2$	${\epsilon }_1=(0-1)\times {10}^{-4}$ molecules${}^{-1}$ s${}^{-1}$
Deactivation of $y_{p2}$	$y_{p2}\stackrel{{\beta }_2}{⟶}{y}_2$		${\beta }_2{y}_{p2}$	${\beta }_2=0.05$ s${}^{-1}$
Activation of $x_3$	$x_3+x_{p2}\stackrel{k_{23x}}{⟶}{x}_{p3}+x_{p2}$		$k_{23x}{x}_{p2}{x}_3$	$k_{23x}=5\times {10}^{-5}$ molecules${}^{-1}$ s${}^{-1}$
Deactivation of $x_{p3}$	$x_{p3}\stackrel{{\alpha }_3}{⟶}{x}_3$		${\alpha }_3{x}_{p3}$	${\alpha }_3=0.05$ s${}^{-1}$
Activation of $y_3$	$y_3+y_{p2}\stackrel{k_{23y}}{⟶}{y}_{p3}+y_{p2}$		$k_{23y}{y}_{p2}{y}_3$	$k_{23y}=5\times {10}^{-5}$ molecules${}^{-1}$ s${}^{-1}$
Deactivation of $y_{p3}$	$y_{p3}\stackrel{{\beta }_3}{⟶}{y}_3$		${\beta }_3{y}_{p3}$	${\beta }_3=0.05$ s${}^{-1}$

and

Table 2. Reactions and corresponding parameter values for the MAPK network motif of S. cerevisiae [6,9,66], related to slower relaxation rate. Other Parameters are $s_x=s_y=10$ molecules/cell, $x_{T1}=x_1+x_{p1}=250$ molecules/cell, $x_{T2}=x_2+x_{p2}=1700$ molecules/cell, $x_{T3}=x_3+x_{p3}=5000$ molecules/cell, $y_{T1}=y_1+y_{p1}=250$ molecules/cell, $y_{T2}=y_2+y_{p2}=1700$ molecules/cell and $y_{T3}=y_3+y_{p3}=5000$ molecules/cell. The kinetic schemes adopted in the present work follows the model of Heinrich et al. [63].

Description	Reaction		Propensity Function	Rate Constant
Activation of $x_1$	$x_1+s_x\stackrel{k_x}{⟶}{x}_{p1}+s_x$		$k_x{s}_x{x}_1$	$k_x={10}^{-4}$ molecules${}^{-1}$ s${}^{-1}$
Deactivation of $x_{p1}$	$x_{p1}\stackrel{{\alpha }_1}{⟶}{x}_1$		${\alpha }_1{x}_{p1}$	${\alpha }_1=0.01$ s${}^{-1}$
Activation of $y_1$	$y_1+s_y\stackrel{k_y}{⟶}{y}_{p1}+s_y$		$k_y{s}_y{y}_1$	$k_y={10}^{-4}$ molecules${}^{-1}$ s${}^{-1}$
Deactivation of $y_{p1}$	$y_{p1}\stackrel{{\beta }_1}{⟶}{y}_1$		${\beta }_1{y}_{p1}$	${\beta }_1=0.01$ s${}^{-1}$
Activation of $x_2$	$x_2+x_{p1}\stackrel{k_{12x}}{⟶}{x}_{p2}+x_{p1}$		$k_{12x}{x}_{p1}{x}_2$	$k_{12x}={10}^{-4}$ molecules${}^{-1}$ s${}^{-1}$
Activation of $x_2$	$x_2+y_{p1}\stackrel{{\epsilon }_2}{⟶}{x}_{p2}+y_{p1}$		${\epsilon }_2{y}_{p1}{x}_2$	${\epsilon }_2=(0-1)\times {10}^{-4}$ molecules${}^{-1}$ s${}^{-1}$
Deactivation of $x_{p2}$	$x_{p2}\stackrel{{\alpha }_2}{⟶}{x}_2$		${\alpha }_2{x}_{p2}$	${\alpha }_2=0.01$ s${}^{-1}$
Activation of $y_2$	$y_2+y_{p1}\stackrel{k_{12y}}{⟶}{y}_{p2}+y_{p1}$		$k_{12y}{y}_{p1}{y}_2$	$k_{12y}={10}^{-4}$ molecules${}^{-1}$ s${}^{-1}$
Activation of $y_2$	$y_2+x_{p1}\stackrel{{\epsilon }_1}{⟶}{y}_{p2}+x_{p1}$		${\epsilon }_1{x}_{p1}{y}_2$	${\epsilon }_1=(0-1)\times {10}^{-4}$ molecules${}^{-1}$ s${}^{-1}$
Deactivation of $y_{p2}$	$y_{p2}\stackrel{{\beta }_2}{⟶}{y}_2$		${\beta }_2{y}_{p2}$	${\beta }_2=0.01$ s${}^{-1}$
Activation of $x_3$	$x_3+x_{p2}\stackrel{k_{23x}}{⟶}{x}_{p3}+x_{p2}$		$k_{23x}{x}_{p2}{x}_3$	$k_{23x}={10}^{-5}$ molecules${}^{-1}$ s${}^{-1}$
Deactivation of $x_{p3}$	$x_{p3}\stackrel{{\alpha }_3}{⟶}{x}_3$		${\alpha }_3{x}_{p3}$	${\alpha }_3=0.01$ s${}^{-1}$
Activation of $y_3$	$y_3+y_{p2}\stackrel{k_{23y}}{⟶}{y}_{p3}+y_{p2}$		$k_{23y}{y}_{p2}{y}_3$	$k_{23y}={10}^{-5}$ molecules${}^{-1}$ s${}^{-1}$
Deactivation of $y_{p3}$	$y_{p3}\stackrel{{\beta }_3}{⟶}{y}_3$		${\beta }_3{y}_{p3}$	${\beta }_3=0.01$ s${}^{-1}$

) while keeping ${\epsilon }_1$ constant. Figure 1A shows that mutual information between $x_{p1}$ and $x_{p2}$ kinases decays with the increment of ${\epsilon }_2$. Augmentation of ${\epsilon }_2$ includes a competition between $x_{p1}$ and $y_{p1}$ to phosphorylate the $x_2$ kinase. During phosphorylation, mutual association is originated, and signal transduction is ensued. Thus, for the low value of ${\epsilon }_2$, maximum level of mutual information is attained due to minimal phosphorylation competition. On the other hand, minimum level of mutual information is propagated at high ${\epsilon }_2$ value due to maximum phosphorylation contribution of $y_{p1}$. In Figure 1B, mutual information between $x_{p1}$ and $y_{p2}$ is plotted, which shows a constant value as a function of ${\epsilon }_2$. This happens as ${\epsilon }_2$ has no influence in the alteration of mutual information. The same logic is applicable to the mutual information between $y_{p1}$ and $y_{p2}$ shown in Figure 1C. In Figure 1D, mutual information between $y_{p1}$ and $x_{p2}$ increases as a function of ${\epsilon }_2$, as ${\epsilon }_2$ is only responsible for establishing the cross-talk between $y_{p1}$ and $x_{p2}$. This result implies that with the enhancement of cross-talk the process of signal integration through $y_{p1}$ increases. The same profiles can be generated as a function of ${\epsilon }_1$, while keeping ${\epsilon }_2$ fixed. These results together indicate that $\mathcal{I}(x_{p1};x_{p2})$ and $\mathcal{I}(y_{p1};x_{p2})$ depend on ${\epsilon }_2$, whereas $\mathcal{I}(x_{p1};y_{p2})$ and $\mathcal{I}(y_{p1};y_{p2})$ depend on ${\epsilon }_1$. In Figure 1A–D, the green lines drawn for slower relaxation rate (see Table 2) always maintains a lower mutual information value compared to the red lines drawn for faster relaxation rate (see Table 1). Relaxation rates of the corresponding kinases i.e., $x_{p1}$, $x_{p2}$, $x_{p3}$, $y_{p1}$, $y_{p2}$ and $y_{p3}$ are $-J_{x1x1}=({\alpha }_1+k_x{s}_x)$, $-J_{x2x2}=({\alpha }_2+k_{12x}\langle x_{p1}\rangle +{\epsilon }_2\langle y_{p1}\rangle )$, $-J_{x3x3}=({\alpha }_3+k_{23x}\langle x_{p2}\rangle )$, $-J_{y1y1}=({\beta }_1+k_y{s}_y)$, $-J_{y2y2}=({\beta }_2+k_{12y}\langle y_{p1}\rangle +{\epsilon }_1\langle x_{p1}\rangle )$ and $-J_{y3y3}=({\beta }_3+k_{23y}\langle y_{p2}\rangle )$, respectively, where the angular bracket $\langle \cdots \rangle$ indicates the deterministic copy number at long time limit (see Appendix). An input signal can reliably flow downstream if relaxation rate (or degradation rate) of a cascade protein is higher than that of its upstream cascade proteins [35]. For red line, we consider higher degradation rate for $x_{p2}$ and $x_{p3}$ ($y_{p2}$ and $y_{p3}$) compared to $x_{p1}$ ($y_{p1}$). Thus, faster relaxation rates are attained under this condition with high information propagation capacity.

Next, we quantify mutual information between two parallel kinases ($x_{pi}$ and $y_{pj}$, with $i=j$) of the two equivalent interacting MAPK pathways. The inter pathway coupling is unidirectional when either ${\epsilon }_1$ or ${\epsilon }_2$ is zero but is bidirectional when both are non-zero. In this situation, both variables ($x_{pi}$ and $y_{pj}$) do not interact with each other but are regulated by a common kinase regulon incorporating coordinated fluctuations into these variables. In other words, quantification of mutual information actually evaluates the extent of cross-correlation between these two variables. We observe zero mutual information value between $x_{p1}$ and $y_{p1}$, as these are uncorrelated. In Figure 1E, we show mutual information between $x_{p2}$ and $y_{p2}$ as a function of ${\epsilon }_2$ keeping ${\epsilon }_1$ fixed. The profile shows an increasing trend as cross-talk parameter ${\epsilon }_2$ increases. Similarly, in Figure 1F, mutual information between $x_{p3}$ and $y_{p3}$ is shown with a similar trend as in Figure 1E. Interestingly, for faster relaxation time scale, mutual information between similar cascade kinases increases while moving from second ($x_{p2}$ and $y_{p2}$) to third ($x_{p3}$ and $y_{p3}$) cascade. On the other hand, an opposite trend is observed for slower relaxation time scale. This characteristic trend is further shown in Figure 1G,H using bar diagram. These results together suggest that fluctuations due to faster relaxation rate transduce correlated fluctuations in a better way compared to the slower one. In Figure 1E, mutual information is high for slower relaxation rate than the faster one, as slower rate parameters yield high level of $x_{p2}$ and $y_{p2}$ which in turn incorporate extra fluctuations that help to increase mutual association. A similar result is also observed in Figure 1F. Identical mutual information profiles of $\mathcal{I}(x_{p2};y_{p2})$ and $\mathcal{I}(x_{p3};y_{p3})$ can be generated as function of ${\epsilon }_1$ keeping ${\epsilon }_2$ fixed. These results suggest that both the cross-talk parameters ${\epsilon }_1$ and ${\epsilon }_2$ contribute equally to the development of an association between two parallel pathways.

Both the mutual information between $x_{p2}$ and $y_{p2}$, $x_{p3}$ and $y_{p3}$ are capable of providing a satisfactory explanation of enhancement of cross-talk with the increment of inter pathway interaction parameters (${\epsilon }_1$ and ${\epsilon }_2$). Under equivalent interactions condition (${\epsilon }_1={\epsilon }_2$), each pathway shares its information with other to an equal extent and is quantified not only by $\mathcal{I}(x_{p2};y_{p2})$ and $\mathcal{I}(x_{p3};y_{p3})$ but also by $\mathcal{I}(x_{p1};y_{p2})$ and $\mathcal{I}(y_{p1};x_{p2})$. However, characterization of the direction of information transduction is difficult under unequal condition (${\epsilon }_1\ne {\epsilon }_2$). Except the equivalent condition (${\epsilon }_1={\epsilon }_2$) where the net information ( $\mathcal{I}(y_{p1};x_{p2})-\mathcal{I}(x_{p1};y_{p2})$ ) flow is zero, it has a definite value with directionality (positive or negative value) at all other conditions. Since the definition of mutual information is symmetric in nature and usage of the same is difficult to provide directionality of information propagation, we define a dimensionless quantity, net information transduction (D) using $\mathcal{I}(x_{p1};y_{p2})$ and $\mathcal{I}(y_{p1};x_{p2})$ as

$$D=\frac{\mathcal{I}(y_{p1};x_{p2})-\mathcal{I}(x_{p1};y_{p2})}{\mathcal{I}(y_{p1};x_{p2})+\mathcal{I}(x_{p1};y_{p2})}.$$

(1)

The above expression implies that it is maximal ($D=1$) when $\mathcal{I}(x_{p1};y_{p2})$ is zero, i.e., no information propagation from $x_{p1}$ to $y_{p2}$ (${\epsilon }_1=0$). It is minimal ($D=-1$) when $\mathcal{I}(y_{p1};x_{p2})$ is zero, which specifies zero information propagation from $y_{p1}$ to $x_{p2}$ (${\epsilon }_2=0$). In Figure 1I, we show the profile of D as a function of ${\epsilon }_2$ while keeping ${\epsilon }_1$ fixed, where the value of D changes from negative to positive as ${\epsilon }_2$ increases. It suggests that at low ${\epsilon }_2$, information flowing from X to Y pathway dominates over the flow from Y to X. In other words, in this regime, the net information flow is accounted for by X → Y, leading to a negative value of D. On the other hand, at high ${\epsilon }_2$, the direction of net information propagation is from Y to X due to reverse situation and generates a positive D value. The opposite scenario can be observed if one generates the profile of D as a function of ${\epsilon }_1$ for fixed ${\epsilon }_2$ In this connection, it is important to mention that both the relaxation time scale limits generate a similar profile of D. As a result, both the profiles of D exhibit data collapse when depicted as a function of ${\epsilon }_2$ for fixed ${\epsilon }_1$ (Figure 1I) or vice versa. This observation indicates that normalised profiles of D are independent of relaxation time scales. In Figure 1J, we also show a 2d-surface plot of D as a function of both ${\epsilon }_1$ and ${\epsilon }_2$ for faster relaxation time scale (Table 1). The surface plot indicates zero (or near to zero) value of D along the diagonal region (${\epsilon }_1\approx {\epsilon }_2$). However, the off diagonal region is positive for ${\epsilon }_1<{\epsilon }_2$ and negative for ${\epsilon }_1>{\epsilon }_2$.

In Figure 2A,B, we show two 2d-surface plots of mutual information between $x_{p2}$ and $y_{p2}$, $x_{p3}$ and $y_{p3}$ kinases, respectively, as a function of two cross-interaction parameters ${\epsilon }_1$ and ${\epsilon }_2$ under faster relaxation time scale (Table 1). Both figures show maximum mutual information at high values of the two parameters. Since, ${\epsilon }_1$ and ${\epsilon }_2$ are equally responsible for developing the cross-correlation between two pathways, one can check the effect of maximisation of mutual information by increasing any one of these two parameters. Although we can quantify the influence of cross-talk with the help of two variable mutual information, $\mathcal{I}(x_{p2};y_{p2})$ and $\mathcal{I}(x_{p3};y_{p3})$, it is difficult to get an insight how the static and dynamic populations of the phosphorylated kinases are correlated. To this end, we have checked such correlation as shown in Figure 2C. Dynamic correlation can be measured using stochastic trajectories of the two variables in a single cell. If sufficient association between two trajectories exist, then correlated fluctuations are observed i.e., one trajectory closely follows the other. Otherwise, uncorrelated fluctuations (trajectories do not follow each other) are seen in the absence of cross-talk. In Figure 2C, we show stochastic time series of different phosphorylated kinases under varied conditions. Here, four different sets of ${\epsilon }_1$ and ${\epsilon }_2$ parameters have been used - mentioned as I, II, III and IV in Figure 2A,B. The stochastic time series exhibit correlated fluctuations at high ${\epsilon }_1$ and ${\epsilon }_2$. However, uncorrelated time series are seen at low ${\epsilon }_1$ and ${\epsilon }_2$. Each of these time series are generated from a single run of stochastic simulation and represents the dynamics in a single cell.

Next, we have executed 1000 independent stochastic simulations and measured the steady state phosphorylated kinase levels from each run to draw a scatter plot. Here, the collection of symbols mimics the behaviour of a hypothetical population of 1000 genetically identical cells. Each symbol in the scatter plots is due to a single stochastic run and represents copies of phosphorylated kinases ($x_{pi}$ and $y_{pj}$, $i,j=1,2,3$) in a single cell expressed at steady state. In other words, each symbol signifies the steady state phosphorylated protein levels in an individual cell within a hypothetical population. Using this concept, we measure static correlation among the pairs of phosphorylated kinases ($(x_{p1},y_{p1})$, etc.) produced in different cells within the population. The motivation behind creating a scatter plot is to measure the static correlation between two phosphorylated protein levels expressed in different cells [57,58,59]. In Figure 2C, we show static correlation between different phosphorylated kinases expressed in different cells. For plots with high ${\epsilon }_1$ and ${\epsilon }_2$ values, most of the symbols are aligned diagonally in a narrow strip but for low ${\epsilon }_1$ and ${\epsilon }_2$, symbols are distributed in a much larger space.

Correlated variation of these kinases in the population is observed along the diagonal direction. However, an uncorrelated variation along the off-diagonal direction reflects an increase in variability between two phosphorylated kinases with respect to each other in a population. Therefore, enhancement of cross-talk decreases variability between two kinases in a population. These results imply that cross-talk develops correlated fluctuations between kinases in a population, thereby assisting in the successful development of a robust adaptation machinery as observed in the bet-hedging program under diverse environmental conditions [22]. Here, an increase in the correlation between two different protein pool (with high Pearson’s correlation coefficient) is observed for maximal cross-talk. On the other hand, a decrease in correlation is seen for minimal cross-interaction. Similar behaviour can be seen under slower relaxation time scale (Table 2) as shown in Figure 3. The primary difference between the nature of correlation between $(x_{p2},y_{p2})$ and $(x_{p3},y_{p3})$ are visible from Figure 2 and Figure 3. In Figure 2, the correlation between $(x_{p3},y_{p3})$ is always higher than $(x_{p2},y_{p2})$ for all four conditions. On the other hand, in Figure 3, it shows an opposite trend. We note that coordinated fluctuations of different MAPK proteins are generated due to cross-talk between two parallel pathways that ultimately lead to similar cellular function with high reliability. Fluctuations in the transcription factor are much faster than the transcript dynamics. On the other hand, fluctuations in the transcription factor (known as extrinsic fluctuations) are comparable with the gene switching (on-off) time scale. Such time scale agreement gets reflected in the dynamics of the transcript [36].

2.2. Three Variable Mutual Information

In the foregoing discussion, we have shown the effect of cross-talk in terms of conventional two variable mutual information. However, as cross-interaction between two pathways develops a complex network, a comprehensive study of three variable mutual information provides an extra insight. In the present study, three variable mutual information is defined as $\mathcal{I}(x_{p1},y_{p1};x_{p2})$ and $\mathcal{I}(x_{p1};y_{p2},x_{p2})$ for two types of branched pathways (see Section 3.2). One of the branched pathways is two inputs ($x_{p1}$ and $y_{p1}$) and one output ($x_{p2}$ or $y_{p2}$) motif where two input signals are integrated into a single output. The other is one input ($x_{p1}$ or $y_{p1}$) and two outputs ($x_{p2}$ and $y_{p2}$) motif where the input signal is bifurcated into two outputs. In this subsection, we investigate the efficacy of such signal integration as well as signal bifurcation. Since marginal and joint distributions of all cascade proteins are considered approximately Gaussian, we adopt multivariate mutual information theory [49,50,51] to analytically estimate three variable mutual information [55]. Each branched motif consists of two signal propagating channels that work together. It is thus interesting to investigate whether these signalling channels perform separately and what significant change arises in the estimation of three variable mutual information. The change in the magnitude of mutual information is defined by net synergy [50] and is evaluated using the theory of partial information decomposition [49,51,55] in terms of the difference between three variable mutual information and two corresponding two variable mutual information. The value of net synergy is either positive or negative; a positive value indicates synergy (extra information) whereas negative value measures redundancy (deficit of information) [50,55].

In Figure 4A, we show mutual information, $\mathcal{I}(x_{p1},y_{p1};x_{p2})$ of two inputs and one output model as a function of ${\epsilon }_2$ for a fixed value of ${\epsilon }_1$. The profile shows a bifunctional behaviour with the increment of ${\epsilon }_2$; initially it decreases up to a certain value of ${\epsilon }_2$, and then it increases. At low ${\epsilon }_2$, a minimal amount of signal is propagated from $y_{p1}$ to $x_{p2}$. Consequently, the motif reduces to a single input-output motif and the motif regains its native form due to the significant contribution of ${\epsilon }_2$. In Figure 1A, we show that two variable mutual information between $x_{p1}$ and $x_{p2}$ decreases with the increment of ${\epsilon }_2$. Similar situation arises in Figure 4A for low value of ${\epsilon }_2$. On the contrary, $\mathcal{I}(y_{p1};x_{p2})$ increases with the increment of ${\epsilon }_2$ (Figure 1D). Thus, two opposing effects work together to generate the convex profile. In Figure 4B, we plot net synergy of the motif as a function of ${\epsilon }_2$ for a fixed value of ${\epsilon }_1$ and it is seen to increase monotonically. It is pertinent to mention here that for this motif, one always gets a positive net synergy value as a function of ${\epsilon }_2$. This result implies that an integrating signalling motif transduces more information compared to the summation of two isolated channels. The extra information, i.e., synergy facilitates fidelity of the output kinase. Intuitively, the sum of the reduction in the uncertainty (cross-correlation) of the output kinase contributed by each input signal is lower than the reduction in the uncertainty of the output provided by both signals together. This phenomenon implicates the aspect of integration of multiple signals in cellular signalling network motif as observed in V. harveyi quorum-sensing circuit [60,61].

In Figure 4C, we show mutual information $\mathcal{I}(x_{p1};x_{p2},y_{p2})$ of one input and two outputs motif with the increment of ${\epsilon }_2$ for a fixed value of ${\epsilon }_1$. The mutual information value decreases with ${\epsilon }_2$ since propagation of information from $x_{p1}$ to $x_{p2}$ is only inhibited by the cross-interaction. However, ${\epsilon }_2$ does not have any influence in information propagation from $x_{p1}$ to $y_{p2}$ and remains unaltered. Thus, three variable mutual information profile follows a decreasing trend. Figure 4D shows decreasing trend of net synergy profile as a function of ${\epsilon }_2$ for a fixed value of ${\epsilon }_1$. Importantly, for this motif negative values of net synergy are observed irrespective of the value of ${\epsilon }_2$. This indicates redundancy in the information transmission in this composite motif compared to the sum of the individual one. Naturally, predictability about the output kinases decreases when two isolated signal propagation channels work together to form a bifurcated signal transduction motif. This result implies that although bifurcated signalling model reduces mutual information, it has a biological significance of the activation of multiple signalling channels in the presence of a single input as identified in the chemotaxis system of E. coli [62]. In all figures (Figure 4A–D), the red lines are plotted for faster relaxation rate constants (Table 1) of $x_{p2},x_{p3},y_{p2}$ and $y_{p3}$ and the green lines are for slower relaxation rate constants (Table 2).

3. Materials and Methods

In Figure 5, we show a schematic diagram of two interacting parallel MAPK pathways (named as X and Y). Each MAPK pathway consists of three kinase components, i.e., $x_1$, $x_2$, $x_3$ (X pathway) and $y_1$, $y_2$, $y_3$ (Y pathway) [7,28,63,64,65]. $x_i$ and $x_{pi}$ represent dephosphorylated and phosphorylated form of a kinase proteins, respectively, and the same applies to $y_j$ and $y_{pj}$ (here $i,j=1,2,3$). The first cascade protein of a MAPK pathway gets phosphorylated with an exposure to the external stimulus. While phosphorylated, it positively regulates the phosphorylation of its own downstream kinase along with the kinase of the other pathway. The phosphorylated intermediate kinase regulates phosphorylation of the last kinase. To maintain the pool of phosphorylated kinase within a cell, a dephosphorylation process is in action with the help of phosphatase molecules. The cross-pathway interactions between two parallel MAPK pathways are denoted by the dashed lines in Figure 5 along with the cross-interaction rate parameters ${\epsilon }_1$ and ${\epsilon }_2$. $S_x$ and $S_y$ are the two extra-cellular signals acting on the X and Y pathway, respectively.

Both pathways get causally correlated through cross-interactions, and a cross-talk develops as a consequence. Causal relationships are frequently examined in various circumstances that are subjected to stochastic fluctuations [23,44,50,55]. In the present manuscript, we quantify the causal relationship in terms of mutual information. Here, the two cross-interaction parameters ${\epsilon }_1$ and ${\epsilon }_2$ play a significant role in establishing different levels of cross-talk. The parameter ${\epsilon }_1$ controls information flow from X to Y pathway ($x_{p1}\to y_{p2}$), but the parameter ${\epsilon }_2$ is responsible for Y to X pathway ($y_{p1}\to x_{p2}$) information flow. In this connection, it is important to mention that during mating process, both pheromone and filamentous growth pathways are activated to a roughly equal extent, whereas during invasive growth process, only filamentous growth pathway is activated [26]. These observations corroborate with our model development. In our calculation, we only consider the post-translationally modified forms of all MAPK proteins. Thus, in the model, the total population of a MAPK protein is the sum of the phosphorylated and the unphosphorylated form of the protein and is considered to be constant ($(x_i+x_{pi})=x_{Ti}=(y_j+y_{pj})=y_{Tj}=$ constant, here $i=j$). In addition, we consider a physiologically relevant parameter set for our calculation [6,9,66].

3.1. Two Variable Mutual Information

Adopting Shannon’s information theory [41,42], we have calculated two variable mutual information between two phosphorylated kinases,

$$\mathcal{I}(x_{pi};y_{pj})=\sum _{x_{pi}}\sum _{y_{pj}}p(x_{pi},y_{pj}){log}_2\left[\frac{p(x_{pi},y_{pj})}{p\left(x_{pi}\right)p\left(y_{pj}\right)}\right].$$

(2)

A generalised index $x_{pi}$ and $y_{pj}$ have been considered to represent the copy number of two different phosphorylated kinases. Similarly, $p\left(x_{pi}\right)$ and $p\left(y_{pj}\right)$ are the marginal and $p(x_{pi},y_{pj})$ is the joint probability distributions associated with the corresponding kinases. For the calculation of mutual information between two kinases of X signalling pathway, we have replaced $y_{pj}$ by $x_{pi}$ (where $i\ne j$) and the reverse replacement has been followed for Y signalling pathway. For the estimation of mutual information between two equivalent kinases ($x_{pi}$ and $y_{pj}$) of the respective pathways, we have used the same formula for $i=j$ condition. Mutual information can also be written in the form of the entropy function. Hence, Equation (2) can be redefined as

$$\mathcal{I}(x_{pi};y_{pj})=H\left(x_{pi}\right)+H\left(y_{pj}\right)-H(x_{pi},y_{pj}).$$

(3)

Here, $H\left(x_{pi}\right)$ and $H\left(y_{pj}\right)$ are individual and $H(x_{pi},y_{pj})$ is total entropy of the respective kinases. In the present study, both probability distribution functions (marginal and joint) are approximately considered to be Gaussian as experimental studies on MAPK pathway show Gaussian dynamics [67]. Thus, using Gaussian channel approximation [41,42,55], Equation (3) takes the reduced form

$$\mathcal{I}(x_{pi};y_{pj})=\frac{1}{2}{log}_2\left[\frac{{\sigma }_{x_{pi}}^2{\sigma }_{y_{pi}}^2}{{\sigma }_{x_{pi}}^2{\sigma }_{y_{pi}}^2-{\sigma }_{x_{pi}{y}_{pi}}^4}\right],$$

(4)

where ${\sigma }_{x_{pi}}^2$ and ${\sigma }_{y_{pj}}^2$ are variances and ${\sigma }_{x_{pi}{y}_{pj}}^2$ is covariance of the corresponding kinases (for detailed calculation see Appendix). At this point it is important to mention that no prior knowledge is required about the nature of probability distribution function for evaluating mutual information using Equation (2). For exact or approximate Gaussian distribution, one can reduce Equation (2) to Equation (4) applying Gaussian channel approximation. However, for systems with non-Gaussian distribution, one can still use Equation (2) with proper analytical expressions of probability distribution functions that may contribute expressions of higher moments in Equation (4).

In the present work, all expressions of two variable mutual information are calculated using Equation (4). The analytical results are then validated by evaluating probability distribution functions (Equation (2)) using exact numerical simulation [43]. The two variable mutual information value is bounded within a scale $0\le \mathcal{I}(x_{pi};y_{pj})\le min(H\left(x_{pi}\right),H\left(y_{pj}\right))$. To quantify the association between two equivalent kinases, we have used Pearson’s correlation coefficient (${\rho }_{ij},i=j$) [68]

$${\rho }_{ij}=\frac{{\sigma }_{x_{pi}{y}_{pj}}^2}{{\sigma }_{x_{pi}}{\sigma }_{y_{pj}}}.$$

(5)

3.2. Three Variable Mutual Information

The three variable mutual information are calculated for both signal integration and signal bifurcation motif. In the first motif, two phosphorylated input kinases interact with one output kinase. Hence the complete description of mutual information is given by

$$\begin{array}{ccc}\hfill \mathcal{I}(x_{p1},y_{p1};x_{p2})& =& \sum _{x_{p1},y_{p1}}\sum _{x_{p2}}p(x_{p1},y_{p1},x_{p2}){log}_2\left[\frac{p(x_{p1},y_{p1},x_{p2})}{p(x_{p1},y_{p1})p\left(x_{p2}\right)}\right],\hfill \end{array}$$

(6)

where $p(x_{p1},y_{p1},x_{p2})$ and $p(x_{p1},y_{p1})$ are the joint distribution functions of the corresponding components. On the other hand, $p\left(x_{p2}\right)$ is the marginal distribution of phosphorylated $x_2$ kinase. One can also write Equation (6) in terms of the respective entropy

$$\mathcal{I}(x_{p1},y_{p1};x_{p2})=H(x_{p1},y_{p1})+H\left(x_{p2}\right)-H(x_{p1},y_{p1},x_{p2}).$$

(7)

Similarly, using Gaussian approximation [41,42,55], one can reduce Equation (7) into the following form

$$\mathcal{I}(x_{p1},y_{p1};x_{p2})=\frac{1}{2}{log}_2\left[\frac{{\sigma }_{x_{p2}}^2({\sigma }_{x_{p1}}^2{\sigma }_{y_{p1}}^2-{\sigma }_{x_{p1}{y}_{p1}}^4)}{|{\Delta }_1|}\right],$$

(8)

with

$$\begin{array}{ccc}\hfill |{\Delta }_1|& =& \left(\begin{array}{ccc}{\sigma }_{x_{p1}}^2& {\sigma }_{x_{p1}{y}_{p1}}^2& {\sigma }_{x_{p1}{x}_{p2}}^2\\ {\sigma }_{y_{p1}{x}_{p1}}^2& {\sigma }_{y_{p1}}^2& {\sigma }_{y_{p1}{x}_{p2}}^2\\ {\sigma }_{x_{p2}{x}_{p1}}^2& {\sigma }_{x_{p2}{y}_{p1}}^2& {\sigma }_{x_{p2}}^2\end{array}\right).\hfill \end{array}$$

Here, the magnitude of three variable mutual information is bounded within a scale $0\le \mathcal{I}(x_{p1},y_{p1};x_{p2})\le min(H(x_{p1},y_{p1}),H\left(x_{p2}\right))$. Using PID formalism, the three variable mutual information can be decomposed into two parts [50,55]. As a result, the net synergy expression becomes

$$\begin{array}{ccc}\hfill \Delta \mathcal{I}(x_{p1},y_{p1};x_{p2})& =& \mathcal{I}(x_{p1},y_{p1};x_{p2})-\mathcal{I}(x_{p1};x_{p2})-\mathcal{I}(y_{p1};x_{p2})\hfill \\ \hfill & =& \frac{1}{2}{log}_2\left[\frac{({\sigma }_{x_{p1}}^2{\sigma }_{y_{p1}}^2-{\sigma }_{x_{p1}{y}_{p1}}^4)({\sigma }_{x_{p1}}^2{\sigma }_{x_{p2}}^2-{\sigma }_{x_{p1}{x}_{p2}}^4)({\sigma }_{y_{p1}}^2{\sigma }_{x_{p2}}^2-{\sigma }_{y_{p1}{x}_{p2}}^4)}{|{\Delta }_1|{\sigma }_{x_{p1}}^2{\sigma }_{y_{p1}}^2{\sigma }_{x_{p2}}^2}\right].\hfill \end{array}$$

(9)

Furthermore, one can calculate mutual information for the signal bifurcating motif with the help of associated distribution functions

$$\begin{array}{ccc}\hfill \mathcal{I}(x_{p1};x_{p2},y_{p2})& =& \sum _{x_{p1}}\sum _{x_{p2},y_{p2}}p(x_{p1},x_{p2},y_{p2}){log}_2\left[\frac{p(x_{p1},x_{p2},y_{p2})}{p\left(x_{p1}\right)p(x_{p2},y_{p2})}\right].\hfill \end{array}$$

(10)

and the entropy representation of Equation (10) is

$$\mathcal{I}(x_{p1};x_{p2},y_{p2})=H\left(x_{p1}\right)+H(x_{p2},y_{p2})-H(x_{p1},x_{p2},y_{p2}).$$

(11)

Using Gaussian approximation [41,42,55] Equation (11) becomes

$$\mathcal{I}(x_{p1};x_{p2},y_{p2})=\frac{1}{2}{log}_2\left[\frac{{\sigma }_{x_{p1}}^2({\sigma }_{x_{p2}}^2{\sigma }_{y_{p2}}^2-{\sigma }_{x_{p2}{y}_{p2}}^4)}{|{\Delta }_2|}\right],$$

(12)

with

$$\begin{array}{ccc}\hfill |{\Delta }_2|& =& \left(\begin{array}{ccc}{\sigma }_{x_{p1}}^2& {\sigma }_{x_{p1}{x}_{p2}}^2& {\sigma }_{x_{p1}{y}_{p2}}^2\\ {\sigma }_{x_{p2}{x}_{p1}}^2& {\sigma }_{x_{p2}}^2& {\sigma }_{x_{p2}{y}_{p2}}^2\\ {\sigma }_{y_{p2}{x}_{p1}}^2& {\sigma }_{y_{p2}{x}_{p2}}^2& {\sigma }_{y_{p2}}^2\end{array}\right).\hfill \end{array}$$

In addition, the three variable mutual information value is bounded within a range $0\le \mathcal{I}(x_{p1};x_{p2},y_{p2})\le min(H\left(x_{p1}\right)$, $H(x_{p2},y_{p2}))$. In this case, one can also use the theory of PID to decompose the three variable mutual information into two parts and calculate the net synergy [50]

$$\begin{array}{ccc}\hfill \Delta \mathcal{I}(x_{p1};x_{p2},y_{p2})& =& \mathcal{I}(x_{p1};x_{p2},y_{p2})-\mathcal{I}(x_{p1};x_{p2})-\mathcal{I}(x_{p1};x_{p2})\hfill \\ \hfill & =& \frac{1}{2}{log}_2\left[\frac{({\sigma }_{x_{p2}}^2{\sigma }_{y_{p2}}^2-{\sigma }_{x_{p2}{y}_{p2}}^4)({\sigma }_{x_{p1}}^2{\sigma }_{x_{p2}}^2-{\sigma }_{x_{p1}{x}_{p2}}^4)({\sigma }_{x_{p1}}^2{\sigma }_{y_{p2}}^2-{\sigma }_{x_{p1}{y}_{p2}}^4)}{|{\Delta }_2|{\sigma }_{x_{p1}}^2{\sigma }_{x_{p2}}^2{\sigma }_{y_{p2}}^2}\right].\hfill \end{array}$$

(13)

For analytical calculation, we have adopted Equations (4), (5), (8), (9), (12) and (13) which contain only variance and covariance expressions, whereas we adopt numerical simulation for evaluation of the expressions given in Equations (2), (6) and (10).

At this point it is important to mention that we validate our analytical calculation by exact stochastic simulation, commonly known as stochastic simulation algorithm or Gillespie algorithm [43]. The validation signifies how much closer the system dynamics with the Gaussian statistics. Corroboration of analytical and simulation results indicate a valid consideration of linear noise approximation. In our numerical simulation, we have used ${10}^7$ trajectories. We note that the dynamics of different MAPK pathways exhibit different temporal evolution [8]. However, the average physiological time scale for the activation of MAPK pathways is ∼ 1 h [15]. Keeping this information in mind, each stochastic simulation was executed for 5000 s (∼ 1.38 h) to ensure that the system dynamics reaches steady state. Using the final value of the variables of each run we carry out the calculation of different variances and covariances to evaluate different expressions of mutual information.

4. Conclusions

To summarize, we have investigated evolutionarily conserved yeast MAPK signalling pathway. In our phenomenological model, we study two parallel MAPK signalling pathways where one signalling pathway in addition to its cognate pathway activates the non-cognate pathway through cross-talk with an emphasis to understand the change in the dynamical behaviour of the system in the presence of cross-talk at the single cell level. The model nonlinear Langevin equations have been solved under the purview of LNA to quantify the variance and the covariance associated with the different phosphorylated kinase. These quantities assist in the evaluation of mutual information (two variable and multivariate) under Gaussian channel approximation. Quantification of mutual information has been carried out with the variation of two cross-talk parameters ${\epsilon }_1$ and ${\epsilon }_2$. The two variable mutual information shows that cross-talk establishes correlation in the signal propagation among the two pathways. To represent a better insight into the directionality of the net information flow, we have defined a new dimensionless parameter (net information transduction D), which varies on a scale of $-1$ to $+1$. Depending on the sign of D, we have deciphered the fidelity of one pathway compared to the other.

We show that cross-talk generates correlated fluctuations at the population level. A minimum and a maximum degree of coordination are observed at the low and high level of cross-talk, respectively. Our analysis thus suggests that coordinated fluctuations are the causal effect of cross-talk in MAPK signalling pathways. Furthermore, we demonstrate the impact of correlated fluctuations in the reduction of variability between two different kinases using scatter plots. At the high degree of cross-talk, scatter plots show high correlation coefficient compared to the lower level of cross talk. These results together imply that cross-talk not only develops synchronisation in a cell but also reduces variability due to the development of correlation between two different phosphorylated kinase levels in a population. Depending on the number of inputs and outputs, we have identified two types of signalling motifs from the composite network. Also, quantification of multivariate mutual information allows us to calculate the net synergy associated with these two different motifs. The signal integration motif (two inputs and one output) reveals high fidelity, whereas the signal bifurcation motif (one input and two outputs) shows redundancy in information propagation.

Based on the aforesaid theoretical discussion, we suggest a satisfactory explanation about the synchronisation in the outputs - a causal effect of cross-talk in parallel MAPK signalling pathways. Nevertheless, one question apparently arises - what is the importance of such synchronisation in cellular physiology? Such functional correlation is possibly required for both the outputs to perform in a combined way to regulate several essential downstream genes. Several experimental results on MAPK cross-talk in S. cerevisiae provide interesting evidence that corroborates with our theoretical analysis. Phosphorylated Fus3 and Kss1 are both responsible for the activation of transcription factor Ste12 that regulates different downstream genes [64]. Additionally, both activated Fus3 and Hog1 assist in arresting the cell cycle in G1 phase temporarily [5,8]. Cross-talk is also highly significant for the eukaryotic cells where the promoter of TATA binding proteins is solely controlled by MAPK signalling pathways [15,69], whereas these binding proteins are essential for the expression of most nuclear genes. Also, they act as a potential vehicle for developing coordination among the multiple disparate classes of genes. Thus, coordinated signalling of MAPK pathways paves the way for TATA binding proteins to establish the association among large-scale nuclear genes. Gene regulation in S. cerevisiae is known to be controlled by more than one transcription factors that bind cooperatively at many promoter sites. This phenomenon suggests that coordinated fluctuations between the outputs of MAPK signalling pathways are necessary to express the gene product in a controlled way. It is also noticed that coordinated fluctuations among gene products are developed through transcriptional as well as translational cross-talk [25,57,58,70]. We propose that it could be more convenient for a cell to establish a functional connection among all intracellular processes if the correlation is initiated in the signalling pathway, not solely in the gene regulation stage. In fact, one interesting signature which was observed in different experiments is that cross-talk is prominent at the low concentration level that is manifested in diverse environmental cues [22,24]. Thus, in these situations, fluctuations in the cellular components are very high, and it is improbable for cells to adopt a constructive decision for survival [71]. Our results indicate that such decision making program becomes easy when correlated fluctuations among the essential proteins are successfully implemented through the bet-hedging program [22].

Overall, we suggest that synchronisation between MAPK signalling pathways is a result of cross-talk. Our analytical calculation supplemented by exact numerical simulation is a general approach and can be applied to other cross-talk pathways to quantify the strength of cross-interactions. In future, we plan to address the influence and physiological relevance of cross-talk in other network motifs. Our theoretical observations in the present work could be verified upon the quantification of phosphorylated kinase protein in a single cell using flow cytometry and time lapse microscopy [72,73,74]. These experimental approaches can be implemented to measure the amount of intra-cellular phosphorylated kinases by treating cells with external stimuli, fixing and permeabilizing cells with appropriate chemicals, and then staining with phospho-specific antibodies for different kinases. After that, one can quantify the intensity of phosphorylated kinases in individual cells of a colony. Using these data, distribution profiles of the concentration of phosphorylated kinases could be developed. These quantifiable distribution profiles could be used to quantify the mutual information.