Information theoretical study of cross-talk mediated signal transduction in MAPK pathways

Information theoretical study of cross-talk mediated signal transduction in MAPK pathways

Alok Kumar Maity amaity@ucla.edu Department of Chemistry, University of Calcutta, 92 A P C Road, Kolkata 700009, India Present address: Department of Chemistry and Biochemistry, University of California, Los Angeles, Los Angeles, CA 90095, United States of America    Pinaki Chaudhury pinakc@rediffmail.com Department of Chemistry, University of Calcutta, 92 A P C Road, Kolkata 700009, India    Suman K Banik skbanik@jcbose.ac.in Department of Chemistry, Bose Institute, 93/1 A P C Road, Kolkata 700009, India
Abstract

Biochemical networks related to 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. Our 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 and scattered plot suggest that the cross-talk generates synchronization within a cell as well as in a cellular population. Depending on the number of input and output, we identify signal integration and signal bifurcation motif that arise due to inter-pathway connectivity in the composite network. Analysis using partial information decomposition quantifies the net synergy in the information propagation through these branched pathways.

I Introduction

The decision making processes at the cellular level are initiated by some specialized signaling networks Tyson et al. (2003); Alon (2006). 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 signaling networks, few signaling networks are incorporated mostly in the eukaryotes due to the evolutionary prospect Lyons and Prasad (2012). One such signaling network is the mitogen-activated protein kinase (MAPK) pathway that plays the central role to attune with extra-cellular signal in eukaryotic cells Bardwell (2006); Ferrell (1996); Huang and Ferrell (1996); Saito (2010); Voliotis et al. (2014). Although different MAPK pathways with diverse inputs and outputs belong to a higher living species, they are sometimes interconnected through overlapping sets of signaling components. Depending on the interconnections, MAPK pathways can be classified into different groups that use one or more than one common signaling components. Moreover, as a result of cross-interaction, a single regulon regulates multiple targets in addition to its own target. Such type of signal association is defined as cross-talk. Cross-coupling in the signaling 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 Laub and Goulian (2007); Posas et al. (1996); Rowland and Deeds (2014); Siryaporn et al. (2010); Trach and Hoch (1993). In the eukaryotic system, cross-talk has been identified in numerous situationsGustin et al. (1998); Kunkel and Brooks (2002); Oeckinghaus et al. (2011). Furthermore, cross-talk and several of its variance have also been identified at different stages of gene regulation Iborra et al. (2004); Mather et al. (2013); Mauri and Klumpp (2014); Riba et al. (2014); Tsimring (2014).

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 despite the continuous development of a large number of evolutionary descent with modifications in the cellular interactions within the surviving trait. This character indicates that crosstalk might have a definite functional potential to build up synchronized cellular regulations by spending the storage energy. If this is true, how a cell balances 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 Saito (2010). Synchronization is necessary to attain natural activity but needs to maintain a threshold value. Otherwise too much synchronization may lead to physiological disorder like epilepsy Palus et al. (2001). Inter-pathway cross-talk becomes prominent due to limitation of common resources, defined as overloaded condition. However, cross-talk effect becomes faint in the underloaded condition, where level of available resources is satisfactory Cookson et al. (2011); Komili and Silver (2008); Mather et al. (2013); Mauri and Klumpp (2014); Tsimring (2014). 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 Tsimring (2014). Cross-talk also facilitates synchronization in different organs such as cardio-respiratory interaction, brain and tissues Palus et al. (2001).

To address the functionality of cross-talk, we undertake a representative cluster of networks with overlapping sets of regulatory components. To be specific, we focus on MAPK pathway, a well-studied eukaryotic signaling machinery, conserved with three kinase cascades. In S. cerevisiae, five MAPK signaling pathways are present among which only three (pheromone response, filamentous growth response and osmostress adaptation) use a common kinase protein Ste11 Bardwell (2006); Saito (2010). 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)gets activated whereas high external osmolarity activates the osmoadaption cascade (Ste11 Pbs2 Hog1) Bardwell (2006); Saito (2010). Due to intercascade correlation among the three signaling pathways, one pathway can be activated by the signal of another pathway in 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 Bardwell (2006); Bardwell et al. (2007); Komarova et al. (2005); McClean et al. (2007); Saito (2010); Ubersax and Ferrell (2007). 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 obvious queries (i) is it possible for an individual signaling pathway to convey its input signal reliably downstream without experiencing any influence from the other pathways signal? (ii) Since the inter-pathway connectivity is known not to allow the uniqueness of transduced signals - what are the physiological advantages of cross-association? (iii) Is there any participation of pathway output in the cooperative regulation of a downstream target in a synchronized manner? (iv) How is it possible for correlated pathways to keep up static as well as dynamic synchronization in a single cell environment that is prevalently stochastic in nature? (v) Does this association have any capability to control the cell-to-cell variability?

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 Bardwell (2006); Saito (2010). Both pathways get stimulated 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 stochastic differential equation. We solve the coupled set of nonlinear Langevin equations using linear noise approximation van Kampen (2011); Gardiner (2009) and calculate the auto variance of each and every kinase and covariance between two different kinases (see Sec. II and Appendix). Recent theoretical development Grima (2015) shows that linear noise approximation is not only limited for high copy number but also exact up to second moments of any chemical species participating 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. In addition, linear noise approximation 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 Elf and Ehrenberg (2003); Gardiner (2009); Maity et al. (2014, 2015); Paulsson (2004, 2005); van Kampen (2011); Grima et al. (2011); Thomas et al. (2013a, b); Grima (2015). To classify the signal transduction efficacy through two pathways in the presence of cross-association, we quantify two as well as three variable mutual information. Distributions of all kinase proteins are approximately considered Gaussian, allowing us to adopt a reduced expression of mutual information Cover and Thomas (2012); Shannon (1948). The reduced equation mainly depends on the auto variance and the covariance of the corresponding kinase. We validate our analytical calculation by exact stochastic simulation Gillespie (1976). 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 synchronization Granger (1969); Palus et al. (2001). 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 very difficult 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. We also verify inter-pathway synchronization with the help of coordinated fluctuations of stochastic trajectories of two parallel kinases. This result implies how two kinases are synchronized within a cell. To understand this phenomenon further, we investigate how much association is developed among the steady state population of these kinases from cell-to-cell. In the second subsection, we quantify three variable mutual information when both the channels of information flow work separately. Applying the theory of partial information decomposition Barrett (2015); Schneidman et al. (2003), we quantify the 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.

Ii The Model

In Fig. 1 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 pathway) and , , (Y pathway) Heinrich et al. (2002); Huang and Ferrell (1996); McClean et al. (2007); Suderman and Deeds (2013); Tănase-Nicola et al. (2006). and represent dephosphorylated and phosphorylated form of a kinase protein, respectively, and the same applies to and (here ). 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 Fig. 1 along with the cross-interaction rate parameters and . and 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 Barrett (2015); Granger (1969); Palus et al. (2001); Schneidman et al. (2003). In the present manuscript, we quantify the causal relationship in terms of mutual information. Here, the two cross-interaction parameters and play a significant role in establishing different levels of cross-talk. The parameter controls information flow from X to Y pathway (), but the parameter is responsible for Y to X pathway () 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 Bardwell et al. (2007). 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 a constant ( constant, here ). In addition, we consider a physiologically relevant parameter set for our calculation Bardwell et al. (1996); Ferrell (1996); Voliotis et al. (2014).

Figure 1: (color online) Schematic diagram of two parallel MAPK (equivalent and identical) signaling pathways (X and Y). Each pathway consists of three successively connected cascade kinases, MAPKKK (red), MAPKK (green) and MAPK (blue). The first activated kinase facilitates the activation of the second one and then the second kinase regulates the activation of the last one. Both signaling pathways are exposed to two different signals ( and ). Cross-talk is developed due to inter-pathway interactions. and are the cross-interaction parameters and the directionality of these interactions are and , respectively.

ii.1 Two variable mutual information

Adopting Shannon’s information theory Cover and Thomas (2012); Shannon (1948), we have calculated two variable mutual information between two phosphorylated kinases,

(1)

A generalised index and have been considered to represent the copy number of two different phosphorylated kinases. Similarly, and are the marginal and is the joint probability distributions associated with the corresponding kinases. For the calculation of mutual information between two kinases of X signaling pathway, we have replaced by (where )and the reverse replacement has been followed for Y signaling pathway. For the estimation of mutual information between two equivalent kinases ( and ) of the respective pathways, we have used the same formula for condition. Mutual information can also be written in the form of the entropy function. Hence, Eq. (1) can be redefined as

(2)

Here, and are individual and is total entropy of the respective kinases. In the present study, both probability distribution functions (marginal as well as joint) are approximately considered to be Gaussian. Thus, using Gaussian channel approximation Barrett (2015); Cover and Thomas (2012); Shannon (1948), Eq. (2) takes the reduced form

(3)

where and are auto variances and 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 Eq. (1). For exact or approximate Gaussian distribution, one can reduce Eq. (1) to Eq. (3) applying Gaussian channel approximation. However, for systems with non-Gaussian distribution, one can still use Eq. (1) with proper analytical expressions of probability distribution functions that may contribute expressions of higher moments in Eq. (3).

In the present work, all expressions of two variable mutual information are calculated using Eq. (3). The analytical results are then validated by evaluating probability distribution functions (Eq. (1)) using exact numerical simulation Gillespie (1976). In our numerical simulation we have used trajectories and smallest bin size of . The two variable mutual information value is bounded with in a scale . To quantify the association between two equivalent kinases, we have used Pearson’s correlation coefficient () Whitlock and Schluter (2009)

(4)

.

ii.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

(5)

where and are the joint distribution functions of the corresponding components. On the other hand, is the marginal distribution of phosphorylated kinase. One can also write Eq. (5) in terms of the respective entropy

(6)

Similarly, using Gaussian approximation Barrett (2015); Cover and Thomas (2012); Shannon (1948), one can reduce Eq. (6) into the following form

(7)

with

Here, the magnitude of three variable mutual information is bounded within a scale . Using partial information decomposition formalism, the three variable mutual information can be decomposed into two parts Barrett (2015); Schneidman et al. (2003). As a result, the net synergy expression becomes

(8)

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

(9)

and the entropy representation of Eq. (9) is

(10)

Using Gaussian approximation Barrett (2015); Cover and Thomas (2012); Shannon (1948) Eq. (10) becomes

(11)

with

In addition, the three variable mutual information value is bounded within a range , . In this case, one can also use the theory of partial information decomposition to decompose the three variable mutual information into two parts and calculate the net synergy Schneidman et al. (2003)

(12)

For analytical calculation, we have adopted Eqs. (3, 4, 7, 8, 11 and 12) which contain only auto variance and covariance expressions, whereas we adopt numerical simulation for evaluation of the expressions given in Eqs. (1, 5 and 9). 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 Gillespie (1976). The validation signifies how much closer the system dynamics with the Gaussian statistics. In the following section, corroboration of analytical and simulation results indicate a valid consideration of Gaussian approximation.

Description Reaction Propensity function Rate constant
Activation of molecules s
Deactivation of s
Activation of molecules s
Deactivation of s
Activation of molecules s
Activation of molecules s
Deactivation of s
Activation of molecules s
Activation of molecules s
Deactivation of s
Activation of molecules s
Deactivation of s
Activation of molecules s
Deactivation of s
Table 1: Reactions and corresponding parameter values for the MAPK network motif of S. cerevisiae Bardwell et al. (1996); Ferrell (1996); Voliotis et al. (2014), related to faster relaxation rate.222Other Parameters are molecules/cell, molecules/cell, molecules/cell, molecules/cell, molecules/cell, molecules/cell and molecules/cell. The kinetic schemes adopted in the present work follows the model of Heinrich et al Heinrich et al. (2002).
Description Reaction Propensity function Rate constant
Activation of molecules s
Deactivation of s
Activation of molecules s
Deactivation of s
Activation of molecules s
Activation of molecules s
Deactivation of s
Activation of molecules s
Activation of molecules s
Deactivation of s
Activation of molecules s
Deactivation of s
Activation of molecules s
Deactivation of s
Table 2: Reactions and corresponding parameter values for the MAPK network motif of S. cerevisiae Bardwell et al. (1996); Ferrell (1996); Voliotis et al. (2014), related to slower relaxation rate.444Other Parameters are molecules/cell, molecules/cell, molecules/cell, molecules/cell, molecules/cell, molecules/cell and molecules/cell. The kinetic schemes adopted in the present work follows the model of Heinrich et al Heinrich et al. (2002).

Iii Results and Discussion

iii.1 Two variable mutual information

The parameters and control the signaling channel X and Y, respectively. In Fig. 2, we show the mutual information profile as a function of for two different sets of parameters while keeping constant. Fig. 2A shows that mutual information between and kinases decays with the increment of . Augmentation of includes a competition between and to phosphorylate the kinase. During phosphorylation, mutual association is originated, and signal transduction is ensued. Thus, for the low value of , 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 value due to maximum phosphorylation contribution of . In Fig. 2B, mutual information between and is plotted, which shows a constant value as a function of . This happens as has no influence in the alteration of mutual information. The same logic is applicable to the mutual information between and shown in Fig. 2C. In Fig. 2D, mutual information between and increases as a function of , as is only responsible for establishing the cross-talk between and . This result implies that with the enhancement of cross-talk the process of signal integration through increases. The same profiles can be generated as a function of , while keeping fixed. These results together indicate that and depend on , whereas and depend on . In Fig. 2A-2D, the dotted lines drawn for slower relaxation rate (see Table II) always maintains a lower mutual information value compared to the solid lines drawn for faster relaxation rate (see Table I). Relaxation rates of the corresponding kinases i.e., , , , , and are , , , , and , respectively, where the angular bracket 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 Maity et al. (2015). For solid line, we consider higher degradation rate for and ( and ) compared to (). Thus faster relaxation rates are attained under this condition with high information propagation capacity.

Figure 2: (color online) Two variable mutual information and net information transduction as a function of cross-talk parameter. A, B, C, D, E and F - Two variable mutual information profiles , , , , and as a function of cross-interaction parameter for a fixed value of . In all figures, solid (with open circle) and dotted (with open diamond) lines are generated using faster (Table I) and slower (Table II) relaxation rate parameters, respectively. The symbols are generated using stochastic simulation algorithm Gillespie (1976) and the lines are due to theoretical calculation. G and H - Bar diagram of two variable mutual information of three parallel cascade kinases under an equivalent cross-talk condition () for faster (Table I) and slower (Table II) relaxation rate parameters, respectively. I - Net information transduction as a function of cross-interaction parameter for a fixed value of . The solid (with open circle) and the dotted (with open diamond) lines are due to faster (Table I) and slower (Table II) relaxation rate parameters, respectively. The figure indicates data collapse for two relaxation rate parameters. The symbols are generated using stochastic simulation algorithm Gillespie (1976) and the lines are obtained from theoretical calculation. J - 2d-surface plot of net information transduction as a function of two cross-talk parameters and for faster (Table I) relaxation rate parameters.

Next, we quantify mutual information between two parallel kinases ( and , with ) of the two equivalent interacting MAPK pathways. The inter pathway coupling is unidirectional when either or is zero but is bidirectional when both are non-zero. In this situation, both variables ( and ) do not interact with each 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 and , as these are uncorrelated. In Fig. 2E, we show mutual information in between and as a function of keeping fixed. The profile shows an increasing trend as cross-talk parameter increases. Similarly, in Fig. 2F, mutual information between and is shown with a similar trend as in Fig. 2E. Interestingly, for faster relaxation time scale, mutual information between similar cascade kinases increases while moving from second ( and ) to third ( and ) cascade. On the other hand, an opposite trend is observed for slower relaxation time scale. This characteristic trend is further shown in Fig. 2G and 2H 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 Fig. 2E, mutual information is high for slower relaxation rate than the faster one, as slower rate parameters yield high level of and which in turn incorporate extra fluctuations that help to increase mutual association. A similar result is also observed in Fig. 2F. Identical mutual information profiles of and can be generated as function of keeping fixed. These results suggest that both the cross-talk parameters and contribute equally to the development of an association between two parallel pathways.

Both the mutual information between and , and are capable of providing a satisfactory explanation of enhancement of cross-talk with the increment of inter pathway interaction parameters ( and ). Under equivalent interactions condition (), each pathway shares its information with other in an equal extent and is quantified not only by and but also by and . However, characterization of the direction of information transduction is difficult under unequal condition (). Except the equivalent condition () where the net information ( ) 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 () using and as

(13)

The above expression implies that it is maximum () when is zero, i.e., no information propagation from to (). It is minimum () when is zero, which specifies zero information propagation from to (). In Fig. 2I, we show the profile of as a function of while keeping fixed, where the value of changes from negative to positive as increases. It suggests that at low , 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 . On the other hand, at high , the direction of net information propagation is from Y to X due to reverse situation and generates a positive value. The opposite scenario can be observed if one generates the profile of as a function of for fixed In this connection, it is important to mention that both the relaxation time scale limits generate a similar profile of . As a result, both the profiles of exhibit data collapse when depicted as a function of for fixed (Fig. 2I) or vice versa This observation indicates that normalised profiles of are independent of relaxation time scales. In Fig. 2J, we also show a 2d-surface plot of as a function of both and for faster relaxation time scale (Table I). The surface plot indicates zero (or near to zero) value of along the diagonal region (). However, the off diagonal region is positive for and negative for .

Figure 3: (color online) 2d-surface plots of two variable mutual information, stochastic time trajectories and scattered plots Gillespie (1976). A and B - 2d-surface plot of two variable mutual information and as a function of two cross-talk parameters and for faster (Table I) relaxation rate parameters. In both figures I, II, III and IV correspond to four different values of and . C - Stochastic time trajectories and steady state population of two parallel kinases for four different sets of and . For CI, CII, CIII and CIV we have used and , and and , respectively. In each scattered plot, represents analytical value of Pearson’s correlation coefficient.

In Fig. 3A and 3B, we show two 2d-surface plots of mutual information between and , and kinases, respectively, as a function of two cross-interaction parameters and under faster relaxation time scale (Table I). Both figures show maximum mutual information at high values of the two parameters. Since, and are equally responsible for developing the cross-correlation between two pathways, one can check the effect of maximization of mutual information by the increment of any of these two parameters. Although we can quantify the cross-talk with the help of two variable mutual information, and , 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 in Fig. 3C. Dynamical correlation is applicable to characterize the stochastic trajectories of two variables in a single cell. If sufficient association between two trajectories exist, then correlated fluctuation is observed i.e., one trajectory closely follows the other. Otherwise, an uncorrelated fluctuations (trajectories do not follow each other) are observed in the absence of cross-talk. In Fig. 3C, we show stochastic time series of different kinases under different conditions (four different sets of and parameters have been used - mentioned as I, II, III and IV in Fig. 3A and 3B). These time series are generated from a single run of stochastic simulation. We show that stochastic time series exhibit correlated fluctuations at high and but are uncorrelated time series are observed at low and . Static correlation implies cell-to-cell population variability between two parallel cascade kinases. In Fig. 3C, we show scattered plots for different kinases under different cross-interaction conditions. Here, each symbol generated from a single realization of stochastic trajectory represents copies of phosphorylated kinase at steady state. For plots with high and values, most of the symbols are disposed diagonally in a narrow strip but for low and , symbols are distributed in a much larger volume of space. Mutual association incorporates correlated variation in populations which is along the diagonal axis. An uncorrelated variation along the off-diagonal axis represents no such association. Therefore, cross-talk enhances correlated variation among populations of the kinases. These results imply that cross-talk not only develops association within a cell but is also capable of generating cell-to-cell association that assists in successful development of a significant robust adaptation machinery as observed in the bet-hedging program under diverse environmental conditions Tsimring (2014). In addition, similar behavior can be observed under slower relaxation time scale (Table II) as shown in Fig. 4. The primary difference between the nature of correlation between and are visible from Fig. 3 and Fig. 4. In Fig. 3, the correlation between is always higher than for all four conditions. On the other hand, in Fig. 4, it shows an opposite trend.

Figure 4: (color online) 2d-surface plots of two variable mutual information, stochastic time trajectories and scattered plots Gillespie (1976). A and B - 2d-surface plot of two variable mutual information and as a function of two cross-talk parameters and for slower (Table II) relaxation rate parameters. In both figures I, II, III and IV correspond to four different values of and . C - Stochastic time trajectories and steady state population of two parallel kinases for four different sets of and . For CI, CII, CIII and CIV we have used and , and and , respectively. In each scattered plot, represents analytical value of Pearson’s correlation coefficient.

iii.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, two types of branched pathways have been considered. One is two inputs ( and ) and one output ( or ) motif where two input signals are integrated into a single output. The other is one input ( or ) and two outputs ( and ) motif where the input signal is bifurcated into two output signals. 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 as approximately Gaussian, we adopt multivariate mutual information theory to analytically estimate three variable mutual information Barrett (2015); Schneidman et al. (2003). Each branched motif consists of two signal propagating channels that work together. It is thus interesting to investigate whether these signaling 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 and is evaluated using the theory of partial information decomposition 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) Barrett (2015); Schneidman et al. (2003).

Figure 5: (color online) Three variable mutual information as a function of cross-talk parameter. A and B - Three variable mutual information (A) and net synergy (B) for signal integration motif. Schematic diagram of signal integration motif in MAPK composite network (see inset in A). C and D - Three variable mutual information (C) and net synergy (D) for signal bifurcation motif. Schematic diagram of signal bifurcation motif in MAPK composite network (see inset in C). All the figures are drawn as a function of cross-interaction parameter for a fixed value of . Here solid (with open circle) and dotted (with open diamond) lines are drawn for faster (Table I) and slower (Table II) relaxation rate parameters, respectively. The symbols are generated using stochastic simulation algorithm Gillespie (1976) and the lines are obtained from theoretical calculation.

In Fig. 5A, we show mutual information, of two inputs and one output model as a function of for a fixed value of . The profile shows a bifunctional behaviour with the increment of ; initially it decreases up to a certain value of , and then it increases. At low , a minimal amount of signal is propagated from to . Consequently, the motif reduces to a single input-output motif and the motif regains its native form due to the significant contribution of . In Fig. 2A, we show that two variable mutual information between and decreases with the increment of . Similar situation arises in Fig. 5A for low value of . On the contrary, increases with the increment of (Fig. 2D). Thus, two opposing effects work together to generate the convex profile. In Fig. 5B, we plot net synergy of the motif as a function of for a fixed value of 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 . This result implies that an integrating signaling 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 signaling network motif as observed in V. harveyi quorum-sensing circuit Waters and Bassler (2005); Mehta et al. (2009)

In Fig. 5C, we show mutual information of one input and two outputs motif with the increment of for a fixed value of . The mutual information value decreases with since propagation of information from to is only inhibited by the cross-interaction. However, does not have any influence in information propagation from to and remains unaltered. Thus, three variable mutual information profile follows a decreasing trend. Fig. 5D shows decreasing trend of net synergy profile as a function of for a fixed value of . Importantly, for this motif negative values of net synergy are observed irrespective of the value of . 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 signaling model reduces mutual information, it has a biological significance of the activation of multiple signaling channels in the presence of a single input as identified in the chemotaxis system of E. coli Armitage (1999). In all figures (Fig. 5A-5D), the solid lines are plotted for faster relaxation rate constants (Table I) of and and the dotted lines are for slower relaxation rate constants (Table II).

Iv Conclusion

To summarize, we have investigated evolutionarily conserved yeast MAPK signaling pathway. In our phenomenological model, we study two parallel MAPK signaling pathways where one signaling pathway in addition to its own activates the other pathway through cross-talk with an emphasis to understand the change in the dynamical behavior 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 linear noise approximation to quantify the auto variance and the covariance associated with the different phosphorylated kinase. These quantities assist in the evaluation of mutual information (two variable and three variable) under Gaussian channel approximation. Quantification of mutual information has been carried out with the variation of two cross-talk parameters and . The two variable mutual information shows that cross-talk establishes an association of 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 ), which varies on a scale of to . Depending on the sign of , 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 low and high level of cross-talk, respectively. Our analysis thus suggests that coordinated fluctuations are the causal effect of cross-talk in MAPK signaling pathways. Furthermore, we have evaluated the impact of correlated association in maintaining cell-to-cell population variability of kinases using scattered plots. At the high degree of cross-talk, scattered plots show high correlation coefficient compared to lower level of cross talk. These results together imply that cross-talk not only develops synchronization in a cell but also among the cellular population. Depending on the number of inputs and outputs, we have identified two types of signaling motifs from the composite network. In addition, quantification of three variable 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 synchronization in the outputs - a causal effect of cross-talk in parallel MAPK signaling pathways. Nevertheless, one question obviously arises - what is the importance of such synchronization 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 corroborate with our theoretical analysis. Phosphorylated Fus3 and Kss1 are both responsible for the activation of transcription factor Ste12 that regulates different downstream genes Suderman and Deeds (2013). Additionally, both activated Fus3 and Hog1 assist in arresting the cell cycle in G1 phase temporarily Bardwell (2006); Saito (2010). Cross-talk is also highly significant for the eukaryotic cells where the promoter of TATA binding proteins is solely controlled by MAPK signaling pathways Gustin et al. (1998); White and Sharrocks (2010), 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 signaling of MAPK pathways paves the way for TATA binding proteins to establish 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 signaling 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 Dunlop et al. (2008); Komili and Silver (2008); Munsky et al. (2012); Stewart-Ornstein et al. (2012). 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 signaling pathway, not solely in the gene regulation stage. In fact, one interesting signature which has been observed in different experiments is that cross-talk is prominent at low concentration level that is manifested in diverse environmental cues Cookson et al. (2011); Tsimring (2014). Thus, in these situations, fluctuations in the cellular components are very high and it is very unlikely for cells to adopt a constructive decision for survival Junttila et al. (2008). 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 Tsimring (2014).

Overall, we suggest that synchronization between MAPK signaling 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 Gao and Stock (2013); Schulz et al. (2012); Yaginuma et al. (2014). 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 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.

Acknowledgements.
We thank Debi Banerjee, Sandip Kar and Jayanta Mukhopadhyay for critical reading and for making constructive suggestions. AKM acknowledges University Grants Commission for a research fellowship (UGC/776/JRF(Sc)). SKB acknowledges financial support from Council of Scientific and Industrial Research, Government of India (01(2771)/14/EMR-II) and Bose Institute (Institutional Programme VI - Development of Systems Biology), Kolkata.

Appendix A Calculation of auto variance and covariance

The MAPK network motif shown in Fig. 1 is explicated through stochastic Langevin equations. Each pathway (X or Y) is activated by the initiation of an extra cellular signal ( or). When the first cascade kinase is activated, it regulates the activation of downstream kinases of the same as well as the parallel pathway through cross-interaction. Once activated, the second kinase regulates the activation of the last kinase. The activated and the deactivated states can be identified in terms of phosphorylated ( and ) and dephosphorylated ( and ) forms of each kinase (), respectively. To construct the theoretical model of the composite MAPK network motif, we have considered the total population (phosphorylated and dephosphorylated form) of all kinases to be a constant ( constant; ). Thus, for X pathway the stochastic differential equations for ,