Original article Biomath 3 (2014), 1404131, 1–11 B f Volume ░, Number ░, 20░░ BIOMATH ISSN 1314-684X Editor–in–Chief: Roumen Anguelov B f BIOMATH h t t p : / / w w w . b i o m a t h f o r u m . o r g / b i o m a t h / i n d e x . p h p / b i o m a t h / Biomath Forum On Nonlinear Dynamics of the STAT5a Signaling Protein Elena Nikolova, Ivan Jordanov, Nikolay Vitanov Institute of Mechanics Bulgarian Academy of Sciences, Sofia, Bulgaria Emails: elena@imbm.bas.bg, i jordanov@email.bg, vitanov@imbm.bas.bg Received: 2 September 2013, accepted: 13 April 2014, published: 29 May 2014 Abstract—In this paper we model dynamics of cross talk between MEK/ERK and JAK/STAT sig- naling pathways by means of nonlinear ordinary differential equations. The considered system of four ordinary differential equations is reduced to one ordinary differential equation, representing dy- namics of the phosphorylated STAT5a signaling protein. We show that the diffusion together with the corresponding biochemical reactions is likely to play a critical role in governing the dynamical behavior of the considered signaling protein. By the modified method of simplest equation to the described reaction-diffusion equation we obtain an analytical solution which explains drop and jump propagation of the STAT5a protein concentration. Keywords-STAT5a signaling protein; PDEs; Modi- fied method of simplest equation; analytical solution; drop and jump propagation; I. Introduction The idea to modeling proteins interactions in the intracellular environment by spatial-temporal systems is based on the following circumstance: While the dynamics of the interacting proteins and their molecular pathways and networks can be described by reaction models, the heterogeneous distributions of the protein concentrations are not taken into consideration. It is proved, however, that similar inhomogeneous protein distributions in the form of cellular jump or drop propagation play an important role in the control of the main processes in the cell. In this way the cellular com- plexity appears to be space-temporal, expressed mathematically by reaction-diffusion models. In fact, the traditional approximation scheme of a well-stirred reactor is most used and somewhat plausible simplification. However, the concentra- tion gradients of cell enzymes that modulate signal transduction refute this simplification [1–4]. The role of diffusion in reaction-diffusion systems of the cell becomes significant when reactions are relatively faster (but not so very) than diffusion rates. Sometimes the term ‘molecular crowding’ is used to denote more specific type of spatial distribution [5, 6]. The biochemical essence of this phenomenon is based on the fact that the state of phosphorylation of target molecules with spatially separated membrane-localized protein ki- nases and cytosolic phosphatases depends essen- tially on diffusion. Then the very high protein density in the intracellular space, or so called molecular crowding, mentioned above can enhance the spatial effect. Consequently, molecular crowd- Citation: Elena Nikolova, Ivan Jordanov, Nikolay Vitanov, On Nonlinear Dynamics of the STAT5a Signaling Protein, Biomath 3 (2014), 1404131, http://dx.doi.org/10.11145/j.biomath.2014.04.131 Page 1 of 11 http://www.biomathforum.org/biomath/index.php/biomath http://dx.doi.org/10.11145/j.biomath.2014.04.131 E Nikolova et al., On Nonlinear Dynamics of the STAT5a Signaling Protein ing can also alter protein activities and break down classical reaction kinetics [5]. Therefore here we introduce reaction-diffusion modeling and its computational tools to concrete example of ERK (Extracellular signal-Regulated Kinase) and STAT (Signal Transducers and Activators of Transcrip- tion) protein interaction. In [7] the cross talk between the MEK/ERK and JAK/STAT signaling pathways has been modeled in a form of a system of nonlinear ordinary differential equations for the protein concentrations (JAK is an abbreviation of Janus Kinase, and MEK is an abbreviation of MAPK/ERK Kinase). Next, in the same paper the spatial modeling of the mentioned interaction has been performed, by introducing an appropriate diffusion-reaction scheme. By the accomplished stability analysis in [7] the authors concluded that in terms of the Turing bifurcation in ERK and STAT dynamical model the mentioned above crowding effect can be interpreted as a process of stabilization of the dissipative structures inherent to the considered intracellular process. II. Themethod of the simplest equation There are many approaches for obtaining an- alytic solutions of nonlinear partial differential equations [8–10]. In this paper we will use the modified method of the simplest equation. The method is a modified version of the method of the simplest equation, created by N. A. Kudryashov [11, 12] that is based on the fact that after ap- plication of an appropriate ansatz a large class of nonlinear PDEs can be reduced to nonlinear ODEs of the kind (P means polynomial) P(F(ξ), dF dξ , d2 F dξ2 , ...) = 0 (1) and for some equations of the kind (1) particular solutions can be obtained which are finite series F(ξ) = n∑ i=0 ai[Φ(ξ)] i (2) constructed by the solution Φ(ξ) of a simpler equation referred to as the simplest equation. The simplest equation can be the equation of Bernoulli, equation of Riccati, etc. The substitution of (2) in (1) leads to the polynomial equation P = σ0 + σ1Φ + σ2Φ 2 + ... + σrΦ r = 0 (3) where the coefficients σi, i = 0, 1, ..., r depend on the parameters of the equation and on the parameters of the solutions. Equating all these coefficients to zero, i.e., by setting σi = 0, i = 0, 1, ..., r (4) one obtains a system of nonlinear algebraic equa- tions. Each solution of this system leads to a solution of kind (2) of (1). In order to obtain a non-trivial solution by the above method we have to ensure that σr contains at least two terms. To do this within the scope of the modified method of the simplest equation we have to balance the highest powers of Φ that are obtained from the different terms of the solved equation of kind (1). As a result of this we obtain an additional equation between some of the parameters of the equation and the solution. This equation is called a balance equation [13,14]. We note that the method of the simplest equation and its modified version are closely connected to the problem for obtaining meromorphic solutions of nonlinear partial dif- ferential equations [15, 16]. By the methodology described in [15, 16] one can obtain other inter- esting classes of solutions of nonlinear PDEs such as rational solutions for example. In addition we stress that by means of the traveling wave ansatz one reduces the nonlinear PDE to a nonlinear ODE and after this if an appropriate simplest ODE exists then a particular solution can be obtained that usually depends on as many parameters of the problem as possible. In many cases such particular solutions are among the few possible exact analytic solutions of the studied nonlinear PDE. III. The interaction between ERK and STAT signaling pathways A. A reaction model and its reduction The specificity of biological responses is often achieved in a combinatorial fashion through the Biomath 3 (2014), 1404131, http://dx.doi.org/10.11145/j.biomath.2014.04.131 Page 2 of 11 http://dx.doi.org/10.11145/j.biomath.2014.04.131 E Nikolova et al., On Nonlinear Dynamics of the STAT5a Signaling Protein concerted interaction of signaling pathways [17]. Here we will discuss the eventual interaction be- tween the units of MEK/ERK and JAK/STAT sig- naling pathways. In cell signaling the pathways are understood as networks of recurrent biochemical reactions, by which the information transmission (in a signal form) is accomplished. The distur- bance of the intracellular signal transmission from the membrane receptors to the nucleus genes is assumed as a general reason for a cancer diseases progress. In [17] the authors discussed the basic features of interaction domains, and suggested that rather simple binary interactions can be used in sophisticated ways to generate complex cellular responses. Moreover, in [18], the protein STAT was found to play important roles in numerous cellular processes including immune responses, cell growth and differentiation, cell survival and apoptosis, and oncogenesis. The STAT target genes include SOCS/CIS, a class of inhibitory proteins that interfere with STAT signaling through several mechanisms. (SOCS is an abbreviation of Sup- pressor Of Cytokine Signaling and CIS means Cytokine Inducible SH2 domain containing). The protein SOCS/CIS can inhibit the STAT phospho- rylation and block the access of STAT to recep- tors of JAKs or both [19]. On the other hand, SOCS-3 can bind to sequester the Ras-GAP proten (The Ras proteins are members of the MEK/ERK- pathway)[20]. However, this is not only the way of interaction between STAT and ERK pathways. In [21–22], the authors suggested that the STAT5 functional capacity of binding DNA could be in- fluenced by the Mitogen-Activated Protein Kinase (MAPK)-pathway. Later on, in [23] the inter- actions between STAT5a and ERK1 (or ERK2) signaling proteins was considered. A simple bio- chemical diagram of this interaction is presented in Fig. 1. According to Fig. 1 in [7] the following system of ordinary differential equations for the kinetics of STAT5a/S phosphorylation and ERK activation was constructed: de1 dt = −k1e1 s1 + k2e2, d s1 dt = −k1e1 s1 + k3 s2 + I de2 dt = k1e1 s1 − k2e2, d s2 dt = k1e1 s1 − k3 s2 − I (5) Fig. 1. Biochemical diagram of the interaction between STAT5a and ERK proteins. Here e1, e2, s1 and s2 are state variables, representing concentrations of ERK-inactive, ERK-active, STAT5a-non-phosphorylated and STAT5a-phosphorylated proteins respectively. The following initial conditions are put on the system (5): e1(0) = 10 −3, s1(0) = 10 −2, e2(0) = 0, s2(0) = 0 They correspond to the initial concentrations of the above–mentioned proteins in mM units [27]. Moreover, k1 is proportional to the frequency of collisions of ERK and STAT5a protein molecules and presents a rate constant of reactions of as- sociations; k2 and k3 are constants of exponen- tial growths and disintegration; I is an inhibitor source respectively. The source I inhibits the phos- phorylation of non-phosphorylated STAT5a. More concrete interpretation of the inhibitor I can be given by considering the role of thee SOCS pro- teins in linking the JAK/STAT pathway. Biological responses elicited by the JAK/STAT pathway are modulated by inhibition of JAK (and respective attenuation of STAT) by a member of the SOCS proteins [24, 25]. In [7] the inhibitor I is presented mathematically in the following manner: Biomath 3 (2014), 1404131, http://dx.doi.org/10.11145/j.biomath.2014.04.131 Page 3 of 11 http://dx.doi.org/10.11145/j.biomath.2014.04.131 E Nikolova et al., On Nonlinear Dynamics of the STAT5a Signaling Protein I = k0Σ (6) where Σ is a constant concentration of SOCS proteins and k0 is a reaction rate constant of inhibition respectively. It is clear that the considered interaction between ERK and STAT pathways can occur only if Σ is sufficiently small, i.e. phosphorylation of the protein STAT5a exists. Therefore we assume the concentration of SOCS proteins to be sufficiently small. In addition in [7] the authors reduced (5) to the following two-dimensional system: de2 dt = k1 E0S 0 − (k1S 0 + k2)e2 − k1 E0 s2 + k1e2 s2 d s2 dt = −k0Σ + k1 E0S 0 − k1S 0e2 (7) −(k1 E0 + k3)s2 + k1e2 s2 taking into account that only two equations of the four ones are independent. This means that between the concentrations e1. e2, s1 and s2 the following relations exist: e1 + e2 = E0, s1 + s2 = S 0 (8) where E0 and S 0 are initial values of the sums of corresponding concentrations of inactive and active ERKs and non-phosphorylated and phos- phorylated STATs. We assume for (7) that the inequality e2 << s2 holds, in view of the fact that the amount of ERK molecules is essentially smaller than the amount of STAT5a ones [21–23]. In this way the inequality e2/s2 = � << 1 will be valid, and (7) can be written as: � de2 dt = k1 E0S 0 − (k1S 0 + k2)e2 − k1 E0 s2 + k1e2 s2 d s2 dt = −k0Σ + k1 E0S 0 − k1S 0e2 (9) −(k1 E0 + k3)s2 + k1e2 s2 Next, in accordance with the QSSA (Quasi- Steady-State Approximation) theorem [26] we consider the first equation of (9) to be linear with respect to e2 and we treat it as an attached system [26], i.e. we take into consideration that e2 is a sufficiently small ”constant“. Further according to the requirements of the QSSA theorem we prove that the attached system has a stable steady state (then well-known Lyapunov definition of stability is satisfied). After replacing the steady state value e02 = k1 E0(S 0 − s2) k1(S 0 − s2) + k2 > 0 (10) in the second equation of (9) (the degenerate system), the quasi-stationary approximation of (7) (or (5)) is obtained in the form d s2 dt = k21 E0 s 2 2 − 2k 2 1 E0S 0 s2 + k 2 1 E0S 2 0 k1 s2 − k1S 0 − k2 − (11) −(k1 E0 + k3)s2 + k1 E0S 0 − k0Σ Next, for our convenience we rewrite (11) in the following manner: d s dt = As2 − Bs + C Ds − E − F s + G (12) where we denote the variable s2 only by s. According to (11) the new coefficients A = k21 E0, B = 2k 2 1 E0S 0, C = k 2 1 E0S 2 0, D = k1, E = k1S 0 + k2, F = k1 E0 + k3, (13) G = k1 E0S 0 − k0Σ are positive because we assumed that the initial concentration of the STAT5a protein is large, the ERK concentration is smaller than the STAT concentration, but bigger than the concentration of the SOCS proteins, and k0, k1, k2 and k3 are positive constants. Biomath 3 (2014), 1404131, http://dx.doi.org/10.11145/j.biomath.2014.04.131 Page 4 of 11 http://dx.doi.org/10.11145/j.biomath.2014.04.131 E Nikolova et al., On Nonlinear Dynamics of the STAT5a Signaling Protein B. Stability analysis of the steady state solution of (12) The fixed points of (12) can be found by solving the nonlinear equation As2 − Bs + C Ds − E − F s + G = 0 (14) They correspond to the following stationary concentrations s0: s01,2 = −(EF + GD − B) 2(A − DF) (1 ∓ (15) ∓ √ 1 − 4(A − F D)(C − EG) (EF + GD − B)2 ) at Ds − E , 0, which are positive, because A − DF < 0, (16) and EF + GD − B > 0, C − EG < 0, (17) if the SOCS concentration Σ is sufficiently smaller than the total ERK (E0) and STAT (S 0) concentrations as we assumed in the previous paragraph. According to the nonlinear character of the first term of (14) we can present it in the following form: As2−Bs+C Ds−E = (Ds−E)(C1 s+C2 ) Ds−E = C1 s + C2 (18) The new coefficients C1 and C2 can be found by the equality As2 − Bs + C = DC1 s 2 + (DC2 − EC1)s − EC2 (19) or C1 = A D , C2 = − C E (20) Moreover, the following relationship among the coefficients of the basic form exists too. B = CD E + AE D (21) We replace the coefficients (20) in (18) and put it in (14). In this way we obtain the following steady state of (12): s03 = D(C − EG) E(A − DF) (22) which is positive in view of (16) and (17). How- ever, it will be valid only if (21) is satisfied. Taking into account (13), the last condition will hold if the coefficients k1 or k2 are sufficiently small (almost approaching zero), i.e. associations between ERK and STAT5a protein molecules or inactivation of active ERK molecules are very slow processes. Let us now analyze the stability of this equi- librium. For the purpose we postulate the substi- tution s = s0 + η, where η is a small variation (perturbation) from the steady state value s0 The corresponding variational equation has the form: dη dt = wη (23) where w = (2As0 − B)(Ds0 − E) − D(As20 − Bs0 + C) (Ds0 − E)2 − F (24) It is easy to show that for s0 = s01,2 the coefficient w will have negative value. This is a sufficient condition for verification of the asymptotically stable character of the corresponding equilibrium states. On contrary, for s0 = s03, w will have positive value, i.e. the corresponding equilibrium state will be unstable. Therefore we proved that the smallest (s01) and biggest (s 0 2) steady states are stable, and the intermediate steady state s03 is unstable, respectively. Biomath 3 (2014), 1404131, http://dx.doi.org/10.11145/j.biomath.2014.04.131 Page 5 of 11 http://dx.doi.org/10.11145/j.biomath.2014.04.131 E Nikolova et al., On Nonlinear Dynamics of the STAT5a Signaling Protein Finally, the results in this section can be summarized as follows: 1) The dynamical behavior of the interaction between ERK and STAT5a proteins (in particular, between MEK/ERK and JAK/STAT signaling pathways) near its quasi-stationary state determines only by the behavior of the STAT5a signaling protein; 2) Stabilization of the considered process in time will observe at smallest (about 10−3mM ) and biggest (over 10−2mM) STAT5a concentrations. Destabilization of the same process will occur at intermediate steady state concentrations of the STAT5a protein, which will lead to initiation of some pathological process (for example, cancerogenesis). IV. Space-temporalmodeling of the STAT5a signaling protein A. A reaction-diffusion model of the STAT5a sig- naling protein. Stability analysis of the inhomoge- neous distribution of the STAT5a concentration. We introduce systems with distributed variables (reaction-diffusion models) when the connections between neighbor points of the space are taken into account by the diffusion law of molecular motion from the higher to lower concentrations. Here we take into account the reaction-diffusion effect, described in the introduction of this article to the reaction model of the STAT5a signaling protein. As a result we obtain the following one-dimensional model with dis- tributed parameters ∂s ∂t = As2 − Bs + C Ds − E − F s + G + Ds ∂2 s ∂r2 (25) where r is the space coordinate from the cell membrane to the nucleus, and Ds is the diffusion coefficient of the STAT5a concentration. In order to analyze qualitatively and solve quantitatively (25), it is necessary to fix initial and boundary conditions for the function s(r, t) in the form: s(r, 0) = 0, s|r=0r=l = 0 where l is the distance between the membrane and nucleus. Next, we consider the equation (25) and search for solutions of the kind of traveling waves: s = s(ξ) = s(r−vt) , where v is the velocity of the wave. In this way, by introducing the traveling- wave coordinate ξ = r − vt, the equation (25) transforms to the following ordinary differential equation of second order: Ds s ′′ + vs′ + As2 − Bs + C Ds − E − F s + G = 0 (26) and next - in the following system of ordinary differential equations of first order: s′ = y (27) y′ = 1 Ds (−vx − As2 − Bs + C Ds − E + F s − G) where ′ denotes d/dξ. By analogy with the previous paragraph the fixed points of (27) can be found by solving the equation (14). Next we will analyze stability of the fixed points of (27) (in particular of (26) or (25)) in the phase plane (s, s’). For the purpose we introduce the substitutions s = s0 + $ and y = y0 + ζ, where $ and ζ are small variations (perturbations) from the inhomogeneous equilibrium (s0, s′0). We substitute the last expressions in (27), and obtain the following variational equations: d$ dξ = ζ dζ dξ = − w Ds $− v Ds ζ where w is a coefficient, presented by (24). The corresponding characteristic equation has the form: µ2 + ςµ + τ = 0 (28) where Biomath 3 (2014), 1404131, http://dx.doi.org/10.11145/j.biomath.2014.04.131 Page 6 of 11 http://dx.doi.org/10.11145/j.biomath.2014.04.131 E Nikolova et al., On Nonlinear Dynamics of the STAT5a Signaling Protein ς = v Ds , τ = w Ds (29) are Routh-Hurwitz coefficients. According to the RouthHurwitz stability criterion they must have positive values to ensure the stable character of the corresponding steady state. Here the diffusion coefficient Ds must be positive and the wave velocity v must be different from 0. Thereby the stability of inhomogeneous steady state of (27) (or (26)) is determined by signs in front of the coefficients v and w. We can observe a stable steady state solution of (27) if v > 0. If v < 0 the corresponding steady state will be always unstable. Thus, at v > 0 the smallest (s01) and biggest (s 0 2) steady states are unstable in view of the negative values of w, and the intermediate steady state s03 is stable because w is positive. Finally, we can conclude that the average steady state concentrations of the STAT5a protein are indicative for the dissipative structures existence, but too low and too high ones are not. In this sense, the well-known crowding effect will hold, and only increase or decrease of the steady state distributions could assure the dissipative structure stability in the cell. B. Traveling wave solutions of (26). Application of the modified method of simplest equation The equation (26) (in particular (25)) will have solutions for such values of s, at which its non- linearity can be reduced to a polynomial non- linearity. The reduction of a non-linearity to a weak non-linearity is a commonly used approach in the nonlinear dynamics. In this way, here we develop the last three terms of (26) in a Taylor series centered in s = 0. We retain only the terms up to cubic power and obtain the following approximation of (26): Ds ∂2 s ∂ξ2 + v ∂s ∂ξ + αs3 + βs2 + χs + δ = 0 (30) where the new coefficients have the form: α = ADE2 − BD2 E3 + CD3 E4 , β = A E − BD E2 + CD2 E3 χ = F − BE + CD E2 , δ = G − C E (31) In this case we have not interested in the con- crete values of the above–given coefficients, be- cause their signs (positive or negative) do not affect our further analytical investigations. In addition, the cubic polynomial approximation means that we accept a weak non–linearity (but not linearization) of (26), i.e. D (the rate constant k1) is sufficiently smaller than E (the rate constant k2) to assure the approximation validity. Indeed, the last inequal- ity follows from the biochemical consideration that the processes of ERK inactivation and STAT dephosphorylation are faster than that of ERK and STAT interaction. The last is of molecular recognition type [21–23]. Next we will consider the wave propagation of the STAT5a density. For the purpose we will apply the methodology from Section I I to (31). In order to solve (31) we constrict a solution as finite series s(ξ) = n∑ i=0 aiφ i, ( dφ dξ )2 = r∑ =0 cφ  (32) where φ is the solution of a simpler equation (referred to as simplest equation according to Section I I), and ai and c j are parameters that we will determine below. We substitute (33) in (31) and obtain an equation that contains powers of φ. Next, we balance the highest power arising from the second derivative in (31) with the highest power arising in the term containing s3 in the same equation. The resulting balance equation is r = 2n + 2, n = 2, 3, ... In the simplest case, if n = 2, then r = 6. Here we will use an equation of kind of Bernoulli as simplest equation. For the purpose we assume that c0 = c1 = c3 = c5 = 0; c2 = p2; c4 = 2pq; c6 = q2 , 0, and search for a solution in the form: s(ξ) = a0 + a1φ + a2φ 2, dφ dξ = pφ + qφ3 (33) Biomath 3 (2014), 1404131, http://dx.doi.org/10.11145/j.biomath.2014.04.131 Page 7 of 11 http://dx.doi.org/10.11145/j.biomath.2014.04.131 E Nikolova et al., On Nonlinear Dynamics of the STAT5a Signaling Protein We note that at another choice of the parameters c j we can obtain another kind of dφ dξ , but we aim to use some simpler equation which has an exact solution according to the modified method of simplest equation. The substitution of (34) in (31) leads to the following system of relationships for the parameters of the solution (the system is of kind (4)): 8Dsa2q2 −αa32 = 0 3αa1a22 + 3Dsa1q 2 = 0 −α((2a0a2 + a21)a2 + 2a 2 1a2 + a 2 2 ∗ a0)+ +2va2q + 12Dsa2 pq −βa22 = 0 4Dsa1 pq + va1q −α(4a0a1a2 + (2a0a2 + a21)a1)− −2βa1a2 = 0 4Dsa2 p2 + 2va2 p −χa2 −α(a20a2 + 2a0a 2 1+ +(2a0a2 + a21)a0 −β(2a0a2 + a 2 1) = 0 (34) −3αa20a1 − 2βa0a1 + va1 p −χa1 + Dsa1 p 2 = 0 −βa20 −χa0 −αa 3 0 −δ = 0 The system (35) implies that a1 = 0 and q is a free parameter. The solution of (35) is: a0 = 1 6 3√ λ α − 2 3 3χα−β2 α 3√ λ − β 3α, a2 = 2 √ 2αDs q α p = − √ 2(−Ds 3√ λ2 + 12Dsχα− (35) −4Dsβ2 + v √ 2αDs 3√ λ)/(12Ds √ αDs 3√ λ) where λ = 36χβα− 108α2 − 8β3 + +12 √ 3(4χ3α−χ2β2 − 18χβαδ + 27δ2α2 + 4δβ3)α The expression for the solitary wave depends on the solution of the differential equation in (34) and it is given by: s(ξ) = a0 + a2 √ p exp [2p(ξ + c)] 1 − q exp [2p(ξ + c)] (36) Fig. 2. Graph of the solution s(ξ) at α = 1; β = 1; χ = 100; δ = 10; k = 1; Ds = 11. for the case p > 0, q < 0 and s(ξ) = a0 + a2 √ p exp [2p(ξ + c)] 1 + q exp [2p(ξ + c)] (37) for the case p < 0, q > 0, where a0, a2, p are given by (36). The values (37) or (38) describe ’kink’ waves, which can be interpreted as fronts of some density propagation of the STAT5a signaling protein in the intracellular space, as we will show in the next paragraph of the paper. Two examples of these kinks are shown in Fig. 2 and Fig. 3. We note, however, that the figures give rather an illustrative idea of the possible spatial–temporal behavior of the considered object due to the lack of experi- mental data for the model parameters. C. Discussions In Section I I I we demonstrated that the lowest and highest steady state values of the STAT5a concentration are realized practically in the ho- mogeneous case corresponding to the absence of appropriate response in the form of ERK or SOCS protein production initiated by the cell signaling. Biomath 3 (2014), 1404131, http://dx.doi.org/10.11145/j.biomath.2014.04.131 Page 8 of 11 http://dx.doi.org/10.11145/j.biomath.2014.04.131 E Nikolova et al., On Nonlinear Dynamics of the STAT5a Signaling Protein Fig. 3. Graph of the solution s(ξ) at α = 1; β = −1; χ = 100; δ = 10; k = −1; Ds = 11. Then stabilization of the STAT5a signaling protein occurs. When, however, similar response (ERK or SOCS proteins increase) is available, the STAT5a concentration dramatically increases or decreases near the nucleus membrane and some inhomoge- neous distribution of the STAT5a molecules takes place in the cytosol. This means that the diffusion involves in the process and it can be presented by (25) (or (26)). In this case the STAT5a concentra- tion will propagate in the intracellular space as a ”kink“ wave, as we demonstrated in the previous paragraph. We choose the wave propagates along the positive direction of axis , i.e. - from the cell membrane to the nucleus one according to the direction of the signal transduction. Equation (25) (in particular (26) or (27) in view of the introduced traveling-wave ansatz) has the same steady state values for s as in the homogeneous case which we considered in Section I I I. However in inho- mogeneous case, presented by (25) or (26), the lowest and the highest steady states are unstable, and the intermediate one is stable. Therefore the steady state values in inhomogeneous case change the type of their character in comparison with the homogeneous one. In this way at additional cell signals, initiating increase of the active ERK molecules (but in absence of additional SOCS protein production), the STAT5a concentration will increase from the lowest destabilized value s01 to stabilized one s03, and density (concentration) jump will be observed. A similar example is shown in Fig.2. On the other hand at additional new signal (ligand), which can inhibit the STAT phosphorylation (For example, the SOCS protein production increases, but there is not ERK pro- tein production), the STAT5a concentration can decrease from its highest destabilized value s02 to stabilized one s03, and density (concentration) drop will be realized in the intracellular space (Fig. 3). Thereby the inhomogeneous effects observed in the STAT5a protein dynamics are predetermined by the levels of the ERK and SOCS proteins. V. Conclusion In this paper, we have re-considered the reaction model of the cross talk between ERK and STAT signaling pathways. By applying a well-known QSSA theorem to the described model we have demonstrated that near the quasi-stationary state of the considered process its dynamics determines by the behavior of the phosphorylated STAT5a protein. Moreover, if the STAT5a concentration is homogeneous distributed in the cell cytoplasm its possible minimum or maximum levels will support stabilization of the cell signaling. But, if some inhomogeneous distribution of the STAT5a protein appears in the intracellular space as a result of external signals destabilization of the above men- tioned levels will occur. In order to ensure the sta- bilization of the cell signaling process, the STAT5a density (concentration) will move in the cell from its unstable levels to stable ones. By applying the modified method of simplest equation to the corre- sponding reaction-diffusion model we have derived an analytic expression (solution) for the spatial STAT5a wave propagation . The ’kink’ form of the considered protein density wave supposes eventual drop and jump concentration effects. These effects are a sign for cell inhomogeneity, which deviates the cell from the homeostatic norms and enables the initiation of pathological processes. Biomath 3 (2014), 1404131, http://dx.doi.org/10.11145/j.biomath.2014.04.131 Page 9 of 11 http://dx.doi.org/10.11145/j.biomath.2014.04.131 E Nikolova et al., On Nonlinear Dynamics of the STAT5a Signaling Protein References [1] S. Khurana, S. Kreydiyyeh, A. Aronzon, et al., Asym- metric signal transduction in polarized ileal Na+- absorbing cells: carbachol activates brush-border but not basolateral-membrane PIP2-PLC and translocates PLC-1 only to the brush border, Biochemical Journal, vol. 313, pp. 509–518, 1996. [2] T. L. Holdaway-Clarke, J. A. Feijo, G. R. Hackett, J. G. Kunkel, P. K. Hepler, Pollen tube growth and the intracellular cytosolic calcium gradient oscillate in phase while extracellular calcium influx is delayed, Plant Cell, vol. 9, pp. 1999–2010, 1997. http://dx.doi.org/10.1105/tpc.9.11.1999 [3] H. Lam, J.-Y. Matroule, C. Jacobs-Wagner,The asym- metric spatial distribution of bacterial signal transduction proteins coordinates cell cycle events, Developmental Cell, vol. 5, pp. 149–159, 2003. http://dx.doi.org/10.1016/S1534-5807(03)00191-6 [4] T. Y. Belenkaya, C. Han, D. Yan, et al., Drosophila Dopmorphogen movement is independent of dynamin- mediated endocytosis but regulated by the glypican mem- bers of heparan sulfate proteoglycans, Cell, vol. 119, pp. 231–244, 2004.http://dx.doi.org/10.1016/j.cell.2004.09.031 [5] U. Bhalla, Signaling in small subcellular volumes. I. Stochastic and diffusion effects on individual pathways, Biophys.J., vol. 87, pp. 733–744, 2004. http://dx.doi.org/10.1529/biophysj.104.040469 [6] S. Schnell, T. Turner, Reaction kinetics in intracellular en- vironments with macromolecular crowding: simulations and rate laws, Prog.Biophys. Mol.Biol., vol. 85, pp. 235– 260, 2004. http://dx.doi.org/10.1016/j.pbiomolbio.2004.01.012 [7] V. Petrov, N. Georgiev, G. Georgiev, Dynamical Model of ERK and STAT Proteins Interaction, Proceedings of the 2006 WSEAS Int. Conf. on Cellular and Molecular Biology, Biophysics and Bioengineering, Athens, Greece, July 14-16, 2006, pp. 147–153. [8] W. Hereman, A. Nuseir, Symbolic methods to construct exact solutions of nonlinear partial differential equations, Math. Comput. Simul., vol. 43, pp. 13–27, 1997. http://dx.doi.org/10.1016/S0378-4754(96)00053-5 [9] Y.-Z. Peng, Exact solutions for some nonlinear partial dif- ferential equations, Phys. Lett. A , vol. 314, pp. 401–408, 2003. http://dx.doi.org/10.1016/S0375-9601(03)00909-5 [10] W. Malfliet, W. Hereman, The tanh method: I. Exact solutions of nonlinear evolution and wave equations, Phys. Scr., vol. 54, pp. 563–568, 1996 . http://dx.doi.org/10.1088/0031-8949/54/6/003 [11] N. A. Kudryashov, Simplest equation method to look for exact solutions of nonlinear differential equations, Chaos Solitons Fractals, vol. 24, pp. 1217–1231, 2005. http://dx.doi.org/10.1016/j.chaos.2004.09.109 [12] N. A. Kudryashov, N. B. Loguinova, Extended sim- plest equation method for nonlinear differential equations, Commun. Nonlinear Sci. Numer. Simul., vol. 14, pp. 3507–3529, 2009. http://dx.doi.org/10.1016/j.cnsns.2009.01.023 [13] N. K. Vitanov, Z. I. Dimitrova, H. Kantz, Modified method of simplest equation and its application to non- linear PDEs, Appl. Math. Comput., vol. 216, pp. 2587– 2595, 2010. http://dx.doi.org/10.1016/j.amc.2010.03.102 [14] N. K. Vitanov, Modified method of simplest equa- tion: powerful tool for obtaining exact and approximate traveling-wave solutions of nonlinear PDEs, Commun. Nonlinear Sci. Numer. Simul., vol. 16, pp. 176–1185, 2011. [15] M. V. Demina, N. A. Kudryashov, Explicit expressions for meromorphic solution of autonomous nonlinear or- dinary differential equations, Commun. Nonlinear Sci. Numer. Simul., vol. 16, pp. 1127–1134, 2011. http://dx.doi.org/10.1016/j.cnsns.2010.06.035 [16] M. V. Demina, N.A. Kudryashov, From Laurent series to exact meromorphic solutions: the Kawahara equation, Phys. Lett. A, vol. 374, pp. 4023–4029, 2010. http://dx.doi.org/10.1016/j.physleta.2010.08.013 [17] T. Pawson, M. Raina, P. Nash, Interaction domains: from simple binding events to complex cellular behavior, FEBS Letters, vol. 513, pp. 2–10, 2002. http://dx.doi.org/10.1016/S0014-5793(01)03292-6 [18] K. Shuai, Modulation of STAT signaling by STAT- interacting proteins, Oncogene, vol.19, pp. 2638–2644, 2000. http://dx.doi.org/10.1038/sj.onc.1203522 [19] W. S. Alexander, Suppressors of cytokine signaling (SOCS) in the immune system, Nat Rev Immunol, vol. 2, pp. 410–416, 2002. [20] N. A. Calcano, D. Sanden, J. A. Johnston, Tyrosine- phosphorilated SOCS-3 inhibits STAT activation but binds to p120 Ras GAP and activates Ras, Nat Cell Biol., vol. 3, pp. 460–465, 2001. http://dx.doi.org/10.1038/35074525 [21] T. J. Wood, D. Silva, P. E. Lobie, T. J. Pircher et al., Mediation of growth hormone-dependent transcrip- tional activation by mammary gland factor/Stat5, J. Biol. Chem., vol. 270, pp. 9448–9453, 1995. http://dx.doi.org/10.1074/jbc.270.16.9448 [22] T. J. Pircher, A. Flores-Morales, A. L. Mui, A. R. Saltiel et al, Mitogen-activated protein kinase kinase inhibition decreases growth hormone stimulated transcription me- diated by STAT5, Mol. Cell. Endocrinol., vol. 133, pp. 169–176, 1997. http://dx.doi.org/10.1016/S0303-7207(97)00164-0 [23] T. J. Pircher, H. Petersen, J. A. Gustafsson, L. A. Hal- dosen, Extracellular signal-regulated kinase (ERK) inter- acts with signal transducer and activator of transcription (STAT)5a, Molecular Endocrinology, vol. 157, pp. 555– 565, 1999. http://dx.doi.org/10.1210/mend.13.4.0263 [24] S. E. Nicholson, T. A. Willson, N. Farley, R. Starr et al., A Mutational analysis of the SOCS proteins suggest a dual domain requirement but distinct mechanisms for in- hibition of LIF and IL-6 signal transduction, The EMBO Journal, vol. 18, pp. 375–385, 1999. http://dx.doi.org/10.1093/emboj/18.2.375 Biomath 3 (2014), 1404131, http://dx.doi.org/10.11145/j.biomath.2014.04.131 Page 10 of 11 http://dx.doi.org/10.1105/tpc.9.11.1999 http://dx.doi.org/10.1016/S1534-5807(03)00191-6 http://dx.doi.org/10.1016/j.cell.2004.09.031 http://dx.doi.org/10.1529/biophysj.104.040469 http://dx.doi.org/10.1016/j.pbiomolbio.2004.01.012 http://dx.doi.org/10.1016/S0378-4754(96)00053-5 http://dx.doi.org/10.1016/S0375-9601(03)00909-5 http://dx.doi.org/10.1088/0031-8949/54/6/003 http://dx.doi.org/10.1016/j.chaos.2004.09.109 http://dx.doi.org/10.1016/j.cnsns.2009.01.023 http://dx.doi.org/10.1016/j.amc.2010.03.102 http://dx.doi.org/10.1016/j.cnsns.2010.06.035 http://dx.doi.org/10.1016/j.physleta.2010.08.013 http://dx.doi.org/10.1016/S0014-5793(01)03292-6 http://dx.doi.org/10.1038/sj.onc.1203522 http://dx.doi.org/10.1038/35074525 http://dx.doi.org/10.1074/jbc.270.16.9448 http://dx.doi.org/10.1016/S0303-7207(97)00164-0 http://dx.doi.org/10.1210/mend.13.4.0263 http://dx.doi.org/10.1093/emboj/18.2.375 http://dx.doi.org/10.11145/j.biomath.2014.04.131 E Nikolova et al., On Nonlinear Dynamics of the STAT5a Signaling Protein [25] T. Tamiya, I. Kashiwagi, R. Takahashi, H. Yasukawa, A. Yoshimura, Suppressors of Cytokine Signaling (SOCS) Proteins and JAK/STAT Pathways: Regulation of T- Cell Inflammation by SOCS1 and SOCS3, Arterioscler Thromb Vasc Biol., vol. 31, pp. 980–985, 2011. http://dx.doi.org/10.1161/ATVBAHA.110.207464 [26] A. N. Tichonov, Systemy differentsialnyh uravneniy, soderjashchie malye parametry pri proizvodnyh, Matem- aticheskiy sbornik, vol. 31, pp. 575-586, 1952 (in Rus- sian). [27] V. Kotev, S. Nikolov, Stability Analysis of Time Delay Model of Crosstalk between ERK and STAT5a Interac- tion, Bioautomation, vol. 7, pp. 90–98, 2007. Biomath 3 (2014), 1404131, http://dx.doi.org/10.11145/j.biomath.2014.04.131 Page 11 of 11 http://dx.doi.org/10.1161/ATVBAHA.110.207464 http://dx.doi.org/10.11145/j.biomath.2014.04.131 Introduction The method of the simplest equation The interaction between ERK and STAT signaling pathways A reaction model and its reduction Stability analysis of the steady state solution of (12) Space-temporal modeling of the STAT5a signaling protein A reaction-diffusion model of the STAT5a signaling protein. Stability analysis of the inhomogeneous distribution of the STAT5a concentration. Traveling wave solutions of (26). Application of the modified method of simplest equation Discussions Conclusion References