New Research In
Physical Sciences
Social Sciences
Featured Portals
Articles by Topic
Biological Sciences
Featured Portals
Articles by Topic
 Agricultural Sciences
 Anthropology
 Applied Biological Sciences
 Biochemistry
 Biophysics and Computational Biology
 Cell Biology
 Developmental Biology
 Ecology
 Environmental Sciences
 Evolution
 Genetics
 Immunology and Inflammation
 Medical Sciences
 Microbiology
 Neuroscience
 Pharmacology
 Physiology
 Plant Biology
 Population Biology
 Psychological and Cognitive Sciences
 Sustainability Science
 Systems Biology
Toward a microscopic model of bidirectional synaptic plasticity

Contributed by Leon N Cooper, May 29, 2009 (received for review January 15, 2009)
Abstract
We show that a 2step phospho/dephosphorylation cycle for the αamino3hydroxy5methyl4isoxazole proprionic acid receptor (AMPAR), as used in in vivo learning experiments to assess longterm potentiation (LTP) induction and establishment, exhibits bistability for a wide range of parameters, consistent with values derived from biological literature. The AMPAR model we propose, hence, is a candidate for memory storage and switching behavior at a molecularmicroscopic level. Furthermore, the stochastic formulation of the deterministic model leads to a mesoscopic interpretation by considering the effect of enzymatic fluctuations on the Michelis–Menten average dynamics. Under suitable hypotheses, this leads to a stochastic dynamical system with multiplicative noise whose probability density evolves according to a Fokker–Planck equation in the Stratonovich sense. In this approach, the probability density associated with each AMPAR phosphorylation state allows one to compute the probability of any concentration value, whereas the Michaelis–Menten equations consider the average concentration dynamics. We show that bistable dynamics are robust for multiplicative stochastic perturbations and that the presence of both noise and bistability simulates LTP and longterm depression (LTD) behavior. Interestingly, the LTP part of this model has been experimentally verified as a result of in vivo, onetrial inhibitory avoidance learning protocol in rats, that produced the same changes in hippocampal AMPARs phosphorylation state as observed with in vitro induction of LTP with highfrequency stimulation (HFS). A consequence of this model is the possibility of characterizing a molecular switch with a defined biochemical set of reactions showing bistability and bidirectionality. Thus, this 3enzymesbased biophysical model can predict LTP as well as LTD and their transition rates. The theoretical results can be, in principle, validated by in vitro and in vivo experiments, such as fluorescence measurements and electrophysiological recordings at multiple scales, from molecules to neurons. A further consequence is that the bistable regime occurs only within certain parametric windows, which may simulate a “historydependent threshold”. This effect might be related to the Bienenstock–Cooper–Munro theory of synaptic plasticity.
Recent experiments (1, 2) have demonstrated the possibility of describing the phosphorylation process of α amino3hydroxy5methyl4isoxazole propionic acid (AMPA) receptors (AMPARs) by using a 2step chain of enzymatic reactions. The enzymatic reactions are 2 phosphorylations and 2 dephosphorylations, catalyzed by 2 specific kinases and 1 nonspecific phosphatase. A peculiarity of these enzymes is that their activities are directly or indirectly controlled by calcium concentration (Ca ^{++}) (3 –5). The intracellular calcium concentration depends on the transmembrane potential and consequently on the stimulus protocol. The relation between lowfrequency stimulations (LFS), highfrequency stimulation (HFS) and AMPA phosphorylation state has been experimentally shown (1), as well as the identification of the involved enzymes: 2 specific kinases, protein kinase cAMPdependent (PKA) and calcium calmodulin kinase of tpe II (CaMKII), and 1 nonspecific phosphatase, the protein phosphatase 1 (PP1) (2). The proposed AMPAR phospho/dephosphorylation cycles were basically of 2 types, one based on 3 phosphorylation states (A, A _{p} and A _{p} ^{p}) and the other based on 4 states: (A, A _{p}, A ^{p}, A _{p} ^{p}), where A, A _{p}, A ^{p}, and A _{p} ^{p} indicate respectively the AMPAR not phosphorylated, phosphorylated at Serine 831, phosphorylated at Serine 845, and phosphorylated at both sites, due to the observation of distinct phosphorylation states in the GluR1 subunit of AMPAR (4, 5). Recently it has been shown that the kinetic scheme based on 3 AMPAR phosphorylation states can account for a real learning paradigm and that the AMPA phosphorylation state is not only related to longterm potentiation (LTP) but also marks synapses in a historydependent way (6). More precisely, 1trial inhibitory avoidance learning in rats produces the same changes in hippocampal glutamate receptors as induction of LTP with HFS and causes a spatially restricted increase in the amplitude of evoked synaptic transmission. It is therefore of great interest to investigate whether this kinetic scheme can exhibit bistable behavior and if this can be mapped onto potentiated and depotentiated states such as LTP and longterm depression (LTD). During the last few years, much effort has been invested in developing bistable models able to model memory and learning phenomena (3, 7 –9) both from a mathematical and biophysical–molecular point of view. This kind of enzymatic chain, as well as the simplest switch based on a single kinase and phosphatase, has also been extensively studied in noisy environments to explain switching phenomena in signaltransduction pathways and hypersensitivity in phenotypic switching (10 –12).
In this paper, we show that a combined deterministic and stochastic approach to the 2step model of enzymatic reactions acting on AMPARs allows us to simulate historydependent bidirectional synaptic plasticity. First we prove the existence of bistable behavior in the deterministic model that could explain an LTP mechanism. Then we show that the presence of fluctuations in the enzyme concentrations allows us to study the model as a stochastic dynamical system, whose diffusion regime simulates the existence of LTP behavior. A more realistic microscopic description of the enzymatic reactions should be based on a single molecule model such as the chemical master equation (13). But the stochastic differential equation model allows a mesoscopic description that generalizes the average approach based on Michaelis–Menten (MM) models. Our approach is consistent with the assumption that the steadystate hypothesis at lowenzyme concentrations should take into account fluctuation effects. Then multiplicative noise is introduced in the reaction velocities, and we associate a probability density with each value of the complexes' and products' concentrations, which can be interpreted as reaction variables. The reaction progression is described as a transition between unbounded and bounded states, as in Markov modeling of ionchannel transitions (14). The stochastic model of the AMPA phosphorylation cycle can be used both to study the average behavior of the system and to look for phase transitions as result of an external perturbation.
AMPA Phosphorylation Dynamics
The kinetic scheme of Fig. 1 is composed of 4 enzymatic MM reactions: E + S ⇋ ES → E + P where k _{f},k _{b} and k _{c} are, respectively, the forward, backward, and catalytic constants. k _{m} = (k _{b} + k _{c})/k _{f} is the MM constant and E, S, ES, and P are the enzyme, the substrate, the enzyme–substrate complex and the product concentrations, respectively. By using the steadystate assumption Ė S = 0 and the total enzyme conservation law (E + ES = E _{T}), the scheme in Fig. 1 can be written in terms of reaction velocities. For the sake of simplicity, we associate the index 1 to the PKA reaction, 2 to the first PP1 dephosphorylation, 3 to the CaMKII and 4 to the second PP1 reaction. With these assumptions, we write the deterministic equations for the AMPARs' phosphorylation/dephosphorylation cycle as follows (see also SI Appendix): where the velocities are (see SI Appendix):
PKA _{T}, PP1_{T}, and CaMKII _{T} respectively denote the total amount of PKA, PP1, and CaMKII. The quantities A, A _{P}, and A _{P} ^{P} are the relative concentrations of AMPA channels in each of the 3 phosphorylation states. Then system 2 has a linear first integral of motion due to the conservation law A + A _{P} + A _{P} ^{P} = 1. The fixed points of system 2 are determined by the condition
that geometrically corresponds to the intersection of 2 conics (see SI Appendix) in the general case. To perform analytic calculations, we assume that the MM constants satisfy the symmetric condition k _{2} ^{m} = k _{4} ^{m}.
Analysis of the Model
In order to study the dynamical properties of system 2, we reduce the parameter number by setting k _{i} ^{c} = 1 for i = 1,2,3,4 in Eq. 2: This corresponds to a time scaling and a redefinition of the enzyme concentrations. Taking into account the symmetric condition and the conservation law A + A _{P} + A _{P} ^{P} = 1, we get the simplified system The A equilibria of system 3 can be easily characterized as the solutions of a cubic function of the form f(A,c _{1},c _{3}) = 0 where we introduce the enzyme ratio concentrations The ratios c _{1,3} are the control parameters for the dynamical properties given the values of the MM constants. By keeping the c _{1} parameter fixed, the ratio c _{3} modulates the strength of the doublephosphorylation reaction by simulating the change in the CaMKII _{T} concentration. In Fig. 2, we plot the equilibrium condition for 3 values of the c _{3} parameter, pointing out the existence of a generic bifurcation phenomenon when the c _{3} parameter (i.e. the CaMKII concentration) increases.
We observe that at low CaMKII concentrations we have only one stable equilibrium for A ≃ .97 and A _{P} ^{P} ≃ .02, whereas at high CaMKII concentrations we have a single equilibrium for A ≃. 03 and A _{P} ^{P} ≃. 92. But there exists an entire interval of c _{3} values, in which we have both the equilibria separated by an unstable equilibria at A ≃ .47 and A _{P} ^{P} ≃ .50. In such a case system 3 undergoes bistable dynamics. This is a general phenomenon that does not depend critically on the chosen parameter values. A numerical study for the existence of bistable dynamics in the parameter space (c _{1},c _{3}) shows that the bistability region has a large area (see Fig. 3). The area of the bistability region depends also on the MM constants and it moves toward the axes as k _{1} ^{m} becomes smaller.
The dynamics of system 3 are dissipative, so that, depending on the initial condition in the (A,A _{P} ^{P}) space, the trajectory will be attracted asymptotically by one of the stable solutions. In Fig. 4, we plot the attraction basin of the stable fixed points in the case k _{1} ^{m} = .1, k _{2} ^{m} = .05: We observe both basins have a comparable area so that the dynamics are analogous to dissipative dynamics with a double well potential. The bistability regime in system 2 can simulate LTP phenomenon by means of a hysteresis mechanism. Consider an AMPA channel population whose states are mostly not phosphorylated. An external stimulus is simulated by assuming a time dependence in the CaMKII _{T} concentration according to in system 3. LTP means that if the external stimulus is sufficiently strong and persistent, most of the channels will transit from the nonphosphorylated equilibrium to the doublephosphorylated state and they will remain in the new equilibrium even when the signal stops. In Fig. 5, we show 3 examples of solutions A(t) of system 3, where we have introduced the time dependence described by Eq. 5. In the simulation, we set CaMKII _{T} ^{0} = 2.2, using the same parameters as in Fig. 4, so that at t = 0 the frozen system has bistable dynamics. Each solution has an initial condition near the A ≃ .97 stable solution, but we used different Δ in the simulation with Ω = 10^{−2} × 2π (in arbitrary time units). Figure 5 shows the existence of a threshold value for Δ, at which the trajectory performs a transition from the attraction basin of the stable nonphosphorylated solution to the attraction basin of the doublephosphorylated stable solution. This behavior can be explained by considering the competition between the adiabatic change of the phasespace structure due to the modulation (see Fig. 5) (we remark that we set Ω ≪ 1 in the simulations) and the relaxation process to stable equilibria. By using an adiabatic approach, we can study the “frozen system” and compute the time during which only a single equilibrium is present in the phase space; if such a time is greater than the relaxation time toward the equilibrium itself, then most of the initial condition will pass from the initial nonphosphorylated equilibrium to the doublephosphorylated one.
However, in order to show the applicability of the model described by Eq. 3 to real situations, one has to study the robustness of the bistability regime when stochastic perturbations are present. This requirement is a key point if we would like to apply the model to describe AMPA channel dynamics at the microscopic level, because fluctuations of the enzyme concentration cannot be averaged. We have to generalize our dynamical model to take account of the effect of multiplicative noise on the bistability regime. This reckoning will be possible by considering system 3 to be a stochastic dynamical system according to the Stratonovich interpretation (see SI Appendix).
Multiplicative Noise Effect
Random perturbations in real systems are unavoidable, and in modeling AMPA channel dynamics the introduction of random fluctuations on the enzyme concentrations is a first step toward microscopic dynamical modeling. We explicitly study the effect of a random perturbation in the CaMKII _{T} concentration because this can be considered a control parameter of the dynamics represented in Eq. 3. As is well known, the stochastic perturbations are mathematical models of physical phenomena whose application requires the realization of specific assumptions. The justification for such assumptions will require further studies using a chemicalmasterequation approach to the modeling of AMPAchannel reactions that is beyond the scope of this paper. By relaxing the steadystate condition used to compute the reaction velocities in Eq. 2, we explicitly consider the effect of fluctuations in the CaMKII concentration (see SI Appendix). Then we substitute the CaMKII _{T} concentration in Eq. 3 with the expression where ξ(t) is a given stochastic process with zero mean value and unitary variance, and ɛ is a small parameter. A key point to note is that the fluctuation correlation time must be fast with respect to a characteristic evolution time scale of the unperturbed system; i.e., the evolution time scale of the MM equations (average dynamics) is much slower than the fluctuation correlation time scale, so that the random variables ξ(t _{1}) ξ(t _{2}) can be considered independent when t _{1}≠t _{2}. It is then possible to approximate ξ(t) by using white noise and system 3 becomes a stochastic dynamical system. But multiplicative noise requires a correct interpretation of a stochastic differential equation and the Stratonovich interpretation turns out to be more suitable for biophysical applications (see SI Appendix) (11). The dynamics of an ensemble of AMPA channels is described in a probabilistic manner by using the distribution function ρ(A, A _{P} ^{P}, t), which satisfies a Fokker–Planck equation where the velocities are defined by the unperturbed system in Eq. 2, and the diffusion coefficient is given by The longterm statistical properties are defined by the stationary solution ρ_{s}(A, A _{P} ^{P}) of Eq. 7. The fluctuations have a 2fold effect. On the one hand, they introduce a diffusion in the (A,A _{P} ^{P}) space, which implies LTD; on the other hand, they introduce an effective potential that could modify the average phase space by means of a stochastic bifurcation phenomenon (15) (see SI Appendix). In Fig. 6, we plot the stationary distribution function ρ_{s}(A, A _{P} ^{P}) on the A space, computed for different noise levels ɛ: We remark that at low noise level ɛ = .5, the distribution has a single peak at A ≪ 1. As the noise increases, the stationary distribution bifurcates into a bimodal distribution, and a new peak at A ≃ 1 appears. Therefore, the statistical properties of the system are described by a bimodal distribution as in the unperturbed case.
Finally, we have studied the fluctuation effect when we consider a modulated time dependence of CaMKII _{T} in according to the expression where we set Δ(t) = Δ_{0} if 0 < t < 2π/Ω and Δ = 0 otherwise. The model seen in Eq. 9 simulates the presence of a stimulating signal for a time interval π/Ω together with continuous fluctuations in the enzyme concentration. The Monte Carlo simulations show the model is able to show both LTP and LTD. We have computed the distribution function ρ(A,t) for Ω ≪ 1 so that an adiabatic approach can be used to describe the dynamics. When the CaMKII _{T} concentration increases, the distribution concentrates at A ≪ 1, indicating that most of the channels are in the doublephosphorylated state (see Fig. 7). Then the distribution remains concentrated around a doublephosphorylated state even when the external stimulus disappears (LTP). But the fluctuations induce a slow diffusion in the phase space that makes the distribution relax around a nonphosphorylated state (see Fig. 8), simulating LTD. The noise level defines the LTD time scale: For small noise most of channels would remain trapped in the doublephosphorylated state for an extremely long time according to Kramer's theory (16).
Discussion and Conclusion
We have demonstrated the existence of a bistable regime and switching properties of a 2 step phosphorylation cycle acting on AMPA channels catalyzed by 3 different enzymes, as described in the biological literature (1, 2) and observed in in vivo learning experiment. The proposed model reproduces experimental results on the induction of LTP. Moreover, we formulate specific predictions concerning the dynamical induction of LTD that can be tested experimentally. The incorporation of multiplicative noise gives a framework for a microscopic formulation of synaptic plasticity based on the socalled chemical master equation. This approach is formulated in terms of state probabilities instead of concentrations and can be used when the number of involved molecules is low, as in the case of AMPAR in a single synapse; it also efficiently describes AMPAR trafficking. The bistability property is a central requirement for memory storage and has been proposed as the key feature in models of synaptic plasticity and LTP induction and persistence (3, 9) based on autophosphorylation of CaMKII. Our model differs from these models in that the autophosphorylation of CaMKII is not required, and LTP is a direct consequence of the bistability character of the 2step enzymatic cycle. We believe that the incorporation of CaMKII autophosphorylation can lead to a further stable state in this model. The stochastic formulation of the deterministic scheme gives a better biophysical framework for this kind of model in the limit of low number of molecules and opens the intriguing possibility of a “constructive” role of noise in enhancing and maintaining memory storage. A challenging question is the link between this model and the Bienenstock–Cooper–Munro (BCM) theory of synaptic plasticity. In mathematical terms, the 2 systems seem to have the same bistable behavior so that the 3state AMPA model is equivalent to a 2input BCM model. A future study will be necessary to compare the 2 models in the case of n + 1 states and n inputs in order to verify the existence of n stable states Experimental tests of this model are feasible a relatively straightforward way in vitro and in cellular system by patchclamp recording and dynamical fluorescence microscopy.
Acknowledgments
This work was partially supported by the Italian Ministry of Research and University and National Institute of Nuclear Physics Single Heat Effects Induced by LowDose Irradiation (SHEILA) experiment grants (to G.C.C. and A.B.); by a Brown/Bologna Exchange Program Grant (to A.B.); and by the Institute for Brain and Neural Systems at Brown University.
Footnotes
 ^{1} To whom correspondence should be addressed. Email: leon_cooper{at}brown.edu

Author contributions: G.C.C., A.B., and L.N.C. designed research, performed research, contributed new reagents/analytic tools, analyzed data and wrote the paper.

The authors declare no conflict of interest.

This article contains supporting information online at www.pnas.org/cgi/content/full/0905988106/DCSupplemental.

Freely available online through the PNAS open access option.
References
 ↵
 ↵
 ↵
 Lisman JA
 ↵
 Castellani GC,
 Quinlan EM,
 Cooper LN,
 Shouval HZ
 ↵
 Castellani GC,
 Quinlan EM,
 Bersani F,
 Cooper LN,
 Shouval HZ
 ↵
 S. Whithlock MF
 ↵
 ↵
 ↵
 Pi HJ,
 Lisman JE
 ↵
 ↵
 Samoilov M,
 Plyasunov S,
 Arkin AP
 ↵
 Brandman O,
 Ferrell JE, Jr,
 Li R,
 Meyer T
 ↵
 Van Kampen NG
 ↵
 Qian H,
 Elson EL
 ↵
 Arnold L
 ↵
Citation Manager Formats
Sign up for Article Alerts
Article Classifications
 Biological Sciences
 Neuroscience