GluN2B-NMDAR subunit contribution on synaptic plasticity: A phenomenological model for CA3-CA1 synapses

Synaptic plasticity is believed to be a key mechanism underlying learning and memory. We developed a phenomenological N-methyl-D-aspartate (NMDA) receptor-based voltage-dependent synaptic plasticity model for synaptic modifications at hippocampal CA3-CA1 synapses on a hippocampal CA1 pyramidal neuron. The model incorporates the GluN2A-NMDA and GluN2B-NMDA receptor subunit-based functions and accounts for the synaptic strength dependence on the postsynaptic NMDA receptor composition and functioning without explicitly modeling the NMDA receptor-mediated intracellular calcium, a local trigger of synaptic plasticity. We embedded the model into a two-compartmental model of a hippocampal CA1 pyramidal cell and validated it against experimental data of spike-timing-dependent synaptic plasticity (STDP), high and low-frequency stimulation. The developed model predicts altered learning rules in synapses formed on the apical dendrites of the detailed compartmental model of CA1 pyramidal neuron in the presence of the GluN2B-NMDA receptor hypofunction and can be used in hippocampal networks to model learning in health and disease.


. Introduction
Long-term synaptic plasticity has been proposed to be the cellular substrate of learning and memory in the brain (Malenka and Nicoll, 1999;Malenka and Bear, 2004). In the hippocampal CA1 area, CA3 Schaffer collateral-CA1 pyramidal neuron synapses can undergo long-term potentiation (LTP) and long-term depression (LTD), triggered by high or low frequency presynaptic stimulation (Collingridge et al., 1983;Dudek and Bear, 1992;Mulkey and Malenka, 1992;Bliss and Collingridge, 1993;Goh and Manahan-Vaughan, 2013;Pousinha et al., 2017). Spike-timing-dependent synaptic plasticity (STDP) is a bidirectional form of synaptic modifications, induced by correlated pre-and postsynaptic neuronal activation, where precise timing of spikes is a major determinant of the direction and magnitude of synaptic modifications (Markram et al., 1997;Bi and Poo, 1998;Debanne et al., 1998;Feldman, 2012). The induction of LTP, LTD, and STDP in excitatory synapses at Schaffer collateral pathway requires activation of N-methyl-D-aspartate receptors .
GluN2B-NMDAR is strongly coupled with calciumcalmodulin-dependent protein kinase II (CaMKII), a key protein that induces downstream signaling cascades mediating LTP expression, learning, and memory. CaMKII leads to phosphorylation of synaptic proteins, increase in the number of active AMPARs or their single-channel conductance (Park et al., 2021;Yasuda et al., 2022). During LTP induction, Ca 2+ influx through GluN2B-NMDAR directly activates CaMKII and leads to synapse strengthening (Shipton and Paulsen, 2014). Experimental data shows that disruption of GluN2B-NMDAR/CaMKII interactions downregulates CaMKII activation and autophosphorylation, prevents LTP in hippocampus, and impairs spatial learning in transgenic mice (Zhou et al., 2007). GluN2B-NMDA type of receptor is crucial for normal learning in hippocampus in vivo . Moreover, GluN2B-NMDAR subunit is implicated in variety of psychiatric disorders like dementia, Alzheimer's disease (Liu et al., 2019), and Schizophrenia (Kocsis, 2012). GluN2B-NMDAR plays a critical role in synaptic plasticity and cognitive impairment in Alzheimer's disease (Parameshwaran et al., 2008;Pousinha et al., 2017Pousinha et al., , 2019. In modeling studies of synaptic plasticity, the challenge is to integrate knowledge at molecular, synaptic, neuronal, and microcircuit levels, to understand the underlying LTP and LTD mechanisms in synapses and transfer the knowledge to large network simulations. Numerous computational studies of synaptic plasticity exist, and the models can be grouped into three main classes that employ phenomenological, optimal, and biophysical approach. Phenomenological models are abstract, and encode data and intuitions about synaptic plasticity taking into account spike timing (Gerstner et al., 1996;Kempter et al., 1999;Kistler and van Hemmen, 2000;. Optimal models use some optimality criterion to deduce the rules of synaptic modifications (Toyoizumi et al., 2005;Pfister et al., 2006). Biophysical models rely on biologically realistic descriptions of the electrophysiological and biochemical processes, usually are based on intracellular calcium dynamics, and involve detailed biochemical reactions to explain synaptic plasticity (Bhalla and Iyengar, 1999;Senn et al., 2001;Shouval et al., 2002;Badoual et al., 2006;Graupner and Brunel, 2007;Pi and Lisman, 2008;Clopath et al., 2010;Urbanczik and Senn, 2014;Migliore et al., 2015;Saudargiene et al., 2015;Sacramento et al., 2018;Ebner et al., 2019;Mäki-Marttunen et al., 2020;Chindemi et al., 2022).
In biophysical models a widely used approach is to investigate the molecular networks underlying synaptic plasticity such as CaMKII and protein phosphatase competition activated by NMDAR-mediated calcium (Graupner and Brunel, 2007;Pi and Lisman, 2008;Saudargiene et al., 2015). A well-known calcium control hypothesis states that low calcium levels in dendritic spine do not evoke any changes, intermediate calcium levels depress the synapse and high calcium transients potentiate the synapse (Shouval et al., 2002). Biophysical models, embedded into detailed compartmental models, account for the factors shaping synaptic plasticity-different membrane mechanisms of the dendritic tree, dendritic integration, morphological features, pattern of pre-and postsynaptic spiking (Poirazi and Papoutsi, 2020). The models include complex biochemical reactions of calcium induced kinase and phosphatase activation that underlie synaptic modifications (Bhalla and Iyengar, 1999;Graupner and Brunel, 2007;Saudargiene et al., 2015;Jędrzejewska-Szmek et al., 2017;Mäki-Marttunen et al., 2020), represent molecular cascades applying simplified functions, dependent on postsynaptic NMDAR-mediated intracellular calcium transients (Shouval et al., 2002;Graupner and Brunel, 2012;Standage et al., 2014;Chindemi et al., 2022), use formulation on the level of postsynaptic voltage (Clopath et al., 2010;Meissner-Bernard et al., 2020), utilize a kinetic model of synapse upregulation and downregulation mediated by NMDAR and based on the precise timing of pre and post spikes (Senn et al., 2001), or describe the weight change in a phenomenological way taking into account spike timing (Gerstner et al., 1996;. Phenomenological models are efficient, but lack biological realism; on the other hand, detailed models, sensitive to NMDAR functioning, are not easily applied in network simulations as they include many complex biochemical reactions, large parameter space, and are computationally expensive. The models that account for the NMDAR subunit properties and are suitable to analyze learning properties in networks are still lacking.
Different forms of LTP and LTD coexist that have different induction and expression mechanisms. In this study we focus on the GluN2B-NMDAR subunit effect on LTP induction in STDP, .
/fnsyn. . high and low frequency stimulation protocols. The aim of this work is to build a phenomenological NMDAR-based synaptic plasticity model that separates the influence of GluN2A-NMDAR and GluN2B-NMDAR subunits, and gain insight into the GluNR2B-NMDAR effect on synaptic modifications of hippocampal CA3-CA1 synapses. We modeled synaptic plasticity induced by a STDP protocol and high and low frequency stimulation, and explored the impact of GluN2B-NMDAR subunit properties on the modification of synaptic strength of the synapses spatially distributed on the apical dendrites of CA1 pyramidal neuron. The novelty of the work is the approach to include the influence of the specific NMDA receptor subunits on synaptic plasticity as the separate mediators of LTP and LTD. We assume that LTP is mainly mediated by GluN2B-NMDAR (Morishita et al., 2007;Andrade-Talavera et al., 2016;Pousinha et al., 2017), and LTD is triggered by other mechanisms, possibly including GluN2A-NMDAR. Experimental studies showed that GluN2A-NMDAR blockade prevented LTD induction in the CA1 region of hippocampal slices (Bartlett et al., 2007;Li et al., 2007). The study of Morishita et al. (2007) also suggested that GluN2A-NMDAR might be responsible for LTD as the application of the GluN2B-NMDAR antagonist ifenprodil did not prevent the induction of LTD by low frequency stimulation.
We analyzed the effect of GluN2B-NMDAR on the properties of learning at the synapses of hippocampal CA1 pyramidal neuron. The modeling results provide insights into the learning rules of hippocampal CA1 pyramidal neuron synapses in healthy and GluN2B-NMDAR hypofunction conditions.

. Methods
We developed a model of synaptic modifications based on the well-established phenomenological models (Clopath et al., 2010;Meissner-Bernard et al., 2020) and integrated the separated influence of postsynaptic NMDAR subunits GluN2A-NMDAR and GluN2B-NMDAR to account for the crucial effect of GluN2B-NMDAR in hippocampal synaptic plasticity. We utilized two computational models of CA1 pyramidal neuron: a modified two-compartmental Pinsky-Rinzel model (Pinsky and Rinzel, 1994;Ferguson and Campbell, 2009) to validate the extended synaptic plasticity model, and a multicompartmental model (Migliore et al., 2018) to study the influence of GluN2B-NMDAR properties on synaptic strength modifications at a cluster of CA3-CA1 synapses distributed randomly onto apical dendrites of CA1 neuron.

. . NMDAR-based voltage-dependent synaptic plasticity model
We extended a voltage-based model of synaptic plasticity (Clopath et al., 2010;Meissner-Bernard et al., 2020) by including the effect of postsynaptic NMDAR subunits GluN2A-NMDAR and GluN2B-NMDAR. The instantaneous weight change d dt (t)w consists of two additive NMDAR-dependent LTD and LTP contributions, where w max and w min set the limits for synaptic weight w.
The LTP component d dt w LTP (t) is expressed as the product of the NMDAR-dependent function φ NMDA + (t) and a low-filtered membrane potential V + (t): where A + is the LTP amplitude parameter. Similarly, the LTD component d dt w LTD (t) is proportional to the product of the NMDAR-dependent function φ NMDA − (t) and a low-filtered membrane potential V − (t): where X(t) is a presynaptic activity variable, and A − is the LTD amplitude parameter. The contribution of the postsynaptic GluN2A-NMDAR and GluN2B-NMDAR gated channel is captured by the newly introduced Hill function φ  , and θ φ [ * ] differ for the LTD and LTP contributions. Function φ NMDA + (t) governs LTP induction and is caused by the filtered NMDAR conductance g NMDA + (t). Function φ NMDA − (t) accounts for LTD, and is triggered by the filtered NMDAR conductance g NMDA − (t). The moving threshold function θ φ [ * ] lowers φ NMDA [ * ] (t) activity and implements competition between LTP and LTD: The moving threshold θ φ − (t) is increasing, if φ NMDA + (t) is strongly activated and LTP is induced, thus vetoing LTD. Threshold θ φ + (t) may also increase if φ NMDA − (t) accumulates, and leads to LTD.
The filtered NMDAR-dependent variables g NMDA [ * ] (t) for LTP and LTD components are described: where τ NMDA [ * ] is a time constant, and g NMDA [ * ] (t) is a conductance of postsynaptic NMDAR that incorporates both GluN2A-NMDAR and GluNB-NMDAR subunits with a different weighting coefficient k 2B[ * ] for LTP and LTD components:  2017), we assume that LTP is mainly governed by GluN2B-NMDAR subunit, and LTD is mediated by GluN2A-NMDAR (or other) subunit. We set the coefficient of GluN2B-NMDAR effect on LTP k 2B+ = 0.8, and coefficient of GluN2B-NMDAR subunit has a slower inactivation time than GluN2A-NMDAR subunit. Kinetic parameters of forward and backward binding rates are adjusted (Cull-Candy et al., 2001).

FIGURE
Schematic diagram of the main components of synaptic plasticity model. Presynaptic action potential activates NMDAR, induces the NMDAR conductance, composed of GluN A-NMDAR and GluN B-NMDAR subunits g NMDA GluN A and g NMDA GluN B , and triggers variables φ NMDA+ and φ NMDA− that account for the LTP and LTD contribution, respectively. The LTP variable φ NMDA+ , once activated, inhibits the LTD variable φ NMDA− to prevent LTD, and vice versa. Postsynaptic local membrane potential V is low-pass filtered, and the resulting LTP and LTD variables V + and V − are multiplied by the corresponding NMDAR-dependent variables φ NMDA+ and φ NMDA− to form the LTP and LTD components of the weight w.
Variables V − (t) and V + (t) are the functions of the filtered membrane potential V(t) at the synapse location, and contribute to the LTD and LTP components: where τ [ * ] is a time constant and θ [ * ] is a threshold for the LTP and LTD components. A presynaptic activity variable X in Equation (3) is calculated as a low pass filter of the presynaptic spike train Schematic representation of synaptic plasticity model is shown in Figure 1. The presynaptic activity triggers NMDAR synaptic conductance g NMDA (t), composed of GluN2A-NMDAR and GluN2B-NMDAR subunits, and induces variables φ NMDA + (t) and φ NMDA − (t). Once activated, the LTP variable φ NMDA + (t) inhibits the LTD variable φ NMDA − (t) preventing LTD induction, and vice versa-φ NMDA − (t) may reduce the activity of φ NMDA + (t). Postsynaptic activity is captured by a local membrane potential V(t) that is thresholded using the thresholds θ + and θ − , and lowpass filtered resulting in V + (t) and V − (t). The LTP component d dt w LTP (t) is obtained by multiplying V + (t) and φ NMDA + (t), and the , and X(t). The weight change of the AMPAR strength w is composed of the scaled LTP and LTD parts d dt w LTP (t) and d dt w LTD (t). The novelty of the model is that it captures the separate specific influence of GluN2A-NMDAR and GluN2B-NMDAR subunits on LTP and LTD induction using the NMDAR subunit-dependent activation functions φ NMDA − (t) and φ NMDA + (t) ( Equation 4). The GluN2A-NMDAR and GluN2B-NMDAR subunitdependent functions (Equations 4, 6, 7) shape the LTP and LTD components (Equations 2, 3). The filtered NMDAR synaptic conductance-dependent variables g NMDA [ * ] (t) (Equation 6) can be interpreted as the intracellular Ca 2+ concentration, and functions φ NMDA [ * ] (t) (Equation 4) reflect the activation of intracellular Ca 2+ -triggered second messenger cascades underlying synaptic plasticity induction. Specifically, φ NMDA + (t) may represent phosphorylation of CaMKII, leading to LTP, and φ NMDA − (t) may indicate dephosphorylation of protein phosphatase 2B (PP2B, calcineurin), responsible for LTD. The description of the signaling pathways is simplified and implemented in a phenomenological manner using the NMDAR-dependent functions. The model does not require the estimation of the intracellular calcium concentration at a postsynaptic site and relies on the NMDAR properties and local membrane potential.
The parameters of the synaptic plasticity model are given in Table 1. The parameters of the GluN2A-NMDAR and GluN2B-NMDAR synaptic conductances are presented in Table 2.

. . Computational models of CA pyramidal neuron
We employed two computational models of CA1 pyramidal neuron: a modified two-compartmental Pinsky-Rinzel model for synaptic plasticity model validation (Pinsky and Rinzel, 1994;Ferguson and Campbell, 2009) and a compartmental detailed model (Migliore et al., 2018) for analysis of GluN2B-NMDARdependent synaptic plasticity properties at CA3-CA1 synapses distributed on the apical dendrites of CA1 pyramidal neuron in the stratum radiatum (SR) region.
. . . Two-compartmental model of CA pyramidal neuron A two-compartmental Pinsky-Rinzel model consisted of a somatic and dendritic compartments connected by the coupling conductance (Pinsky and Rinzel, 1994;Ferguson and Campbell, 2009). The somatic compartment had five ionic current channels: inward Na + current I Na,s , outward delayed rectifier K + current I KDR,s , inward Ca 2+ current I Ca,s , outward short-duration voltage and Ca 2+ -dependent K + current I KCa,s , and outward long-duration Ca 2+ -dependent after hyperpolarization (AHP) K + current I KAHP,s . The dendritic compartment had three ionic current channels: inward Ca 2+ current I Ca,d , outward short-duration voltage and Ca 2+ -dependent K + current I KCa,d , outward long-duration Ca 2+ -dependent AHP potassium current I KAHP,d . Both compartments had leak current I leak,s , I leak,d . A single synapse containing AMPAR and GluN2A-NMDAR/GluN2B-NMDAR was formed on the dendritic compartment. The two-compartment model and synaptic currents are described in Supplementary material. The partial blockade of GluN2B-NMDAR-gated channel was simulated by reducing the conductance g NMDA GluN2B (Equation 7).

. . . Multicompartmental model of CA pyramidal neuron
A multicompartmental model of a CA1 pyramidal cell oh140807_A0_idA (Migliore et al., 2018) consisting of 175 compartments was used, and it included 11 ionic current channels and a leak current. The ionic currents were the following: inward Na + current I Na ; four types of K + currents: outward delayed rectifier K + current I KDR , transient A-type K + current I KA , currents I KM , I KD ; three types of inward Ca 2+ currents: N-type current I CaN , L-type current I CaL , T-type current I CaT ; two types of Ca 2+dependent K + currents: outward short-duration voltage and Ca 2+dependent K + current I KCa and I Cagk current; and the nonspecific I h current. Ionic channels were uniformly distributed in all dendritic compartments except I KA and I h , which increased with distance from the soma. Channels were described using a conventional Hodgkin-Huxley formalism, and peak conductances of each channel were optimized for soma, axon, basal, and apical dendrite compartments and validated against experimental data. Intracellular calcium concentration was described by a simple Ca 2+ extrusion mechanism with a single exponential decay. The multicompartmental model and synaptic currents are described in Migliore et al. (2018).
A cluster of 50 AMPARs and GluN2A/GluN2B-NMDARs containing synapses, distributed randomly on the apical dendrites of the neuron in the SR region at 140 µm from the soma with a synaptic density of 0.8 synapses/µm of dendrite (Gasparini et al., 2004;Bezaire et al., 2016) was formed to model synaptic modifications.
The ratio of AMPAR/NMDAR gated channel currents was replicated using the experimental protocol used in (Pousinha et al., 2017). The voltage was clamped at -65 mV for AMPAR gated channel current estimation, and at +40 mV for NMDAR gated channel assessment. The peak AMPAR current was compared with the NMDAR current 60 ms after the onset of stimulus. The maximal conductances of the AMPAR and NMDAR-gated channels were set to ensure this ratio to be equal to 4 as in Pousinha et al. (2017). The partial blockade of GluN2B-NMDAR-gated channel was simulated by lowering the conductance g NMDA GluN2B (Equation 7).

. . Stimulation protocols for synaptic plasticity induction at CA -CA synapses
Synapses were stimulated using the activation patterns applied in the following electrophysiological studies of synaptic plasticity: • STDP induction protocol (Wittenberg and Wang, 2006;Inglebert et al., 2020). Presynaptic input was paired with a doublet of postsynaptic action potentials 60 times at 5 Hz; 5 times at 5 Hz frequency; and 30 times at 1 Hz. Temporal difference T was measured between pre-and a second postsynaptic spike. In addition, a presynaptic spike was paired with a single postsynaptic action potential 60 times at 5 Hz. Temporal difference T was measured between a pre-and a postsynaptic spike (Wittenberg and Wang, 2006). Pairing frequency was increased from 1 to 50 Hz for spike pairs with a temporal difference T = 10 ms between a pre-and a single
• Frequency-dependent synaptic plasticity induction protocol (Pousinha et al., 2017). Presynaptic input was stimulated at 100 Hz for 1 s (LTP protocol) or at 1 Hz for 100 s (LTD protocol). To estimate the change in the excitatory postsynaptic potential (EPSP), a presynaptic stimulus was delivered before and after the conditioning stimulation, and the resulting ratio between the maximal values of the resulting EPSPs was calculated.
Simulations were performed in the Python and NEURON simulation environment (version 8.0.0) (Hines and Carnevale, 1997
We used a two-compartment neuron model with a single synapse on its dendrite, and applied the STDP induction protocol by pairing presynaptic activity with a doublet of postsynaptic action potentials (Wittenberg and Wang, 2006) (Figure 2). For pre-post stimulation protocol and T = 10 ms, the presynaptic activation precedes a second action potential ( Figure 2A1, blue triangle and black line, respectively) and results in opening of GluN2A-NMDAR g NMDA GluN2A (red line) and GluN2B-NMDAR g NMDA GluN2B (green line) ( Figure 2A2); filtered NMDAR-dependent variables g NMDA − (red) and g NMDA + (green) for LTD and LTP induction ( Figure 2A3) favor activation of the LTP function φ NMDA + ( Figure 2A4, green line). The product of φ NMDA + and V + forms d dt w LTP (t) and leads to the increased weight w ( Figure 2A5). For post-pre stimulation protocol and T = −10 ms, the presynaptic activation follows a second somatic action potential ( Figure 2B1, blue triangle and black line, respectively), NMDAR activation is weaker ( Figure 2B2), failing to sufficiently activate φ NMDA + , but strong enough to induce φ NMDA − ( Figure 2B4, green and red lines, respectively). Functions φ NMDA − , V − , and X(t) combine into d dt w LTD (t) and result in the decreased weight w ( Figure 2B5).
The STDP curve of the weight change w induced by pairing the presynaptic input with a doublet of postsynaptic action potential 60 times at 5 Hz frequency with temporal difference T ∈ [-100; 100 ms] is presented in Figure 3A. For the positive T window from 0 ms up to 40 ms a synapse undergoes LTP, and LTD is obtained for anti-causal pairings and causal pairings within the T interval [40; 100 ms]. Shorter stimulation of five pairings at 5 Hz results in a potentiation-only plasticity rule ( Figure 3B) as the duration is not sufficient for the accumulation of the LTD mechanisms activity (φ NMDA − in our model). Pairings at low frequency of 1 Hz results in LTD only ( Figure 3C). A single postsynaptic action potential paired with a presynaptic action potential 60 times at 5 Hz evokes LTD, as the activation of LTP variable is too weak (φ NMDA + in our model, Figure 3D). The modeled STDP weight modifications replicate the experimental data (Wittenberg and Wang, 2006).
We investigated the weight change dependence on frequency of postsynaptic pairings and a number of postsynaptic spikes. When a presynaptic action potential was paired with a single postsynaptic action potential at T = 10 ms, the synapse underwent LTD for very low repetition frequencies and switched to LTP for the increasing frequency above 5 Hz ( Figure 4A). A single postsynaptic spike, paired with the input activity at T = 10 ms, induced LTD, while two, three, and four postsynaptic spikes led to LTP ( Figure 4B). The results qualitatively reproduces the experimental observations on LTP recovery with increasing pairing frequency and postsynaptic spike number (Inglebert et al., 2020).
Next, we applied the frequency-dependent stimulation protocol using the same two-compartment model of CA1 pyramidal neuron with a single synapse ( Figure 5). In subthreshold regime ( Figures 5A1-A5), stimulation of presynaptic input at 100 Hz resulted in opening of Glu2NA-NMDAR and Glu2NB-NMDAR channels ( Figure 5A2), activation of φ NMDA + and inhibition of φ NMDA − (Figures 5A3, A4, green and red lines respectively), and increase in weight w ( Figure 5A5). In suprathreshold regime ( Figures 5B1-B5), the same protocol led to the generation of somatic action potentials ( Figure 5B1, black line), high activity of NMDAR channels, strong increase in φ NMDA + , and inhibition of φ NMDA − ( Figure 5B4) causing strong LTP ( Figure 5B5). Low frequency stimulation at 1 Hz only slightly opened Glu2NA-NMDAR and Glu2NB-NMDAR channels ( Figure 5C2) that was sufficient to activate φ NMDA − , but not φ NMDA + ( Figure 5C4, red line) and induce LTD ( Figure 5C5). The results indicate that the model is suitable to account for the synaptic changes using frequency dependent LTD and LTP protocols in subthreshold and suprathreshold regimes.
The model validation analysis showed that the NMDARdependent synaptic plasticity model was capable of reproducing the experimental STDP weight change curves for different frequencies and postsynaptic patterns, high-frequency, and low-frequency stimulation protocols.

. . Synaptic plasticity depends on GluN B-NMDAR properties in a synapse cluster on a CA pyramidal neuron
We employed the developed model of synaptic plasticity and a compartmental detailed model of a CA1 pyramidal neuron (Migliore et al., 2018) to analyze the dependence of synaptic modifications on the GluN2B-NMDAR functioning using STDP and frequency-dependent stimulation protocols. We modeled weight modifications at clustered synapses on the apical branches of CA1 pyramidal neuron and measured EPSP change in soma before and after the stimulation protocol. Experimental data and computational modeling studies suggest that synapses tend to form tightly-packed groups or clusters on the dendrites of neurons [for review see Kastellakis and Poirazi (2019) and Miry et al. (2021)]. Such nearly-synchronous activated inputs . /fnsyn. . carry similar information onto the same dendrite and enable emerging of memory engrams. Thus, we formed a cluster of 50 effectively activated AMPAR and GluN2A-NMDAR/GluN2B-NMDAR containing synapses on the apical dendrites of CA1 pyramidal neuron. During stimulation, each synapse on the dendritic branch developed its weight depending on the NMDAR-gated synaptic conductance function (Equation 4) that sensed local depolarization and presynaptic glutamate release. Figure 6 shows the evolvement of 50 synaptic weights in a cluster during LTP and LTD stimulation protocols.
High frequency presynaptic stimulation at 100 Hz for 1 s depolarized membrane potential ( Figure 6A1; blue linemembrane potential V d at a location of a randomly chosen synapse; black line-membrane potential in soma V s ) and activated GluN2A-NMDAR and GluN2B-NMDAR channels g NMDA GluN2A and g NMDA GluN2B ( Figure 6A2, black and blue lines, respectively) that led to the increase of synaptic weights ( Figure 6A3). The histograms of the synaptic weights ( Figure 6A4) shows that the weights distributed in the interval from 1 up to the predefined maximum value w max = 2.5. Some synapses were only slightly potentiated due to the low local membrane potential and weakly activated NMDAR. The EPSP increased by 190% if compared to the EPSP before the conditioning stimulation ( Figure 6A5, green line vs. black line). Low-frequency stimulation at 1 Hz for 100 s induced a small membrane depolarization at the synapse location ( Figure 6B1, red line), weak NMDAR activation ( Figure 6B2), and resulted in the weakening of synaptic strength of all synapses (Figures 6B3,  B4). After the stimulation, the somatic EPSP decreased to 61% ( Figure 6B5). The synaptic plasticity model embedded into a detailed model of a CA1 pyramidal cell qualitatively reproduced the experimental results. Experimental data showed that 500 pulses at 1 Hz induced 57% LTD and 100 pulses at 100 Hz led to 191% LTP in hippocampal CA1 pyramidal neurons in rats (Pousinha et al., 2017).
We investigated the influence of the partial and full blockade of GluN2B-NMDAR on synaptic plasticity outcome using frequencydependent stimulation and STDP protocols (Figure 7). First, presynaptic input was stimulated at 100 Hz for 1 s (LTP protocol), and the normalized EPSP to the pre-LTP baseline value was estimated. The blockade of the GluN2B-NMDAR synaptic conductance, leaving 0.3 fraction its active baseline value, resulted in the decrease of LTP from 190% to 164%, while the full blockade of the GluN2B-NMDAR led to LTD, the decrease to 90% of somatic EPSP. The impairment of GluN2B-NMDAR did not affect LTD leaving it to 63% ( Figure 7A). The model of synaptic plasticity qualitatively reproduced experimental data of GluN2B-NMDAR inhibitor ifenprodil effect on LTP (Pousinha et al., 2017). It was shown that ifenprodil dose-dependently inhibited LTP, evoked by high frequency stimulation in CA3-CA1 synapses. The maximal inhibition efficacy was observed with the increasing ifenprodil concentration of 5 µm that converted 190% LTP to 75% LTD, but it did not affect LTD. The results show that the synaptic plasticity model can capture the influence of GluN2B-NMDAR properties at a single synapse by decreasing the synaptic weight changes in response to the impaired GluN2B-NMDAR functioning, and quantitatively follows the experimental data (Pousinha et al., 2017).
Next, we applied the STDP stimulation protocol with the temporal difference T = ± 20 ms between a presynaptic activity and a second postsynaptic spike ( Figure 7B).
GluN2B-NMDAR blockade led to the LTP switch into LTD for pre-post pairings and left LTD intact.
The results illustrate that synaptic plasticity is strongly affected by Glu2NB-NMDAR subunit properties. Hypofunction of GluN2B-NMDAR abolishes LTP induction and does not affect LTD for STDP protocol. For high frequency stimulation, LTP switches to LTD, and leaves LTD intact for low frequency stimulation. The results quantitatively align well with the experimental evidence on Glu2NB-NMDAR importance in shaping LTP at hippocampal synapse (Morishita et al., 2007;Andrade-Talavera et al., 2016;Pousinha et al., 2017). GluN2B-NMDAR may act as an additional modulatory mechanism of synaptic plasticity, and further experimental and computational studies are needed to understand the importance of the NMDAR subunit composition, and its effect on synaptic plasticity.

. Discussion
We developed an NMDAR subunit-dependent voltage-based synaptic plasticity model of synaptic weight modifications at hippocampal CA3-CA1 synapses. We extended the computational model of STDP (Clopath et al., 2010;Meissner-Bernard et al., 2020) by simultaneously incorporating the GluN2A-NMDAR and GluN2B-NMDAR components to account for the specific functions of NMDAR subunits in synaptic learning. The model of synaptic plasticity was validated against the experimental data (Wittenberg and Wang, 2006;Pousinha et al., 2017;Inglebert et al., 2020) and reproduced STDP and frequency dependent LTP and LTD. Furthermore, the results show that this synaptic plasticity model is able to account for the impairment of LTP in GluN2B-NMDAR hypofunction conditions as in experimental studies (Pousinha et al., 2017), demonstrating that synaptic plasticity depends on GluN2B-NMDAR properties, and dysfunction of GluN2B-NMDAR leads to LTP impairment and its transformation to LTD.
The developed model can be interpreted as a phenomenological model, standing at the intersection with the class of biophysical models of synaptic plasticity. The model captures the influence of the GluN2B-NMDAR subunit on synaptic modifications. Usually, the calcium-based models of synaptic plasticity do not distinguish between the NMDAR subunits as mediators of calcium influx. Our approach enriches the model with the new features of specific NMDAR effect on synaptic plasticity.
Our model uses the functions of the NMDAR subunit dynamics that in an abstract form accounts for the CaMKII and phosphatase activation, does not require to model dendritic spines and estimate intracellular calcium concentration, a main trigger of synaptic plasticity. The formalism proposed captures synapse-specific mechanisms that define synaptic plasticity-the local non-linear activation of NMDAR, its subunit composition and functioning.
We chose the voltage-based approach as the voltage traces are conventionally recorded in the experimental setups, well-described .
/fnsyn. . by the mathematical formalism, enabling the model to be usable for future network level simulations. The level of modeling may also be more detailed and rely on intracellular calcium dynamics, as in e.g., Shouval et al. (2002) and Graupner and Brunel (2012), while focusing on the NMDAR subunit effects on LTP and LTD induction.
We explored the impact of the NMDAR properties on the somatic EPSP changes using a biologically realistic multicompartmental CA1 pyramidal neuron and taking into account the influence of the spatial distribution of the synapses.
The results indicate that GluN2B-NMDAR regulates the amount of synaptic strength on the dendritic tree and the resulting EPSP changes in soma. The hypofunction of GluN2B-NMDAR leads to the impairment of LTP and gradual switch to LTD. The study extends the experimental observations of Pousinha et al. (2017) and predicts the pattern of the GluN2B-NMDAR .
/fnsyn. .  functioning-mediated synaptic plasticity measured as changes in somatic EPSP for specific synapse cluster for STDP induction protocol. In detailed biophysical modeling studies, long-term synaptic plasticity depends on intracellular calcium influx, but the sources of calcium is usually not taken into consideration.
Here, we discuss the importance of considering the mediating role of Glu2NB-NMDAR to study LTP in STDP and frequency dependent synaptic plasticity. As NMDARs can undergo activitydependent long-term plasticity (Hunt and Castillo, 2012), this work shows the importance to consider the state of NMDARs, not only in the modeling studies of learning and memory, but also in physiological experiments. The model offers a possibility to include GluN2B-NMDAR effective contribution to the synaptic weight modifications. The incorporation of the separated influence of the GluN2A-NMDAR and GluN2B-NMDAR functioning makes the model a candidate to explore learning in the diseased brain, as the Glu2NB-NMDAR normal functioning is crucial for healthy CA3-CA1 synapses, and its dysfunction is observed in cognitive deficits in neurological diseases (Kocsis, 2012;Pousinha et al., 2017Pousinha et al., , 2019Adell, 2020). The model of synaptic plasticity is principally based on the critical role of postsynaptic NMDAR in LTP and LTD induction in adult CA3-CA1 synapses. GluN2A-NMDAR and .
/fnsyn. . solution compositions, and other variables linked to the different experimental settings. In general, it is hypothesized that the GluN2A-to-GluN2B ratio defines the magnitude and sign of frequency-induced synaptic plasticity and shifts the LTP and LTD threshold. Higher GluN2A-to-GluN2B ratio requires stronger stimulation to induce LTP, confirming the critical role of GluN2B in LTP (Paoletti et al., 2013).
Recent experimental evidence shows that synaptic plasticity has different induction mechanisms depending on the NMDAR position (pre-or post-synaptic) and subunit composition, .

FIGURE
Blockade of GluN B-NMDAR prevents LTP for high-frequency and STDP stimulation protocols. (A) Presynaptic input was stimulated at Hz for s (LTP protocol, green bars) and at Hz for s (LTD protocol, red bars). LTP was impaired, and LTD remained intact. (B) Presynaptic activity was paired with a doublet of postsynaptic action potentials times at Hz frequency with a temporal di erence T = + ms (LTP protocol, green bars) and T = -ms (LTD protocol, red bars). LTP was abolished, and LTD was preserved.
developmental stage of animal, or experimental settings including the type of protocol used. For example, presynaptic NMDARs at the CA3-CA1 synapse can mediate a pre-synaptic form of STDP LTD (t-LTD) in young mice (P13-021) (Andrade-Talavera et al., 2016). This t-LTD is lost in adult synapses when applying the same post-pre protocol and even converts LTP in adult animals not requiring NMDARs anymore (Pérez-Rodríguez et al., 2019;Falcón-Moya et al., 2020). Such diversity in mechanisms of synapse plasticity, even within one type of synapse, shows the limitation of the synaptic plasticity model proposed and indicates the need to extend the study by including other mechanisms such as the involvement of presynaptic NMDAR, group I metabotropic glutamate receptor (mGluR), astrocytic signaling.
The limitation of this study is the phenomenological nature of the synaptic plasticity model that relies on the NMDAR subunit-dependent synaptic plasticity functions, but does not include detailed molecular pathways of the possible LTP and LTD induction mechanisms such as protein kinase A (PKA), CaMKII, protein phosphatase 2A and 2B (PP2A, PP2B) activation and competition. The potential direction of synaptic plasticity studies is the extension of the detailed biophysical models to account for the influence of the postsynaptic NMDAR subunit effects on the biochemical pathways of CaMKII, PKA, PP2A/2B activation in LTP and LTD induction. On the other hand, the model allows reducing the complexity of the description of the molecular network underlying LTP and LTD. We were also confronted with sometimes seemingly contradicting experimental data of synaptic outcomes that probably depend on the developmental stage of the animal, extracellular or intracellular solution compositions, and other variables linked to the different experimental settings. For example, pre-post pairings of synaptic activity lead to LTD for 5 Hz stimulation (Wittenberg and Wang, 2006), but to LTP in Inglebert et al. (2020), and it might be explained by the extracellular calcium concentration (Inglebert et al., 2020).
The model could in future studies be extended to include more complexity and account not only for postsynaptic Glu2NA NMDAR and Glu2NB-NMDAR, but also for the influence of mGluR activation, endocannabinoind, astrocytic signaling, presynaptic Glu2NC/D-NMDAR mediating effect on synaptic plasticity (Andrade- Talavera  Current trend in neuroscience is shifting the focus toward studying the relationships between different levels and scales of brain organization to better understand the nervous systems and, ultimately, human behavior. Strong emphasis in the field is to connect these different levels and scales by multiscale techniques, to allow better exploration of the information flow between the cellular-, molecular-, and network/circuit-level phenomena, and the cognitive processes and behavior. Therefore, the presented model of synaptic plasticity that enables linking synapse-specific NMDAR function to the cell and network behavior is an important step toward understanding learning in hippocampal networks.

Data availability statement
The original contributions presented in the study are included in the article/Supplementary material, further inquiries can be directed to the corresponding author.