Accuracy of position determination in Ca 2+ signaling
AAccuracy of position determination in Ca signaling Vaibhav H. Wasnik,
1, 2
Peter Lipp, and Karsten Kruse NCCR Chemical Biology, Departments of Biochemistry and Theoretical Physics,University of Geneva, 1211 Geneva, Switzerland Indian Institute of Technology Goa, Ponda-403401, India Institute for Molecular Cell Biology,Research Centre for Molecular Imaging and Screening,Center for Molecular Signaling (PZMS), Medical Faculty,Saarland University, Homburg/Saar, Germany (Dated: July 22, 2019)
Abstract
A living cell senses its environment and responds to external signals. In this work, we studytheoretically, the precision at which cells can determine the position of a spatially localized transientextracellular signal. To this end, we focus on the case, where the stimulus is converted into therelease of a small molecule that acts as a second messenger, for example, Ca , and activates kinasesthat change the activity of enzymes by phosphorylating them. We analyze the spatial distributionof phosphorylation events using stochastic simulations as well as a mean-field approach. Kinasesthat need to bind to the cell membrane for getting activated provide more accurate estimates thancytosolic kinases. Our results could explain why the rate of Ca detachment from the membrane-binding conventional Protein Kinase C α is larger than its phosphorylation rate. a r X i v : . [ q - b i o . CB ] J u l . INTRODUCTION Living cells respond to external stimuli. Often these stimuli consist in binding a ligand to acell surface receptor [1], but cells also react to mechanical forces that are applied to them [2],to changes in ambient light or temperature. In this way, cells can migrate towards favourableand away from unfavourable environmental conditions [3], sense the density of cells of theirown kind [4], initiate or inhibit cell division, or trigger developmental programs [5]. Theability to sense external stimuli is thus of paramount importance for the success of a single-or multi-celled species to survive: The better a cell can read out signals, the better it willdo.This has led to investigations of the physical limits of cellular signalling. Berg and Purcelldetermined the conditions under which an organism can optimally determine the concentra-tion of a molecule in its environment in case the cell uses independent receptors [6]. In thebacterium
Escherichia coli some kind of chemoreceptor clusters at one cell end [7], and phys-ical studies revealed how this clustering can enhance sensitivity to external signals [8–13].Subsequently, these investigations have been generalized to account also for the energeticsof ligand-receptor binding [14, 15]. Whereas
E. coli is too small to directly sense spatialgradients in the concentration of chemoattractants and thus relies to this end on detectingtemporal concentration changes while moving, eukaryotic cells like neutrophils and leuko-cytes of the human immune system, the budding yeast
Saccharomyces cerevisiae , or the slimemold
Dictyostelium discoideum can directly sense spatial gradients [16]. Physical constraintson the accuracy of sensing spatial gradients by single cells have been established [17–19]. Inthis context, spatial aspects of the intracellular signalling network have been considered [20].In addition to sensing spatial gradients, eukaryotic cells can also respond to spatiallylocalized signals. For example, neurons reinforce or weaken synapses in response to transientstimuli [21, 22]. Another example is provided by cells from the immune system, whichpolarize upon making contact with an antigen presenting cell [23].In this work, we address the question of how the spatial resolution at which a cell candetect a stimulus depends on the intracellular mechanism of reading it out. Stimuli of anykind are typically transformed into a cell response by activating or inactivating proteinsthrough adding or removing phosphate groups (phosphorylation and dephosphorylation)to and from certain amino acids. Proteins leading to phosphorylation are called kinases,2hereas phosphatases carry out dephosphorylation. Eventually, changing the phosphory-lation of proteins can lead to changes of the cytoskeletal organization on short time scalesand to modifications of the cellular gene expression profile on long time scales to give buttwo examples. Commonly, external stimuli do not directly (de-)phosphorylate proteins, butfirst elicit the release of a second messenger, for example, cyclic Adenosine monophosphate(cAMP), inositol triphosphate (IP ), diacylglycerol (DAG), or Ca . In this way, the samemachinery can be used to respond to a variety of different stimuli. At the same time it raisesthe question of how to obtain a specific response to a given signal.The second messenger Ca is read out by different proteins. Notably, a number of kinasesare activated by Calmodulin, a peptide diffusing in the cytoplasm that changes conformationafter binding to a Ca ion. In contrast, the ubiquitously expressed conventional ProteinKinase C α (PKC α ) directly binds Ca . For activation it requires also binding to DAG inthe plasma membrane of a cell. Strikingly, the detachment rate of Ca from PKC α hasbeen measured to be of the order of 50ms, whereas the phosphorylation of a target proteinby PKC α takes about 500ms [24]. This suggests that PKC α is an inefficient kinase andraises the question of the evolutionary benefit of a high Ca detachment rate.In this work, we develop a framework for studying the spatial distribution of phospho-rylation events inside a cell. We start by studying a toy model that allows us to introducesome notation and to present the tools we will use for analyzing the stochastic processescorresponding to phosphorylation by a cytosolic and by a membrane-binding kinase. Wewill then study the spatial distribution of phosphorylation events in response to a singleCa for a cytosolic and a membrane-binding kinase, where we find that the latter yields amore accurate read-out than the former. We will then discuss responses to Ca puffs andinvestigate the influence of background phosphorylation on the precision of estimating thesite of Ca entry. A short account of some of the results presented in this work has beengiven in Ref. [Letter]. II. A TOY MODEL
In order to introduce some quantities as well as some methods that we will use later toanalyze the localization of Ca influx into the system, we will study a toy model in thissection. It can be interpreted as describing kinase-dependent phosphorylation following the3 a target ν p P i a) ν l b)c) d) xz FIG. 1. (color online) Toy model for position determination of a Ca entry site. a) Illustrationof the system. The Ca ion enters through the membrane at x = 0 and immediately activatesa kinase, which phosphorylates target proteins at rate ν p . The Ca ion is lost immediatelyafter detaching from the kinase at rate ν l . b) Distribution P of estimated Ca entry sites ˆ x obtained from 10 stochastic simulations for ν p = 1. c) Estimation error (cid:96) as a function of thephosphorylation rate ν p . Circles are from stochastic simulations. The yellow full line represents theerror for ν p (cid:28)
1, Eq. (22), the red dashed line the error for ν p → ∞ , Eq. (30). d) Distribution P ofestimated Ca entry sites ˆ x obtained from 10 stochastic simulations for ν p = 0 .
01. In (b, d), thered line presents a Gaussian fit to the data, whereas the yellow line is the normalized distribution P for ν p (cid:28) entry of a single Ca ion, see Fig. 1a. The kinase is activated by immediately bindingthe Ca ion at the latter’s entry site. After inactivation of the kinase following Ca detachment, the Ca ion is immediately lost. A. The Master equation
Consider activation of the kinase at x = 0. We will restrict attention to the dynamicsalong the x -direction and assume that the system has no boundaries in that direction.The activated kinase phosphorylates at a constant rate ν p , while diffusing with a diffusionconstant D . The kinase is inactivated at rate ν l . Let n denote the spatial distribution ofphosphorylation events. The state of the system is given by the probability P a [ n ( ξ ); x, t ]for having an active kinase at x at time t and a distribution of phosphorylation events4 ( ξ ). Similarly, P i [ n ( ξ ); t ] is the corresponding distribution when the kinase is inactive. Inthis case its position is irrelevant, because the kinase cannot be activated again such that nofurther phosphorylation events can cccur. The Master equation governing the time evolutionof these distributions is given by ∂ t P i [ n ( ξ ); t ] = ν l (cid:90) d x P a [ n ( ξ ); x, t ] (1) ∂ t P a [ n ( ξ ); x, t ] = D∂ x P a [ n ( ξ ); x, t ] − ν l P a [ n ( ξ ); x, t ]+ ν p {P a [ n ( ξ ) − δ ( ξ − x ); x, t ] − P a [ n ( ξ ); x, t ] } , (2)where δ denotes the Dirac distribution. The initial condition at t = 0 is given by P a [ n ( ξ ) = 0; x, t = 0] = δ ( x ) (3)and all other probabilities equal to zero.For the analysis of the Master equation, we will scale time by ν l and space by (cid:112) D/ν l .The only remaining dimensionless parameter is then ν p /ν l . We will keep the notation ν p forthe dimensionless phosphorylation rate. B. Number of phosphorylation events
It is instructive to first neglect the spatial aspects of the phosphorylation dynamics and todetermine the distribution of the number of phosphorylation events. Consider the probabilitydistribution P a ( n, t ) that n phosphorylation events have taken place at time t and with thekinase being active. The corresponding distribution when the kinase is inactive is P i ( n, t ).The Master equation for P a and P i readsdd t P i ( n ) = P a ( n ) (4)dd t P a (0) = − P a (0) − ν p P a (0) (5)dd t P a ( n ) = − P a ( n ) + ν p ( P a ( n − − P a ( n )) , (6)where n ≥ n ≥ P a (0 , t = 0) = 1, P a ( n, t = 0) = 0 for all n >
0, and P i ( n, t = 0) = 0 for all n ≥
0. Thesolution to these equations is P i ( n, t ) = (1 − r ) r n − (1 − r ) (cid:34) n (cid:88) k =0 ( ν p t ) k k ! r n − k (cid:35) e − (1+ ν p ) t (7) P a ( n, t ) = ( ν p t ) n n ! e − (1+ ν p ) t (8)with r = ν p ν p (9)for all n ≥
0. For t → ∞ we have P i ( n ) = (1 − r ) r n (10) P a ( n ) = 0 , (11)which yields the average number N p ≡(cid:104) n (cid:105) = (cid:80) ∞ n =0 n ( P i ( n, t ) + P a ( n, t )) of phosphoryla-tion events, namely N p = ν p . The variance around this value is var( n ) ≡(cid:104) n (cid:105) − (cid:104) n (cid:105) = ν p (1 + ν p ) = (1 + N p ) N p . C. Estimated position of Ca release and estimation error For a given distribution of phosphorylation events n ( ξ ), we estimate the position of Ca release by computing the distribution average:ˆ ξ n ( ξ ) = (cid:90) d ξ ξn ( ξ ) / (cid:90) d ξ n ( ξ ) , (12)which obviously only exists, if (cid:82) d ξ n ( ξ ) (cid:54) = 0. If n ( ξ ) = 0 for all ξ , the position cannot beestimated. In this way, we obtain a distribution P of estimated positions. Explicitly, thisdistribution is given by P ( ˆ ξ, t ) = (cid:90) D n ( ξ ) δ ( ˆ ξ n ( ξ ) − ˆ ξ ) P [ n ( ξ ); t ] , (13)where ˆ ξ n ( ξ ) is the average of the distribution n ( ξ ) according to Eq. (12) and P [ n ( ξ ); t ] ≡P i [ n ( ξ ); t ] + (cid:82) d x P a [ n ( ξ ); x, t ] is the probability functional of phosphorylation distributions.6e will use the variance (cid:96) of the distribution P as a measure of the error in localizing thesite of Ca release. Given that the mean of the estimated positions is zero, we obtain (cid:96) = (cid:90) d ˆ ξ ˆ ξ P ( ˆ ξ, t ) . (14)Combining Eqs. (12)-(14), we get (cid:96) = (cid:90) D n ( ξ ) (cid:18) (cid:82) d ξ ξn ( ξ ) (cid:82) d ξ n ( ξ ) (cid:19) P [ n ( ξ ); t ] . (15) D. Solution of the Master equation
In general, it is not possible to solve the Master equation (1) and (2) and to determinethe estimation error analytically. Let us thus start our analysis by a stochastic simulation.To this end, we use a variant of the Gillespie algorithm: we define the total rate of possibleevents as ν tot ≡ ν p and draw the time ∆ t until occurrence of the next event from thedistribution ν tot exp ( − ν tot t ). Then, the kind of event is determined by drawing a numberfrom a uniform distribution on the interval [0 , ν p ]. If the number is smaller than 1, thenthe kinase is inactivated and the simulation stops. In the opposite case, the position of thekinase is changed by drawing a random number from a Gaussian distribution with zero meanand variance 2∆ t and adding this value to the current position of the kinase. This positionis then recorded as the site of phosphorylation and the simulation continues with drawingthe time to the next event. From the distribution of positions of phosphorylation eventsobtained in this way, we calculate the mean position, which we take to be the estimatedposition at which the Ca ion entered the system. If there was no phosphorylation event,then the position of Ca entry is not estimated.In Figure 1b, we show the distribution of estimated Ca entry sites for ν p = 1 obtainedform 10 simulation runs. The dependence of the error (cid:96) on the phosphorylation rate ν p is given in Fig. 1c, where for each value at least 10 simulation runs have been performed.The error is about 2 for ν p → ν p an increasing number of phosphorylations and thusposition measurements occur on average for a Ca ion. In the limit ν p → ∞ , the error is (cid:96) ≈ .We can make some analytical progress by considering the limiting cases of very small andlarge phosphorylation rates. 7 . The limit ν p (cid:28) For very small phosphorylation rates ν p (cid:28)
1, the probability of having a trajectory withtwo or more phosphorylation events is negligible. In that case, the Master equation (1) and(2) reduces to ∂ t P a = ∂ x P a − (1 + ν p ) P a (16) ∂ t P a = ∂ x P a + ν p δ ( x − ξ ) P a − P a (17) ∂ t P i = (cid:90) d x P a ≡ ¯ P a , (18)where P a ( x, t ) is the probability of having, at time t , an active kinase at x without anyphosphorylation, whereas P a ( ξ ; x, t ) is the probability of the active kinase being at x andwhere a phosphorylation event had occurred at position ξ . Similarly, P i ( ξ ; t ) denotes thecorresponding probability after inactivation of the kinase. We do not show the dynamicequation for the probability of having lost the kinase before it phosphorylated, because it isirrelevant for estimating the position of the Ca release site. In the last equation, we haveintroduced the marginal distribution of the phosphorylation position after integrating outthe position of the kinase, ¯ P a ( ξ, t ) = (cid:82) d x P a ( ξ ; x, t ). Integrating Eq. (17) with respect tothe kinase position x , we obtain its dynamic equation ∂ t ¯ P a ( ξ, t ) = ν p P a ( ξ, t ) − ¯ P a ( ξ, t ) . (19)The solution to Eq. (16) is P a ( x, t ) = 1 √ πt exp (cid:26) − (1 + ν p ) t − x t (cid:27) (20)The distribution of the phosphorylation position P ( ξ ) ≡ ¯ P a ( ξ ) + P i ( ξ ) is obtained from theexpression for P a through P ( ξ, t ) = ν p (cid:82) t d t (cid:48) P a ( ξ, t (cid:48) ). In the limit t → ∞ , we get P ( ξ ) = ν p (cid:112) ν p exp (cid:8) − (cid:112) ν p | ξ | (cid:9) , (21)such that (cid:96) = 21 + ν p . (22)For small enough values of ν p , the distribution Eq. (21) agrees well with the distribution ofestimated entry sites obtained from stochastic simulations, see Fig. 1d. The expression ofthe error as given by Eq. (22) gives a good approximation for ν p (cid:46) . . Continuous phosphorylation We now consider the limit of a very large phosphorylation rate, such that phosphorylationoccurs at any point of the kinase’s trajectory. In that case, the estimated position of releaseˆ x of a single Ca as defined in Eq. (12) is obtained by the average position of the kinasealong its trajectory x ( t ), that is, ˆ x T = 1 T (cid:90) T d t x ( t ) , (23)where T is the time span between entry of the Ca and loss of the kinase. The error isthen given by the variance (cid:96) = (cid:90) ∞ d T (cid:10) ˆ x T (cid:11) e − T , (24)where the exponential factor accounts for the probability of finding a trajectory of duration T . To calculate the the expectation value (cid:104) ˆ x T (cid:105) , we note that the trajectories x ( t ) are solutionsto a Langevin equation ˙ x ( t ) = ζ ( t ) (25)for 0 ≤ t ≤ T , where ζ is a fluctuating “force” that obeys a Gaussian distribution at eachtime t with zero mean, (cid:104) ζ ( t ) (cid:105) = 0, and (cid:104) ζ ( t ) ζ ( t (cid:48) ) (cid:105) = 2 δ ( t − t (cid:48) ). The solution to this equationis x ( t ) = (cid:82) t d t (cid:48) ζ ( t (cid:48) ), such that (cid:10) ˆ x T (cid:11) = 1 T (cid:90) T d t (cid:90) T d¯ t (cid:104) x ( t ) x (¯ t ) (cid:105) (26)= 1 T (cid:90) T d t (cid:90) T d¯ t (cid:90) t d t (cid:48) (cid:90) ¯ t d t (cid:48)(cid:48) (cid:104) ζ ( t (cid:48) ) ζ ( t (cid:48)(cid:48) ) (cid:105) (27)= 1 T (cid:90) T d t (cid:90) T d¯ t | t − ¯ t | (28)= 23 T. (29)From this we obtain using Eq. (24) (cid:96) = 23 , (30)which is in good agreement with the numerical results presented in Fig. 1c.9 . The mean-field limit A natural approximation is to use a mean-field ansatz and to pose that the rate ofphosphorylation at position x is proportional to the probability of finding the kinase atthis position. Let p a denote the probability of finding an active kinase at x , p a ( x, t ) = (cid:82) D n ( ξ ) P a [ n ( ξ ); x, t ], where the functional integral extends over all possible phosphoryla-tion distributions n ( ξ ) with n ( ξ ) ≥ ξ . Furthermore, let p i denote the probability ofhaving no kinase, p i ( t ) = (cid:82) D n ( ξ ) P i [ n ( ξ ); t ]. The time-evolution of p a obeys ∂ t p a = ∂ x p a − p a . (31)Once p a is known, then p i = 1 − (cid:82) d x p a ( x ). Furthermore, the probability of having n phosphorylation events at position x at time t , P ( n, x, t ), is determined by˙ P (0 , x, t ) = − ν p p a ( x, t ) P (0 , x, t ) (32)˙ P ( n, x, t ) = ν p p a ( x, t ) ( P ( n − , x, t ) − P ( n, x, t )) (33)for n ≥
1. The normalization conditions read (cid:80) ∞ n =0 P ( n, x, t ) = 1 for all x and all t and (cid:82) ∞−∞ d x p a ( x, t ) + p i ( t ) = 1 for all t . The initial condition is P (0 , x, t = 0) = 1 for all x and p a ( x, t = 0) = δ ( x ).The solution to these equations is P ( n, x, t ) = ν np n ! ¯ p a ( x, t ) n exp {− ν p ¯ p a ( x, t ) } (34)with ¯ p a ( x, t ) = (cid:90) t d t (cid:48) p a ( x, t (cid:48) ) (35)and p a ( x, t ) = (cid:114) π t exp (cid:26) − t − x t (cid:27) . (36)We are interested in the distribution as t → ∞ . In that case,¯ p a ( x ) = 12 exp {−| x |} , (37)where the bar indicates the distribution for t → ∞ .10n the spirit of the mean-field ansatz , we replace the expression for the estimation error,see Eq. (15), by (cid:96) = (cid:82) d x x ˆ n ( x ) (cid:82) d x ˆ n ( x ) . (38)Here, ˆ n is the mean number of phosphorylation events at x . As we will see below inSect. IV B, the expression for the estimation error as defined in Eq. (15) differs from theabove expression if the probability distribution is given by Eq. (34). However, the two ex-pressions for the error are the same in the limit of a small number of phosphorylation events, (cid:104) n (cid:105) (cid:28) P ( n, x, t ) in the limit t → ∞ , we get for the mean numberof phosphorylation events at x ˆ n ( x ) = ν p ¯ p a ( x ) and thus (cid:96) = 1 . (39)In the rescaled units, this is just the diffusion length of the activated kinase, (cid:112) D/ν l in theoriginal units. Note, that it holds for arbitrary values of ν p >
0, which is different from theresults of the stochastic simulations. Indeed, the mean-field equations can only be expectedto work well in cases, when the number of phosphorylation events is small, that is, when twoor more phosphorylation events for a single Ca ion are rare. In fact, if there is only onephosphorylation event, then the spatial distribution of the active kinase determines directlythe distribution of phosphorylation events. This is not true for subsequent phosphorylations.In spite of this obvious failure of the mean-field approximation to determine the dependenceof the error on ν p , it does give, however, the correct order of magnitude of the error, whichvaries between 2/3 and 2, see Fig. 1d. III. PHOSPHORYLATION DYNAMICS IN RESPONSE TO A SINGLE CA ION
We will now study the phosphorylation dynamics and the ensuing spatial distributionsof phosphorylation events in two scenarios that reflect essential properties of Ca activatedkinases in response to a single Ca ion entering the system, see Fig. 2. A single Ca ion entering a real cell is unlikely to elicit a response. However, when treating the case of aCa puff below, we will assume that all Ca ions are independent of each other. Therefore,studying the response to a single Ca ion is appropriate. As we will see in Sect. IV, however,11 a target ν a ν d ν p P i Ca PKC ν a ν d ν p CaM P i target ν b ν u a) b) ν l ν l xz FIG. 2. (color online) Schematic of the two scenarios of Ca dependent phosphorylation. a)Cytosolic kinase: A cytosolic kinase is activated directly by Ca . The Ca attaches at rate ν a and detaches at rate ν d . The activated kinase phosphorylates target proteins at rate ν p and Ca is lost at rate ν l . b) Membrane-binding kinase: In this case, after attaching Ca , the kinase stillneeds to bind to the membrane before it is active. Binding and unbinding occur at rates ν b and ν u . the passage from one Ca ion to a Ca puff is not trivial, because not all Ca ions leadto a phosphorylation event.The Ca ion enters the system at x = z = 0 at t = 0, where the z -direction is thedirection perpendicular to the membrane, which is located at z = 0. In the first scenario,the Ca ion diffuses in the cytoplasm until it either binds to and thereby activates a kinaseor gets lost from the system. The Ca ion can detach and reattach to the kinase and canonly be lost from the system, when not being attached to the kinase. In the second scenario,after attaching a Ca ion, the kinase still needs to bind to the membrane to be activated.As in the toy model, we will consider the average position of phosphorylation events inresponse to a Ca ion as an estimate of the Ca entry site. This is adequate if the rateof dephosphorylation is much smaller than the rate at which the Ca ion is eventually lostform the system. It also requires that the phosphorylated target proteins remain immobileif the estimate obtained in this way should give a good proxy for the localization of the cellresponse. The dynamics of the target proteins, however, depends on the protein at handand can meaningfully be studied only with a specific cell response in mind. In contrast, theresults we obtain are independent of the specific response and provide a general lower boundfor the estimation error. 12 . The Master equations In this section, we present the Master equations for the two different scenarios of phos-phorylation dynamics. We will consider only the situation in which kinases are abundantand uniformly distributed, such that the attachment rate of Ca to a kinase is constant.
1. Cytosolic kinase
Consider a kinase that is active immediately after attaching a Ca ion. In fact, there areno cytosolic kinases known that are directly activated by Ca binding. Instead, Ca typi-cally binds to Calmodulin (CaM), which then binds to and thereby activates a Calmodulin-dependent kinase [25]. In our analysis, we neglect the second step. If anything, directactivation by Ca will only increase the accuracy of position determination.For a cytosolic kinase that is directly activated by binding a Ca ion, the system stateis determined by the spatial distribution of phosphorylation events n ( ξ ) and, if the Ca has not been lost from the system, by its position x , where we have to distinguish betweenthe possibility that the Ca ion is either free or attached to a kinase. Let P i be thecorresponding probability distribution, when Ca is present, but not attached to the kinase, P a the probability distribution in the case that Ca is attached to a kinase, and P theprobability distribution, when the Ca ion is lost from the system. Then the Masterequation is given by ∂ t P i [ n ( ξ ); x, t ] = D C ∂ x P i [ n ( ξ ); x, t ] + ν d P a [ n ( ξ ); x, t ] − ( ν a + ν l ) P i [ n ( ξ ); x, t ] (40) ∂ t P a [ n ( ξ ); x, t ] = D K ∂ x P a [ n ( ξ ); x, t ] − ν d P a [ n ( ξ ); x, t ] + ν a P i [ n ( ξ ); x, t ]+ ν p {P a [ n ( ξ ) − δ ( ξ − x ); x, t ] − P a [ n ( ξ ); x, t ] } (41) ∂ t P [ n ( ξ ); t ] = ν l (cid:90) d x P i [ n ( ξ ); x, t ] , (42)where D C and D K are, respectively, the diffusion constants of Ca and the kinase, ν a and ν d are the respective rates of Ca attachment to and detachment from the kinase, ν l is therate at which Ca is lost from the system, and ν p again the rate at which an active kinasephosphorylates. Note, that the dynamics in the z -direction is irrelevant for this scenario.13he initial condition is P i [ n ( ξ ) ≡ x, t = 0] = δ ( x ) (43)and all other probabilities zero reflecting that the Ca enters the system at x = 0. Theprobability distributions obey the normalization condition (cid:82) D n ( ξ ) (cid:8) P + (cid:82) d x [ P i + P a ] (cid:9) =1 for all times t . In contrast to the toy model, we will scale time by ν p and length by (cid:112) D C /ν p .We keep the same notation for the dimensionless parameters.
2. Membrane-binding kinase
The dynamics in the case of a kinase that needs to attach to a membrane for activationfollows the same reasoning. The membrane is assumed to be localized at z = 0 and toextend infinitely into the x -direction. Calcium enters at x = z = 0 into the domain z ≥ P i , when Ca is present, but not attached to the kinase, P a , when the kinase has Ca attached to it, but is not bound to the membrane and thusinactive, and P , when Ca is lost from the system, there is the distribution P b , when thekinase is active, that is, it has Ca attached to it and is bound to the membrane. Theyobey the following Master equation: ∂ t P i [ n ( ξ ); x, z, t ] = D C (cid:0) ∂ x + ∂ z (cid:1) P i [ n ( ξ ); x, z, t ]+ ν d P a [ n ( ξ ); x, z, t ] − ( ν a + ν l ) P i [ n ( ξ ); x, z, t ] (44) ∂ t P a [ n ( ξ ); x, z, t ] = D K (cid:0) ∂ x + ∂ z (cid:1) P a [ n ( ξ ); x, z, t ] − ν d P a [ n ( ξ ); x, z, t ] + ν a P i [ n ( ξ ); x, z, t ] (45) ∂ t P b [ n ( ξ ); x, t ] = ν b P a [ n ( ξ ); x, z = 0 , t ] − ν u P b [ n ( ξ ); x, t ]+ ν p {P b [ n ( ξ ) − δ ( ξ − x ); x, t ] − P b [ n ( ξ ); x, t ] } (46) ∂ t P [ n ( ξ ); t ] = ν l (cid:90) d x (cid:90) ∞ d z P i [ n ( ξ ); x, z, t ] , (47)where ν b and ν u are the respective rates of kinase binding to and unbinding from the mem-brane. The bulk equations are complemented by no-flux boundary conditions for the Ca ,that is, ∂ z P i | z =0 = 0 . (48)14he boundary condition for the kinase accounts for its binding to and unbinding from themembrane − D K ∂ z P a | z =0 = − ν b P a | z =0 + ν u P b . (49)Again, we scale time by ν − p and length by (cid:112) D C /ν p and keep the notation for the nowdimensionless parameters. B. Number of phosphorylation events
Let us first neglect the spatial degrees of freedom and consider only the distribution ofthe number of phosphorylation events.
1. Cytosolic kinase
In the case, we neglect the spatial degrees of freedom, the Master equation (40)-(42) canbe written as ˙ C n = − ( ν a + ν l ) C n + ν d K n (50)˙ K = ν a C − ν d K − K (51)˙ K n = ν a C n − ν d K n − K n + K n − (52)˙ P n = ν l C n , (53)where Eqs. (50) and (53) hold for all n ≥
0, whereas Eq. (52) is valid for n >
0. Here, C n , K n , and P n denote the respective probabilities of having n phosphorylation events, whenthe Ca ion is free, attached to the kinase, or lost from the system.In the limit t → ∞ , we have C n = K n = 0 for all n ≥
0. The distribution of phos-phorylation events is thus entirely determined by P ∞ n ≡ lim t →∞ P n ( t ). From Eq. (53)we have P ∞ n = ν l (cid:82) ∞ d t C n ≡ ν l ¯ C n . Integrating the dynamic equations (50)-(52) withrespect to time from t = 0 to ∞ and using the initial condition C ( t = 0) = 1 and15 n +1 ( t = 0) = K n ( t = 0) = 0 for all n ≥ ν a + ν l ) ¯ C − ν d ¯ K = 1 (54)( ν a + ν l ) ¯ C n − ν d ¯ K n = 0 (55) ν a ¯ C − ν d ¯ K − ¯ K = 0 (56) ν a ¯ C n − ν d ¯ K n − ¯ K n + ¯ K n − = 0 (57)for all n >
0. The bars indicate that the corresponding quantities have been integratedfrom t = 0 to ∞ . Solving these equations, we obtain for the distribution of the number ofphosphorylation events P ∞ = 1 − ν a ν a + ν l + ν d ν l (58) P ∞ n = ν a ν d ν l ( ν a + ν l ) n − ( ν a + ν l + ν d ν l ) n +1 . (59)We present an example for the distribution of the number of phosphorylation events inFig. 3a.For the average number of phosphorylation events N p , CaM and its variance var
CaM we get N p , CaM = ν a ν d ν l (60)var CaM = (cid:18) ν d + N p , CaM (cid:19) N p , CaM . (61)Compared to the toy model, the average N p , CaM contains an extra factor ν a /ν d and thevariance has an additional term 2 /ν d . On Figure 3b, we display the average average numberof phosphorylation events as a function of the detachment rate ν d and for several values of ν l .
2. Membrane-binding kinase
The calculation in the case of a membrane-binding kinase proceeds along the same linesas for the case of a cytosolic kinase. To simplify the task, we will consider the case, whenneither Ca reattaches to a kinase after detaching nor the kinase rebinds after unbinding16
10 20 3000.050.1 -3 -2 -1 -3 -2 -1 a) d) -2 -5 -3 -1 b) e) -2 -5 c) -2 -2 -2 -1 f) FIG. 3. (color online) Number of phosphorylation events. a) The distribution P ∞ n for a cytosolickinase with ν a = ν d = 1 and ν l = 0 . N p , CaM (b) and variance var p , CaM (c) of phosphorylation events for a cytosolic kinase from 10 simulations for each parameter set and according to Eqs. (60) and (61) as a function of ν d for ν a = 1 and ν l = 0 .
01 ( ∗ , yellow), 0 . (cid:52) , black), 1 ( (cid:5) , green), 10 ( (cid:3) , red), and 100 ( ◦ , blue). d)The distribution P ∞ n for a membrane-binding kinase with ν a = 10, ν d = 1, ν b = 10, ν u = 0 . ν l = 1, and D K = 1 from stochastic simulations (blue bars) and according to Eqs. (70) and (71)(red stars). Inset: semilogarithmic plot of the same data. e,f) Mean number N p , PKC and variancevar p , PKC of phosphorylation events as a function of ν d for ν a = 10, ν b = ν u = 1, D K = 0 .
01, and ν l = 0 . (cid:52) , black), 1 ( (cid:5) , green), 10 ( (cid:3) , red), and 100 ( ◦ , blue). Lines are according to Eq. (75). ∂ t C = ∂ z C − ( ν a + ν l ) C (62) ∂ t K = D K ∂ z K + ν a C − ν d K (63)˙ k = ν b K ( z = 0) − ν u k − k (64)˙ k n = − ν u k n − k n + k n − (65)˙ P = (cid:90) ∞ ( ν l C + ν d K ) d z + ν u k (66)˙ P n = ν u k n , (67)where Eqs. (65) and (67) hold for all n ≥
1. In these equations, C and K , respectively,denote the probabilities of the free Ca ion and the Ca ion attached to the kinase diffusingthe cytoplasm. Under the above conditions, phosphorylation cannot have occurred in thesestates. With k n and P n we denote the respective probabilities of having n phosphorylationevents with the kinase bound to the membrane and after the Ca ion is lost from thesystem. These equations are complemented by boundary condition for Eqs. (62) and (63).Explicitly, ∂ z C | z =0 = 0 (68) − D K ∂ z K | z =0 = − ν b K ( z = 0) + ν u k . (69)Initially, the Ca ion is localized at z = 0.For the probability distribution P ∞ n of the number of phosphorylation events we obtain P ∞ = 1 − ¯ k (70) P ∞ n = ν u (1 + ν u ) − n ¯ k (71)for n ≥ k = ν a ν b ν u √ ν a + ν l √ D K (cid:112) D K ( ν a + ν l ) + √ ν d √ D K ν d + ν b . (72)Figure 3d shows an example of the distribution.From these expressions we get for the average number of phosphorylation events and thecorresponding variance N p , PKC = 1 + ν u ν u ¯ k (73)var PKC = (cid:18) ν u − N p , PKC (cid:19) N p , PKC . (74)18imilarly, the distribution of the number of phosphorylation events can be calculated forthe full Master equation (40)-(42), see Appendix A. For the mean value, we find N p , PKC = ν b ν u (cid:32)(cid:114) ν l + (cid:114) D K ν d (cid:33) + D K ν a ν d ν l − / ν a ν d ν l (75)The expression for the variance is very lengthy and not illuminating. The mean value andvariance are shown in Figure 3e,f as a function of ν d and for various values of ν l , where thevariance has been obtained from a numerical solution of the Master equation. C. Spatial distribution of phosphorylation events by a cytosolic kinase
We now turn to the spatial distribution of phosphorylation events for a cytosolic kinase,for which we need to consider the full Master equation presented in Sect. III A. We firstsolve it numerically and then present results of a mean-field analysis.
1. Stochastic simulations
The numerical analysis of the Master equation (40)-(42) is done through simulations asdescribed in Sect. II D with appropriate modifications. In Figure 4a, b, we present examplesof the distribution P of estimated positions. It is non-Gaussian and has exponential tails.In Figure 4c, we present the estimation error as a function of the detachment rate ν d .Initially, it decreases as 1 /ν d . For ν d ≈ ν p , the dependence changes. For large enoughvalues of the loss rate ν l , the error apparently saturates as a function of ν d . Below a criticalloss rate, the error first increases before saturating. As a consequence, there is an optimalvalue of ν d for which the error is minimal. Note, however, that this minimum is not veryprominent. The dependence of the error on the loss rate ν l is similar, see Fig. 4d: initially, itdecreases as 1 /ν l and then saturates. Saturation occurs for ν l > ν p for values of ν d (cid:46) ν d , saturation is observed for increasing values of ν l . In contrastto the dependence on ν d , our simulations do not indicate the existence of an optimum lossrate that would minimize the error. 19 a) -20 -10 0 1000.020.040.06 b) -2 -2 -2 -2 c) d) FIG. 4. (color online) Spatial distribution of phosphorylation events for a cytosolic kinase. a, b)Distributions P (ˆ x ) of estimated Ca entry site from 10 numeric simulations for ν a = 1, ν d = 1(a) and ν a = 10, ν d = 0 . −| x | /λ ) / λ ,yellow) fit. Other parameter values: ν l = 1 and D K = 1. c) Estimation error as a function of ν d . d) Estimation error as a function of ν l . In (c,d) lines represent the mean-field result Eq. (83).Parameter values are ν a = 10, D K = 0 .
01 and ν l = 100 ( ◦ , blue), 10 ( (cid:3) , red), 1 ( (cid:5) , green), 0 . (cid:52) ,black) (c) and ν d = 100 ( ◦ , blue), 10 ( (cid:3) , red), 1 ( (cid:5) , green), 0 . (cid:52) , black) (d).
2. Mean-field analysis
We now perform a mean-field analysis similar to Sect. II E. Let p C denote the probabilityof finding a free Ca , that is, p C ( x, t ) = (cid:82) D n ( ξ ) P i [ n ( ξ ); x, t ]. Analogously, p K ( x, t )denotes the probability of finding a kinase with the Ca bound at position x at time t and p ( t ) the probability that there is no Ca in the system at time t . These quantities obey ∂ t p C = ∂ x p C + ν d p K − ( ν a + ν l ) p C (76) ∂ t p K = D K ∂ x p K − ν d p K + ν a p C (77)˙ p = ν l (cid:90) d x p C . (78)The normalization condition reads (cid:82) ( p C + p K ) d x + p = 1 and the initial condition is p C ( x, t = 0) = δ ( x ).Let us focus on the limit t → ∞ . Similar to the mean-field analysis of the toy model, thedistribution P of having n phosphorylation events at x is given by P ( n, x ) = 1 n ! ¯ p K ( x ) n exp {− ¯ p K ( x ) } , (79)20here the barred quantities indicate as above time-integrated quantities, for example,¯ p K ( x ) = (cid:82) ∞ p K ( x, t ) d t . To obtain the distribution of expected phosphorylation events, weare then left with solving − ∂ x ¯ p C + ν d ¯ p K − ( ν a + ν l ) ¯ p C (80)0 = D K ∂ x ¯ p K − ν d ¯ p K + ν a ¯ p C . (81)In Fourier space, we obtain¯ p K,q = ν a (cid:8) D K q + [ D K ( ν a + ν l ) + ν d ] q + ν l ν d (cid:9) − . (82)We do not need the expression for ¯ p C,q , as the expected ˆ n number of phosphorylation eventsin the limit t → ∞ is ˆ n ( x ) = ¯ p K ( x ). Using (cid:82) ˆ n ( x ) dx = ¯ p K, and (cid:82) x ˆ n ( x ) dx = − ¯ p (cid:48)(cid:48) K, aswell as expression (38) for determining the estimation error, we finally obtain for the errorin the mean-field limit (cid:96) = 2 (cid:26) (cid:96) C + (cid:96) K (cid:18) ν a ν l (cid:19)(cid:27) , (83)where (cid:96) C = 1 /ν l is the diffusion length of Ca and (cid:96) K = D K /ν d (we recall that in therescaled units used here, ν p = D C = 1).The mean-field result reproduces some of the features presented by the simulation results,see Fig. 4c,d: As for the stochastic simulations, the error decays inversely proportional withthe loss- and the detachment rates if ν l , ν d (cid:28)
1. Furthermore, the error saturates if theserates are large with (cid:96) → D K /ν d for ν l → ∞ and (cid:96) → /ν l for ν d → ∞ . However,neither as a function of ν d nor of ν l does the mean-field calculation indicate optimal ratesthat would minimize the error. The mean-field estimation error agrees quantitatively withthe simulation result in the limits of small loss rates ν l and large detachment rates ν d . D. Spatial distribution of phosphorylation events by a membrane-binding kinase
1. Stochastic simulations
In the case of a membrane-binding kinase, the presence of a boundary at z = 0 requiresspecial attention in the stochastic simulation. Whenever the system is in a state, wherethe Ca ion is diffusing in the cytoplasm, or when the kinase is bound to the membrane21nd active, we use a Gillespie-like algorithm as explained for the toy model. If the kinase isbound to the membrane at x and thus active, it can either phosphorylate or unbind from themembrane. In the first case, we record the position of the phosphorylation event, otherwisethe system state is changed and the new coordinates of the now unbound kinase are ( x, ion, then the new position ( x new , z new ) isdetermined as in Sec. II D. Should z new <
0, which is outside the considered domain, thenthe particle is assumed to have been reflected and the z coordinate of the Ca is set to − z new .In case, the Ca ion is attached to the kinase, which itself is residing in the cytoplasm,then binding to the membrane needs to be considered. To do so, we use in this state ascheme with continuous space and discrete time steps of length ∆ t and employ the methodspresented in Ref. [26]: During each time step, we first determine the new position ( x new , z new )of the particle at t + ∆ t , as explained in Sec. II D. If z new <
0, then the kinase has “crossed”the membrane and one has to determine, whether it bound to the membrane during thisprocess. To this end, a new random number between 0 and 1 is drawn. If it is smaller than1 − ν b √ π ∆ t/ (2 √ D K ), then the kinase has not bound to the membrane, but instead wasreflected and the new position is ( x new , − z new ). We then determine if the Ca has detachedand change the state if necessary. In the opposite case, the kinase binds to the membrane at( x new ,
0) and the system state is changed accordingly. Even if z new >
0, the kinase might stillhave bound to the membrane. To determine, whether this happened a random number be-tween 0 and 1 is drawn. If it is smaller than exp {− ( z new z old ) / ( D K ∗ ∆ t ) } ν b √ π ∆ t/ (2 √ D K ),then the kinase bound to the membrane at ( x new ,
0) [27]. If the kinase has not bound to themembrane, we check whether the Ca detached and change the state if necessary. The sizeof the time step dt is chosen to be ∆ t = 0 . / max { ν b , ν d } .In Figure 5a, b two examples of the distribution of estimated positions are shown. Asin the previous cases, the distributions are not Gaussian, but instead have an exponentialtail. The dependence of the estimation error on the detachment rate ν d and the loss rate ν l are shown in Fig. 5c,d. Overall, the behavior is similar to the case of a cytosolic kinase:After an initial decrease of the error with ν d and ν l , the error saturates. As a function of ν d , saturation occurs around ν d ≈ ν p . In contrast to the cytosolic kinase, a clear minimumof the error as a function of ν d cannot be detected even for small loss rates. Finally, let usnote that the estimation error is independent of the membrane binding and unbinding rates22 ) b) -2 -2 -2 -2 c) d) -5 0 500.10.20.3 -5 0 500.10.20.3 -1 FIG. 5. (color online) Spatial distribution of phosphorylation events for a membrane-bindingkinase. a, b) Distributions P (ˆ x ) of estimated Ca entry site from 10 numeric simulations for ν a = 1 (a) and ν a = 10 (b). Lines show a Gaussian (red) and an exponential (exp( −| x | /λ ) / λ ,yellow) fit. Other parameter values: ν d = 1, ν l = 1, ν b = 1, ν u = 1, and D K = 1. c) Estimationerror as a function of ν d . d) Estimation error as a function of ν l . Inset: estimation error as afunction of ν b /ν u . In (c,d) lines represent the mean-field result Eq. (83). Parameter values are ν a = 10, D K = 0 .
01 and ν l = 100 ( ◦ , blue), 10 ( (cid:3) , red), 1 ( (cid:5) , green), 0 . (cid:52) , black) (c) and ν d = 100 ( ◦ , blue), 10 ( (cid:3) , red), 1 ( (cid:5) , green), 0 . (cid:52) , black) (d). ν b and ν u , as long as they have finite values, see Fig. 5d, inset.
2. Mean-field analysis
The mean-field analysis proceeds along the same lines as for the cytosolic kinase. Asabove, let p C denote the probability of finding a free Ca , p C ( x, z, t ) = (cid:82) D n ( ξ ) P i [ n ( ξ ); x, z, t ].Analogously, p K denotes the probability distribution for a Ca bound kinase in the cytosol, p k the one for the membrane-bound kinase, and p ( t ) the probability that the Ca is lostfrom the system at time t . These quantities obey ∂ t p C = (cid:0) ∂ x + ∂ z (cid:1) p C + ν d p K − ( ν a + ν l ) p C (84) ∂ t p K = D K (cid:0) ∂ x + ∂ z (cid:1) p K − ν d p K + ν a p C (85) ∂ t p k = ν b p K ( z = 0) − ν u p k (86)˙ p = ν l (cid:90) d x (cid:90) d z p C . (87)23hese equations are complemented by boundary conditions at the membrane. Explicitly, ∂ z p C | z =0 = 0 and − D K ∂ z p K | z =0 = − ν b p K ( z = 0) + ν u p k . The normalization conditionreads (cid:82) d x (cid:82) d z ( p C + p K ) + (cid:82) p k d x + p = 1 and the initial condition is p C ( x, z, t = 0) = δ ( x ) δ ( z ). Finally, the distribution of phosphorylation events in space for t → ∞ is given byˆ n ( x ) = (cid:82) ∞ p k ( x, t ) d t .To obtain the latter, we first integrate Eqs. (84)-(86) with respect to time, which yields − δ ( x ) δ ( z ) = (cid:0) ∂ x + ∂ z (cid:1) ¯ p C + ν d ¯ p K − ( ν a + ν l ) ¯ p C (88)0 = D K (cid:0) ∂ x + ∂ z (cid:1) ¯ p K − ν d ¯ p K + ν a ¯ p C (89)0 = ν b ¯ p K ( z = 0) − ν u ¯ p k , (90)where the bars indicate the time-integrated quantities as above. From Eq. (90) the boundarycondition for ¯ p K at z = 0 is seen to be ∂ z ¯ p K | z =0 = 0. Furthermore, it shows that ˆ n ( x ) = ν b ¯ p K ( x, z = 0) /ν u . The solution for ¯ p C and ¯ p K is easiest obtained after performing a Fouriertransform with respect to x and a cosine transform with respect to z . It yields − − (cid:0) q + k (cid:1) ¯ p C,qk + ν d ¯ p K,qk − ( ν a + ν l ) ¯ p C,qk (91)0 = − D K (cid:0) q + k (cid:1) ¯ p K,qk − ν d ¯ p K,qk + ν a ¯ p C,qk , (92)where the indices q and k denote the wavenumbers in x - and z -direction, respectively. Thesolution for the time-integrated distribution of the cytosolic kinase bound to Ca is¯ p K,qk = ν a (cid:8)(cid:2) D K (cid:0) q + k (cid:1) + ν d (cid:3) (cid:0) q + k + ν a + ν l (cid:1) − ν d ν a (cid:9) − . (93)From this expression, we eventually get for the error (cid:96) = 12 (cid:2) (cid:96) + (cid:96) C (cid:96) K (cid:3) . (94)Using this expression for the estimation error, we can write the average number of phospho-rylation events, Eq. (75) as N p , PKC ≡ (cid:90) ˆ n ( x ) d x = ν b ν u (cid:2) (cid:96) + (cid:96) C (cid:96) K (cid:3) − / ν a ν d ν l . (95) E. Comparison between the two scenarios
The simulation results show that the estimation error of the measured position for thecytosolic kinase is always larger than for the membrane-binding kinase, if we compare simula-tions with the same parameter values, see Fig. 6. This result also obtained by the mean-field24 -2 FIG. 6. (color online) Ratio of the estimation error (cid:96) of a cytosolic and (cid:96) of a membrane-binding kinase as a function of ν d . Parameter values are ν a = 10, ν b = 1, ν u = 1, D K = 0 .
01 and ν l = 0 . (cid:52) , black), 1 ( (cid:5) , green), 10 ( (cid:3) , red), 100 ( ◦ , blue). expressions for (cid:96) CaM and (cid:96)
PKC . However, the mean-field result for the estimation error ratiodoes not represent the functional dependence of the ratio on ν d well. Instead, it is ratherconstant with a value that is close to the maximal value of the ratio in the interval of ν d displayed on Fig. 6.Since for a membrane-binding kinase the estimation error is independent of the membranebinding and unbinding rates, ν b and ν u , we cannot meaningfully compare the number ofphosphorylation events for the two scenarios. Note, however, that by increasing the ratio ν b /ν u , the number of phosphorylation events can be increased without affecting the accuracyof position estimate in the case of a membrane-binding kinase. This is another advantageof membrane-binding kinases over cytosolic kinases: signals will be transmitted with higherfidelity in case more target proteins are phosphorylated. IV. RESPONSES TO CA PUFFS
We will now turn to situations, in which more than one Ca is present in the system.Since all particles are independent of each other, one might expect that the error (cid:96) simplyscales as N − / if N Ca is the number of Ca ions. However, since only those events arecounted in which a phosphorylation took place, this expectation is not met. Following thepresentation of simulation results, we will apply a mean-field ansatz to express the error for25 ) b)c) d) FIG. 7. (color online) Estimation error for a Ca ion puff. a, b) Estimation error as a functionof the number N Ca of Ca ions in a puff from simulations (blue circles) and the mean-field result(106) (blue line) for a cytosolic (a) and a membrane-binding kinase (b). The green dashed linein (b) is a fit of Eq. (106) to the simulation data with rescaled parameters (cid:96) and (cid:104) n (cid:105) . c, d)Estimation error for a puff of N Ca = 1000 Ca ions as a function of the loss rate from simulations(blue circles) and Eq. (106) (red full line) for a cytosolic (c) and a membrane-binding kinase (d).Parameters are D K = 0 . ν a = 1 (a,d), ν a = 10 (b), and ν a = 0 .
01 (c), ν d = 100 (a) and ν d = 1(b-d), ν l = 10 (a,b), ν b = 1 (b,d), and ν u = 1 (b,d). a puff in terms of the average phosphorylation profile ˆ n for a single Ca ion. A. Stochastic simulations
Simulations are done as described in Sects. III C 1 and III D 1. For a puff of N Ca Ca ions, we ran N simulations and recorded the positions of all phosphorylation events duringthese simulations. We then obtained the estimated position by calculating their average.For each data point we performed at least 2 · simulations.In Figure 7a,b, we present the estimation error for a puff as a function of the number N Ca of Ca ions per puff. For both kinds of kinases, the estimation error decreases monotonicallywith increasing N Ca . The data points fall onto a sigmoidal curve: initially the accuracy ofthe estimate increases less than for larger values of N Ca . Note, that with increasing numberof Ca ions, the gap between the estimation error for the membrane-binding kinase andthe cytosolic kinase gets wider. 26n the case of a cytosolic kinase and for a Calcium puff with N Ca = 1000, the dependenceof the estimation error as a function of ν l changes qualitatively with respect to the caseof a single Ca ion: In the latter case the error decreases monotonically whereas for thepuff a minimum and a maximum exist, respectively, around ν l = 1 and 2, see Fig. 7c. Incontrast, for a membrane-binding kinase the dependence remains monotonically decreasing,see Fig. 7d. However, the estimation error saturates for smaller values of ν l compared to thecase of a single Ca ion. B. Mean-field analysis
We will now apply the mean-field ansatz introduced in Sect. II E to the case of Ca puffs, where we focus directly on the limit t → ∞ . For the general expression (15) forthe estimation error, we need the probability P [ n ( ξ )] for a particular realization of thephosphorylation profile n ( ξ ). In the mean-field approximation, the probability P ( n, ξ ) ofhaving n phosphorylation events at position ξ is given by Eq. (34), where in the generalcase ¯ p a has to be replaced by the time-integrated probability to find the particle in thephosphorylating state at position ξ . Explicitly, for a cytosolic kinase it is ¯ p K , whereas forthe membrane-binding kinase it is ¯ p k . In contrast to the case of a single Ca ion, we cannotuse expression (38) for the error, because it does not depend on the average number ofphosphorylation events (cid:104) n (cid:105) , which obviously increases with the number of Ca ions N Ca .Therefore, we will now calculate the exact expression of the estimation error (15) in themean-field approximation.From P ( n, ξ ), we obtain the mean-field probability P mf [ n ( ξ )] of a particular realization n ( ξ ) through P mf [ n ( ξ )] = N (cid:89) ξ P ( n, ξ ) (96)for all n ( ξ ) with (cid:82) d ξ n ( ξ ) (cid:54) = 0. The normalization factor N assures that (cid:82) D (cid:48) n ( ξ ) P [ n ( ξ )] =1, where the prime indicates summation over all distributions n ( ξ ) except for n ( ξ ) ≡ P mf [ n ( ξ )] = e − (cid:82) d ξ (cid:48) ˆ n ( ξ (cid:48) ) − e − (cid:82) d ξ (cid:48) ˆ n ( ξ (cid:48) ) (cid:89) x n ( ξ )! (ˆ n ( ξ )) n ( ξ ) , (97)27here ˆ n ( ξ ) is the average phosphorylation profile, such that P mf [ n ( ξ )] = e − N p − e − N p (cid:89) ξ n ( ξ )! (ˆ n ( ξ )) n ( ξ ) . (98)Using this probability distribution in expression (15) for the estimation error, we get (cid:96) = (cid:90) D (cid:48) n ( ξ ) (cid:82) d ξ (cid:82) d ξ ξ ξ n ( ξ ) n ( ξ ) (cid:0)(cid:82) d ξ (cid:48) n ( ξ (cid:48) ) (cid:1) P [ n ( ξ )] (99)= e − N p − e − N p (cid:89) ξ ∞ (cid:88) n ( ξ )=0 (cid:82) d ξ (cid:82) d ξ ξ ξ n ( ξ ) n ( ξ ) (cid:0)(cid:82) d ξ (cid:48) n ( ξ (cid:48) ) (cid:1) n ( ξ )! ˆ n ( ξ ) n ( ξ ) (100)= e − N p − e − N p (cid:89) ξ ∞ (cid:88) n ( ξ )=0 (cid:90) d ξ (cid:90) d ξ ξ ξ ˆ n ( ξ ) δδ ˆ n ( ξ ) ˆ n ( ξ ) δδ ˆ n ( ξ ) ˆ n ( ξ ) n ( ξ ) n ( ξ )! 1 (cid:0)(cid:82) d ξ (cid:48) n ( ξ (cid:48) ) (cid:1) (101)= e − N p − e − N p (cid:90) d ξ (cid:90) d ξ ξ ξ ˆ n ( ξ ) δδ ˆ n ( ξ ) ˆ n ( ξ ) δδ ˆ n ( ξ ) (cid:89) ξ ∞ (cid:88) n ( ξ )=0 ˆ n ( ξ ) n ( ξ ) n ( ξ )! 1 (cid:0)(cid:82) d ξ (cid:48) n ( ξ (cid:48) ) (cid:1) (102)= e − N p − e − N p (cid:90) d ξ ξ ˆ n ( ξ ) δδ ˆ n ( ξ ) (cid:89) ξ ∞ (cid:88) n ( ξ )=0 ˆ n ( ξ ) n ( ξ ) n ( ξ )! 1 (cid:0)(cid:82) d ξ (cid:48) n ( ξ (cid:48) ) (cid:1) (103)= e − N p − e − N p (cid:90) d ξ ξ ˆ n ( ξ ) δδ ˆ n ( ξ ) ∞ (cid:88) N =1 (cid:0)(cid:82) d ξ ˆ n ( ξ ) (cid:1) N N ! N (104)= e − N p − e − N p (cid:82) d ξ ξ ˆ n ( ξ ) (cid:82) d ξ ˆ n ( ξ ) ∞ (cid:88) N =1 (cid:0)(cid:82) d ξ ˆ n ( ξ ) (cid:1) N N ! N (105)= (cid:96) e − N p − e − N p ∞ (cid:88) N =1 N N p N ! N . (106)In the final expression (cid:96) is the estimation error in the mean-field approximation, see Eq. (38).For N p (cid:28)
1, we have (cid:96) = (cid:96) as announced in Sect. II E. In the case of a Ca puff with N Ca Ca ions, N p = N Ca (cid:104) n (cid:105) , where (cid:104) n (cid:105) is the average number of phosphorylation eventsfor one Ca ion, because we assume that all Ca ions are independent of each other andthat there is an excess of kinases. As the final expression for (cid:96) shows, the dependenceon N Ca is more complicated than the usual 1 /N Ca -dependence of the variance in case onehas N Ca independent measurements. The reason is that not all Ca ions produce at leastone phosphorylation event, such that the position of Ca influx cannot be estimated for allCa ions. However, for N Ca (cid:29)
1, we find again (cid:96) ∼ (cid:96) N Ca .28e will now use this general expression in combination with the results for ˆ n ( x ) obtainedfrom the mean-field analysis. For the case of a cytosolic kinase, we obtain (cid:96) , puff = 2 (cid:26) (cid:96) C + (cid:96) K (cid:18) ν a ν l (cid:19)(cid:27) e − N Ca νaνdνl − e − N Ca νaνdνl (cid:88) N N ! N N N Ca ν Na ν Nd ν Nl . (107)Application of the general expression (106) to the case of a membrane-binding kinase istricky. This expression depends on N p , PKC , which in turn can be changed by changing either N Ca or ν b /ν u . In the simulations, however, only changing N Ca affects the estimation error,whereas it is independent of the ration ν b /ν u as we had seen above for the case of a singleCa ion. Since these two effects are not separated in the mean-field expression one cannotexpect it to describe the dependence of the error on the number of Ca ions. We thusrefrain from giving the mean-field result for the case of a membrane-binding kinase.In Figure 7, we present the estimation error for a puff obtained from the mean-fieldtreatment as a function of N Ca . In case of a cytosolic kinase, where the error for a singleCa ion is given by the mean-field result, the dependence on the number of Ca matchesthe simulation results perfectly. From Equation (107), we see that the error only decreasessignificantly, when N Ca ∼ ν d ν l /ν a . For the parameters chosen in Fig. 7, we get ν d ν l /ν a =1000, which matches well the simulation data. As we have argued before, we cannot expectthe mean-field error for puffs to describe simulation results for a membrane-binding kinase.For the parameters chosen in Fig. 7b, it is indeed off. However, by appropriately rescaling (cid:96) for N Ca = 1 and (cid:104) n (cid:105) , the expression (106) provides a fit to the data. In Figure 7d, we seethat for increasing values of ν l , we obtain agreement between Eq. (106) and the simulationresults. These results show that the mean-field expression does capture important aspectsof the estimation error even in the case of a membrane-binding kinase. V. ESTIMATING THE SITE OF CA RELEASE IN PRESENCE OF BACK-GROUND PHOSPHORYLATION
Living cells have a cytosolic Ca concentration of roughly 100 nM [28]. Consequently,a fraction of calmodulin and PKC α are active even in absence of an external signal. Howdoes the corresponding background phosphorylation affect the accuracy of the estimatedposition of the Ca release site? On general grounds, cells might be expected to suppressthe influence of the background by employing a threshold mechanism: a cellular response29s only elicited if the number of phosphorylation events exceeds a certain value. Still, it isinteresting to account explicitly for background phosphorylation in our analysis.In the presence of background phosphorylation, our theoretic approach has to be modifiedto some extent. Above, we considered the distribution of all phosphorylation sites that weregenerated by a Ca ion or -puff, independently of when they occurred. If we appliedthe same approach in presence of background phosphorylation, then the background wouldalways outcompete the signal. We thus introduce a rate ν dp of dephosphorylation of targetproteins. In this way, a phosphorylated protein contributes only during a time 1 /ν dp tothe cellular response, which we still take to be given by the spatial distribution of thephosphorylation events. In general, it is now a time-dependent quantity, whereas before,we considered the accumulated distribution of all phosphorylation events following a signal.We will assume that phosphorylated proteins do not move. Let us denote the number ofphosphorylated proteins resulting from the signal by N s and those from the background by N bg . Let us note again that N s depends on time and N bg does not. To arrive at a singlenumber for the error in estimating the position of Ca release, we consider the time pointat which N s is maximal.We calculate the error of the position estimate by a weighted mean of the error from thephosphorylated proteins resulting from the signal and those from the background. Since,background phosphorylation is independent of phosphorylation in response to the signal,the corresponding variances and thus errors simply add up. The error (cid:96) s associated with theresponse to the signal is calculated as before. For the error resulting from the background,we assume that the corresponding phosphorylation events are uniformly distributed in thecell, such that (cid:96) is given by the size of the cell. The total error then is (cid:96) = N s, max N s, max + N bg (cid:96) s + N bg N s, max + N bg (cid:96) . (108)In Figure 8, we present the error for estimating the Ca release site for a cytosolic anda membrane-binding kinase in presence of background phosphorylation. As expected, if thebackground phosphorylation exceeds a certain threshold, the signal is masked and the errorequals the size of the cell, such that any information about the site of Ca release is lost.30 ) b) FIG. 8. (color online) Estimation error in presence of background phosphorylation as a function ofthe number N bg of phosphorylated kinases due to background phsophorylation for (a) a cytosolickinase with ν d = 0 . ν d = 1 (blue circles), 10 (red squares), and 100 (green stars). Other parametervalues are D K = 0 . ν a = 0 . ν l = 10, ν dp = 0 . ν b = 1, and ν u = 0 . VI. DISCUSSION
In this work, we have presented a framework for studying cellular responses to localizedsignals. For concreteness, we have considered the signal to be given by a localized Ca influx and the response to be represented by the spatial distribution of phosphorylationevents of either a cytosolic or a membrane-binding kinase. We very much simplified thecellular response. For example, we considered direct activation of the cytosolic kinase bybinding a Ca ion, although often Ca activates Calmodulin, which in turn activates thekinase. Furthermore, we assumed that the membrane-binding kinase is activated directlyby binding to the membrane, whereas the Protein Kinase C α , for example, requires bindingto Diacylglycerol (DAG) in the membrane for activation. In spite of these simplifications,we expect our general results to be valid also in more realistic situations. This holds, no-tably, for the tendency of the estimation error to decrease with increasing Ca detachment.Furthermore, membrane-binding kinases should provide a better spatial localization of thesignal than cytosolic kinases. Beyond the specific question of how cells can localize externalsignals, our framework can also be applied in various other situations, in which a stochasticbirth process is coupled to diffusion.We defined the error in determining the position at which the Ca ions entered the sys-tem. To this end, we considered the variance of the distribution of estimated positions fromthe generated phosphorylation events. The cell response does not end with the phosphory-lation of target proteins. Ultimately, the Ca signal should lead to endo- or exocytosis or31o a change in the strength of a synapse. However, it is difficult to see how the outcomes ofthese processes could spatially be more precise than the phosphorylation ’signal’ into whichthe Ca signal has been transformed.Let us estimate the error in the light of measured parameter values for the Protein KinaseC α . The diffusion constant of free cytoplasmic Ca is about 500 µ m /s [29], whereascytoplasmic PKC α has a diffusion constant on the order of 10 µ m /s [30]. With ν p ≈ ν d ≈ ν l ≈ /s [31], the estimation error for a single Ca ion is about 7 µ m.This justifies neglecting any boundaries in the lateral direction as a typical cell diameteris 50 µ m. For a Ca puff this value decreases further with an increasing number of Ca ions. With regards to our assumption that membrane-binding kinases are immobile, notethat the diffusion length on the membrane (cid:112) D K, mem /ν u ≈ . µ m. In this estimate, wetook the diffusion constant D K, mem of PKC α to be 100 times smaller than its cytoplasmicdiffusion constant [32] and ν ≈ α it should be interesting to consider its dynam-ics in more detail. In addition to adding the effect of binding DAG to PKC α , it might beinteresting to consider the full cascade of a signal activating Phospholipase C that breaksPhosphatidylinositol-biphosphate (PIP ) into DAG and Inositol-triphosphate (IP ). Thelatter diffuses and can open nearby internal Ca stores, making the activation of PKC α much more involved than considered here. In addition, PKC α forms clusters on the mem-brane that affects the lifetime of its activated state [33, 34]. Studying the dynamics of PKC α might yield insights into how its localized activation can help cells to obtain a specific re-sponse to external signals even though they get transformed into general purpose secondmessengers [35, 36]. Appendix A: Full Master equation for the number of phosphorylation events for amembrane-binding kinase
In case, Ca can reattach to the kinase after detachment and the kinase can rebindto the membrane after unbinding, the Master equation for the number of phosphorylation32vents is ∂ t C n = ∂ z C n + ν d K n − ( ν a + ν l ) C n (A1) ∂ t K n = D K ∂ z K n − ν d K n + ν a C n (A2)˙ k = ν b K ( z = 0) − ν u k − k (A3)˙ k n = ν b K n ( z = 0) − ν u k n − k n + k n − (A4)˙ P n = ν l (cid:90) ∞ C n d z (A5)where Eqs. (A1), (A2), and (A5) hold for all n ≥
0, whereas Eq. (A4) holds for all n ≥ C n and K n , respectively, denote the probabilities of having n phospho-rylation events with a free Ca ion and with it bound to a cytosolic kinase, whereas k n and P n denote the respective probabilities with the kinase bound to the membrane and after theCa ion is lost from the system. These equations are complemented by boundary conditionfor Eqs. (A1) and (A2). Explicitly, ∂ z C n | z =0 = 0 (A6) − D K ∂ z K n | z =0 = − ν b K n ( z = 0) + ν u k n . (A7)Finally, initially, the Ca is located at z = 0.In the limit t → ∞ , the distribution P ∞ n of the number of phosphorylation events is givenby P ∞ n = ν l (cid:82) ∞ ¯ C n d z , where the bar indicates integration of C n from 0 to ∞ with respectto time. For the barred quantities, the equations read ∂ z ¯ C + ν d ¯ K − ( ν a + ν l ) ¯ C = − δ ( z ) (A8) ∂ z ¯ C n + ν d ¯ K n − ( ν a + ν l ) ¯ C n = 0 (A9) D K ∂ z ¯ K − ν d ¯ K + ν a ¯ C = 0 (A10) D K ∂ z ¯ K n − ν d ¯ K n + ν a ¯ C n = 0 (A11) ν b ¯ K ( z = 0) − ν u ¯ k − ¯ k = 0 (A12) ν b ¯ K n ( z = 0) − ν u ¯ k n − ¯ k n + ¯ k n − = 0 , (A13)33here n ≥
1. Their solution can be written as¯ C ( z ) = ˜ I ( z ) + ˆ C (cid:0) λ e − λ z − λ e − λ z (cid:1) (A14)¯ K ( z ) = ν a (cid:20) I ( z ) − ˆ C (cid:18) λ D K λ − ν d e − λ z − λ D K λ − ν d e − λ z (cid:19)(cid:21) (A15)¯ C n ( z ) = ˆ C n (cid:0) λ e − λ z − λ e − λ z (cid:1) (A16)¯ K n ( z ) = − ν a ˆ C n (cid:18) λ D K λ − ν d e − λ z + λ D K λ − ν d e − λ z (cid:19) (A17)with I ( z ) = (cid:90) ∞ cos( kz )( D K k + ν d ) ( k + ν a + ν l ) − ν a ν d d k (A18)˜ I ( z ) = (cid:90) ∞ ( D K k + ν d ) cos( kz )( D K k + ν d ) ( k + ν a + ν l ) − ν a ν d d k. (A19)Using Equations (A3) and (A4) as well as the boundary conditions (A7), we findˆ C = ν b I (0)( λ − λ ) [(1 + ν u ) A + B ] (A20)ˆ C = − ν b ν u A (1 + ν u ) A + B ˆ C (A21)ˆ C n = (cid:20) A + B (1 + ν u ) A + B (cid:21) n − ˆ C (A22)with A = D K λ λ ( λ + λ ) (A23) B = ν b ν a ν d ( λ λ + ν a + ν b ) . (A24)This eventually leads to the probability distribution P ∞ = 1 − E (A25) P ∞ n = ν u A (1 + ν u ) A + B (cid:20) A + B (1 + ν u ) A + B (cid:21) n − E (A26)for n ≥ E = ν l ν b (1 + ν u ) A + B λ + λ λ λ I (0) . (A27)34 CKNOWLEDGMENTS
We acknowledge funding through SFB 1027 by Deutsche Forschungsgemeinschaft. Thecomputations were performed at University of Geneva on the Baobab cluster. [1] B. Alberts, A. Johnson, J. Lewis, M. Raff, K. Roberts, and P. Walter,
Molecular Biology ofthe Cell , 5th ed., edited by B. Alberts (Garland Science, 2008).[2] P. A. Janmey and C. A. McCulloch, Annu. Rev. Biomed. Eng. , 1 (2007).[3] K. F. Swaney, C.-H. Huang, and P. N. Devreotes, Annu. Rev. Biophys. , 265 (2010).[4] W.-L. Ng and B. L. Bassler, Annu. Rev. Genet. , 197 (2009).[5] A. J. Engler, S. Sen, H. L. Sweeney, and D. E. Discher, Cell , 677 (2006).[6] H. C. Berg and E. M. Purcell, Biophys. J. , 193 (1977).[7] J. R. Maddock and L. Shapiro, Science , 1717 (1993).[8] D. Bray, M. D. Levin, and C. J. Morton-Firth, Nature , 85 (1998).[9] T. A. Duke and D. Bray, Proc. Natl. Acad. Sci. USA , 10104 (1999).[10] B. A. Mello and Y. H. Tu, Proc. Natl. Acad. Sci. USA , 8223 (2003).[11] B. A. Mello, L. Shaw, and Y. H. Tu, Biophys. J. , 1578 (2004).[12] B. A. Mello and Y. H. Tu, Proc. Natl. Acad. Sci. USA , 17354 (2005).[13] R. G. Endres and N. S. Wingreen, Proc. Natl. Acad. Sci. USA , 13040 (2006).[14] W. Bialek and S. Setayeshgar, Proc. Natl. Acad. Sci. USA , 10040 (2005).[15] W. Bialek and S. Setayeshgar, Phys. Rev. Lett. , 258101 (2008).[16] P. J. M. Van Haastert and P. N. Devreotes, Nat. Rev. Mol. Cell Bio. , 626 (2004).[17] R. G. Endres and N. S. Wingreen, Proc. Natl. Acad. Sci. USA , 15749 (2008).[18] W.-J. Rappel and H. Levine, Phys. Rev. Lett. , 228101 (2008).[19] W.-J. Rappel and H. Levine, Proc. Natl. Acad. Sci. USA , 19270 (2008).[20] Y. Xiong, C.-H. Huang, P. A. Iglesias, and P. N. Devreotes, Proc. Natl. Acad. Sci. USA ,17079 (2010).[21] K. Deisseroth, H. Bito, and R. W. Tsien, Neuron , 89 (1996).[22] D. G. Wheeler, C. F. Barrett, R. D. Groth, P. Safa, and R. W. Tsien, J. Cell Biol. , 849(2008).
23] M. L. Kapsenberg, Nat. Rev. Immunol. , 984 (2003).[24] E. A. Nalefski and A. C. Newton, Biochemistry , 13216 (2001).[25] There are also Calmodulin phosphatases, which dephosphorylate target proteins. Cum granosalis the results developed in this work also apply to phosphatases.[26] R. Erban and S. J. Chapman, Phys. Biol. , 16 (2007).[27] S. S. Andrews and D. Bray, Phys. Biol. , 137 (2004).[28] R. Milo and R. Phillips, Cell Biology by the Numbers , 1st ed. (Garland Science, 2015).[29] B. S. Donahue and R. F. Abercrombie, Cell Calcium , 437 (1987).[30] M. Schaefer, N. Albrecht, T. Hofmann, T. Gudermann, and G. Schultz, FASEB J. , 1634(2001).[31] G. D. Smith, J. E. Keizer, M. D. Stern, W. J. Lederer, and H. Cheng, Biophys. J. , 15(1998).[32] J. Lippincott-Schwartz, E. Snapp, and A. Kenworthy, Nat. Rev. Mol. Cell Bio. , 444 (2001).[33] M. Bonny, X. Hui, J. Schweizer, L. Kaestner, A. Zeug, K. Kruse, and P. Lipp, Sci. Rep. ,36028 (2016).[34] C. J. Swanson, R. F. Sommese, K. J. Petersen, M. Ritt, J. Karslake, D. D. Thomas, andS. Sivaramakrishnan, PLoS ONE , e0162331 (2016).[35] J. H. Horne and T. Meyer, Science , 1690 (1997).[36] C. Maasch, S. Wagner, C. Lindschau, G. Alexander, K. Buchner, M. Gollasch, F. C. Luft,and H. Haller, FASEB J. , 1653 (2000)., 1653 (2000).