Biochemically detailed model of LTP and LTD in a cortical spine (Maki-Marttunen et al 2020)

 Download zip file 
Help downloading and running models
"Signalling pathways leading to post-synaptic plasticity have been examined in many types of experimental studies, but a unified picture on how multiple biochemical pathways collectively shape neocortical plasticity is missing. We built a biochemically detailed model of post-synaptic plasticity describing CaMKII, PKA, and PKC pathways and their contribution to synaptic potentiation or depression. We developed a statistical AMPA-receptor-tetramer model, which permits the estimation of the AMPA-receptor-mediated maximal synaptic conductance based on numbers of GluR1s and GluR2s predicted by the biochemical signalling model. We show that our model reproduces neuromodulator-gated spike-timing-dependent plasticity as observed in the visual cortex and can be fit to data from many cortical areas, uncovering the biochemical contributions of the pathways pinpointed by the underlying experimental studies. Our model explains the dependence of different forms of plasticity on the availability of different proteins and can be used for the study of mental disorder-associated impairments of cortical plasticity."
1 . Mäki-Marttunen T, Iannella N, Edwards AG, Einevoll GT, Blackwell KT (2020) A unified computational model for cortical post-synaptic plasticity. Elife [PubMed]
Model Information (Click on a link to find other models with that property)
Model Type: Synapse;
Brain Region(s)/Organism: Neocortex;
Cell Type(s): Neocortex spiking regular (RS) neuron;
Channel(s): I Calcium;
Gap Junctions:
Transmitter(s): Glutamate; Norephinephrine; Acetylcholine;
Simulation Environment: NEURON; NeuroRD;
Model Concept(s): Long-term Synaptic Plasticity;
Implementer(s): Maki-Marttunen, Tuomo [tuomomm at];
Search NeuronDB for information about:  I Calcium; Acetylcholine; Norephinephrine; Glutamate;
Ca_HVA.mod *
Ca_LVAst.mod *
CaDynamics_E2.mod *
epsp.mod *
Ih.mod *
Im.mod *
K_Pst.mod *
K_Tst.mod *
Nap_Et2.mod *
NaTa_t.mod *
NaTs2_t.mod *
ProbUDFsyn2.mod *
ProbUDFsyn2group.mod *
ProbUDFsyn2groupdet.mod *
SK_E2.mod *
SKv3_1.mod *
TITLE AMPA and NMDA receptor with presynaptic short-term plasticity 

AMPA and NMDA receptor conductance using a dual-exponential profile
presynaptic short-term plasticity based on Fuhrmann et al. 2002
Implemented by Srikanth Ramaswamy, Blue Brain Project, July 2009
Etay: changed weight to be equal for NMDA and AMPA, gmax accessible in Neuron



        RANGE tau_r_AMPA, tau_d_AMPA, tau_r_NMDA, tau_d_NMDA, Nsyns
        RANGE Use, u, Dep, Fac, u0, weight_factor_NMDA
        RANGE i, i_AMPA, i_NMDA, g_AMPA, g_NMDA, e, gmax


        tau_r_AMPA = 0.2   (ms)  : dual-exponential conductance profile
        tau_d_AMPA = 1.7    (ms)  : IMPORTANT: tau_r < tau_d
	tau_r_NMDA = 0.29   (ms) : dual-exponential conductance profile
        tau_d_NMDA = 43     (ms) : IMPORTANT: tau_r < tau_d
        Use = 1.0   (1)   : Utilization of synaptic efficacy (just initial values! Use, Dep and Fac are overwritten by BlueBuilder assigned values) 
        Dep = 100   (ms)  : relaxation time constant from depression
        Fac = 10   (ms)  :  relaxation time constant from facilitation
        e = 0     (mV)  : AMPA and NMDA reversal potential
	mg = 1   (mM)  : initial concentration of mg2+
    	gmax = .001 (uS) : weight conversion factor (from nS to uS)
    	u0 = 0 :initial value of u, which is the running value of Use
	Nsyns = 10 : How many synapses are there actually
        weight_factor_NMDA = 1

The Verbatim block is needed to generate random nos. from a uniform distribution between 0 and 1 
for comparison with Pr to decide whether to activate the synapse or not


double nrn_random_pick(void* r);
void* nrn_random_arg(int argpos);

extern int ifarg(int iarg);
extern int vector_capacity(void* vv);
extern void* vector_arg(int iarg);


        v (mV)
        i (nA)
	i_AMPA (nA)
	i_NMDA (nA)
        g_AMPA (uS)
	g_NMDA (uS)
 	space       : A pointer to the vector containing the synapse times. Note that the underlying vector should not be touched after initialization by setVec().


        A_AMPA       : AMPA state variable to construct the dual-exponential profile - decays with conductance tau_r_AMPA
        B_AMPA       : AMPA state variable to construct the dual-exponential profile - decays with conductance tau_d_AMPA
	A_NMDA       : NMDA state variable to construct the dual-exponential profile - decays with conductance tau_r_NMDA
        B_NMDA       : NMDA state variable to construct the dual-exponential profile - decays with conductance tau_d_NMDA


        LOCAL tp_AMPA, tp_NMDA
	A_AMPA = 0
        B_AMPA = 0
	A_NMDA = 0
	B_NMDA = 0
	tp_AMPA = (tau_r_AMPA*tau_d_AMPA)/(tau_d_AMPA-tau_r_AMPA)*log(tau_d_AMPA/tau_r_AMPA) :time to peak of the conductance
	tp_NMDA = (tau_r_NMDA*tau_d_NMDA)/(tau_d_NMDA-tau_r_NMDA)*log(tau_d_NMDA/tau_r_NMDA) :time to peak of the conductance
	factor_AMPA = -exp(-tp_AMPA/tau_r_AMPA)+exp(-tp_AMPA/tau_d_AMPA) :AMPA Normalization factor - so that when t = tp_AMPA, gsyn = gpeak
        factor_AMPA = 1/factor_AMPA
	factor_NMDA = -exp(-tp_NMDA/tau_r_NMDA)+exp(-tp_NMDA/tau_d_NMDA) :NMDA Normalization factor - so that when t = tp_NMDA, gsyn = gpeak
        factor_NMDA = 1/factor_NMDA



        SOLVE state METHOD cnexp
        mggate = 1 / (1 + (mg/4.1 (mM))*exp(0.063 (/mV)*(-v))) :mggate kinetics - Spruston et al. 1995
        g_AMPA = gmax*(B_AMPA-A_AMPA) :compute time varying conductance as the difference of state variables B_AMPA and A_AMPA
	g_NMDA = gmax*(B_NMDA-A_NMDA) * mggate :compute time varying conductance as the difference of state variables B_NMDA and A_NMDA and mggate kinetics
        i_AMPA = g_AMPA*(v-e) :compute the AMPA driving force based on the time varying conductance, membrane potential, and AMPA reversal
	i_NMDA = g_NMDA*(v-e) :compute the NMDA driving force based on the time varying conductance, membrane potential, and NMDA reversal
	i = i_AMPA + i_NMDA


        A_AMPA' = -A_AMPA/tau_r_AMPA
        B_AMPA' = -B_AMPA/tau_d_AMPA
	A_NMDA' = -A_NMDA/tau_r_NMDA
        B_NMDA' = -B_NMDA/tau_d_NMDA

NET_RECEIVE (weight, Pv, Pr, u, myInd, tsyn (ms), Pv_tmp){
	:printf("NMDA weight = %g\n", weight_NMDA)


        :Randomize which of the synapses is activated. Note that an additional random number is generated by rand() - this may interfere with the random number order in parallel simulations.
          void** vv = (void**)(&space);
          double *x;
          int nx = vector_instance_px(*vv, &x);
          int myInd = rand()%((int)Nsyns);
          _args[4] = myInd;
          _args[5] = x[myInd];                //tsyn
          _args[1] = x[myInd+(int)Nsyns];     //Pv
          _args[3] = x[myInd+2*((int)Nsyns)]; //u
	::printf("NET_RECEIVE_beg: Pv = %g, Pr = %g, u = %g, myInd = %g, tsyn = %g, t = %g\n", Pv, Pr, u, myInd, tsyn, t)
	:printf("NET_RECEIVE_beg:  myInd = %g/%g, Pv = %g, u = %g, tsyn = %g, t = %g. ", myInd, Nsyns, Pv, u, tsyn, t)

        : calc u at event-
        if (Fac > 0) {
              u = u*exp(-(t - tsyn)/Fac) :update facilitation variable if Fac>0 Eq. 2 in Fuhrmann et al.
        } else {
              u = Use  
        if(Fac > 0){
              u = u + Use*(1-u) :update facilitation variable if Fac>0 Eq. 2 in Fuhrmann et al.

        Pv_tmp  = 1 - (1-Pv) * exp(-(t-tsyn)/Dep) :Probability Pv for a vesicle to be available for release, analogous to the pool of synaptic
                                                  :resources available for release in the deterministic model. Eq. 3 in Fuhrmann et al.
        Pr  = u * Pv_tmp                          :Pr is calculated as Pv * u (running value of Use)
        Pv_tmp  = Pv_tmp - u * Pv_tmp             :update Pv as per Eq. 3 in Fuhrmann et al.
        :printf("Pv = %g\n", Pv)
        :printf("Pr = %g\n", Pr)

	if (erand() < Pr){
            tsyn = t
	    Pv = Pv_tmp
            A_AMPA = A_AMPA + weight*factor_AMPA
            B_AMPA = B_AMPA + weight*factor_AMPA
            A_NMDA = A_NMDA + weight*weight_factor_NMDA*factor_NMDA
            B_NMDA = B_NMDA + weight*weight_factor_NMDA*factor_NMDA
            :printf ( "R! Pr = %g\n" , Pr )
          } else {
            ::printf("Not released! value = %g, Pr = %g\n", erand(), Pr )
            :printf ( "NR! Pr = %g\n" , Pr )
	:printf("NET_RECEIVE_end: Pv = %g, Pr = %g, u = %g, myInd = %g, tsyn = %g, t = %g\n", Pv, Pr, u, myInd, tsyn, t)

          x[myInd] = _args[5];
          x[myInd+(int)Nsyns] = _args[1];
	  x[myInd+2*((int)Nsyns)] = _args[3];

         * This function takes a NEURON Random object declared in hoc and makes it usable by this mod file.
         * Note that this method is taken from Brett paper as used by netstim.hoc and netstim.mod
         * which points out that the Random must be in negexp(1) mode
        void** pv = (void**)(&_p_rng);
        if( ifarg(1)) {
            *pv = nrn_random_arg(1);
        } else {
            *pv = (void*)0;

FUNCTION erand() {
	    //FILE *fi;
        double value;
        if (_p_rng) {
                :Supports separate independent but reproducible streams for
                : each instance. However, the corresponding hoc Random
                : distribution MUST be set to Random.negexp(1)
                value = nrn_random_pick(_p_rng);
		        //fi = fopen("RandomStreamMCellRan4.txt", "w");
                //fprintf(fi,"random stream for this simulation = %lf\n",value);
                //printf("random stream for this simulation = %lf\n",value);
                return value;
                : the old standby. Cannot use if reproducible parallel sim
                : independent of nhost or which host this instance is on
                : is desired, since each instance on this cpu draws from
                : the same stream
                erand = exprand(1)
        :erand = value :This line must have been a mistake in Hay et al.'s code, it would basically set the return value to a non-initialized double value.
                       :The reason it sometimes works could be that the memory allocated for the non-initialized happened to contain the random value
                       :previously generated (or if _p_rng is always a null pointer). However, here we commented this line out.

PROCEDURE setVec() {    : Sets the times of firing of each synapse. This should be done only once for each ProbAMPANMDA2group,
                        : before the running of the simulation, and the underlying vector should be untouched after that.
  void** vv;
  vv = (void**)(&space);
  *vv = (void*)0;
  if (ifarg(1)) {
    *vv = vector_arg(1);
    Nsyns = vector_capacity(*vv)/3;

PROCEDURE printVec() { : Prints the previous times of firing of each synapse.
    void** vv = (void**)(&space);
    double *x;
    int nx = vector_instance_px(*vv, &x);
    int i1;
    for (i1=0; i1<Nsyns;i1++) {
      printf("tsyns[%i] = %g, Pv[%i] = %g, u[%i] = %g\n", i1, x[i1], i1, x[i1+(nx/3)], i1, x[i1+2*(nx/3)]);

Loading data, please wait...