Hélène Marie

Romain Veltz
email: romain.veltz@inria.fr

Eduardo R Zimmer

Evandro Manica

Elie L Bienenstock

Leon N Cooper

Paul W Munro

66CHAPTER 5. A STOCHASTIC MODEL OF HIPPOCAMPAL SYNAPTIC PLASTICITY WITH GEOMETRICAL READOUT OF ENZYMES -ARTICLE

Keywords: computational neuroscience, synaptic plasticity, learning rules, memory, hippocampus, enzymes modelling, experimental conditions

How learning occurs has been a long-standing question in neuroscience. Since the first demonstration that the strength wiring up neurons can be persistent, the study of neuronal connections, the synapses, became a path to understanding memory formation. In the 70s, the first electrophysiology methods to modify the synaptic strength were discovered, leading to the evidence of how synapses subjected to stimulation are plastic. Such a form of synaptic plasticity was predicted two decades before by a theoretical synaptic rule coined by the neuropsychologist Donald Hebb. The possibility that a devised rule could explain memory motivated the birth of theories providing new questions and mechanistic representations of the brain's functioning.

The diversification of techniques allowed researchers to investigate in depth the nature of synaptic rules. However, the heterogeneity of experimental conditions adopted by different laboratories implicated that the same stimulation pattern could produce different synaptic modifications. The observed heterogeneity in the methods and outcomes have hindered the formalization of a coherent view on how synaptic plasticity works. To fill this gap, this thesis developed a stochastic computational model of the rat CA3-CA1 glutamatergic synapse to explain and gain insights into how experimental conditions affect plasticity outcomes. I uncovered a new plasticity rule that accounts for methodological differences such as developmental aspects, extracellular medium and temperature influences on synaptic plasticity outcome. The model relies on an expanded version of the previous methods to predict synaptic plasticity, modified to handle combined dynamics. That is achieved by introducing a geometrical readout to interpret the dynamics of two calcium-binding enzymes controlling plasticity induction. In this way, the model covers classical and recent stimulation paradigms (e.g. STDP, FDP) using a single rule parameter set. Finally, the model's robustness is tested for in vivo-like spike time irregularity showing how different protocols converge iii to the same outcome when regularity is altered. This model allows one to obtain testable predictions since it links the simulated variables to the specificity needed to describe a plasticity protocol. Although the model is specific to a single CA3-CA1 synapse, the study's insights may be generalized to other types, enabling a deeper understanding of the rules of synaptic plasticity and learning.

Résumé

La façon dont l'apprentissage se produit est une question de longue date à laquelle la recherche en neuroscience tente de répondre. Depuis la première démonstration que la force de câblage des neurones peut être persistante, l'étude des connexions neuronales, les synapses, est devenue une voie de compréhension pour la formation de la mémoire. Dans les années 70, les premières méthodes d'électrophysiologie pour modifier la force synaptique a été découverte, menant à la preuve que les synapses stimulées ont des propriétés plastiques. Une telle forme de plasticité synaptique avait été prédite deux décennies auparavant par une règle synaptique théorique inventée par le neuropsychologue Donald Hebb. La possibilité qu'une règle élaborée puisse expliquer la mémoire a motivé la naissance de théories apportant de nouvelles questions et des représentations mécanistiques du fonctionnement du cerveau.

La diversification des techniques expérimentales a permis aux chercheurs d'enquêter en profondeur sur la nature des règles synaptiques. Cependant, l'hétérogénéité des conditions expérimentales adoptées par différents laboratoires impliquait que le même schéma de stimulation pourrait produire différentes modifications synaptiques. L'hétérogénéité observée dans les méthodes et les résultats ont entravé la formalisation d'une vision cohérente du fonctionnement de la plasticité synaptique. Pour combler cette lacune, pendant cette thèse, j'ai développé un modèle stochastique neurocomputationel de la synapse CA3-CA1 glutamatergique de rat pour expliquer et obtenir des informations sur la manière dont les conditions expérimentales affectent les résultats de cette plasticité. J'ai découvert une nouvelle règle de plasticité qui tient compte des différences méthodologiques telles que les aspects développementaux, l'influence du milieu extracellulaire et de la température sur l'issue de la plasticité synaptique. Le modèle repose sur une version étendue des méthodes précédentes pour prédire la plasticité synaptique, modifiée pour gérer la dynamique combinée. Cela est possible en introduisant v une lecture géométrique afin d'interpréter la dynamique de deux enzymes de liaison au calcium contrôlant l'induction de la plasticité. De cette façon, le modèle couvre les paradigmes de stimulation classiques et récents (par exemple STDP, FDP, BTSP) en utilisant un seul jeu de paramètres. Enfin, la robustesse du modèle est testée dans un contexte mimant l'irrégularité des décharges de potentiels d'actions in vivo montrant comment différents protocoles convergent vers la même issue de plasticité synaptique lorsque la régularité est altérée. Ce modèle permet d'obtenir des prédictions testables de façon expérimentale car il relie les variables simulées à la spécificité requise pour décrire un protocole de plasticité. Bien que le modèle soit spécifique à la synapse CA3-CA1, les résultats de cette étude peuvent être généralisés à d'autres types de synapses, permettant une meilleure compréhension des règles de la plasticité synaptique et de l'apprentissage.

Mots clés : neurosciences computationnelles, plasticité synaptique, règles d'apprentissage, mémoire, hippocampe, modélisation des enzymes, conditions expérimentales vi vii I dedicate this thesis to my family, partners, friends and teachers with a special mention to the person that taught me about love and perseverance, Serenita (aka Nena).

viii

List of abbreviations

Introduction

We live in a world of constant change, in which the environment forces organisms to adapt. Such adaptability is achieved through learning and memory, which are properties beyond simply sensing the environment. Learning is rooted to experience and can modify how information is processed, and therefore, how thoughts, actions and behaviours occur. The learning property is so fundamental that it is even observed in brainless fungi [START_REF] Nakagaki | Maze-solving by an amoeboid organism[END_REF], gifted to machines (Friedman et al. 2001) and imprinted in materials (Tian et al. 2017). For brained animals, neurons can learn through a process called plasticity, which is currently understood as the induction of persistent changes in the neuronal function, constitution, and interactions with other neurons. Such long-lasting changes [START_REF] Abraham | How long will long-term potentiation last?[END_REF] are triggered by the activity of brain cells that can connect in various ways. A perspective on this diversity is represented by the number of possible connections among the neurons of a worm (c. elegans) which can host a large number of configurations stemming from only 9 of its 302 neurons dedicated to movement [START_REF] Rakowski | Optimal synaptic signaling connectome for locomotory behavior in Caenorhabditis elegans: Design minimizing energy cost[END_REF]. Furthermore, each neuron connection has a different strength which can change under sufficient stimulation. The hotspot of these changes is called the synapse, the locus where neurons "touch" each others (συναψιζ, syn συν, ensemble and haptein αψιζ, touching -Sherrington 1910). Owning to understand learning, neuroscientists have investigated the conditions ruling synaptic efficacy. However, experiments have shown that the same activity patterns do not always produce the same plasticity outcome. This thesis presents a new plasticity rule to account for the critical conditions controlling how neurons wire up together. The synaptic rule proposed in this thesis considers biophysical conditions to explain the diversity of plasticity outcomes. Also, it builds over the original ideas of John Lisman and other names of neuroscience which the historical contributions to plasticity are drafted next.

From Cajal to Lisman

Historically, theoretical and experimental works have alternated for the primacy of the progress of synaptic plasticity understanding as shown in Figure 1.1. A mini-historical review can start by referring to Ramón y Cajal 1894 and his description of neuronal structures using the Golgi staining process. This unveiled the tiny neuronal structures called dendritic spikes (espinas) sprouting from dendrites with unprecedented details. At the beginning of the twentieth century, the mathematical formalization of neuronal theories was in its infancy. For instance, the McCulloch-Pitts model [START_REF] Mcculloch | A logical calculus of the ideas immanent in nervous activity[END_REF] laid the concepts still in use in artificial neuronal networks nowadays. Soon after, the first plasticity rule describing how a source neuron takes part in changing the excitability of its downstream neurons was stated in Hebb 1949. Hebb's rule precedes the first electrophysiology methods to induce plasticity in the hippocampus developed in the 70s (Bliss and Lømo 1973;Lynch et al. 1977). Simultaneously to the experimental confirmation of plasticity induction in 1966-77, Sejnowski 1977 anticipated the "bidirectional learning rules" based on the precise spiking times of connected neurons, currently known as spike time-dependent plasticity (STDP). Experimental evidence for this rule would be found years after by [START_REF] Markram | Action potentials propagating back into dendrites trigger changes in efficacy of single-axon synapses between layer V pyramidal neurons[END_REF] later by Bi and Poo 1998. The stages of the STDP rule understanding can be seen from different levels of description and usage, 1) an abstract implementation level brought by [START_REF] Sejnowski | Storing covariance with nonlinearly interacting neurons[END_REF], without mechanistic understanding 2) an experimental confirmation by Markram and Sakmann 1995 and a more popular definition given by Bi and Poo 1998, and 3) a current widespread use of STDP in neuromorphic computers and bio-inspired materials [START_REF] Prezioso | Spike-timing-dependent plasticity learning of coincidence detection with passively integrated memristive circuits[END_REF].

The different levels of implementation of the BCM (Bienenstock, Cooper, and Munro) rule (Bienenstock et al. 1982) shared a similar path to the STDP rule. The BCM rule is based on the spike frequency (instead of spike times) and applied to visual cortex data (Wiesel and Hubel 1963), and such frequency-dependent rule was initially uncovered by Dunwiddie and Lynch 1978. Nowadays, the BCM and STDP rules are taken as canonical in bio-inspired applications without proper consideration of their physiological relevance (Tigaret et al. 2016). The Bienenstock et al. 1982 rule uses a threshold to separate opposite types of plasticity. The idea of threshold in biology can be traced back to the first neuronal model McCulloch-Pitts model. This thesis will expand the threshold concept to predict plasticity.

With the possibility of visualizing the calcium dynamics inside the neurons thanks to the green fluorescent protein (GFP) discovery, more specifically cal-the plasticity bidirectionality is controlled by two-thresholds over either high or moderate calcium levels. His original work also suggests an enzymatic maintenance role for the kinase and phosphatase interactions.

Bi Guo-Qiang (毕国强) and Mu-ming Poo (蒲慕明) using cultured hippocampal neurons helped to popularize the spike time-dependent plasticity which currently is largely applied in bioinspired applications.

Figure 1.1:

A short roadmap of plasticity experiment and theory hallmarks. A short roadmap of plasticity experiment and theory hallmarks through the 20th century. Note that the only woman cited was Brenda Milner, who has continuously contributed to the understanding of memory (even now at 102 years old). Women's contributions to neuroscience are not widely available or mentioned; fortunately, a project to highlight these has been pushed forward in a timely and needed rescue of the scientific heritage.

cium indicators with a GFP domain, plasticity rules branched and diversified. For instance, parallel to the synaptic rules based on neurons spiking times, Lisman 1989 proposed a new model hypothesizing how calcium-binding molecules induce and maintain plasticity at the molecular level motivated by precursor studies on protein phosphorylation pathways from [START_REF] Nestler | Protein phosphorylation in the brain[END_REF]. This thesis will also expand such original idea.

Despite the advances brought by the plasticity models, the growing number of plasticity experiments has caused few theories to expire since they were not sufficiently general or had central assumptions proved false. As a benefit, these improved the synergy between theory and experimental investigation, motivating new questions on how plasticity occurs. However, mostly overlooking experimental heterogeneity. The problem addressed in this thesis is how to incorporate the experimental diversity in a unified theory capable of providing testable predictions.

Reductionist representations of synaptic mechanisms

Synaptic models aiming to elucidate plasticity induction have been built on simplified representations of the biological entities and interactions in different granularities [START_REF] Brette | Frankenstein models and the double-edged sword of modelling tools[END_REF][START_REF] Marr | Vision: A computational investigation into the human representation and processing of visual information[END_REF]. Similarly, non-mechanistic approaches (e.g. deep learning) can make accurate predictions. However, they do not attend a biological description level (e.g. components of the synaptic nano/microdomains, calcium sources) to explain the plasticity mechanisms (Kievit et al. 2011). In this thesis, representations of synaptic mechanisms, or reduction of real biological entities Kievit et al. 2011), will be used to gain insights on the process of plasticity induction. Unavoidably, the simulations of mechanisms activating neurons are described only partially by models. For example, the equations describing ion channels kinetics and their collective behaviour simplify the different states, localizations and compositions of these units. An intermediate approach is to model individual ion channel contributions to the neuron electrophysiology of using a stochastic approach. Here, synaptic mechanisms are implemented using a partially deterministic Markov process (see Annexe 1). The postsynaptic mechanisms are modelled with more biophysical detail, including different stochastic sources (Ribrault et al. 2011), and signals originated in the soma are modelled phenomenologically.

Unified theories of plasticity

The growing experimental evidence solidified that precise temporal spike patterns alone are insufficient to explain plasticity. The characterization of mechanisms propelled the emergence of unified theories for synaptic plasticity. These theories expanded the mechanistic understanding of synaptic plasticity, gaining further insight into the induction process by encompassing multiple plasticityrelated modifiers: spikes, calcium, neuron morphology, calcium-binding molecules, temperature, extracellular and extracellular ion concentration, astrocytes, etc. Such a growing number of interactions poses an analysis challenge since the interpretation of multiple dynamics did not advance simultaneously with the higher number of mechanisms included. This thesis will address such methodological shortcomings by proposing a novel prediction strategy more suitable to multidimensional simulations.

Research problem

It is hard to mimic nature. Laboratories tend to use their recipe variations (Tebaykin et al. 2018) without justifications. For instance, several ex vivo studies in which the brain slices are extracted and prepared for electrophysiology do not recreate the neurons' conditions. For example, non-physiological temperatures are used since the first description of LTP induction (Bliss and Lømo 1973). The mentioned problem is not only restricted to ex vivo experimentation. In vivo recordings are also affected by the sampling process. Electrodes or optogenetic probes inserted in living brains can cause inflammation or heat damage during stimulating [START_REF] Misra | Preventing neuronal damage and inflammation in vivo during cortical microelectrode implantation through the use of Poloxamer P-188[END_REF]. These mixed experimental approaches can cause a reproducibility problem, and as mentioned by Wittenberg et al. 2006, laboratory differences can be a source of divergent plasticity outcomes even using a similar stimulation. This thesis intends to provide a mechanistic explanation of how plasticity outcomes, sometimes qualitatively conflicting and sampled from different experimental conditions, are part of the same phenomenon.

The thesis is structured as follows: chapter 2 introduces plasticity-related terms; chapter 3 discusses unified theories, their strategies and mechanisms; chapter 4 formalizes the gap in these unified theories and presents the conceptual approach developed in the results; chapter 5 states the context in which this research was done and how it led to the article containing the main results, methods and supplemental results; chapter 6, adds additional results; and, chapter 7 discusses the results and comments on future perspectives. Ribrault, Claire, Ken Sekimoto, and Antoine Triller (2011). "From the stochasticity of molecular processes to the variability of synaptic transmission". In: Nature Reviews Neuroscience 12.7, pp. 375-387. Sejnowski, Terrence J (1977). "Storing covariance with nonlinearly interacting neurons". In: Journal of mathematical biology 4.4, pp. 303-321. Sherrington, Charles Scott (1910). "Flexion-reflex of the limb, crossed extension-reflex, and reflex stepping and standing". In: The Journal of physiology 40.1-2, pp. 28-121. Tebaykin, Dmitry, Shreejoy J Tripathy, Nathalie Binnion, Brenna Li, Richard C Gerkin, and Paul Pavlidis (2018). "Modeling sources of interlaboratory variability in electrophysiological properties of mammalian neurons". In: Journal of neurophysiology 119.4, pp. 1329-1339.

Tian, He et al. (2017). "Emulating bilingual synaptic response using a junction-based artificial synaptic device". In: ACS nano 11.7, pp. 7156-7163.

Tigaret, Cezar M, Valeria Olivo, Josef HLP Sadowski, Michael C Ashby, and Jack R Mellor (2016).

"Coordinated activation of distinct Ca 2+ sources and metabotropic glutamate receptors encodes Hebbian synaptic plasticity". In: Nature communications 7, p. 10289.

Wiesel, Torsten N and David H Hubel (1963). "Single-cell responses in striate cortex of kittens deprived of vision in one eye". In: Journal of neurophysiology 26.6, pp. 1003-1017.

Chapter 2

Synaptic plasticity and its mechanisms

Synaptic plasticity regulates the transmission efficacy between pre and postsynaptic neurons. This process starts with a Ca 2+ influx across the neuronal membrane following the activation of postsynaptic contents such as receptors, ion channels, enzymes, and vesicle-filled neurotransmitters. The activation of synaptic contents triggers structural changes, provoking the addition, removal or remodelling of structural proteins. Also, the synaptic content undergoes functional changes, causing mechanisms to alter their biophysical properties, such as their conductance (charge movement), permeability (ion flux) or molecular conformation. The combination of functional and structural changes may lead to a potentiated or depressed synaptic efficacy. However, despite the accumulated knowledge about synaptic plasticity, the phenomenon is not fully understood due to the heterogeneity of electrophysiological approaches (Tebaykin et al. 2018) and the lack of biophysical characterization of synaptic contents [START_REF] Heil | Analysis of proteins in computational models of synaptic plasticity[END_REF]. This chapter sums up plasticity forms, synaptic components, mechanisms, the synaptic sensibility to extracellular ion concentrations and my main hypothesis concerning plasticity induction.

The timescales of plasticity

Pre and postsynaptic mechanisms such as ion channels, receptors, enzymes are regularly activated, adapt at different timescales according to their kinetics and available resources. Thanks to the interaction of these mechanisms, neurons can communicate and induce plasticity. Figure 2.1 exemplifies the short-and longterm timescales of different plasticity forms discussed next.

Figure 2.1: Short-and long-term plasticity timescales. A) AMPAr desensitization occurs when glutamate binds to the AMPAr, and the pore does not open-depending on how deep is the desensitization, the full recovery of AMPA current can take from 100 to 300 ms given the AMPAr subunits composition (Robert and Howe 2003). B) Short-term facilitation is thought to be caused by Ca 2+ influx through VGCCs which accumulates in the presynaptic terminal, enhancing the probability of vesicle release (The Residual Ca 2+ Hypothesis) (Debanne et al. 1996). However, other mechanisms can cause facilitation, such as buffers saturation and spike broadening [START_REF] Jackman | The mechanisms and functions of synaptic facilitation[END_REF]. C) Under constant firing the vesicle pools become empty causing short-term depression (Fernández-Alfonso and Ryan 2004; Rizzoli and Betz 2005). D) Long-term plasticity phases showing structural changes occurring in dendritic spines [START_REF] Bosch | Structural and molecular remodeling of dendritic spine substructures during long-term potentiation[END_REF]).

Short-term plasticity

Short-term plasticity (STP) consists of adaptations at fast timescales producing immediate changes in the excitability or resources availability (Deperrois and Graupner 2020;Tsodyks and Markram 1997). In the literature [START_REF] Zucker | Short-term synaptic plasticity[END_REF]Rizzoli and Betz 2005), STP refers to a phenomena, such as receptor desensitization (∼1 to 25 ms) and ion-channels inactivation, facilitation of vesicle release (∼50 to 200 ms), or depletion of neurotransmitter filled vesicles in the presynaptic terminals (∼50 to 2000 ms). The listed STP examples adapt in milliseconds to seconds, and their effect exceeds such range (Deperrois and Graupner 2020). Also, STP can be characterized as an all-or-none process, such as the vesicle release and its consequent short-term depression or facilitation evaluated in multiple trials. Although short-term plasticity refers to millisecond range adaptations, faster processes (fraction of milliseconds) are also prone to adaptations (Ribrault et al. 2011;Dobrunz et al. 1997).

Long-term plasticity

Two forms of long-term plasticity, long-term potentiation (LTP) and long-term depression (LTD), reflects the long-lasting effects on synaptic efficacy. LTP was first characterized by Bliss and Lømo 1973 with local field potential in the CA3-CA1 pathway using a high-frequency stimulation (HFS). Years after, Lynch et al. 1977 would demonstrate LTD in the same pathway using a slow-frequency stimulation (LFS). Both forms of plasticity and their importance for understanding memory were debated in the context of hippocampal lesion studies [START_REF] Douglas | The hippocampus and behavior[END_REF]. LTP or LTD have protein reorganization phases. In these phases, the dendritic spines malleable body are reshaped minutes after or during the induction protocol [START_REF] Bosch | Structural and molecular remodeling of dendritic spine substructures during long-term potentiation[END_REF]. Initially, LTP and LTD undergo structural reorganization and stabilization (early long-term plasticity), consolidation, and later new protein synthesis (late long-term plasticity), see Figure 2.1D. Fluorescence markers on multiple structural proteins identified the rhythms of LTP induction from the early phase, occurring until ∼1 hour after the stimulation, to late potentiation phase, occurring after ∼ 1 hour [START_REF] Bosch | Structural and molecular remodeling of dendritic spine substructures during long-term potentiation[END_REF]). Usually, structural protein studies investigate the different plasticity phases with a time resolution of minutes or dozens of seconds to follow long term changes [START_REF] Bosch | Structural and molecular remodeling of dendritic spine substructures during long-term potentiation[END_REF][START_REF] Thomazeau | Dissociation of functional and structural plasticity of dendritic spines during NMDAR and mGluR-dependent long-term synaptic depression in wild-type and fragile X model mice[END_REF]. Few other studies have measured long-term plasticity phases at shorter time scales [START_REF] Frost | Singlemolecule discrimination of discrete perisynaptic and distributed sites of actin filament assembly within dendritic spines[END_REF]. Also, [START_REF] Bosch | Structural and molecular remodeling of dendritic spine substructures during long-term potentiation[END_REF] characterization of early and late-LTP was performed at room temperature, in which enzyme activity is expected to react (under equivalent stimulation protocol) slower than physiological temperatures, potentially modifying the timescales of structural changes [START_REF] Thomazeau | Dissociation of functional and structural plasticity of dendritic spines during NMDAR and mGluR-dependent long-term synaptic depression in wild-type and fragile X model mice[END_REF]. Next, receptors and ion channels, the primary Ca 2+ sources provoking structural and functional changes, will be discussed.

Receptors and ion-channels

NMDA receptor

N-methyl-D-aspartate receptors (NMDAr) are the primary source of synaptically induced Ca 2+ influx in dendritic spines (Yuste et al. 1999); they have a tetramer arrangement composed of subunits GluN1, GluN2 (A-D) and GluN3 [START_REF] Salussolia | Arrangement of subunits in functional NMDA receptors[END_REF]. The NMDAr's subtype expression is age-dependent, producing receptors with different kinetics and functions within the brain. Both glutamate and glycine activate NMDArs leading to ionotropic (K + , Na + and Ca 2+ ), and metabotropic activity depending on NMDAr composition. A concept to understand the role of NMDAr is the coincidence detection in which the NMDAr-Mg 2+ unblocking occurs when dendritic spines are sufficiently depolarized [START_REF] Collingridge | Excitatory amino acids in synaptic transmission in the Schaffer collateral-commissural pathway of the rat hippocampus[END_REF]. That enhances NMDAr conductance since the NMDAr blocking is alleviated by the BaP current reaching dendritic spines. Another form of NMDAr conductance modulation occurs by changing extracellular Ca 2+ as shown in Maki and Popescu 2014.

AMPA receptor

α-amino-3-hydroxy-5-methyl-4-isoxazolepropionic receptor (AMPAr) are tetrameters composed of GluA1-A4 subunits. These glutamatergic receptors are majorly expressed in the dendritic spines, and their central role is to relay excitatory neuronal impulses [START_REF] Lu | Subunit composition of synaptic AMPA receptors revealed by a single-cell genetic approach[END_REF]. It is of particular interest for plasticity since the inclusion and removal of AMPArs from the postsynaptic density (PSD) reflects the plasticity outcomes illustrated in Figure 2.2. Due to this, AMPArs density in den-dritic spines correlates with PSD size. Unlike NMDAr, given the synapse type, AMPArs have a specific distance-dependent expression to compensate EPSP attenuation towards the soma [START_REF] Nicholson | Distance-dependent differences in synapse number and AMPA receptor expression in hippocampal CA1 pyramidal neurons[END_REF]. Note also that AMPAr current can have its magnitude changed even without the inclusion of new units due to the phosphorylation of its subunits [START_REF] Incontro | The CaMKII/NMDA receptor complex controls hippocampal synaptic transmission by kinase-dependent and independent mechanisms[END_REF]).

Voltage gated Ca 2+ channels

An important Ca 2+ source in plasticity induction is the voltage-gated Ca 2+ channels (VGCCs) which are a class of ion channels permeable to Ca 2+ ions, activated by changes in the membrane potential. VGCCs (L, N, P/Q, R and T) have specific distribution and variations within the brain and non-neuronal cells (Magee and Johnston 1995). A classification for the VGCCs subtypes attributes their activation to a particular membrane depolarization level. For instance, low-, intermediaryand high-voltage activated channels, such as T-type VGCCs, activated at resting membrane conditions and L-type VGCCs activated at high voltage (Magee and Johnston 1995). Although VGCCs provides less Ca 2+ (synaptic induced Yuste et al. 1999) compared to NMDAr in dendritic spines there are VGCC-dependent forms plasticity [START_REF] Adermark | Combined activation of L-type Ca2+ channels and synaptic transmission is sufficient to induce striatal long-term depression[END_REF]. However, the VGCCs major contribution is the presynaptic vesicle release triggered by responding to a Ca 2+ influx caused by presynaptic action potentials [START_REF] He | New insights into interactions of presynaptic Calcium Channel subtypes and SNARE proteins in neurotransmitter release[END_REF]). This effect is the major stochastic component in the signal transduction since its ineffective opening causes vesicle release failures (Yuste et al. 1999).

GABA receptors

So far, we have reviewed only excitatory mechanisms (VGCC, AMPAr and NM-DAr). However, plasticity requires also the participation of inhibitory mechanisms (Meredith et al. 2003). For instance, the gamma-aminobutyric acid receptor (GABAr) [START_REF] Rowlett | Different GABAA receptor subtypes mediate the anxiolytic, abuse-related, and motor effects of benzodiazepine-like drugs in primates[END_REF]) is activated by the presence of GABA neurotransmitters which opens for the passage for Cl -and K + , leading to inhibitory currents to cause cell polarization [START_REF] Pelkey | Hippocampal GABAergic inhibitory interneurons[END_REF]. GABAr controls neuronal excitability and the EPSPs integration since inhibitory currents act in different parts of the neuron, such as soma, dendrites and spines [START_REF] Somogyi | Defined types of cortical interneurone structure space and spike timing in the hippocampus[END_REF]. By interacting with Ca 2+ sources, GABAr can modify plasticity induction outcomes as shown by blocking studies (Meredith et al. 2003;Buchanan and Mellor 2007).

Small Potassium ion channel

A class of potassium channel is activated in a [Ca 2+ ] i -dependent manner, hyperpolarizing the membrane potential and reducing the NMDAr and VGCCs enhancement caused by membrane depolarization. For instance, the small potassium (SK) channel which interacts with Ca 2+ sources (Griffith et al. 2016) to reduce Ca 2+ influx, therefore affecting plasticity outcomes (Tigaret et al. 2016). This class of channels limits the Ca 2+ influx and the Ca 2+ -induced excitotoxicity [START_REF] Trombetta-Lima | Calcium-activated potassium channels: Implications for aging and age-related neurodegeneration[END_REF]).

The Ca 2+ dynamics

Plasticity of dendritic spines depends on intracellular Ca 2+ elevations and it is mainly controlled by Ca 2+ sources such as [START_REF] Yuste | Dendritic spines[END_REF]): NMDAr, Ca 2+ permeable AMPAr, VGCCs, Kainate receptors, Ca 2+ storages from endoplasmatic reticulum (ER). However, Ca 2+ dynamics impact other neuronal functions such as neuronal transmission, plasticity, neurite growth, apoptosis, ATP production, protein expression, to cite a few (Knot et al. 2005). Due to the multiple functions of Ca 2+ , its influx and extrusion are tightly controlled in dendritic spines. Therefore, dendritic spines actively compartmentalise Ca 2+ with its high extrusion efficiency (outside the cell), which is supposed to exceed Ca 2+ diffusion through the spine neck (Sabatini and Svoboda 2000). Such efficiency at extruding Ca 2+ from the tiny dendritic spine's cytoplasm is granted by Ca 2+ pumps (e.g. PMCA), mobile buffers (e.g. Calmodulin) or also Ca 2+ ER storages (e.g. SERCA pumps) which can prevent depotentiation [START_REF] Perez-Alvarez | Endoplasmic reticulum visits highly active spines and prevents runaway potentiation of synapses[END_REF]. Despite this tight control, internal Ca 2+ fluctuates due to tonic T-type VGCCs flickering, keeping a basal Ca 2+ concentration of around 50 nM (Maravall et al. 2000). The main method to indirectly probe the Ca 2+ dynamics is through fluorescent Ca 2+ -binding proteins. However, this form of indirect visualization has some drawbacks such as laserinduced phototoxicity. The laser exciting the fluorescent proteins in the synaptic domains can modify the receptors and ion channels' kinetic parameters since the temperature increases and it can cause photodamage (Podgorski and Ranganathan 2016;Schmidt and Oheim 2018) given the laser parameters.

The Ca 2+ hypothesis for plasticity

Ca 2+ is the primary signalling ion in the brain (Knot et al. 2005) and since the discovery that intracellular levels of Ca 2+ are fundamental to induce both LTP and LTD [START_REF] Lynch | Intracellular injections of EGTA block induction of hippocampal long-term potentiation[END_REF]Mulkey and Malenka 1992) neuroscientists have searched for a relationship between plasticity outcomes and Ca 2+ . Lisman 1989 coined a hypothesis linking the Ca 2+ amplitude to plasticity directionality at the same theoretical work that proposes a molecular mechanism for memory storage based on enzyme activity. Lisman statement can be rewritten in the following way: high levels of Ca 2+ influx lead to potentiation, while moderate levels lead to depression. Such a hypothesis, illustrated in Figure 2.3 (Left), finds support in a frequencydependent plasticity experiment done by [START_REF] Hansel | Different threshold levels of postsynaptic [Ca2+] i have to be reached to induce LTP and LTD in neocortical pyramidal cells[END_REF] using Ca 2+ fluorescence imaging. Despite that, direct observation of Ca 2+ without fluorescent dye prevents a proper evaluation of the amplitude hypothesis. Moreover, evidence against the Ca 2+ amplitude hypothesis is given by Tigaret et al. 2016, which shows a high Ca 2+ amplitude (measured with Ca 2+ dye) induced not LTP, while a moderate level induced LTP does. Also by Nevian and Sakmann 2006 that shows the same Ca 2+ amplitude can induce both LTP and LTD. The two last examples make such a claim without measuring the whole Ca 2+ plasticity protocols due to phototoxicity issues. Another variation of it the Ca 2+ hypothesis, the Ca 2+ duration hypothesis suggests the duration of Ca 2+ intracellular elevations dictates plasticity directionality (Mizuno et al. 2001). That can be stated as brief and high Ca 2+ levels induce LTP, while sustained and moderate Ca 2+ levels induce LTD, as depicted in Figure 2.3 (Middle). Despite such hypothesis has received favourable support by Yang et al. 1999;Mizuno et al. 2001, precise control of Ca 2+ is hard to be achieved using only spikes (Yang et al. 1999). Therefore, part of the support of the duration hypothesis comes from computational studies [START_REF] Gamble | The dynamics of free calcium in dendritic spines in response to repetitive synaptic input[END_REF]. Other ideas are proposed, such as the Ca 2+ location hypothesis (Figure 2.3,Right) that states depending on the source of Ca 2+ (e.g. NMDAr or VGCCs), a different directionality would be selected. However, it falls in similar technical limitations from the previous hypothesis when controlling Ca 2+ entry (Evans and Blackwell 2015). Given the shortcomings, it is interesting to point out a simple enquiring by [START_REF] Malenka | Long-term depression: not so depressing after all[END_REF] It is well established that a rise in postsynaptic Ca 2+ is necessary for LTP induction. How can a rise in postsynaptic Ca2+ also be responsible for LTD? [...]. Ca 2+ is extruded very fast from the spines. Therefore it has a relatively fast decay time (50 -500 ms) when compared to enzymes. This thesis will try to interpret the dynamics of "what comes after Ca 2+ " since new measurements are available for the immediate Ca 2+ -cascade, and this can lead to a new interpretation of Malenka's observation and hopefully make it more clear how different plasticity experiments, with different experimental conditions, are part of the same phenomenon.

Calmodulin

Calmodulin (CaM) is a calcium-binding enzyme highly expressed in dendritic spines interacting with LTP and LTD pathways (see Figure 2.4), and other calcium-related regulatory mechanisms (e.g. pain [START_REF] Hasan | Calmodulin is responsible for Ca 2+-dependent regulation of TRPA1 Channels[END_REF], vesicle release Pang et al. 2010). It is a mobile buffer that efficiently detects the Ca 2+ influx due to its high affinity and fast binding [START_REF] Faas | Calmodulin as a direct detector of Ca 2+ signals[END_REF]. CaM dependent enzymes, such as CaMKII and CaN, which inhibition and blocking can disrupt plasticity induction (O'Connor et al. 2005) are discussed next.

CaMKII

The Ca 2+ /Calmodulin-dependent kinase II, or CaMKII, has four isoforms (α, β, γ, and δ, Gaertner et al. 2004) associated with the induction of long term plasticity (LTP and LTD) (Coultrap and Bayer 2012). This enzyme is composed of 12 subunits, forming a holoenzyme that can phosphorylate (activate by adding a phosphate group) its neighbouring units allowing the enzyme to be activated longer in the absence of Ca 2+ . Such autophosphorylation property was hypothesized to have a maintenance role in dendritic spines (Lisman 1989). When CaMKII is activated, it would reach a stable state due to its autophosphorylation, which would surpass protein turnover by having new included CaMKII becoming self-phosphorylated at Ca 2+ basal levels and saturation of dephosphorylation pathways (Lisman 1985;Miller et al. 2005). Despite that, a maintenance role through autophosphorylation was not confirmed in brain slices by Förster resonance energy transfer (FRET)-sensor experiments recording CaMKII activity (Lee et al. 2009). Instead, experiments have shown that CaMKII activity returns to basal levels instead of staying in an activated (self-sustained) state (Chang et al. 2017;Chang et al. 2019). Other evidence for the lack of maintenance role was obtained by applying CaMKII inhibitors after an LTP protocol, producing no changes in the LTP expression [START_REF] Chen | Is persistent activity of calcium/calmodulin-dependent kinase required for the maintenance of LTP?[END_REF]. Despite that, CaMKII was found to modify the NMDAr function by biding on its GluN2B subunit and forming an NMDAr-CaMKII complex necessary for LTP induction and possibly maintenance [START_REF] Incontro | The CaMKII/NMDA receptor complex controls hippocampal synaptic transmission by kinase-dependent and independent mechanisms[END_REF]. The maintenance role of NMDAr-CaMKII was suggested by test-tube essays (Urakubo et al. 2014). To date, a self-sustained CaMKII activity to justify a maintenance role has not been confirmed for LTP (Michalski 2014;Coultrap and Bayer 2012). More studies are necessary to evaluate the NMDAr-CaMKII complex maintenance role in more physiological conditions (Urakubo et al. 2014). Rather than LTP maintenance, the CaMKII induction role was demonstrated by knocking-out CaMKIIα, the most abundant isoform in excitatory synapses together with CaMKIIβ, causing impairment of learning and memory capabilities [START_REF] Silva | Impaired spatial learning in alpha-calcium-calmodulin kinase II mutant mice[END_REF][START_REF] Incontro | The CaMKII/NMDA receptor complex controls hippocampal synaptic transmission by kinase-dependent and independent mechanisms[END_REF]). Still, CaMKII may have a maintenance role for LTD (Coultrap and Bayer 2012). Another aspect of Lisman's hypothesis is that the phosphatase activity causes CaMKII to dephosphorylate. However, by blocking a series of synaptic phosphatases, Otmakhov et al. 2015 found CaMKII dephosphorylation unaltered.

CaMKII has structural and functional roles. The formation of the NMDAr-CaMKII complex increases the Ca 2+ influx and CaMKII phosphorylated AMPAr has increased conductance. Because CaMKII can increase Ca 2+ influx, CaMKII excitability path-ways have been highlighted as a drug target to prevent excitotoxicity [START_REF] Coultrap | CaMKII in cerebral ischemia[END_REF]. Structurally, CaMKII isoforms differentially interact with dendritic spines structural proteins [START_REF] Okamoto | The roles of CaMKII and F-actin in the structural plasticity of dendritic spines: a potential molecular identity of a synaptic tag?[END_REF], for instance, by bundling together with F-actin filaments that compose dendritic spine body and stimulating neurite growth [START_REF] Okamoto | The roles of CaMKII and F-actin in the structural plasticity of dendritic spines: a potential molecular identity of a synaptic tag?[END_REF]). The roles associated with this enzyme have been the focus of theoretical studies owning to predict plasticity and better understand memory. However, the complete CaMKII enzyme is a modelling challenge due to the staggering number of state combinations its 12 subunits can form (Pharris et al. 2019).

Calcineurin

PP2, or Calcineurin (CaN) and its isoforms (α, β, γ), are a Calmodulin activated phosphatases mainly associated with LTD [START_REF] Creamer | Calcineurin[END_REF]O'Connor et al. 2005). It stands in Lisman's hypothesis as a mediator of LTD and reversal pathway of LTP. Indeed, CaN is relevant to LTD since inhibiting the phosphatase leads to lower LTD amplitude. As CaMKIIβ, CaN has structural and functional plasticity roles such as the dephosphorylation of AMPAr leading to LTD [START_REF] Woolfrey | Coordination of protein phosphorylation and dephosphorylation in synaptic plasticity[END_REF] and NMDAr desensitization by interacting with the GluN2A subunit [START_REF] Rycroft | Inhibitory interactions of calcineurin (phosphatase 2B) and calmodulin on rat hippocampal NMDA receptors[END_REF]. Structural plasticity effects include regulation microtubulesstabilizing proteins such as tau (axons) and MAP2 (dendrites) [START_REF] Hoffman | The protein serine/threonine phosphatases PP2A, PP1 and calcineurin: A triple threat in the regulation of the neuronal cytoskeleton[END_REF]. CaN-tau regulation is a relevant path to study tauopathies and other neurodegenerative diseases such as Alzheimer's. Excessive tau phosphorylation may lead to protein bundle accumulation with a cascade of deleterious effects from slower synaptic transmission (L. Zhou et al. 2017) to neuronal death [START_REF] Fath | Tau-mediated cytotoxicity in a pseudohyperphosphorylation model of Alzheimer's disease[END_REF]. Still, CaN apoptotic pathways can be triggered by excessive Ca 2+ influx (F. [START_REF] Shibasaki | Calcineurin functions in Ca (2+)-activated cell death in mammalian cells[END_REF] and be explored in neuroprotective strategies [START_REF] Erin | In vitro hypoxia and excitotoxicity in human brain induce calcineurin-Bcl-2 interactions[END_REF]. Fujii et al. 2013 suggested that CaN and CaMKII enzymatic activity could represent a pulse number and frequency decoder in neuron somata. However, in dendritic spines, CaMKII measurements (Chang et al. 2017) showed a faster saturation profile given frequency (from 0.49 to 7.8 Hz) and pulse number (more than eight pulses to saturation in room temperature). Given these contradictory evidences, CaN-CaMKII as an enzymatic decoder needs to be investigated in dendritic spines where these enzymes are primarily expressed instead of neuron somata. Another experiment on CaN-CaMKII dynamics reinforces that CaN and CaMKII respectively induces LTD and LTP in a frequency-dependent manner (O'Connor et al. 2005).

CaN and CaMKII joint activity

In O'Connor et al. 2005 experiment, by inhibiting CaMKII, LTP is abolished, and the same is valid for CaN and LTD. Despite this association, CaMKII and CaN are not exclusively related to potentiation or depression [START_REF] Woolfrey | Coordination of protein phosphorylation and dephosphorylation in synaptic plasticity[END_REF]. Interestingly, one of the first transgenic animal models was a CaMKIIα knockout designed to study the learning and memory dependence of this enzyme [START_REF] Silva | Impaired spatial learning in alpha-calcium-calmodulin kinase II mutant mice[END_REF]. Note that both enzymes are related to plasticity and apoptotic pathways, disrupted CaMKII-NMDA complex induced excitotoxicity (Coultrap and Bayer 2012) and CaN-induced apoptosis (F. [START_REF] Shibasaki | Calcineurin functions in Ca (2+)-activated cell death in mammalian cells[END_REF]. Despite that, the thresholds separating the healthy and pathological functions of these enzymes are unknown.

Stimulation protocols

Experimentalists can artificially manipulate Ca 2+ influx through dendritic spines using spiking times (or the firing patterns) to investigate how neuronal activity induces plasticity. The stimulation protocol may be induced by external electrodes in brain slices exposing (or not) neuronal circuitry (Bliss and Lømo 1973); optogenetics which requires the expression of light-sensitive ion channels (e.g. channelrhodopsins) (Anisimova et al. 2019); glutamate uncaging (Ellis-Davies 2019) which releases caged glutamate from a light-sensitive molecular matrix, allowing the precise stimulation of individual spines; or also, delivering electrical current using a patch through (or at) the neuronal membrane surface [START_REF] Neher | The extracellular patch clamp: a method for resolving currents through individual open channels in biological membranes[END_REF], with various clamping variations such as current, voltage and dynamical [START_REF] Berecki | Dynamic clamp as a tool to study the functional effects of individual membrane currents[END_REF]. Also, plasticity can be induced chemically, without manipulating neuronal spiking activity [START_REF] Santschi | Chemically induced, activityindependent LTD elicited by simultaneous activation of PKG and inhibition of PKA[END_REF]. A general description of plasticity protocols highlights pre-and postsynaptic spiking times or their rate, duration, brain region, aCSF temperature, preparation steps of the neural tissue, developmental age, etc. Usually, stimulation protocols are heterogeneous like the electrophysiology methods shown in Tebaykin et al. 2018. However, most plasticity protocols use firing structures with highly correlated regular spikes (not observed in vivo measurements (Cui et al. 2018)) and non-physiological experimental conditions (Inglebert et al. 2020;Bi and Poo 1998).

The first electrophysiology measurements of LTP were recorded and induced through local field potentials (LFP) (Bliss and Lømo 1973). It consists of estimating the slope of the electromagnetic field decay produced by the movement of the charges in response to a neural pathway stimulation. This measure is preferred when recording in the extracellular medium since the signal amplitude depends on the source distance from the electrode (e.g. metal or micropipette) [START_REF] Buzsáki | The origin of extracellular fields and currents-EEG, ECoG, LFP and spikes[END_REF]. Other techniques such as patch-clamp variations and intracellular recordings are invasive but allow one to measure the neuron's current and voltage, with a precise estimation of the EPSP amplitude change [START_REF] Covey | Basic Electrophysiological Methods[END_REF]. Electrophysiological techniques have several intricacies. For instance, in in vivo recording, the correct electrode placement has to be verified by surgery.

It is not unusual to have experiments discarded due to technical issues and artefacts (Ahnaou et al. 2020). Automatized electrophysiology is still not widely available and still limited to few types of experiments [START_REF] Annecchino | Progress in automating patch clamp cellular physiology[END_REF]. Other forms of measure of the plasticity outcomes are by imaging, targeting specific proteins linked to the plasticity expression [START_REF] Bosch | Structural and molecular remodeling of dendritic spine substructures during long-term potentiation[END_REF], or imaging genetically encoded indicators (e.g. voltage and calcium) to evaluate the EPSP or its effects before being filtered by the spine neck [START_REF] Jackman | The mechanisms and functions of synaptic facilitation[END_REF]Kwon et al. 2017).

FDP and STDP

The first stimulation protocol descends from the Bliss and Lomo (Bliss and Lømo 1973) which proposes a method to induce LTP using an external stimulation in a rabbit CA3-CA1 pathway [START_REF] Nicoll | A brief history of long-term potentiation[END_REF]. Their technique continued to be used due to its practicality (Dudek and Bear 1992;Dudek and Bear 1993) even after new methods allowing control over pre and postsynaptic spiking times [START_REF] Neher | The extracellular patch clamp: a method for resolving currents through individual open channels in biological membranes[END_REF]) have been described. An example is the frequency-dependent stimulation (FDP) paradigm presented by Dunwiddie and Lynch 1978 that uses regular spikes. However, note that later experiments have used jittered spikes, such as the ratedependent in-duction of plasticity by Sjöström et al. 2001. In the 70s, when LTP and LTD were experimentally confirmed in the hippocampus, simultaneously the Bienstook Cooper Monro (BCM) theory formalized the experiments of visual cortex directional selectivity using a nonlinear function of synaptic efficacy and frequency (Lynch et al. 1977). The BCM popularized the idea that a threshold marks a transition between LTD and LTP over the frequency. For instance, the BCM-like rule found in the hippocampus [START_REF] Wang | Priming-induced shift in synaptic plasticity in the rat hippocampus[END_REF][START_REF] Holland | Primed facilitation of homosynaptic long-term depression and depotentiation in rat hippocampus[END_REF].

The experiment that established a canonical bidirectional plasticity rule was popularized by Bi and Poo 1998 and made possible due to the patch-clamp techniques (Markram et al. 2011). Bi and Poo 1998 discovered that precise timing of pre and postsynaptic spikes and its causality can describe plasticity outcomes observed in cultured neurons, with the causal pairing inducing LTP, and the anticausal one inducing LTD. The bidirectionality of STDP, experimentally demonstrated by Bi and Poo 1998, was a theoretical construct previously used by Sejnowski 1977. Currently, experiments have demonstrated that the synaptic rules STDP have multiple shapes influenced by: age (Meredith et al. 2003), temperature Krelstein et al. 1990, frequency (Sjöström et al. 2001), distance from the soma (Letzkus et al. 2006), aCSF composition (Inglebert et al. 2020), duration (Mizuno et al. 2001). To date, there is no mechanistic justification for how the experimental conditions affect the STDP outcomes. For instance, as Wittenberg and Wang 2006 interestingly points out, a similar STDP causal stimulation protocol can induce LTP [START_REF] Nishiyama | Calcium stores regulate the polarity and input specificity of synaptic modification[END_REF], LTD [START_REF] Christie | The role of dendritic action potentials and Ca2+ influx in the induction of homosynaptic long-term depression in hippocampal CA1 pyramidal neurons[END_REF]Wittenberg and Wang 2006) and no change [START_REF] Pike | Postsynaptic bursting is essential for 'Hebbian'induction of associative long-term potentiation at excitatory synapses in rat hippocampus[END_REF], without telling why such differences arise.

Sensitivity to experimental conditions

Synapses are sensitive, and experimental choices affect plasticity, as will be exemplified further. Despite that, plasticity experiments have a wide range of nonphysiological settings selected due to technical limitations in mimicking the neuron environment. Figure 2.5 shows histograms of the experimental parameters used in hallmark plasticity experimental publications described in Table 3.3 and 3.1. The data was extracted from hallmark plasticity studies reproduced by theoretical models (Graupner and Brunel 2012, Ebner et al. 2019, Jedrzejewska-Szmek et al. 2017, Inglebert et al. 2020). For articles which uses a interval, such as P14-21, the whole interval was counted. This histogram is based on the ones described in Tebaykin et al. 2018.

Preparation of the neural tissue

Electrophysiological recordings are usually done in brain slices and cultured neurons. Brain slicing preparations have variations, and the effects of these are not completely understood. A basic protocol to brain slicing goes as follows [START_REF] Papouin | Obtaining acute brain slices[END_REF]: anaesthetized animals are sacrificed for brain extraction and further isolation of the region of interest while submerged in aCSF; the region of interest is isolated under a specific aCSF slicing solution designed to cause minimal damage to neurons while having their connections severed, then slices are left to recover under a specific aCSF solution until further use in a recording chamber with recirculated aCSF (designed for recording). Each step of this process has several variables that can modify the plasticity outcomes. However, the effects of preparation methods are not widely discussed in the papers using them [START_REF] Papouin | Obtaining acute brain slices[END_REF]. For instance, high Mg 2+ during slicing prevents epileptiform activity in slices, and it is a stabilizing factor for further recording. However, a variation of 0.3 mM in Mg 2+ concentration can change whether or not LTP is induced [START_REF] Capron | The characteristics of LTP induced in hippocampal slices are dependent on slice-recovery conditions[END_REF]. Moreover, given the slicing temperature, dendrite spines can significantly shrink when slicing is done in ice-cold temperatures (Q. [START_REF] Zhou | Shrinkage of dendritic spines associated with long-term depression of hippocampal synapses[END_REF]). An alternative to slices is cultured neurons in which the brain is extracted during embryonic development and processed to differentiate and form connections in vitro [START_REF] Gordon | General overview of neuronal cell culture[END_REF].

Given the diversity of electrophysiology recording done in non-physiological conditions, slicing preparation and the associated recipes usually do not recreate the ideal neuronal environment. Tebaykin et al. 2018 highlighted that electrophysiological choices between laboratories are variable, and aCSF recipes are passed from supervisor to student in a conservative way. Given that synapses are very sensitive to their environment, the lack of uniformity on how close the experiments are to reproduce the extracellular environment does little to produce replicate and comparable data. Despite the lack of uniformity in data, the cost to re-investigate the hallmarks of plasticity neuroscience under uniform experimental conditions is not negligible, and new hypothesis from theoretical plasticity models may justify such effort.

Temperature

Previously, it was mentioned that synapses are sensitive to experimental conditions. An important factor affecting biochemical reactions is the temperature since it defines the average kinetic energy of molecules and ions, therefore receptor activation and plasticity timescales [START_REF] Roelandse | Hypothermia-associated loss of dendritic spines[END_REF]. Krelstein et al. 1990 showed that changing temperature 3-5 • C can induce qualitative changes in the plasticity induced in hibernator hamsters. Also, Wittenberg and Wang 2006 compare different temperature effects on STDP protocols for hippocampal slices in rats without conclusive evidence. In terms of temperature-induced changes, ionic channels such as VGCC can have a kinetic constant 50x different when the temperature passes from physiological to non-physiological range [START_REF] Shibasaki | Effects of body temperature on neural activity in the hippocampus: regulation of resting membrane potentials by transient receptor potential vanilloid 4[END_REF]. Enzymes critical to plasticity induction also have their kinetics modified. For instance, the CaMKII saturation point decreases for higher temperatures under the same stimuli (Chang et al. 2017). Therefore, it is unsafe to presume that an accurate understanding of how plasticity works can be achieved using non-physiological electrophysiology. Still, room temperature recordings are still used since physiological temperatures may reduce neural tissue viability [START_REF] Buskila | Extending the viability of acute brain slices[END_REF].

Extracellular Ca 2+ and Mg 2+

Ca 2+ , as previously mentioned, dictates the plasticity outcomes. Extracellular Ca 2+ affects Ca 2+ dynamics and the Ca 2+ -influx through the modification of Ca 2+reversal potential according to the diffusion laws formalized by the GHK (Goldman-Hodgkin-Katz) equation. Receptors also are susceptible to such variations; for instance, NMDAr conductance reduces in response to high extracellular Ca 2+ (Maki and Popescu 2014). Despite Ca 2+ importance, plasticity studies use a variety of Ca 2+ concentrations in their recipes, which are often non-physiological (Inglebert et al. 2020). Different concentrations can modify the release probability of presynaptic vesicles containing neurotransmitters which is sufficiently heterogeneous even under controlled conditions, within the same brain region (King et al. 2001;Dobrunz et al. 1997). Note that during a plasticity experiment, aCSF is recirculated with a fixed concentration. However, in an intracellular environment, the external Ca 2+ is dynamically subjected to Ca 2+ fluctuations caused by astrocytes, which regulate various extracellular concentrations [START_REF] Bazargani | Astrocyte calcium signaling: the third wave[END_REF] and it is modulated by the astrocyte network [START_REF] Wang | Astrocytes modulate neural network activity by Ca2+-dependent uptake of extracellular K+[END_REF].

Mg 2+ is also a largely available divalent ion (4th major ion in the brain Romani 2011) which interacts with Ca 2+ sources controlling plasticity. It affects NMDAr conductance (Jahr and Stevens 1990) and has its internal concentration tightly regulated [START_REF] Romani | Cellular magnesium homeostasis[END_REF]. During brain slice preparations, Mg 2+ is added to control spontaneous and epileptiform activity. Despite its importance, plasticity protocols have variable choices of Mg 2+ concentrations (see Figure 2.5). Its effects can cause experiments in free-Mg 2+ to induce plasticity with only a single pulse (Mizuno et al. 2001). Tebaykin et al. 2018 have mapped various ions from electrophysiology experiments (not directly to plasticity) and noted that Mg 2+ has the greatest variability (1, 1.3 and 2 mM) among ions used in the aCSF recipes, meanwhile Ca 2+ is well centered around 2 mM despite physiological external Ca 2+ in rats is usually inferior to 1.8 mM (Inglebert et al. 2020). The model developed in this thesis predicts the effects of Mg 2+ variations for different plasticity protocols.

Development

The plasticity induction requirements are not constant during development since structural and biochemical alterations occur within brain cells. For instance, during puberty, there is a high density of dendritic spines and an up-regulated pruning which decreases towards adulthood [START_REF] Boivin | Adolescent pruning and stabilization of dendritic spines on cortical layer 5 pyramidal neurons do not depend on gonadal hormones[END_REF]. Such shift process depends on inhibitory GABAr (>P9-21 in rats [START_REF] Rivera | The K+/Cl-co-transporter KCC2 renders GABA hyperpolarizing during neuronal maturation[END_REF][START_REF] Afroz | Synaptic pruning in the female hippocampus is triggered at puberty by extrasynaptic GABAA receptors on dendritic spines[END_REF]) tonic current, which the deletion of which can cause cognitive impairment. GABAr also has its own shift passing from excitatory to inhibitory (Rinetti-Vargas et al. 2017) and can affect STDP induction (Meredith et al. 2003). Other synaptic receptors and ion channels are also subjected to composition modifications and expression shifts. During development, the NMDAr exchange the preferential expression of GluN2B subunits to the GluN2A, producing NMDAr with a smaller coincidence detection window (Iacobucci and Popescu 2018;[START_REF] Iacobucci | NMDA receptors: linking physiological output to biophysical operation[END_REF]. Furthermore, the second component of the coincidence detection, the BaP, also has a dynamical shift caused by maturation and expression of ion channels (Buchanan and Mellor 2007). Plasticity components in the developmental and adult brain are not the same, and in correlation, the plasticity induction can differ in an age-dependent way (Dudek and Bear 1993;Cao and Harris 2012;Meredith et al. 2003). Despite the various changes that occur in the animal models' first days of life, it is usual to observe studies using a rat cohort with a wide age gap (see Table 3.3 and 3.1). The model presented in this thesis analyses the effect of age in different plasticity protocols.

Next, a mini-review will discuss how different plasticity models have contributed to propose testable predictions.

Chapter 3

Computational models of synaptic plasticity

Computational science fields such as AI or computer vision generally have benchmark datasets to rank models' performance. The same practice is uncommon in computational neuroscience [START_REF] Jolivet | The quantitative single-neuron modeling competition[END_REF] since the specificity of experiments imposes a less flexible approach when the mechanistic understanding is the goal. The plasticity experimentation hallmarks such as STDP from Bi and Poo 1998 or the BCM-rule from [START_REF] Kirkwood | Experience-dependent modification of synaptic plasticity in visual cortex[END_REF] cover only a fraction of the staggering heterogeneity seen in the plasticity experimentation literature. Therefore, they cannot be used as the sole quantitative or qualitative ground truth to judge if a plasticity model performs well. Since the plasticity phenomenon is still being explored, an ideal model can support this exploration by reproducing a comprehensive set of experiments to provide testable predictions and sharper questions [START_REF] Karlin | The eleventh ra fisher memorial lecture-kin selection and altruism[END_REF].

Unified theories for synaptic -modeling approaches

Due to the diversity of plasticity experiments, conceptual advances are necessary to address a large body of experiments. Unified theories (Kievit et al. 2011) have thus been proposed to tackle this issue. These unified theories, exemplified by models, combine different concepts to describe a larger domain of evidence related to synaptic plasticity. In this way, they open new questions. Although the theory status could sound premature compared to long-standing theories in physics (e.g. as stated in Rohrlich 1989 pg. 27), the term is justified by how a model achieves a broader understanding by combining previously established principles, concepts or evidence. For instance, in [START_REF] Costa | Unified pre-and postsynaptic long-term plasticity enables reliable and flexible learning[END_REF] the descrip-tion of both pre-and postsynaptic long-term plasticity are unified to give a new vision on how long-term plasticity is induced. Another example, Ebner et al. 2019 models the evidence that plasticity is deferentially induced given the neuron's morphologic characteristics together with the concept of the enzymatic control of pre-and postsynaptic plasticity. Other examples of unified models approaches are given by Mäki-Marttunen et al. 2020, Shouval et al. 2002, De Pittà and Brunel 2016[START_REF] Singh | Presynaptic endoplasmic reticulum regulates short-term plasticity in hippocampal synapses[END_REF]. Integrating mechanisms into a single model can increase its complexity due to the resulting number of interactions but allows one to incorporate the diverse mechanistic contributions to synaptic plasticity. However, a unified theory is not necessarily built over complexity or the detailed description of all known interactions. More abstract models (Graupner and Brunel 2012;Clopath and Gerstner 2010) can also unify different concepts. Note that a complex model does not imply its quality when reproducing a phenomenon, as shown by Brette 2015, neither, reproducing experiments should be the single goal of modelling. Thus, we will identify how models unifying different concepts contribute to pose testable hypotheses.

Complexity of mechanisms on plasticity models

To understand how mechanistic models differ and propose testable predictions, we need to know how models are assembled and what are their limits. Testable predictions can be quoted as what an experimentalist needs in order to test the model predictions. For instance, an electrophysiology experiment has many details and steps, with few exceptions, overlooked by current models. However, they can affect learning rules as shown in Figure 3.1. Therefore, it is important to define the scope in which the hypothesis drawn from the model are testable. Table 3.4 categorizes the level mechanistic details to specify model capabilities as shown in Table 3.5. Such categorization is done not to undermine other models efforts but to identify the gap in which different model could contribute. To this, four attributes were considered: experimental conditions (temperature, development, aCSF), firing structures (FDP, STDP, irregular firing), localization (presynaptic, morphology, glia, sampling methods), and fast/slow variables (membrane potential, neurotransmitters, Ca 2+ sources, molecular dynamics). With these four attributes, five models are compared as shown in Table 3.5. Next, the four attributes and their relevance for testable predictions are discussed.

Experimental conditions

Modelling plasticity requires qualitative or quantitative data and the data type (e.g. neuroimage, EEG, electrophysiology). Because sampling methods normally vary between labs (Tebaykin et al. 2018) and methodological differences can modify the experimental outcomes (as in Figure 3.1), knowing how the data used in a model is sampled can help to determine the model's scope. Electrophysiology measurements are the most common type of data used in plasticity modelling. Despite that, experimental conditions can vary as shown in the Tables 3.3, 3.1 and3.2 which summarize ex vivo experiments which are reproduced by the five models in the Table 3.5 (except Tigaret et al. 2016;Mizuno et al. 2001;Dudek and Bear 1993). Sampling aspects known to affect plasticity outcomes are: sex [START_REF] Dachtler | Do cortical plasticity mechanisms differ between males and females?[END_REF], age (Dudek and Bear 1993;[START_REF] Kirkwood | Co-regulation of long-term potentiation and experience-dependent synaptic plasticity in visual cortex by age and experience[END_REF], temperature [START_REF] Alonso | Temperature compensation in a small rhythmic circuit[END_REF]Krelstein et al. 1990), aCSF recipe (Inglebert et al. 2020), slicing procedures and brain regions. However, this list of plasticity modifiers is overlooked by most models which are unlikely to propose clear instructions linking experimental conditions and their simulated variables. This limitation of current models will be discussed in more detail later when we describe how models implement experimental conditions. Figure 3.2 depicts an example of implementation of experimental conditions that is classified as complex in Table 3.4. Inglebert et al. 2020). On the Right, plasticity outcomes from TBS affected by temperature (see Krelstein et al. 1990).

Temperature

Most temperatures used in plasticity experiments do not closely follow physiological conditions due to neural tissue viability [START_REF] Buskila | Extending the life span of acute neuronal tissue for imaging and electrophysiological studies[END_REF]) and other exper-imental constraints. Temperature fluctuations can occur during experimentation, induced by the environment, or sampling techniques such as laser or electrodeinduced heating, both poorly controlled. By observing the categorized plasticity models in Table 3.5, most of them overlook the temperature differences in the experiments they reproduce (Tables 3.3 and3.1). For instance, the three STDP experiments reproduced by Jedrzejewska-Szmek et al. 2017's model were recorded at different temperatures. However, the implemented NDMAr and AMPAr are not adapted to the temperature differences, although they are known to be affected by it (Korinek et al. 2010;Postlethwaite et al. 2007). By stating which temperature the model is designed for can help to define its validity. That could improve the bottleneck on the data availability since synaptic or neuronal mechanisms are not experimentally characterized in a sufficiently broad range to support modelling. Unlike plasticity models, few other computational neuroscience models have temperature adaptations, even when using heterogeneous data sources [START_REF] Bartol | Computational reconstitution of spine calcium transients from individual proteins[END_REF][START_REF] Alonso | Temperature compensation in a small rhythmic circuit[END_REF]. For instance, [START_REF] Bartol | Computational reconstitution of spine calcium transients from individual proteins[END_REF] reconstitutes Ca 2+ dynamics in physiological temperatures, adjusting their model to 34 • C with the reference parameters for the model sampled at different temperatures. Also, Alonso and Marder 2020 characterizes network dynamics using data from crabs (Cancer Borealis) in different temperatures. They show that temperature can modulate network rhythms [START_REF] Tang | Robustness of a rhythmic circuit to short-and long-term temperature changes[END_REF], this modulation is expected to be acute in hot-blooded animals, for which temperature range is narrow. Temperature oscillations occur during circadian cycles or pathological circumstances (e.g. fever, hyper/hypothermia). Thus, the implementation of temperature is a relevant factor impacting how we can model memory, learning and normal functioning of the brain.

Developmental age of specimens

Animals undergo constant change during their lifetimes. The genetic variability analysis of the developing synaptosome (Cizeron et al. 2020) pointed out that the synaptic genetic diversity reaches its peak during adulthood, whereas old and young individuals have it less differentiated. An example of developmental impacts over plasticity is shown by Dudek and Bear 1993. They use local field stimulation on the CA3-CA1 pathway to induce LTD at different developmental ages, showing that LTD is more robust for younger than adult rats. Despite the LTP and LTD requirements being not the same during development (Pinar et al. 2017;Dudek and Bear 1993), most models for plasticity overlook age despite recognizing its relevance Graupner and Brunel 2012. Usually, computational mod-els depict age as a class rather than a continuum. For instance, in McKiernan and Marrone 2017 review, age is represented as age classes with different firing properties. The range of animal ages used in plasticity experiments is also variable, usually within few days of difference Nevian and Sakmann 2006 (13-15) or a month as in [START_REF] Weber | Location-dependent synaptic plasticity rules by dendritic spine cooperativity[END_REF]. Despite Ca 2+ sources undergo fast expression changes during development, such as the NMDAr developmental shift [START_REF] Sinclair | Effects of sex and DTNBP1 (dysbindin) null gene mutation on the developmental GluN2B-GluN2A switch in the mouse cortex and hippocampus[END_REF]), these aspects are not explicitly described in plasticity model (Shouval et al. 2002;[START_REF] Philpot | Visual experience and deprivation bidirectionally modify the composition and function of NMDA receptors in visual cortex[END_REF]. Therefore, it is relevant for models to capture these differences since they induce qualitative alteration in electrophysiology and Ca 2+ dynamics within few days (Meredith et al. 2003;Dudek and Bear 1993;Buchanan and Mellor 2007;[START_REF] Kirkwood | Co-regulation of long-term potentiation and experience-dependent synaptic plasticity in visual cortex by age and experience[END_REF]. Table 3.5 shows that developmental aspects are the most overlooked among the compared models for plasticity.

aCSF composition

In ex vivo experiments, an artificial cerebrospinal fluid (aCSF) is used to support the neural tissue viability. The aCSF composition influences neuronal electrophysiology (Tebaykin et al. 2018) by changing the spontaneous neuronal activity [START_REF] Barreto | Ion concentration dynamics as a mechanism for neuronal bursting[END_REF][START_REF] Barreto | Ion concentration dynamics as a mechanism for neuronal bursting[END_REF] or modifying receptor kinetics (Maki and Popescu 2014). An experimental-theoretical work (Inglebert et al. 2020) showed how the external Ca 2+ concentration, [Ca 2+ ] o , causes qualitative changes on the STDP curves obtained from hippocampal rat slices. They modelled the collected data using a nonlinear version of Graupner and Brunel 2012 Ca 2+ -based rule since the linear version was not sufficient (see Table 3.5). They managed to reproduce STDP with different delays, frequencies, and burst modes. However, despite it is well known that extracellular ion concentrations (e.g. [Ca 2+ ] o , [Mg 2+ ] o ) modulates various molecular mechanisms implicated in synaptic plasticity, such as NMDAr conductance (Maki and Popescu 2014) and VGCCs. Despite that the implementation of a biophysically detailed external ion influence has not been done for plasticity models (Higgins et al. 2014;Inglebert et al. 2020). Note that some biophysical models in computational neuroscience use GHK formalism to cover the diversity of extracellular ionic concentrations (Huang et al. 2015)

Firing patterns

Firing patterns are the plasticity experimentation aspect that received the most attention from modellers despite not being the only cause for plasticity induction.

Next, let's describe which protocols are mostly targeted by plasticity models.

Spike-timing-dependent plasticty

STDP protocols control the Ca 2+ influx in dendritic spines through pre and postsynaptic stimulation. In the literature, STDP firing pattern variations include bursts, various frequencies and pulse numbers, inducing different plasticity outcomes shown by Wittenberg and Wang 2006;Tigaret et al. 2016. These variations are also affected by experimental conditions (Inglebert et al. 2020), and by dendritic location as implemented by Ebner et al. 2019. Despite neurocomputational studies have given hold important the canonical STDP curve estimated by Bi and Poo 1998, it was not found to be reproducible in low Ca 2+ settings (Froemke and Dan 2002;Inglebert et al. 2020). Yet, it became a benchmark for some models (De Pittà and Brunel 2016) and neuro-inspired technologies (Zhuravleva et al. 1997). Table 3.5 shows how different plasticity models cover a wide range of STDP experiments, while other models focus on the canonical STDP with little physiological relevance.

Frequency dependent platicity (FDP) Bliss and Lømo 1973 were the first to formalize a method to induce LTP. Further, Lynch et al. 1977 did the same with LTD, leading to the experimental discovery of the frequency-dependent plasticity. The BCM-like and the STDP curves are milestones for plasticity models (Graupner and Brunel 2012;[START_REF] Izhikevich | Relating stdp to bcm[END_REF]. Despite the various FDP experiments, spiking times of EPSP-induced APs are not recorded (Mayr and Partzsch 2010). For instance, in Dudek and Bear 1992 that only mentions the presynaptic frequency and does not record postsynaptic spikes. Therefore, depending on the experimental conditions, superthreshold summated EPSPs can produce AP, affecting the Ca2+ influx. From the modelling perspective, one should be aware that an FDP experiments can have EPSP-induced-AP, thus it is important to have a mechanism to generate APs. Excluding EPSP-induced-AP can hinder the participation of NMDAr coincidence detection (Jahr and Stevens 1990). The models in Table 3.5 do not include short term depression, which is known to limit the high-frequency stimulation due to the fast depletion of neurotransmitter resources (Fernández-Alfonso and Ryan 2004). Other plasticity models (Deperrois and Graupner 2020;[START_REF] Costa | Unified pre-and postsynaptic long-term plasticity enables reliable and flexible learning[END_REF]) implement these presynaptic adaptations, which are relevant to grasp the unreliable nature of neuronal transmission.

Another experimental example is the sliding threshold [START_REF] Kirkwood | Experience-dependent modification of synaptic plasticity in visual cortex[END_REF]) implemented by Shouval et al. 2002 through the NMDAr decay times, however, without explicitly modelling what grants the NMDAr slower or faster kinetics such as the NMDAr shift [START_REF] Philpot | Visual experience and deprivation bidirectionally modify the composition and function of NMDA receptors in visual cortex[END_REF]. Also, other models have suggested the hyperpolarizing-activated currents as the factors controlling the sliding threshold (Narayanan and Johnston 2010), however, without a direct data comparison.

Irregular Firing

The most used experimental data in plasticity models are the regular firing structures with correlated spiking times. However, neurons in vivo are subjected to background noise and irregular spiking times incompatible with the usual experimental approach. As previously noted, in vivo plasticity cannot be explained only by delay times as in STDP (Froemke and Dan 2002), thus suggesting rules that use a single frequency or a single delay time characterize plasticity is of limited interest. Unified theories can be applied to natural firing patterns as shown in previous plasticity models (Graupner et al. 2016;Cui et al. 2018;Bittner et al. 2017). Table 3.4 collects studies of irregular firing in plasticity protocols which are found in various degrees in the models listed in Table 3.5.

Brain regions and structures (localization)

Synaptic plasticity induction can be different between brain region since neurons exhibit various morphologies (e.g. Purkinje and CA1 pyramidal cell), input patterns dependent on the incoming neurocircuitry, genetic expression profiles (Cizeron et al. 2020) and are supported for different densities of non-neuronal cells (e.g. non-uniform distribution astrocytes, [START_REF] Keller | Cell densities in the mouse brain: a systematic review[END_REF]. A formulation of the non-neurocentric plasticity has been approached by unified theories which aim to generalize concepts (De Pittà and Brunel 2016). Depending on which sampling method is used, plasticity outcomes from ex vivo and in vivo can differ due to the background activity and others aspects (Ahnaou et al. 2020). Table 3.4 describes key aspects regulating plasticity in different localization contexts which are used to classify models in Table 3.5.

Presynaptic plasticity pathways

The presynaptic release is pointed out as the largest source of variability in the postsynaptic calcium influx due to its all-or-none behaviour (Yuste et al. 1999). Unified theories have included presynaptic plasticity pathways such as the endocannabinoid system Heifets and Castillo 2009 to predict pre and postsynaptic plastic changes [START_REF] Costa | Unified pre-and postsynaptic long-term plasticity enables reliable and flexible learning[END_REF]Ebner et al. 2019;Deperrois and Graupner 2020;[START_REF] Costa | Synaptic transmission optimization predicts expression loci of long-term plasticity[END_REF]Cui et al. 2018). Also, different forms of vesicle release are implemented by the models as shown in Table 3.5, most usually, a discrete release with or without short-term plasticity (STP) is used, or an averaged STP form in which the effects are illustrated in Figure 6.4A. However, vesicle pool dynamics are generally ignored by models (Jedrzejewska-Szmek et al. 2017). Instead, models (De [START_REF] Pittà | Modulation of synaptic plasticity by glutamatergic gliotransmission: A modeling study[END_REF]Ebner et al. 2019) have used an averaged-form of vesicle release or no adaptation, i.e. does not mimic the depletion of vesicles during sustained inputs. An example of discrete and averaged vesicle release forms is given in this thesis's additional results in Figure 6.4.

Morphology

Different implementations of spine [START_REF] Bartol | Computational reconstitution of spine calcium transients from individual proteins[END_REF] and neuron morphologies (Herz et al. 2006) are found in models to represent time and position of variables relevant to plasticity (Chindemi et al. 2020;Ebner et al. 2019;Jedrzejewska-Szmek et al. 2017). The relevance of spatial characterization intersects with electrophysiology and diffusion properties of enzymes (Yasuda 2017) which has been approached by more detailed models (Blackwell et al. 2019). Furthermore, dendritic spikes contribution to plasticity prediction has been implemented by Ebner et al. 2019. Another example, Chindemi et al. 2020 includes multiple simplified Ca 2+ -dependent rules based on Graupner and Brunel 2012 to estimate plasticity in cortical neurons for various synapses. Table 3.4 shows different levels of neuron morphology implementations and how plasticity models adopt them (Table 3.5).

Non-neuronal cells contribution

The tripartite synapse (pre, post and astrocyte, [START_REF] Tewari | A mathematical model of the tripartite synapse: astrocyte-induced synaptic plasticity[END_REF]Nadkarni and Jung 2007) or quadripartite (including microglia [START_REF] Schafer | The "quad-partite" synapse: Microglia-synapse interactions in the developing and mature CNS[END_REF] are non-neurocentric frameworks to understand plasticity through the synaptic micro and macrodomains. These non-neuron cells regulate the extracellular environment, excitability, neuronal death and development [START_REF] Singer | Noninvasive 40-Hz light flicker to recruit microglia and reduce amyloid beta load[END_REF]. For instance, astrocytes prevent glutamate spillover to nearby spines through glutamate transporters (GLT-1) which recently has been receiving attention in models González et al. 2020. Despite De Pittà and Berry 2019 initiatives to model the participation glia cells on plasticity there is still a gap to be filled as suggested by Table 3.5.

Sampling methods

The rules to induce plasticity differ between brain regions due to their specialization, specific morphology, expression, and connectivity rules. Therefore, models need to be adaptable to capture aspect in different regions, such as done by Graupner and Brunel 2012. Theoretically, a model can evaluate if the experiment done for one region is valid for another one as shown by Ebner et al. 2019 which made a cortex model and applied it to reproduce cortical and hippocampal data. The benefit of having a less abstract parameter space is that a model can be repurposed using available biophysical measures rather than decode how a given abstraction maps an experiment specificity. The meaning of abstraction here is how a model represents physical units and biological entities. Table 3.5 shows how generic plasticity models are when including different brain regions and sampling methods.

Fast/slow variables

Within the dendritic spines, ion channels can open and close in milliseconds, and enzymes can stay activated for minutes. Such duality can bridge short and long time scales Antunes et al. 2016a. For instance, the hippocampus's sliding threshold phenomena was observed using first a slow priming stimulation with a rest period followed by the frequency-dependent protocol [START_REF] Wang | Priming-induced shift in synaptic plasticity in the rat hippocampus[END_REF]. Another experiment requiring slow variables is the Bittner et al. 2017 behavioural spiking time plasticity. They replayed a burst-like-STDP recorded in vivo with long delay times (seconds), finding a bell shape curve that requires a slow variable to bridge long timescales. To include these, models have used slow variables (e.g. based on enzymes) to integrate fast timescales, such as done by Graupner and Brunel 2007 inspired by their CaMKII-CaN model (Graupner and Brunel 2007). A similar strategy is done by Clopath and Gerstner 2010, which integrates fast membrane potential timescales using a low pass filter. Table 3.5 shows how models integrate these different aspects.

Membrane potential

Membrane potential is usually a fast variable that is more accessible to measure than Ca 2+ . Although some models have not included electrophysiological aspects directly (Graupner and Brunel 2012;Inglebert et al. 2020), its characterization was shown to be relevant to understand how NMDAr coincidence contributes to plasticity. BaPs, which are the postsynaptic responses dendritic spines interact with, attenuates under constant use, and it is subjected to the differential expression of ion channels (Ebner et al. 2019). The BaP attenuation controls the Ca 2+ entry in neurons through NMDAr coincidence. The membrane potential is closely related to neuron morphology since it defines the location of excitatory, inhibitory and modulatory inputs (Gerstner et al. 2018), and the dendritic tree integration properties.

Neurotransmitters & receptors

The plasticity phenomena have necessary conditions other than entry of Ca 2+ influx, for instance mGluR1 (Tigaret et al. 2016) and CaMKII activation (Chang et al. 2017). Beyond the classically implemented receptor and ion channels, neuromodulators modify excitability and requirements to plasticity induction [START_REF] He | Modelling the dynamics of CaMKII-NMDAR complex related to memory formation in synapses: The possible roles of threonine 286 autophosphorylation of CaMKII in long term potentiation[END_REF]. Models have included neuromodulators contribution (Blackwell et al. 2019). Alternatively, they have modelled regions that are markedly related to neuromodulation, such as dopaminergic synapses in the striatum (Jedrzejewska-Szmek et al. 2017). Also, this has been done in a more phenomenological fashion Gerstner et al. 2018. As suggested by Magee and Grienberger 2020, neuromodulators participate in vivo plasticity (Isaac et al. 2009;Bittner et al. 2017). Therefore it would be valuable if they were more represented in models.

Ca 2+ sources

Models have introduced both detailed ( [START_REF] Bartol | Computational reconstitution of spine calcium transients from individual proteins[END_REF]Blackwell et al. 2019) and simplified (Shouval et al. 2002;Graupner and Brunel 2012) Ca 2+ dynamics. Usually these simulations where validated with dye fluorescence data which reflects the buffered Ca 2+ . There is a diversity of plasticity models focusing on the Ca 2+ sources, such as endoplasmatic reticulum [START_REF] Singh | Presynaptic endoplasmic reticulum regulates short-term plasticity in hippocampal synapses[END_REF], Ca 2+ diffusion (Jedrzejewska-Szmek et al. 2017), and NMDAr decay time (Shouval et al. 2002).

Due to the diversity of mechanisms, a Ca 2+ linear model was shown to be insufficient to reproduce STDP rules sampled at different [Ca 2+ ] o (Inglebert et al. 2020).

Molecular dynamics

Molecular interactions are linked to different forms of plasticity The data output from molecular simulations tend to be larger than those focusing on only voltage or calcium, and it has complicated interpretation due to the possible unknown interactions which could escape analysis. An example of multiple threshold systems is given by Ebner et al. 2019. Alternatively, as in Blackwell et al. 2019, an average over key molecules during a fixed interval was used to predict how alcohol disturbs plasticity. Simulation with various molecules [START_REF] Kotaleski | Modelling the molecular mechanisms of synaptic plasticity using systems biology approaches[END_REF][START_REF] Bhalla | Synaptic input sequence discrimination on behavioral timescales mediated by reaction-diffusion chemistry in dendrites[END_REF], with redundant pathways, could indicate the molecular compensations in the disruption pathways implicated in plasticity [START_REF] Mizusaki | Neural circuit function redundancy in brain disorders[END_REF]. Therefore, methodological advances are necessary to analyse these systems.

Determinism, stochasticity and averages

A competition in which neuron model had to reproduce the voltage traces as a result of a current injection raised questions on which kind of modelling (biophysical or abstract) could better predict the spiking events of cortical neurons [START_REF] Jolivet | The quantitative single-neuron modeling competition[END_REF]. The collective behaviour of hundreds of stochastic ion channels shaping a neuron's membrane potential has been estimated by deterministic modelling when parameters are optimized to this. However, when scaling down towards smaller compartments, such as dendritic spines, an averaged approach does not represent the dynamics of a few stochastic units (Ribrault et al. 2011). The number of receptors and ion channels falls sharply in dendritic spines (Nimchinsky et al. 2004), such as ∼3 VGCCs to ∼15 NMDAr (Ribrault et al. 2011) demanding a stochastic modelling approach to encompass the inherent stochasticity of such small compartments (Blackwell et al. 2019). An alternative is to include an external global noise, however without a defined source (Graupner and Brunel 2012). Another stochastic event relevant to understand plasticity is the all-ornone release of neurotransmitter vesicles. Despite these events being discrete and depending on the availability of vesicles, plasticity models usually characterize vesicle release using averages or unitary input without short-term plasticity. That is, when using an averaged input a vesicle release failure may never occur, it is plausible that strategies without short-term plasticity introduce more glutamate than expected (see Figure 6.4A). Stochastic contributions of individual receptors and ion-channels are essential to characterize the highly miniaturized dendritic spines due to its small number of components. Furthermore, action potential backpropagation experienced by dendritic spines also can be stochastic, with random invasion of distant dendrites (Short et al. 2017). The richness of stochastic process in the neurons and how they contribute to learning still need to be understood better since multiple descriptions of the stochastic phenomena does not attribute a function to what is measured (Ribrault et al. 2011;Rusakov et al. 2020).

Thresholds strategies

Mechanisms with short and long timescales affect how plasticity is expressed. Thus implementing all known players [START_REF] Meunier | Modulation of synaptic plasticity in the cortex needs to understand all the players[END_REF]) would be out of reach for a human to do without mistakes. Instead of having a cascade of events linking stimulation to the inclusion of AMPAr in the PSD, one could map a variable into the desired outcome without necessarily simulating the intermediary inter-actions, also known as the phenomenological model of cascade induction. For instance, Graupner and Brunel 2012 thresholds Ca 2+ , showing that it is possible to generate various learning rules seen in the literature; to this, they use a Heaviside function detecting when Ca 2+ is over the LTD and/or LTP threshold. When Ca 2+ crosses the threshold, it modifies the synaptic efficacy with a given LTD and LTP associated rate. Furthermore, Graupner and Brunel 2012 model was modified to unify other plasticity mechanisms (Inglebert et al. 2020;Graupner et al. 2016;De Pittà and Brunel 2016;Chindemi et al. 2020). However, since Graupner and Brunel 2012 model extension is also abstract, they could not be merged in a modular way and are restricted by the specificity of the experiment they reproduce. Thresholds are usually applied to fast time scales (e.g. Ca 2+ in Graupner and Brunel 2012 and voltage in Clopath and Gerstner 2010), translated to integrators, filters or indicator functions. A general limitation of the threshold approach is the definition of an upper limit, meaning that a variable contributes to induce plasticity regardless of how high it is above the threshold, disregarding abnormally influxes. This lack of upper limit for excessively high Ca 2+ influx ignores that apoptotic pathways can be triggered for example. Another thresholds issue is that variables contribute to plasticity from the first spike in a spike train, which does not assume the possibility that some spikes are used to communicate only instead of inducing plasticity. For instance, the dendritic spine structure is unchanged for a single EPSP however, after eight EPSPs spine structure reorganization occurs (Chang et al. 2017).

A threshold strategy following original Lisman 1989's Ca 2+ amplitude hypothesis, also discussed in Shouval et al. 2002, uses two thresholds over Ca 2+ to describe synaptic weight change, using the expression below from Graupner and Brunel 2012:

τ ρ = -ρ(1 -ρ)(ρ ⋆ -ρ) + γ p (1 -ρ)Θ(c(t) -θ p ) -γ d ρΘ[c(t) -θ d ] + N oise(t).
Here, c(t) is the Ca 2+ concentration, ρ is the synaptic weight, γ p and γ d are respectively the rate of change for LTP and LTD. This expression represents plasticity induction with a Heaviside function Θ using a LTP and LTD threshold over calcium, respectively, θ p and θ d . Note that the noise is added externally, representing a global measure of noise, such strategy is refined by De Pittà and Brunel 2016 which specifies more precisely the noise sources before adding it externally. Other strategies, instead, apply a similar two-threshold system in more than one variable using a low-pass filter as in Ebner et al. 2019. In this thesis, a new threshold system will proposed that instead of being applied individually to each Fino, Elodie, Vincent Paille, Yihui Cui, Teresa Morera-Herreras, Jean-Michel Deniau, and Laurent Venance (2010). "Distinct coincidence detectors govern the corticostriatal spike timing-dependent plasticity". In: The Journal of physiology 588. 16, pp. 3045-3062. Froemke, Robert C and Yang Dan (2002). "Spike-timing-dependent synaptic modification induced by natural spike trains". In: Nature 416.6879, pp. 433-438.

Froemke, Robert C, Ishan A Tsay, Mohamad Raad, John D Long, and Yang Dan (2006). "Contribution of individual spikes in burst-induced long-term synaptic modification". In: Journal of neurophysiology. Gerstner, Wulfram, Marco Lehmann, Vasiliki Liakoni, Dane Corneil, and Johanni Brea (2018). "Eligibility traces and plasticity on behavioral time scales: experimental support of neohebbian three-factor learning rules". In: Frontiers in neural circuits 12, p. 53.

González, Janneth, Andrés Pinzón, Andrea Angarita-Rodrıéguez, Andrés Felipe Aristizabal, George E Barreto, and Cynthia Martıén-Jiménez (2020). "Advances in Astrocyte Computational Models:

From Metabolic Reconstructions to Multi-omic Approaches". In: Frontiers in neuroinformatics 14, p. 35. Graupner, Michael and Nicolas Brunel (2012). "Calcium-based plasticity model explains sensitivity of synaptic changes to spike pattern, rate, and dendritic location". In: Proceedings of the National Academy of Sciences 109. 10, pp. 3991-3996. -(2007). "STDP in a bistable synapse model based on CaMKII and associated signaling pathways". In: PLoS computational biology 3.11, e221. Graupner, Michael, Pascal Wallisch, and Srdjan Ostojic (2016). "Natural firing patterns imply low sensitivity of synaptic plasticity to spike timing compared with firing rate". In: Journal of Neuroscience 36.44, pp. 11238-11258. Hay, Etay, Sean Hill, Felix Schürmann, Henry Markram, and Idan Segev (2011). "Models of neocortical layer 5b pyramidal cells capturing a wide range of dendritic and perisomatic active properties". In: PLoS Comput Biol 7.7, e1002107.

He, Kaiwen, Marco Huertas, Su Z Hong, XiaoXiu Tie, Johannes W Hell, Harel Shouval, and Alfredo Kirkwood (2015). "Distinct eligibility traces for LTP and LTD in cortical synapses". In: Neuron 88.3, pp. 528-538. Heifets, Boris D and Pablo E Castillo (2009). "Endocannabinoid signaling and long-term synaptic plasticity". In: Annual review of physiology 71, pp. 283-306. Herz, Andreas VM, Tim Gollisch, Christian K Machens, and Dieter Jaeger (2006). "Modeling singleneuron dynamics and computations: a balance of detail and abstraction". In: science 314.5796, pp. 80-85. Higgins, David, Michael Graupner, and Nicolas Brunel (2014). "Memory maintenance in synapses with calcium-based plasticity in the presence of background activity". In: PLoS Comput Biol 10.10, e1003834. Huang, Shiwei, Sungho Hong, and Erik De Schutter (2015). "Non-linear leak currents affect mammalian neuron physiology". In: Frontiers in cellular neuroscience 9, p. 432. Inglebert, Yanis, Johnatan Aljadeff, Nicolas Brunel, and Dominique Debanne (2020). "Synaptic plasticity rules with physiological calcium levels". In: Proceedings of the National Academy of Sciences. issn: 0027-8424. doi: 10.1073/pnas.2013663117.

Chapter 4

The gap in the computational models of plasticity Computational models of plasticity have the potential to provide testable hypothesis to understand learning and memory better. As in the previous chapter, unifying concepts enlarge the domain of questions one can ask. If an ideal model can provide sharper questions, a complex model would grant broader coverage. However, drawbacks from overlaying complexity may not ensure a precise representation of the real phenomenon [START_REF] Brette | What is the most realistic single-compartment model of spike initiation?[END_REF]. For instance, models using data from heterogeneous conditions [START_REF] Bartol | Computational reconstitution of spine calcium transients from individual proteins[END_REF]) can fall on the Frankenstein effect discussed by Brette 2016. That is, there is a risk of false modularity when assembling mechanisms obtained from different sources such as repositories of neuronal mechanisms. Also, there is an analytical challenge when simulating many variables (Ebner et al. 2019;Blackwell et al. 2019). Therefore, a balanced modelling approach encompassed with data validation is a more suitable goal.

As previously discussed in chapter three, most models overlook experimental conditions of their datasets, limiting how one can compare experiment and simulation. More contributions are also required for in vivo-like firing, non-neural cells plasticity modulation, how calcium sources and molecular dynamics are modelled. One can attribute such arrested development to a data shortage. However, similarly to the precursor plasticity models, they could predict experimental outcomes for years to come, even with scarce evidence. Therefore, in this thesis, with the available data, a biophysical model was developed to cover experimental conditions and how they interact with the key calcium sources and the immediate calcium-binding enzymes. The model phenomenologically includes contributions from the more established advances of synaptic plasticity (FDP, STDP, membrane potential, morphology, vesicle release) to balance complexity. As previously discussed, there is an analytical bottleneck when interpreting molecular 57 simulations; next, a new plasticity prediction strategy is defined that expands the interpretability of established unidimensional thresholds and averages.

Limits of unidimensional thresholds

A challenge in modelling is to decide whether adding a mechanism is relevant for a given task and how abstract or complex it should be. Owning to represent synaptic mechanisms, models can assume that relevant factors are hidden (not explicitly described), allowing them to predict plasticity without necessarily simulating all cascade of events linked to plasticity outcomes. In the literature, strategies to bridge cascade of events without describing intermediary steps are done by thresholds (Graupner and Brunel 2012;Shouval et al. 2002), filters (Clopath and Gerstner 2010;Gerstner et al. 2018) or averages over the simulated variables (Blackwell and Jedrzejewska-Szmek 2013) Bear et al. 1987) suggests the existence of a threshold mechanism over calcium and voltage interactions to induce plasticity modifications. Although it is long known that the interaction over two or more variables should be taken into account in a combined way, current models still use unidimensional thresholds, as mentioned in chapter 3. Therefore, new techniques to deal with combined activity are required. In this thesis, we developed a polygonal threshold, named here geometrical readout, which can be used to capture the dynamical patterns of two or more variables simultaneously (Figure 4.2C), instead of individually (Figure 4.2A-B). Instead of relying on a unidimensional threshold (Graupner and Brunel 2012;Clopath and Gerstner 2010), one could define a geometrical shape or volume to readout dynamical pattern inducing plasticity. Combining variables into orbits can potentially improve the separability of plasticity protocol since when the number of dimensions increases, the distance between points also does for a suitable metric (Friedman et al. 2001 pages 22-23). Conversely, this is an issue for methods for data classification (Friedman et al. 2001) since more data is required to fill up the sparsity associated with the increase of dimensions, also known as the "curse of dimensionality". The advantage of using it in the dynamical systems is that such extra space gained with more dimensions can separate better plasticity protocols which can be too similar in lower dimensions. This is illustrated by Guyon and Elisseeff 2003 which shows that even adding a second dimension with noise can increase separability in scatter plots. We depict this increase of separability with a toy model example in Figure 4.1.

Limits are important for biological systems. As mentioned in chapter 2, exces- sive calcium activation triggers apoptotic pathways rather than plasticity (Zhivotovsky and Orrenius 2011; Orrenius et al. 2003). Unidimensional thresholds usually do not consider upper limits, Figure 4.2D-I shows how the geometrical and unidimensional thresholds are saturated by two slow variables with increasingly higher increment size (var1 and var2). Using these two aspects, orbits and a geometrical readout, we developed a model which can predict a wide range of plasticity experiments. The following section proposes a model to cover the exposed gaps. The main question it tries to answer is twofold, how a geometrical threshold can decode the enzymatic activity to predict plasticity and how it can be used to understand the heterogeneity of plasticity outcomes. (Pousinha et al. 2019), decreasing cell excitability and neuron integration properties which are believed to influence cognitive decay during the early stage of Alzheimer's disease. The research goal related to the pathology modelling could not start as a satisfactory synapse model was not readily available to cover the heterogeneity of experimental conditions. The thesis project, therefore, focused on creating an excitatory synapse model. In parallel, Romain Veltz and Cian O'Donnell (University of Bristol) had initiated a collaboration to build this type of model in healthy conditions providing the backbone on which the thesis project started.

The project was funded through the ComputaBrain project granted by Université Côte d'Azur, which intends to model and unify different brain scales to under-stand healthy and pathological aspects of memory processing. This PhD thesis was funded to cover the synaptic scale of this project. The first challenge was to explain the Tigaret et al. 2016's STDP experiments claiming that the canonical STDP rules were not reproducible in brain slices. Tigaret et al. 2016 also used a broad range of stimulation patterns, including burst-STDP, different delays, presynaptic stimulation with burst and frequencies. This article had an apparent paradox caused by the Ca 2+ dye, with similar fluorescence curves for plasticity protocols yielding different outcomes (discussed in Figure 6.5).

Initially, we used Fujii et al. 2013's discovery attributing to CaMKII and CaN somatic amplitude and integral a filter property the ability to decode the firing patterns (number of pulses and frequency). Our first readout predicted plasticity using the simulated enzymatic integral and amplitude (Fujii et al. 2013) M11). It was important to understand how the adaptation of CaMKII occur in Chang et al. 2017since Tigaret et al. 2016 (supp. files) also had a strong adaptation of Ca 2+ -sources (visualized through Ca 2+ dye) which caused the fluorescence to vanish in 20 s after the stimulation (no Ca 2+ influx). Due to this, the hypothesis lead to the investigation and inclusion of sources of adaptation in the model, such as BaP attenuation (Golding et al. 2001), AMPA desensitisation (Robert and Howe 2003), the SK channels (Griffith et al. 2016) and short-term depression and facilitation. From all possibilities, the one which best matches Chang et al. 2017 Dudek and Bear 1992;Dudek and Bear 1993. Later, the model included age-related BaP and GABA changes since, as mentioned by Buchanan andMellor 2007 andMeredith et al. 2003, these are relevant to understand age-related changes in STDP with BaP being more efficient for an adult than younger rats, with implications for the NMDAr coincidence detection.

The aCSF conditions were slightly different between Dudek andBear 1992 andTigaret et al. 2016. However, these differences were even greater considering Meredith et al. 2003 andWittenberg andWang 2006, done at room temperature (see Table 3.3). Wittenberg andWang 2006 andKrelstein et al. 1990 suggested how plasticity rules can be affected by temperature. Due to this, it was decided to extend the temperature control initially over NMDAr using Korinek et al. 2010's NMDAr-dependent EPSP measurements. However, the temperature adjustments from other receptors and ion channels were done later. During NeuroFrance 2019 conference, we discovered through Yanis Inglebert's poster how external Ca 2+ and Mg 2+ modulated STDP. That was important to understand since experiments used to validate the model had different aCSF recipes. After the inclusion of a better Ca 2+ control through GHK formalism, we found that Ca 2+ also affected the vesicle release probability [START_REF] Hardingham | Extracellular calcium regulates postsynaptic efficacy through group 1 metabotropic glutamate receptors[END_REF]King et al. 2001) and the

Article preamble

Plasticity experiments are heterogeneous in their methods, and conditions are thought to represent physiological learning and memory. Therefore, it is vital to understand how plasticity is induced in such extreme cases (non-physiological conditions) since they can also induce long-lasting changes in synaptic efficiency. The article presented next describes a data-driven model to gain insight into the relationship between plasticity outcomes in heterogeneous experimental conditions. It shows a new synaptic rule which reproduces several plasticity experiments using the concepts of geometrical readout and orbits mentioned previously. Also, it offers a hypothesis on how irregular firing from in vivo plasticity occurs and new strategies to induce LTP and LTP based on plasticity protocols. The model also formalises strategies to induce LTP or LTD when less favourable experimental conditions are not met (e.g. LTD and adult rats). For instance, the model predicts that presynaptic burst can rescue LTD induced by LFS in adult rats to levels similar to rat pups using Dudek and Bear 1992;Dudek and Bear 1993's experimental conditions. It also predicts several plasticity maps indexed by the experimental conditions or firing structure variations (in vivo-like and regular spikes, single interspike stimulation interval). This theoretical result extends the reproducibility factor of hallmarks experiments. And, it provides testable predictions using parameter modifications such as temperature, developmental age, aCSF (Mg and Ca), and distance from the soma.

The synapse model which will presented next is formulated as a Piecewise deterministic Markov process, see Annexe 1 for details.

A stochastic model of hippocampal synaptic plasticity with geometrical readout of enzyme dynamics

Yuri E. Rodrigues 1,2,3,* , Cezar M. Tigaret 4 , H él ène Marie 1,3,6 , Cian O'Donnell 5,6 , and Romain Veltz 2,6,* Discovering the rules of synaptic plasticity is an important step for understanding brain learning. Existing plasticity models are either 1) top-down and interpretable, but not flexible enough to account for experimental data, or 2) bottom-up and biologically realistic, but too intricate to interpret and hard to fit to data. We fill the gap between these approaches by uncovering a new plasticity rule based on a geometrical readout mechanism that flexibly maps synaptic enzyme dynamics to plasticity outcomes. We apply this readout to a multi-timescale model of hippocampal synaptic plasticity induction that includes electrical dynamics, calcium, CaMKII and calcineurin, and accurate representation of intrinsic noise sources. Using a single set of model parameters, we demonstrate the robustness of this plasticity rule by reproducing nine published ex vivo experiments covering various spike-timing and frequency-dependent plasticity induction protocols, animal ages, and experimental conditions. The model also predicts that in vivo-like spike timing irregularity strongly shapes plasticity outcome. This geometrical readout modelling approach can be readily applied to other excitatory or inhibitory synapses to discover their synaptic plasticity rules.

To understand how brains learn, we need to identify the rules governing how synapses change their strength in neural circuits. What determines whether each synapse strengthens, weakens, or stays the same? The dominant principle at the basis of current models of synaptic plasticity is the Hebb postulate 1 which states that neurons with correlated electrical activity strengthen their synaptic connections, while neurons active at different times weaken their connections. In particular, spike-timingdependent plasticity (STDP) models 2-4 were formulated based on experimental observations that precise timing of pre-and postsynaptic spiking determines whether synapses are strengthened or weakened 5-8 . However, experiments also found that plasticity induction depends on the rate and number of stimuli presented to the synapse 9,10 , and the level of dendritic spine depolarisation 11-15 . The lack of satisfactory plasticity models based solely on neural spiking prompted researchers to consider more elaborate models based on synapse biochemistry 16 . Following a proposed role for postsynaptic calcium (Ca 2+ ) signalling in synaptic plasticity 17 , previous models assumed that the amplitude of postsynaptic calcium controls long-term alterations in synaptic strength, with moderate levels of calcium causing long-term depression (LTD) and high calcium causing long-term potentiation (LTP) 18,19 . Recent experimental data suggests that calcium dynamics is also important 20-24 . As a result, subsequent phenomenological models of plasticity incorporated slow variables that integrate the fast synaptic input signals, loosely modelling calcium and its downstream effectors 25-32 .

However, even these models do not account for data showing that plasticity is highly sensitive to physiological conditions such as the developmental age of the animal 33-36 , extracellular calcium and magnesium concentrations 37,38 and tissue temperature 39-41 . The fundamental issue is that the components of these phenomenological models do not directly map to biological components of synapses, so they cannot automatically model alterations due to physiological and experimental conditions. This absence limits the predictive power of existing plasticity models.

To tackle this problem, we devised a new plasticity rule based on a bottom-up, data-driven approach by building a biologicallygrounded model of plasticity induction at a single rat hippocampal CA3-CA1 synapse. We focus on this synapse type because of the abundant published experimental data that can be used to quantitatively constrain the model parameters. Compared to previous models in the literature, we aimed for an intermediate level of detail: enough biophysical components to capture the key dynamical processes underlying plasticity induction, but not the full molecular cascade underlying plasticity expression; much of which is poorly quantified 42 . Our model centred on dendritic spine electrical dynamics, calcium signalling and immediate downstream molecules, which we then mapped to synaptic strength change via an conceptually new dynamical, geometric readout mechanism. Crucially, the model also captures intrinsic noise based on the stochastic switching of synaptic receptors and ion channels 43,44 . We found that the model can account for published data from spike-timing and frequency-dependent plasticity experiments, and variations in physiological parameters influencing plasticity outcomes. We also tested how the model responded to in vivo-like spike timing jitter and spike failures, and found that the plasticity rules were highly sensitive to these subtle input alterations. (Methods). The same vesicle depletion rule is applied to both glutamate and GABA. Right: BaP efficiency as function of time. BaP efficiency controls the axial resistance between soma and dendrite in order to phenomenologically capture the distance-dependent BaP attenuation 46,47 (Methods). e, Activated enzyme concentration for CaM, CaN and CaMKII, as function of time for the stimulus 1Pre1Post10, 30 pairing repetitions delivered at 5 Hz. Note that the vertical grey bar is the duration of the stimuli, 6 s.

Results

A multi-timescale model of synaptic plasticity induction. We built a computational model of plasticity induction at a single CA3-CA1 rat glutamatergic synapse (Fig. 1). Our goal was to reproduce results on synaptic plasticity that explored the effects of several experimental parameters: fine timing differences between pre and postsynaptic spiking (Fig. 23); stimulation frequency (Fig. 4); animal age (Fig. 5); external calcium and magnesium (Fig. 6); stochasticity in the firing structure (Fig. 7), temperature and experimental conditions variations (Supplemental Information). Where possible, we set parameters to values previously estimated from experiments, and tuned the remainder within physiologically plausible ranges to reproduce our target plasticity experiments (Methods).

The model components are schematized in Fig. 1a (full details in Methods). For glutamate release, we used a twopool vesicle depletion and recycling system, which accounts for short-term presynaptic depression and facilitation. When glutamate is released from vesicles, it can bind to the postsynaptic α-amino-3-hydroxy-5-methyl-4-isoxazolepropionic acid and N-methyl-D-aspartate receptors (AMPArs and NMDArs, respectively). When dendritic spine voltage depolarises, it activates voltage-gated calcium channels (VGCCs) and removes magne-sium (Mg 2+ ) block from NMDArs. Backpropagating action potentials (BaP) can also cause spine depolarisation. As an inhibitory component, we modelled a gamma-aminobutyric acid receptor (GABAr) synapse on the dendrite shaft 48 . Calcium ions influxing through VGCCs and NMDArs can activate hyperpolarising SK potassium channels 49,50 , bind to calmodulin (CaM) or to a generic calcium buffer. Calcium-bound calmodulin activates two major signalling molecules immediately downstream of Ca/CaM enzymes 51 : Ca 2+ /calmodulin-dependent protein kinase II (CaMKII) or calcineurin (CaN) phosphatase, also known as PP2B 52 . We included these two proteins because CaMKII activation is necessary for Schaffer-collateral LTP 53,54 , while CaN activation is necessary for LTD 55,56 . Later, we show how we map the joint activity of CaMKII and CaN to LTP and LTD.

Synaptic receptors and ion channels have an inherent random behavior, stochastically switching between open and closed states 44 . If the number of receptors or channels is large, then the variability of the total population activity becomes negligible relative to the mean 57 . However individual hippocampal synapses contain only small numbers of receptors and ion channels, for example ∼10 NMDA receptors and <15 VGCCs 58-60 , making their total activation highly stochastic. Therefore, we modelled AMPAr, NMDAr, VGCCs and GABAr as stochastic processes. Presynaptic vesicle release events were also stochastic: glutamate release was an all-or-nothing event, and the amplitude of each glutamate pulse was drawn randomly, modelling heterogeneity in vesicle size 61 . The inclusion of stochastic processes to account for an intrinsic noise in synaptic activation 62 contrasts with most previous models in the literature, which either represent all variables as continuous and deterministic or add an external generic noise source 63-65 .

The synapse model showed nonlinear dynamics across multiple timescales. For illustration, we stimulated the synapse with single simultaneous glutamate and GABA vesicle releases (Fig. 1b). AMPArs and VGCCs open rapidly but close again within a few milliseconds. The dendritic GABAr closes more slowly, on a timescale of ∼10 ms. NMDArs, the major calcium source, closes on timescales of ∼50 ms and ∼250 ms for the GluN2A and GluN2B subtypes, respectively.

To show the typical responses of the spine head voltage and Ca 2+ , we stimulated the synapse with a single presynaptic pulse (EPSP) paired 10 ms later with a single BaP (1Pre1Post10) (Fig. 1c, left). For this pairing, when BaP is triggered immediately after an EPSP, it leads to a large Ca 2+ transient aligned with the BaP due to the NMDArs first being bound by glutamate then unblocked by the BaP depolarisation (Fig. 1c,right).

Single pre or postsynaptic stimulation pulses did not cause depletion of vesicle reserves or substantial activation of the enzymes. To illustrate these slower-timescale processes, we stimulated the synapse with a prolonged protocol: one presynaptic pulse followed by one postsynaptic pulse 10 ms later, repeated 30 times at 5 Hz (Fig. 1d-e). The number of vesicles in both the docked and reserve pools decreased substantially over the course of the stimulation train (Fig. 1d left), which in turn causes decreased vesicle release probability. Similarly, by the 30th pulse, the dendritic BaP amplitude had attenuated to ∼85% (∼70% BaP efficiency; Fig. 1d right) of its initial amplitude, modelling the effects of slow dendritic sodium channel inactivation 47,66 . CaM concentration rose rapidly in response to calcium transients but also decayed back to baseline on a timescale of ∼500 ms (Fig. 1e top). In contrast, the concentration of active CaMKII and CaN accumulated over a timescale of seconds, reaching a sustained peak during the stimu-3/46 5.3. ARTICLE lation train, then decayed back to baseline on a timescale of ∼10 and ∼120 s respectively, in line with experimental data 51,54,67 (Fig. 1e).

The effects of the stochastic variables can be seen in Fig. 1b-d. The synaptic receptors and ion channels open and close randomly (Fig. 1b). Even though spine voltage, calcium, and downstream molecules were modelled as continuous and deterministic, they inherited some randomness from the upstream stochastic variables. As a result, there was substantial trial-totrial variability in the voltage and calcium responses to identical pre and postsynaptic spike trains (grey traces in Fig. 1d). This variability was also passed on to the downstream enzymes CaM, CaMKII and CaN, but was filtered and therefore attenuated by the slow dynamics of CaMKII and CaN. In summary, the model contained stochastic nonlinear variables acting over five different orders of magnitude of timescale, from ∼1 ms to ∼1 min, making it sensitive to both fast and slow components of input signals.

Distinguishing between stimulation protocols using the CaMKII and CaN joint response. It has proven difficult for simple models of synaptic plasticity to capture the underlying rules and explain why some stimulation protocols induce plasticity while others do not. We tested the model's sensitivity by simulating its response to a set of protocols used by 24 in a recent ex vivo experimental study on adult (P50-55) rat hippocampus with blocked GABAr. We focused on three pairs of protocols (three rows in Fig. 2). In each case in 24 's experiments, one of the pair induced LTP or LTD, while the other subtly different protocol caused no change (NC) in synaptic strength. We asked if the model's joint CaMKII-CaN activity could distinguish between each pair of protocols.

The first pair of protocols differed in intensity. A protocol which caused no plasticity consisted of 1 presynaptic spike followed 10 ms later by one postsynaptic spike repeated at 5 Hz for one minute (1Pre1Post10, 300 at 5Hz). The other protocol induced LTP, but differed only in that it included a postsynaptic doublet instead of a single spike (1Pre2Post10, 300 at 5Hz), implying a slightly stronger BaP amplitude initially. For the plots in Fig. 2a, it was not possible to set a single concentration threshold on either CaMKII or CaN that would discriminate between the protocols.

To achieve better separability, we combined the activity of the two enzymes, plotting the joint CaMKII and CaN responses against each other on a 2D plane (Fig. 2b). In this geometric plot, the two protocol's trajectories can be seen to overlap for the initial part of the transient, but then diverge. To quantify trial to trial variability, we also calculated contour maps showing the mean fraction of time the trajectories spent in each part of the plane during the stimulation (Fig. 2c). Importantly, both the trajectories and contour maps were substantially non-overlapping between the two protocols, implying that they can be separated based on the joint CaN-CaMKII activity. We found that the 1Pre2Post10 protocol leads to a weaker response in both CaMKII and CaN, corresponding to the lower blue traces Fig. 2b. The decreased response to the doublet protocol was due to the enhanced attenuation of dendritic BaP amplitude over the course of the simulation 47 , leading to less calcium influx through NMDArs and VGCCs (data not shown).

The second pair of protocols we explored differed in sequenc-ing. We stimulated the synapse model with one causal (EPSP-BaP) protocol involving a single presynaptic spike followed 50 ms later by a doublet of postsynaptic spikes (1Pre2Post50, 300 at 5Hz), repeated at 5 Hz for one minute, which 24 found caused LTP.

The other anticausal protocol involved the same total number of pre and postsynaptic spikes, but with the pre-post order reversed (2Post1Pre50, 300 at 5Hz). Experimentally this anticausal (BaP-EPSP) protocol did not induce plasticity. Notably, the only difference was the sequencing of whether the pre or postsynaptic neuron fired first, over a short time gap of 50 ms. Despite the activations being apparently difficult to distinguish (Fig. 2d), we found that the LTP-inducing protocol caused greater CaN activation than the protocol that did not trigger plasticity. Indeed, this translated to a horizontal offset in both the trajectory and contour map (Fig. 2e-f), demonstrating that another pair of protocols can be separated in the joint CaN-CaMKII plane.

The third pair of protocols differed in both duration and intensity. In line with previous studies, 24 found that a train of doublets of presynaptic spikes separated by 50 ms repeated at a low frequency of 3 Hz for 5 minutes (2Pre50, 900 at 3Hz) induced LTD, while a slightly more intense but shorter duration protocol of presynaptic spike doublets separated by 10 ms repeated at 5 Hz for one minute (2Pre10, 300 at 5Hz) did not cause plasticity. When we simulated both protocols in the model (Fig. 2g-i), both caused similar initial responses in CaMKII and CaN. In the shorter protocol, this activation decayed to baseline within 100 s of the end of the stimulation. However the slower and longer-duration 2Pre50 3Hz 900p protocol caused an additional sustained, stochastically fluctuating, plateau of activation of both enzymes (Fig. 2g). This resulted in the LTD-inducing protocol having a downward and leftward-shifted CaN-CaMKII trajectory and contour plot, relative to the other protocol (Fig. 2h-i). These results again showed that the joint CaN-CaMKII activity may be useful to predict plasticity changes.

A geometrical readout mapping joint enzymatic activity to plasticity outcomes. We found that the simulated CaN-CaMKII trajectories from the two LTP-inducing protocols (1Pre2Post10 and 1Pre2Post50, at Fig. 2a and D respectively) spent a large fraction of time near ∼ 20 µM CaMKII and 7-10 µM CaN. In contrast, protocols that failed to trigger LTP had either lower (2Post1Pre50 and 2Pre10, Fig. 2d and 2g respectively), or higher CaMKII and CaN activation (1Pre1Post10, Fig. 2a). The LTDinducing protocol, by comparison, spent a longer period in a region of sustained but lower ∼ 12µM CaMKII and ∼ 2µM CaN and activation. The plots in Fig. 2c,f andg, show contour maps of histograms of the joint CaMKII-CaN activity, indicating where in the plane the trajectories spent most time. Fig. c andf indicate this measure can be used to predict plasticity, because the NC and LTP protocol histograms are largely non-overlapping. In Fig. 2c, the NC protocol response 'overshoots' the LTP protocol response, whereas in Fig. 2f the NC protocol response 'undershoots' the LTP protocol response. In contrast, when we compared the response histograms for LTD and NC protocols, we found a greater overlap (Fig. 2i). This suggested that in this case the histogram alone was not sufficient to separate the protocols, and that protocol duration is also important. LTD induction (2Pre50) required a more prolonged activation than NC (2Pre10).

To design a geometrical readout mechanism to map the The upper square bumps are caused by the protocol crossing the LTP region, the lower square bumps when the protocol crosses the LTD region (as in panel d). f, Synaptic weight (%) as function of time for each protocol. The weight change is defined as the number (out of 100) of states in the LTP state minus the number of states in the LTD state (panel c). The trajectories correspond to the median of the simulations in panel g. g, Synaptic weight change (%) predicted by the model compared to data (EPSC amplitudes) from 24 (100 samples for each protocol, also for panel h and i). The data (grey dots) was provided by 24 (note an 230% outlier as the red asterisk), red bands indicate data means. h, Predicted mean synaptic weight change (%) as function of delay (ms) and number of pairing repetitions (pulses) for the protocol 1Pre2Post(delay), for delay between -100 and 100. LTD is induced by 2Post1Pre50 after at least 500 pulses. The mean weight change along each dashed line is reported in the STDP curves in panel i. i, Synaptic weight change (%) as function of pre-post delay. Each plot corresponds to a different pairing repetition number (legend). The solid line shows the mean, and the ribbons are the 2nd and 4th quantiles. The red dots are the data means estimated in 24 , also shown in panel g. enzyme activity to plasticity outcomes, we first drew nonoverlapping boxes of LTP and LTD "plasticity regions" in the CaN-CaMKII plane (Fig. 3a). We positioned these regions over the parts of the phase space where the enzyme activities corresponding to the LTP-and LTD-inducing protocols were most different, as shown by trajectories in Fig. 2. We then fixed these regions for all subsequent parts of this study. When a trajectory enters in one of these plasticity regions, it activates LTD or LTP indicator variables (Methods) which encode the joint enzyme activities (trajectories in the phase plots) transitions across the LTP and LTD regions over time (Fig. 3b). These indicator variable drove transition rates in a plasticity Markov chain used to predict LTP or LTD (Fig. 3c), see Methods. The LTD transition rates were slower than the LTP transition rates, to reflect studies showing that LTD requires sustained synaptic stimulation 20-22 . The parameters for this plasticity Markov chain (Methods) were fit to the plasticity induction outcomes from different protocols (Table M1). In the beginning of the simulation, the plasticity Markov chain starts with 100 processes 48 in the state NC, with each variable representing 1% weight change, an abstract measure of synaptic strength that can be either EPSP, EPSC, or field EPSP slope depending on the experiment. Each process can transit stochastically between NC, LTP and LTD states. At the end of the protocol, the plasticity outcome is given by the difference between the number of processes in the LTP and the LTD states (Methods). This readout mechanism acts as a parsimonious model of the complex signalling cascade linking CaMKII and CaN activation to expression of synaptic plasticity 68 . It can be considered as a two-dimensional extension of previous computational studies that applied analogous 1D threshold functions to dendritic spine calcium concentration 18,19,30,31 . Our model is scalable, as it gives M1). Data are represented as normal distributions with the mean and variance extracted from 9 . Note that data from Dudek and Bear are given in field EPSP slope change. f, Prediction for the mean weight change (%) varying the stimulation frequency and pulse number (24x38x100 data points, respectively pulse x frequency x samples). The red dots show the 9 protocol parameters, the corresponding results are shown in panel e.

the possibility for the readout to be extended to dynamics of n different molecules, using n-dimensional closed regions.

In Fig. 3d, we plot the model's responses to seven different plasticity protocols used by 24 by overlaying example CaMKII-CaN trajectories for each protocol with the LTP and LTD regions. The corresponding region occupancies are plotted as function of time in E, and long-term alterations in the synaptic strength are plotted as function of time in F. The three protocols that induced LTP in the 24 experiments spent substantial time in the LTP region, and so triggered potentiation. In contrast, the 1Pre1Post10 (yellow trace) overshoots both regions, crossing them only briefly on its return to baseline, and so resulted in little weight change. The protocol that induced LTD (2Pre50, purple trace) is five times longer than other protocols, spending sufficient time inside the LTD region (Fig. 3f). In contrast, two other protocols that spent time in the same LTD region of the CaN-CaMKII plane (2Post1Pre50 and 2Pre10) were too brief to induce LTD. These protocols were also not strong enough to reach the LTP region, so resulted in no net plasticity, again consistent with 24 's experiments.

We observed run-to-run variability in the amplitude of the predicted plasticity, due to the inherent stochasticity in the model. In Fig. 3g, we plot the distribution of the predicted plasticity from each protocol (colours) alongside the data from 24 's study, finding a good correspondence.

Experimentally, LTP can be induced by few pulses while LTD usually requires longer-duration stimulation 20-22 . We incorporated this effect into the readout model by letting LTP have faster transition rates than LTD (Fig. 3c). 24 found that 300 repetitions of anticausal post-before-pre pairings did not cause LTD, in contrast to the canonical spike-timing-dependent plasticity curve 7 . We hypothesized that LTD might indeed appear with the anticausal 24 protocol (Table M1) if stimulation duration was increased. To explore this possibility in the model, we systematically varied the number of paired repetitions from 100 to 1200, and also co-varied the pre-post delay from -100 to 100 ms. Fig. 3h shows a contour plot of the predicted mean synaptic strength change across for the 1Pre2Post(delay) stimulation protocol for different numbers of pairing repetitions. A LTD window appears after ∼500 pairing repetitions for some anticausal pairings, in line with our hypothesis. The magnitude of LTP also increases with pulse number, for causal positive pairings. For either 100 or 300 pairing repetitions, only LTP or NC is induced (Fig. 3i). The model also made other plasticity predictions by varying 24 's experimental conditions (Fig. S1). In summary, our model readout reveals that the direction and magnitude of the change in synaptic strength can be predicted from the joint CaMKII-CaN activity in the LTP and LTD regions.

Frequency-dependent plasticity. The stimulation protocols used by 24 explored how subtle variations in pre and postsynaptic spike timing influenced the direction and magnitude of plasticity (see Table M1 for experimental differences). In contrast, traditional synaptic plasticity protocols exploring the role of presynaptic stimulation frequency did not measure the timing of co-occurring postsynaptic spikes 9,69,70 . These studies found that long-duration low-frequency stimulation induces LTD, whereas short-duration high-frequency stimulation induces LTP, with a cross-over point of zero change at intermediate stimulation frequencies. In addition to allowing us to explore frequencydependent plasticity (FDP), this stimulation paradigm also gives us further constraints for LTD in the model since in 24 , only one LTD case was available. For FDP, we focused on modelling the 6/46 70CHAPTER 5. A STOCHASTIC MODEL OF HIPPOCAMPAL SYNAPTIC PLASTICITY WITH GEOMETRICAL READOUT OF ENZYMES -ARTICLE experiments from 9 , who stimulated Schaffer collateral projections to pyramidal CA1 neurons with 900 pulses in frequencies ranging from 1 Hz to 50 Hz. In addition to presynaptic stimulation patterns, the experimental conditions differed from 24 in two other aspects: animal age and control of postsynaptic spiking activity (see Table M1 legend). We incorporated both age-dependence and EPSP-evoked-BaPs (Methods). Importantly, the read-out mechanism mapping joint CaMKII-CaN activity to plasticity is the same for all experiments in this work.

Fig. 4a shows the joint CaMKII-CaN activity when we stimulated the model with 900 presynaptic spikes at 1, 3, 5, 10 and 50 Hz 9 . Higher stimulation frequencies drove stronger responses in both CaN and CaMKII activities (Fig. 4a). Fig. 4b,c show the corresponding plasticity region indicator for the LTP/LTD region threshold crossings and the synaptic strength change. From this set of five protocols, only the 50 Hz stimulation drove a response strong enough to reach the LTP region of the plane (Fig. 4a a,d). Although the remaining four protocols drove responses primarily in the LTD region, only the 3 and 5 Hz stimulations resulted in substantial LTD. The 1 Hz and 10 Hz stimulations resulted in negligible LTD, but for two distinct reasons. Although the 10 Hz protocol's joint CaMKII-CaN activity passed through the LTD region of the plane (Fig. 4a,d), it was too brief to activate the slow LTD mechanism built into the readout (Methods). The 1 Hz stimulation, on the other hand, was prolonged, but its response was mostly too weak to reach the LTD region, crossing the threshold only intermittently (Fig. 4b, bottom trace). Overall the model matched well the mean plasticity response found by Dudek and Bear (Fig. 4e), following a classic BCM-like curve as function of stimulation frequency 71,72 .

We then used the model to explore the stimulation space in more detail by varying the stimulation frequency from 0.5 Hz to 50 Hz, and varying the number of presynaptic pulses from 50 to 1200. Fig. 4f shows a contour map of the mean synaptic strength change (%) in this 2D frequency-pulse number space. Under 9 's experimental conditions, we found that LTD induction required at least ∼300 pulses, at frequencies between 1Hz and 3Hz. In contrast, LTP could be induced using ∼50 pulses at ∼20Hz or greater. The contour map also showed that increasing the number of pulses (vertical axis in Fig. 4e) increases the magnitude of both LTP and LTD. This was paralleled by a widening of the LTD frequency range, whereas the LTP frequency threshold remained around ∼20Hz, independent of pulse number. The pulse dependence amplitude increase predicted in Fig. 4 is also valid for 24 experiment shown in Fig. S1f.

Ex vivo experiments in 9 were done at 35°C. However, lower temperatures are more widely used ex vivo because they extend brain slice viability. We performed further simulations testing temperature modifications for 9 's experiment, finding that it had a strong effect on plasticity outcomes (Fig. S2d-f).

Variations in plasticity induction with developmental age.

The rules for induction of LTP and LTD change during development 33,35 , so a given plasticity protocol can produce different outcomes when delivered to synapses from young animals versus mature animals. For example, when 33 tested the effects of low-frequency stimulation (1 Hz) on CA3-CA1 synapses from rats of different ages, they found that the magnitude of LTD decreases steeply with age from P7 until becoming minimal in mature animals >P35 (Fig. 5a, circles). Across the same age range, they found that a theta-burst stimulation protocol induced progressively greater LTP magnitude with developmental age (Fig. 5b, circles). Paralleling this, multiple properties of neurons change during development: the NMDAr switches its dominant subunit expression from GluN2B to GluN2A 73-75 , the reversal potential of the receptor (GABAr) switches from depolarising to hyperpolarizing 34,76,77 , and the action potential backpropagates more efficiently with age 46 . These mechanisms have been proposed to underlie the developmental changes in synaptic plasticity rules because they are key regulators of synaptic calcium signalling 34,46 . However, their sufficiency and individual contributions to the age-related plasticity changes are unclear. We incorporated these mechanisms in the model (Methods) by parameterizing each of the three components to vary with the animal's postnatal age, to test if they could account for the age-dependent plasticity data.

We found that elaborating the model with age-dependent changes in NMDAr composition, GABAr reversal potential, and BaP efficiency, while keeping the same plasticity readout parameters, was sufficient to account for the developmental changes in LTD and LTP observed by 33 (Fig. 5a,b). We then explored the model's response to protocols of various stimulation frequencies, from 0.5 to 50 Hz, across ages from P5 to P80 (Fig. 5c,e). Fig. 5c shows the synaptic strength change as function of stimulation frequency for ages P15, P25, P35 and P45. The magnitude of LTD decreases with age, while the magnitude of LTP increases with age. Fig. 5e shows a contour plot of the same result, covering the age-frequency space.

The 1Hz presynaptic stimulation protocol in 33 did not induce LTD in adult animals 9 . We found that the joint CaN-CaMKII activity trajectories for this stimulation protocol underwent an age-dependent leftward shift beyond the LTD region (Fig. 5f). This implies that LTD is not induced in mature animals by this conventional LFS protocol due to insufficient activation of enzymes. In contrast, 24 and 78 were able to induce LTD in adult rat tissue by combining LFS with presynaptic spike pairs repeated 900 times at 3 Hz. Given these empirical findings and our modelling results, we hypothesized that LTD induction in adult animals requires that the stimulation protocol: 1) causes CaMKII and CaN activity to stay more in the LTD region than the LTP region, and 2) is sufficiently long to activate the LTD readout mechanism. With experimental parameters used by 33 , this may be as short as 300 pulses when multi-spike presynaptic protocols are used since the joint CaMKII-CaN activity can reach the LTD region more quickly than for single spike protocols. We simulated two such potential protocols as predictions: doublet and quadruplet spike groups delivered 300 times at 1 Hz, with 50 ms between each pair of spikes in the group (Fig. 5d). The model predicted that both these protocols induce LTD in adults (green and blue curves), whereas as shown above, the single pulse protocol did not cause LTD (yellow curve). These findings suggest that the temporal requirements for inducing LTD may not be as prolonged as previously assumed, since they can be reduced by varying stimulation intensity. See Fig. S3 for frequency versus age maps for presynaptic bursts. 33 also performed theta-burst stimulation (TBS, Table M1) at different developmental ages, and found that LTP is not easily induced in young rats (Fig. 5b), see also 35 . The model qualitatively 7/46 5.3. ARTICLE matches this trend, and also predicts that TBS induces maximal LTP around P21, before declining further during development (Fig. 5b, green curve). Similarly, we found that high-frequency stimulation induces LTP only for ages >P15, peaks at P35, then gradually declines at older ages (Fig. 5e). Note that in Fig. 5b, we used 6 epochs instead of 4 used by 33 to increase LTP outcome which is known to washout after one hour for young rats 35 . Fig. 5. | Age-dependent plasticity, Dudek and Bear 1993 33 dataset a, Synaptic weight change for 1Pre, 900 at 1 Hz as in 33 . The solid line is the mean and the ribbons are the 2nd and 4th quantiles predicted by our model (same for panel b, c and f). b, Synaptic weight change for Theta Burst Stimulation (TBS -4Pre at 100 Hz repeated 10 times at 5Hz given in 6 epochs at 0.1Hz (Table M1) and 33 . c, Synaptic weight change as function of frequency for different ages. BCM-like curves showing that, during adulthood, the same LTD protocol becomes less efficient. It also shows that high-frequencies are inefficient at inducing LTP before P15. d, Synaptic weight change as function of age. Proposed protocol using presynaptic bursts to recover LTD at ≥ P35 with less pulses, 300 instead of the original 900 from 33 . This effect is more pronounced for young rats. Fig. S3 shows a 900 pulses comparison. e, Mean synaptic strength change (%) as function of frequency and age for 1Pre 900 pulses (32x38x100, respectively, for frequency, age and samples). The protocols in 33 (panel a) are marked with the yellow vertical line. The horizontal lines represent the experimental conditions of panel c. Note the P35 was used for 9 experiment in Fig. 4f. f, Mean time spent for the 1Pre 1Hz 900 pulses protocol showing how the trajectories are left-shifted as rat age increases.

In contrast to 33 's findings, other studies have found that LTP can be induced in hippocampus in young animals (<P15) with STDP. For example, 34 found that at room temperature, 1Pre1Post10 induces LTP in young rats, whereas 1Pre2Post10 induces NC. This relationship was inverted for adults, with 1Pre1Post inducing no plasticity and 1Pre2Post10 inducing LTP (Fig. S5).

Together, these results suggest that not only do the requirements for LTP/LTD change with age, but also that these agedependencies are different for different stimulation patterns. Finally, we explore which mechanisms are responsible for plasticity induction changes across development in the FDP protocol (Fig. S3) by fixing each parameter to young or adult values for the FDP paradigm. Our model analysis suggests that the NMDAr switch iacobucci2017 is a dominant factor affecting LTD induction, but the maturation of BaP 46 is the dominant factor affecting LTP induction, with GABAr shift having only a weak influence on LTD induction for 33 's FDP.

Plasticity requirements during development do not necessarily follow the profile in 33 as shown by 34 's STDP experiment. Our model shows that multiple developmental profiles are possible when experimental conditions vary within the same stimulation paradigm. This is illustrated in Fig. S5a-c 50,81 .

Fig. 6a-c compares our model to 38 's STDP data at different [Ca 2+ ] o and [Mg 2+ ] o . Note that 38 used 150 (100) pairing repetitions for the anti-causal (causal) both delivered at 0.3 Hz. At [Ca 2+ ] o =1.3 mM, Fig. 6a shows the STDP rule induced weak LTD for brief causal delays. At [Ca 2+ ] o = 1.8 mM, in Fig. 6b, the model predicted a full-LTD window. At [Ca 2+ ] o = 3 mM, in Fig. 6c, it predicted a bidirectional rule with a second LTD window for long pre-before-post pairings, previously theorized by 26 . Fig. 6d illustrates the time spent by the joint CaN-CaMKII activity for 1Pre1Post10 using 38 's experimental conditions. Each density plot corresponds to a different panel in Fig. 6 with the respective Ca/Mg. The response under low [Ca 2+ ] o spent most time inside the LTD region, but high [Ca 2+ ] o shifts the trajectory to the LTP region. Fig. S4a presents density plots for the anticausal post-before-pre protocols. 38 fixed the Ca/Mg ratio at 1.5, although aCSF formulations in the literature differ (see Table M1). Fig. S4d shows that varying Ca/Mg ratio and [Ca 2+ ] o for 38 's experiments restrict LTP to Fig. 6. Effects of extracellular calcium and magnesium concentrations on plasticity a, Synaptic weight (%) for a STDP rule with [Ca 2+ ] o = 1.3 mM (fixed ratio, Ca/Mg=1.5). According to the data extracted from 38 , the number of pairing repetitions for positive (negative) delays is 100 (150), both delivered at 0.3 Hz 38 . The solid line is the mean, and the ribbons are the 2nd and 4th quantiles predicted by our model (all panels use 100 samples). b, Same as A, but for [Ca 2+ ] o = 1.8 mM (Ca/Mg ratio = 1.5). c, Same as A, but for [Ca 2+ ] o = 3 mM (Ca/Mg ratio = 1.5). d, Mean time spent for causal pairing, 1Pre1Post10, at different Ca/Mg concentrations. The contour plots are associated with the panels a, b and c. e, STDP and extracellular Ca/Mg. Synaptic weight change (%) for causal (1Pre1Post10, 100 at 0.3 Hz) and anticausal (1Post1Pre10, 150 at 0.3 Hz) pairings varying extracellular Ca from 1.0 to 3 mM (Ca/Mg ratio = 1.5). The dashed lines represent the experiments in the panel a, b and c. We used 21x22x100 data points, respectively calcium x delay x samples. f, Varying frequency and extracellular Ca/Mg for the causal pairing 1Pre1Post10, 100 at 0.3 Hz. Contour plot showing the mean synaptic weight (%) for a single causal pairing protocol (1Pre1Post10, 100 samples) varying frequency from 0.1 to 10 Hz and [Ca 2+ ] o from 1.0 to 3 mM (Ca/Mg ratio = 1.5). We used 21x18x100 data points, respectively calcium x frequency x samples. Ca/Mg>1.5 and [Ca 2+ ] o >1.8 mM.

Our model can also identify the transitions between LTD and LTP depending on Ca/Mg. Fig. 6e shows a map of plasticity as function of pre-post delay and Ca/Mg concentrations and the parameters where LTP is induced for the 1Pre1Post10 protocol. Since plasticity rises steeply at around [Ca 2+ ] o = 2.2 mM, small fluctuations in [Ca 2+ ] o near this boundary could cause qualitative transitions in plasticity outcomes. For anticausal pairings, increasing [Ca 2+ ] o increases the magnitude of LTD (Fig. S4b illustrates this with 38 's data).

38 also found that increasing the pairing frequency to 5 or 10 Hz results in a transition from LTD to LTP for 1Pre1Post10 at [Ca 2+ ] o = 1.8 mM (Fig. S4c), similar frequency-STDP behaviour has been reported in cortex 10 . In Fig. 6f, we varied both the pairing frequencies and [Ca 2+ ] o and we observe similar transitions to 38 . However, the model's transition for [Ca 2+ ] o = 1.8 mM was centred around 0.5 Hz, which was untested by 38 . The model predicted no plasticity at higher frequencies, unlike the data, that shows scattered LTP and LTD (see Fig. S4c). Fig. S1d and S4e shows that 24 's burst-STDP and 38 's STDP share a similar transition structure, but not 9 's FDP.

In contrast to 38 's results, we found that setting low [Ca 2+ ] o for 24 's burst-STDP abolishes LTP, and does not induce strong LTD (Fig. S1d). For 9 's experiment, Fig. S2a [Mg 2+ ] o controls a sliding threshold between LTD and LTP but not [Ca 2+ ] o (Fig. S2b). For another direct stimulation experiment, Fig. S4c shows that in an Mg-free medium, LTP expression requires fewer pulses 21 .

Despite exploring physiological [Ca 2+ ] o and [Mg 2+ ] o 38 use a non-physiological temperature (30°C) which extends T-type VGCC closing times and modifies the CaN-CaMKII baseline (Fig. S5i). Fig. S5g,h shows comparable simulations for physiological temperatures. Overall our model predicts that temperature can change STDP rules in a similar fashion to [Ca 2+ ] o (Fig. S4a,b). In summary, plasticity is highly sensitive to variations in extracellular calcium, magnesium, and temperature (Fig. S1a, Fig. S5d-f; 40 ).

In vivo-like spike variability affects plasticity. In the above sections, we used highly regular and stereotypical stimulation protocols to replicate typical ex vivo plasticity experiments. In contrast, neural spiking in hippocampus in vivo is irregular and variable 78,82 . It is unclear how natural firing variability affects the rules of plasticity induction 27,83-85 . We explored this question using model simulations by adding three distinct types of variability: 1) spike time jitter, 2) failures induced by dropping spikes, 3) independent pre and postsynaptic Poisson spike trains 83 .

We introduced spike timing jitter by adding zero-mean Gaussian noise (s.d. σ ) to pre and postsynaptic spikes, changing spike pairs inter-stimulus interval (ISI). In Fig. 7a, we plot the LTP magnitude as function of jitter magnitude (controlled by σ ) for protocols taken from 24 . With no jitter, σ = 0, these protocols have different LTP magnitudes (corresponding to Fig. 3) and become similar once σ increases. The three protocols with a postsynaptic spike doublet gave identical plasticity for σ = 50 To understand the effects of jittering, we plotted the trajectories of CaN-CaMKII activity (Fig. 7c). 2Post1Pre50 which "undershoots" the LTP region shifted into the LTP region for jitter σ = 50 ms. In contrast, 1Pre1Post10 which "overshoots" the LTP region shifted to the opposite direction towards the LTP region.

Why does jitter cause different spike timing protocols to yield similar plasticity magnitudes? Increasing jitter causes a fraction of pairings to invert causality. Therefore, the jittered protocols became a mixture of causal and anticausal pairings (Fig. 7c). This situation occurs for all paired protocols. So any protocol with the same number spikes will become mixed if the jitter is large enough. Note that despite noise the mean frequency was conserved at 5 ± 13.5 Hz (see Fig. 7e).

Next, we studied the effect of spike removal. In the previous sections, synaptic release probability was ∼60% (for [Ca 2+ ] o = 2.5 mM) or lower, depending on the availability of docked vesicles (Methods). However, baseline presynaptic vesicle release probability is heterogeneous across CA3-CA1 synapses, ranging from ∼ 10 -90% 86,87 and likely lower on average in vivo 88,89 . BaPs are also heterogeneous with random attenuation profiles 47 and spike failures 90 . To test the effects of pre and postsynaptic failures on plasticity induction, we performed simulations where we randomly removed spikes, altering the previously regular attenuation in 24 's protocols.

In Fig. 7b we plot the plasticity magnitude as function of sparsity (percentage of removed spikes). The sparsity had different specific effects for each protocol. 1Pre2Post10 and 1Pre2Post50 which originally produced substantial LTP were robust to spike removal until ∼60% sparsity. In contrast, the plasticity magnitude from both 1Pre1Post10 and 2Post1Pre50 showed a nonmonotonic dependence on sparsity, first increasing then decreasing, with maximal LTP at ∼40% sparsity.

To understand how sparsity causes this non-monotonic effect on plasticity magnitude, we plotted the histograms of time spent in the CaN-CaMKII plane for 2Post1Pre50 for three levels of sparsity: 0%, 30% and 80% (Fig. 7d). For 0% sparsity, the activation spent most time at the border between the LTP and LTD regions, resulting in no plasticity. Increasing sparsity to 30% caused the activation to shift rightward into the LTP region because there was less attenuation of pre and postsynaptic resources. In contrast, at 80% sparsity, the activation moved into the LTD region because there were not enough events to substantially activate CaMKII and CaN. Since LTD is a slow process and the protocol duration is short (60s), there was no net plasticity. Therefore for this protocol, high and low sparsity caused no plasticity for distinct reasons, whereas intermediate sparsity enabled LTP by balancing resource depletion with enzyme activation.

Next we tested the interaction of jitter and spike removal. Fig. 7f shows a contour map of weight change as function of jitter and sparsity for the 2Post1Pre50 protocol, which originally induced no plasticity (Fig. 2). Increasing spike jitter enlarged the range of sparsity inducing LTP. In summary, these simulations (Fig. 7a,b,f andh) show that different STDP protocols have different degrees of sensitivity to noise in the firing structure, suggesting that simple plasticity rules derived from regular ex vivo experiments may not predict plasticity in vivo.

How does random spike timing affect rate-dependent plasticity? We stimulated the model with pre and postsynaptic Poisson Fig. 7. | Jitter and spike dropping effects on STDP; and Poisson spike trains. a, Mean weight (%) for the jittered STDP protocols. The solid line is the mean, and the ribbons are the 2nd and 4th quantiles predicted by our model estimated using 100 samples (same for all panels). b, Mean weight (%) for the same 24 protocols used at panel a subjected to random spike removal (sparsity %). c, Mean time spent (s) varying jittering. Contour plot shows 2Post1Pre50 and 1Pre1Post10 (300 at 5 Hz) with (grey contour plot) and with jittering (coloured contour plot). The circles and squares correspond to the marks in panel a. d, Mean time spent (s) varying sparsity. Contour plot in grey showing 0% sparsity for 2Post1Pre50 300 at 5Hz (see Fig. 2f). The contour plots show the protocol with spike removal sparsities at 0% (NC) 30% (LTP), and 80% (NC). The triangles correspond to the same marks in panel a. e, Distribution of the 50 ms jittering applied to the causal protocol 1Pre1Post10, 300 at 5 Hz in which nearly half of the pairs turned into anticausal. The mean frequency is 5 ± 13.5 Hz. The protocol 2Post1Pre50 will have nearly half of the pairings turning into causal, making them have a similar firing structure and position inside the LTP region. spike trains for 10s, under 9 's experimental conditions. We systematically varied both the pre and postsynaptic rates (Fig. 7h). The 10s stimulation protocols induced only LTP, since LTD requires a prolonged stimulation 21 . LTP magnitude monotonically increased with the presynaptic rate (Fig. 7g,h). In contrast, LTP magnitude varied non-monotonically as function of postsynaptic rate, initially increasing until a peak at 10 Hz, then decreasing with higher stimulation frequencies. This non-monotonic dependence on post-synaptic rate is inconsistent with classic rate-based models of Hebbian plasticity. We also investigated how this plasticity dependence on pre-and post-synaptic Poisson firing rates varies with developmental age (Fig. S2g-i). We found that at P5 no plasticity is induced, at P15 a LTP region appears at around 1 Hz postsynaptic rate, and at P20 plasticity becomes similar to the mature age, with a peak in LTP magnitude at 10 Hz postsynaptic rate.

Discussion

We built a model of a rat CA3-CA1 hippocampal synapse, including key electrical and biochemical components underlying synaptic plasticity induction (Fig. 1). We used a novel geometric readout of CaN-CaMKII dynamics (Fig. 234) to predict the outcomes from a range plasticity experiments with heterogeneous conditions: animal developmental age (Fig. 5), aCSF composition (Fig. 6), temperature (Supplemental files), and in vivo-like firing variability (Fig. 7).

Our model included critical components for plasticity induction at CA3-CA1 synapses: those affecting dendritic spine voltage, calcium signalling, and enzymatic activation. We were able to use model to make quantitative predictions, because its variables and parameters corresponded to biological components. This property allowed us to incorporate the model components' dependence on developmental age, external Ca/Mg levels, and temperature to replicate datasets across a range of experimental conditions. The model is relatively fast to simulate, taking ∼1 minute of CPU time to run 1 minute of biological time. These practical benefits should enable future studies to make experimental predictions on dendritic integration of multiple synaptic inputs 42,91,92 and on the effects of synaptic molecular alterations in pathological conditions. In contrast, abstract models based on spike timing 28,93,94 or simplified calcium dynamics 18,30 must rely on ad hoc adjustment of parameters with less biological interpretability.

The model was built based the new concept that the full temporal activity of CaN-CaMKII over the stimulus duration 51 , and not their instantaneous activity levels 18,19 , is responsible for plasticity changes. We instantiated this concept by analyzing the joint CaN-CaMKII activity in the two-dimensional plane and designing polygonal plasticity readout regions (Fig. 3a). In doing so, we generalised previous work with plasticity induction based on single threshold and a slow variable 26,28,30,95 . Here, we used a two-dimensional readout, but anticipate a straightforward generalisation to higher-dimensions for different cellular processes in neuroscience but also in systems biology more broadly. The central discovery is that these trajectories, despite being stochastic, can be separated in the plane as function of the stimulus (Fig. 3). This is the basis of our new synaptic plasticity rule.

Let us describe the intuition behind our model more concisely.

First, we abstracted away the sophisticated cascade of plasticity expression. Second, the plasticity regions, crossed by the trajectories, are described with a minimal set of parameters, their tuning is quite straightforward and done once and for all, even when the joint activity is stochastic. The tuning of the model is possible thanks to the decoupling of the plasticity process from the spine biophysics which acts as a feedforward input to the plasticity Markov chain and from the distributions of the different trajectories, which are well separated. It is expected that one can find other model versions (parameters or conceptual) instantiating our concept that also match the data well.

In our model, some CaMKII-CaN trajectories overshot the plasticity regions (e.g. Fig. 3d). Although abnormally high and prolonged calcium influx to cells can trigger cell death 96 , the effects of high calcium concentrations at single synapses are poorly understood. Notably, a few studies have reported evidence consistent with an overshoot, where strong synaptic calcium influx does not induce LTP 20,24,97 .

Intrinsic noise is an essential component of the model. How can the synapse reliably express plasticity but be noisy at the same time 43,44 ? Noise can be reduced either by redundancy or by averaging across time, also called ergodicity 98 . However redundancy requires manufacturing and maintaining more components, and therefore costs energy. We propose that, instead, plasticity induction is robust due to temporal averaging by slow-timescale signalling and adaptation processes. These slow variables reduce noise by averaging the faster timescale stochastic variables. This may be a reason why CaMKII uses auto-phosphorylation to sustain its activity and slow its decay time 54,99 . In summary, this suggests that the temporal averaging by slow variables, combined with the separability afforded by the multidimensional readout, allows synapses to tolerate noise while remaining energy-efficient.

We identified some limitations of the model. First, we modelled only a single postsynaptic spine attached to twocompartment neuron for soma and dendrite. Second, the model abstracted the complicated process of synaptic plasticity expression, and even if this replicated the "early" phase of LTP/LTD expression in the first 30-60 minutes after induction, slower protein-synthesis-dependent processes, maintenance processes, and synaptic pruning proceed at later timescales 100 . Third, like most biophysical models, ours contained many parameters (Methods). Although we set these to physiologically plausible values and then tuned to match the plasticity data, other combinations of parameters may fit the data equally well 101,102 due to the ubiquitous phenomenon of redundancy in biochemical and neural systems 103,104 . Indeed synapses are quite heterogeneous in receptor and ion channel counts 58-60,105 , protein abundances 106,107 , and spine morphologies 65,108 , even within the subpopulation of CA1 pyramidal neuron synapses that we modelled here. It remains to be discovered how neurons tune their synaptic properties in this vast parameter space to achieve functional plasticity rules, or implement meta-plasticity 109-111 .

Several predictions follow from our results. Since the model respected the stochasticity of vesicle release 112,113 , NMDAr 60,74,75,114 , and VGCC opening 59,115,116 , the magnitude of plasticity varied from simulation trial to trial (Methods, Fig. 3g,4e). This suggests that the rules of plasticity are inherently stochastic 63,117 and that the variability observed in experiments 9,21,24,33,34,38,40 is not just due to heterogeneity in synapse properties. By running extensive simulations over the space of protocols beyond those tested experimentally (Fig. 3h,i; 4f; 5c,e and f; 6e,f), we made testable predictions for plasticity outcomes. For example, 24 did not find LTD when using classic post-before-pre stimulation protocols, but the model predicted that LTD could be induced if the number of pairing repetitions was extended (Fig. h,i). The model also predicts that the lack of LTD induced by FDP in adults can be recovered using doublets or quadruplet spike protocols (Fig. 5d). We tested the model's sensitivity to spike time jitter and spike failure in the stimulation protocols (Fig. 7). Our simulations predicted that this firing variability can alter the rules of plasticity, in the sense that it is possible to add noise to cause LTP for protocols that did not otherwise induce plasticity.

What do these results imply about the rules of plasticity in vivo? First, we noticed that successful LTP or LTD induction required a balance between two types of slow variables: those that attenuate, such as presynaptic vesicle pools and dendritic BaP, versus those that accumulate, such as slow enzymatic integration 62,118,119 . This balance is reflected in the inverted-U shaped magnitude of LTP seen as a function of post-synaptic firing rate (Fig. 7h). Second, although spike timing on millisecond timescales can in certain circumstances affect the direction and magnitude of plasticity (Fig. 3), in order to drive sufficient activity of synaptic enzymes, these patterns would need to be repeated for several seconds. However, if these repetitions are subject to jitter or failures, as observed in hippocampal spike trains in vivo 82,120 , then the millisecond-timescale information will be destroyed as it gets averaged out across repetitions by the slow integration processes of CaMKII and CaN (Fig. 7a-d). The net implication is that millisecond-timescale structure of individual spike pairs is unlikely to play an important role in determining hippocampal synaptic plasticity in vivo 83,88,121 .

In summary, we presented a new type of biophysical model for plasticity induction at the rat CA3-CA1 glutamatergic synapse. Although the model itself is specific to this synapse type, the study's insights may generalise to other synapse types, enabling a deeper understanding of the rules of synaptic plasticity and brain learning.

Online Methods

Data and code availability

All simulations were performed in the Julia programming language (version 1.4.2). Simulating the synapse model is equivalent to sampling a piecewise deterministic Markov process, and this relies on the Julia package PiecewiseDeterministicMarkovProcesses.jl. These simulations are event-based, and no approximation is made beyond the ones required to integrate the ordinary differential equations by the method LSODA (Livermore Solver for Ordinary Differential Equations). We run the parallel simulations in the Nef cluster operated by Inria. EXPERIMENT

0

Table M1. Table with the parameters extracted from the respective publications. To fit the data associated to publications displaying a parameter interval (e.g. 70 or 100) we used a value within the provided limits. Otherwise, we depict in red the value used to fit to the data. For complete data structure on these publications and the ones used for method validation see github code. We allowed the AP to be evoked by EPSPs for these protocols: 5-7 . Note that 1 used GABA(A)r blocker modelled by turning GABAr conductance to zero. Also, 7 LTD protocol used partial NMDA blocker modelled by reducing NMDA conductance by 97 %.

Modelling procedures Notations

We write 1 A the indicator of a set A meaning that 1 A (x) = 1 if x belongs to A and zero otherwise.

Vesicle release and recycling

Vesicle-filled neurotransmitters from the presynaptic terminals stimulate the postsynaptic side when successfully released. We derived a vesicle release Markov chain model based on a deterministic approach described in 8 on page 183. We denote by (t 1 , • • • ,t n ) the time arrivals of the presynaptic spikes.

Vesicles can be in two states, either belonging to the docked pool (with cardinal D) with fast emptying, or to the reserve pool (with cardinal R) which replenishes D 9 . The docked pool loses one vesicle each time a release 10 occurs with transition D -→ D -1 (Fig. M1C). The reserve pool replenishes the docked one with transition (R, D) → (R -1, D + 1). Finally, the reserve pool is replenished with rate

(R 0 -R)/τ re f D with the transition (R, D) -→ (R + 1, D). TRANSITION RATE INITIAL CONDITION (R, D) → (R -1, D + 1) (D 0 -D) • R/τ D D(0) = D 0 (R, D) → (R + 1, D -1) (R 0 -R) • D/τ R R(0) = R 0 (R, D) -→ (R + 1, D) (R 0 -R)/τ re f R Table M2
. Stochastic transitions used in the pool dynamics. Note that the rates depend on the pool's cardinal 11 .

In addition to the stochastic dynamics in Table M2, each spike t i triggers a vesicle release D -→ D -1 with probability p rel :

p rel (Ca pre , [Ca 2+ ] o , D) = (Ca pre ) s (Ca pre ) s + h([Ca 2+ ] o ) s 1 D>0 , h([Ca 2+ ] o ) = 0.654 + 1.349 1 + e 4•([Ca 2+ ] o -1.708 mM)
which is a function of presynaptic calcium Ca pre and extracellular calcium concentration

[Ca 2+ ] o through the threshold h([Ca 2+ ] o ).
To decide whether a vesicle is released for a presynaptic spike t i , we use a phenomenological model of Ca pre (see Fig. M1A) based on a resource-use function 12 :

Ċa pre = - Ca pre τ pre Ca pre (0) = 0 Ċa jump = 1-Ca jump τ rec -δ decay •Ca jump •Ca pre Ca jump (0) = 1. (1)
Upon arrival of the presynaptic spikes, t ∈ (t 1 , • • • ,t n ), we update Ca pre according to the deterministic jump:

Ca pre -→ Ca pre +Ca jump .

Finally, after Ca pre has been updated, a vesicle is released with probability p rel (Fig. M1B).

Parameters for the vesicle release model is given in Table M3. The experimental constraints to devise a release probability model are given by 13 and 1 . Because [Ca 2+ ] o modifies the release probability dynamics 14 , we fixed an initial release probability of 68 % for [Ca 2+ ] o = 2.5 mM as reported by 1 (initial value in Fig. M1B andD). Additionally, 13 report a 38% reduction in the initial release probability when changing [Ca 2+ ] o from 2.5 mM to 1 mM. Taking these into account, the decreasing sigmoid function in the Fig. M1E depicts our [Ca 2+ ] o -dependent release probability model (p rel ).

Fig. M1E shows that our p rel function is in good agreement with a previous analytical model 14 suggesting that p rel ([

Ca 2+ ] o ) ∝ ([Ca 2+ ] o ) 2 mM -2 .
Our model also qualitatively reproduces the vanishing of calcium dye fluorescence levels after 20 s of theta trains from 1 (in their Supplementary Materials). We interpret their fluorescence measurements as an effect of short-term depression (see Fig. M1B).

Despite our model agreeing with previous works, it is a simplified presynaptic model that does not encompass the vesicle release's highly heterogeneous nature. Vesicle release dynamics are known to be sensitivity to various experimental conditions such as temperature 15 , the age for some brain regions 10 or magnesium concentration 13 . Furthermore, since our model of vesicle dynamics is simple, τ rec in Equation ( 1) has two roles: to delay the p rel recovery caused by Ca pre inactivation (enforced by δ Ca in Equation ( 1)) and to prevent vesicle release after HFS induced depression 9,14 . Later, we incorporate a higher number of experimental parameters Table M3. | Parameter values used in the presynaptic model. Our model does not implement a larger pool called "resting pool" containing ∼ 180 vesicles (CA3-CA1 hippocampus) 16 . Terminology note: In other works, the larger pool with ∼180 vesicles can be found with different nomenclatures such as "reserve pool" 19 or "resting pool" 16 . Furthemore, the nomenclature used in our model for the reserve pool, can also be found as "recycling pool" in 9,16 . ] o = 0.5 mM which can also be difficult to experimentally measure given the rarity of events 13 . Table M4. Parameters for the neuron electrical properties. * The membrane leak conductance in the spine is small since the spine resistance is so high that is considered infinite (> 10 6 MΩ) 21 , therefore the current mostly leaks through the neck. Additionally, the dendrite leak conductance is equally small in order to control the distance-dependent attenuation by the axial resistance term g adapt BaP in Equations 3 and 4.

NAME

20/46

Membrane potential and currents

Our model is built over three compartments, a spherical dendritic spine linked by the neck to a cylindrical dendrite connected to a spherical soma. The membrane potential of these compartments satisfy the equations below (parameters in Table M4). The different currents are described in the following sections.

C sp • Vsp = g neck • (V dend -V sp ) + g sp L • (E rev -V sp ) + I T + I L + I R + I NMDA + I AMPA + I SK (2)
C dend • Vdend = g adapt BaP • (V soma -V dend ) + g neck • (V sp -V dend ) + g dend L • (E rev -V dend ) + I GABA (3)
C soma • Vsoma = g adapt BaP • (V dend -V soma ) + g soma L • (E rev -V soma ) + β age • (I BaP + I Na ) + I K (4)
Action-potential backpropagation (BaP)

Postsynaptic currents

The postsynaptic currents are generated in the soma, backpropagated to the dendritic spine and filtered by a passive dendrite. The soma generates BaPs using a version of the Na+ and K+ channel models developed by 33 . The related parameters are described in Table M5 (the voltage unit is mV).

Sodium channel

Potassium channel

α m (V soma ) = 0.4 • V soma + 30 1 -e -Vsoma+30 7.2 α n (V soma ) = e -0.11•(V soma -13) β m (V soma ) = 0.124 • V soma + 30 e Vsoma+30 7.2 -1 β n (V soma ) = e -0.08•(V soma -13) m inf (V soma ) = α m (V soma ) α m (V soma ) + β m (V soma ) n inf (V soma ) = 1 1 + α n (V soma ) m τ (V soma ) = 1 α m (V soma ) + β m (V soma ) n τ (V soma ) = max 50 • β n (V soma ) 1 + α n (V soma ) ; 2 α h (V soma ) = 0.01 • V soma + 45 e Vsoma+45 1.5 -1 ṅ(V soma ) = n inf -n n τ β h (V soma ) = 0.03 • V soma + 45 1 -e -Vsoma+45
1.5

I K = γ K • n • (Erev K -V soma ) ḣ(V soma ) = α h (V soma ) • (1 -h) -β h (V soma ) • h ṁ(V soma ) = m inf -m m τ I Na = γ Na • m 3 • h • (Erev Na -V soma ).
To trigger a BaP, an external current I BaP is injected in the soma at times t ∈ {t 1 , ...,t n } (postsynaptic input times) for a chosen duration δ in j with amplitude I amp (nA):

I BaP = n ∑ i=1 H(t i ) • (1 -H(t i + δ in j )) • I amp .
The current injected in the soma is filtered in a distance-dependent manner by the dendrite before it reaches the dendritic spine. The distant-dependent BaP amplitude attenuation changes the axial resistance g adapt BaP (see equations 3 and 4) between the dendrite and the soma as follows (Fig. M2C top):

g adapt BaP = β • g di f f • β soma , β soma (d soma ) = 0.1 + 1.4 1 + e 0.02•(d soma -230.3µm) (5)
where d soma is the distance of the spine to the soma and where the factor β is dynamically regulated based on a resource-use equation 12 with a dampening factor β aux changing the size of the attenuation step δ decay : Golding measurements have neurons that are strongly attenuated and weakly attenuated (dichotomy mark represented by the dashed line). However, in this work, we consider only strongly attenuated neurons. b, Attenuation of somatic action potential from 37 and model in response to five postsynaptic spikes delivered at 100 Hz. The value showed for the model is the spine voltage with distance from the soma set to zero (scale 25 ms, 20 mV).c, Top panel shows the β soma used in Equation ( 5) to modify the axial conductance between the soma and dendrite. Bottom panel shows the age-dependent changes in the step of the resource-use equation, in Equation ( 6) that accelerates the BaP attenuation and decreases the sodium currents in the Equation ( 4). d, Probability of evoking an AP multiplied by the successfully evoked AP (p AP (V evoked ) • 1(evoked) for the protocol 1Pre, 300 at 5 Hz (2.5 mM Ca). d, Two-pool dynamics with the stimulation than panel D showing the vesicle release, the reserve and docked pools, and the evoked AP. e, Probability of evoking an AP for the protocol 1Pre 300 pulses at different frequencies (3 and 5 Hz have the same probability).

β = 1 -β τ rec -δ decay • β -1 aux • β • I BaP (t) βaux = 1 -β aux τ aux rec -δ aux • β aux • I BaP (t).
The BaP attenuation model is based on 24 data for strongly attenuating neurons. Therefore, the second type of attenuation (weakly attenuating) in neurons is not considered (dichotomy in Fig. M2A). Fig. M2A compares Golding data to our model and illustrates the effect of BaP attenuation in the upper panels of Fig. M2A andB.

Table M5 shows the BaP attenuation parameters. The plasticity outcomes as function of the dendritic spine distance from the soma are shown in Fig. S2C and Fig. S1E.

Age-dependent BaP adaptation

Age-dependent BaP attenuation modifies the neuronal bursting properties through the maturation and expression of potassium and sodium channels 34 , therefore changing the interaction of polarizing and depolarizing currents (see Fig. M2B) 35,36 . We reproduce 37 somatic attenuation profiles (Fig. M2B) with our model by including an age-dependent BaP amplitude attenuation factor. We define the attenuation factor β age (Fig. M2C

In Equation ( 4), the age effects are introduced by multiplying the attenuation factor β age by the sodium I Na and the external I BaP currents.

AP evoked by EPSP

A presynaptic stimulation triggers a BaP if sufficient depolarization is caused by the EPSPs reaching the soma 27 chapter 13. We included an option to choose whether an EPSP can evoke an AP using an event generator resembling the previous release probability model (p rel ). Like the p rel , the BaPs evoked by EPSPs are estimated before the postsynaptic simulation. To this, we use a variable V evoke which is incremented by 1 at each presynaptic time t ∈ (t 1 , ...,t n ) and has exponential decay:

Vevoke = -V evoke τ v V evoke (0) = 0 V evoke -→ V evoke + 1. ( 7
)
22/46 M5. The Na+ and K+ conductances intentionally do not match the reference because models with passive dendrite need higher current input to initiate action potentials 42 . Therefore we set it to achieve the desired amplitude on the dendrite and the dendritic spine according to the predictions of 24 and 43 .

Since the BaPs evoked by EPSP are triggered by the afferent synapses and are limited by their respective docked pools (D), we use the previous p rel to define the probability of an AP to occur. We test the ratio of successful releases from 25 synapses to decide if a BaP is evoked by an EPSP, setting a test threshold of 80%. Therefore, we express the probability of evoking an AP, p AP (V evoke ), with the following test:

∑ 25 1(rand < p rel (V evoked , [Ca 2+ ] o , D)) 25 > 80 %.
The EPSP summation dynamics on the soma and dendrites depend on the complex neuron morphology 38,39 which was not implemented by our model. Therefore, our "AP evoked by EPSP test" intends to give a simplified way to produce BaPs similar to an integrate-and-fire model 8 .

Previous work suggests that BaPs can be evoked with a ∼5 % probability for low-frequencies 40 in the Dudek and Bear 1992 experiment ([Ca 2+ ] o = 2.5 mM). Our model covers this estimation, but the chance to elicit an AP increases with the frequency 38 . This is captured by the V evoke (in a integrate-and-fire fashion 27 ) as shown in Fig. M2F. The Fig. s M2D and E show how a 5 Hz stimulation evokes APs. The delay between the EPSP and the evoked AP is set to δ delay-AP = 15ms, similar to the EPSP-spike latency reported for CA1 neurons 41 .

AMPAr

Markov chain

The AMPAr is modeled with the Markov chain (Fig. M3) described by 44,45 and adapted to temperature changes according to 46 . Here, we introduce the additional parameters ρ AMPA f , ρ AMPA b to cover AMPAr temperature-sensitive kinetics 46 . The corresponding parameters are given in Table M6. The AMPAr current is the sum of the subcurrents associated to the occupancy of the three subconductance states O2, O3 and O4 of the Markov chain in the Fig. M3 and described as follows:

C0 C1 C2 C3 C4 D0 D1 D2 D3 D4 O2 O3 O4 D 2 2 D 2 3 D 2 4 4k 1 • [Glu] • ρ AMPA f k -1 • ρ AMPA b 3k 1 • [Glu] • ρ AMPA f 2k -1 • ρ AMPA b 2k 1 • [Glu] • ρ AMPA f 3k -1 • ρ AMPA b k 1 • [Glu] • ρ AMPA f 4k -1 • ρ AMPA b 3k 1 • [Glu] • ρ AMPA f k -2 • ρ AMPA b 3k 1 • [Glu] • ρ AMPA f k -1 • ρ AMPA b 2k 1 • [Glu] • ρ AMPA f 2 • k -1 • ρ AMPA b k 1 • [Glu] • ρ AMPA f 3k -1 • ρ AMPA b 2k 1 • [Glu] • ρ AMPA f k -1 • ρ AMPA b k 1 • [Glu] • ρ AMPA f 2k -1 • ρ AMPA b 4δ 0 γ 0 δ 1 γ 1 2δ 1 γ 1 3δ 1 γ 1 4δ 1 γ 1 δ 2 γ 2 2δ 2 γ 2 3δ 2 γ 2 α 2β α 3β α 4β
I AMPA = (Erev AMPA -V sp ) • (γ A2 • O2 + γ A3 • O3 + γ A4 • O4).
The adaptation of the Markov chain from 44 The effects of temperature change on AMPAr dynamics are presented in Fig. M4, which also shows that the desensitisation is unaltered between temperature changes (Fig. M4B andC). The recovery time from desensitisation is the same as at room temperature 44 . Desensitisation measurements are required to account for a temperature-dependent change in the rates of the "vertical" transitions in Fig. M3, see 46 . This can be relevant for presynaptic bursts. ) and the decay time at different temperatures in response to 1 mM glutamate during 1 ms (standard pulse). 46 data (our model) suggests that AMPAr decay time at 35°C is ∼ 0.5 ms (∼ 0.6 ms) and at 25°C is ∼ 0.65 ms (∼ 0.95 ms), this shows a closer match towards more physiological temperatures. b, Desensitisation profile of AMPAr at 35°C showing how many AMPAr are open in response to a glutamate saturating pulse (5 mM Glu during 20 ms) separated by an interval (x-axis). c, Same as in panel B but for 25°C.

24/46

= 1.6 • 10 7 M -1 s -1 44 glu unbinding 1 k -1 = 7400 s -1 44 glu unbinding 2 k -2 = 0.41 s -1 44 closing α = 2600 s -1 44 opening β = 9600 s -1 44 desensitisation 1 δ 1 = 1500 s -1 44 desensitisation 2 δ 2 = 170 s -1 44 desensitisation 3 δ 0 = 0.003 s -1 44 re-desensitisation 1 γ 1 = 9.1 s -1 44 re-desensitisation 2 γ 2 = 42 s -1 44 re-desensitisation 3 γ 0 = 0.83 s -1 44
Table M6. Parameter values for the AMPAr Markov chain and glutamate release affecting NMDAr, AMPAr. Properties of GABA release are the same as those for glutamate.

Postsynaptic Ca 2+ influx

The effects of experimental conditions in the calcium dynamics are due to receptors, ion channels and enzymes. A leaky term models the calcium resting concentration in the Equation (8). The calcium fluxes from NMDAr and VGCCs (T, R, L types) are given in Equation ( 9). The diffusion term through the spine neck is expressed in Equation (10). Finally, the buffer, the optional dye and the enzymatic reactions are given in Equation ( 11) (parameter values given at the Table M7):

Ċa = Ca ∞ -Ca τ Ca + (8)
Ca NMDA + I T + I R + I L 2 • F • A sp + (9) max(Ca ∞ ,Ca/3) -Ca τ CaDi f f - (10)
Ḃu

f f Ca -Ḋ ye + enzymes. (11)
Despite the driving force to the resting concentration, Ca ∞ = 50 nM, the tonic opening of T-type channels causes calcium to fluctuate making its mean value dependent on temperature, extracellular calcium and voltage. The effects of this tonic opening in various experimental conditions are shown in Fig. S4C. To avoid modelling dendritic calcium sources, we use a dampening term as one-third of the calcium level since calcium imaging comparing dendrite and spine fluorescence have shown this trend 49 . The Equation (10) implements the diffusion of calcium from the spine to the dendrite through the neck. The time constant for the diffusion coefficient τ CaDi f f , is estimated as described in 50 . The calcium buffer and the optional dye are described as a two states reaction system 51 :

Ḃu f f Ca = k Bu f f on • (Bu f f con -Bu f f Ca ) •Ca -k Bu f f o f f • Bu f f Ca Ḋ ye = k Fluo5 on • (Fluo5 f con -Dye) •Ca -k Fluo5 o f f • Dye.
We estimated the calcium reversal potential for the calcium fluxes using the Goldman-Hodgkin-Katz (GHK) flux equation described in 52 . The calcium ion permeability, P Ca , was used as a free parameter adjusting a single EPSP to produce a calcium 25/46 5.4. ARTICLE METHODS amplitude of ∼ 3 µM as reported in 53 .

φ = z Ca •V sp • F/R • (T + 273.15K) Φ Ca = -P Ca • z Ca • F • φ • [Ca] i -[Ca] o • e -φ 1 -e -φ (12)
Φ Ca is used to determine the calcium influx through NMDAr and VGCC in the Equations 13,14,15

NMDAr -GluN2A and GluN2B

Markov chain

In the hippocampus, the NMDAr are principally heteromers composed of the obligatory subunit GluN1 and either the GluN2A or GluN2B subunits. These N2 subunits guide the activation kinetics of these receptors with the GluN1/GLUN2B heteromers displaying slow kinetics (∼ 250ms) and the GluN1/GluN2A heteromers displaying faster kinetics (∼ 50ms). We modeled both NMDA subtypes. The NMDAr containing GluN2A is modeled with the following Markov chain 54 where we introduce the additional parameters

ρ NMDA f , ρ NMDA b : A 0 A 1 A 2 A 3 A 4 A O1 A O2 k a • [Glu] • ρ NMDA f k -a • ρ NMDA b k b • [Glu] • ρ NMDA f k -b • ρ NMDA b k c • ρ NMDA f k -c • ρ NMDA b k d • ρ NMDA f k -d • ρ NMDA b k e • ρ NMDA f k -e • ρ NMDA b k f • ρ NMDA f k -f • ρ NMDA b
The NMDAr containing GluN2B is modeled with a Markov chain based on the above GluN2A scheme. We decreased the rates by ∼75% in order to match the GluN2B decay at 25°C as published in 55 .

B 0 B 1 B 2 B 3 B 4 B O1 B O2 s a • [Glu] • ρ NMDA f s -a • ρ NMDA b s b • [Glu] • ρ NMDA f s -b • ρ NMDA b s c • ρ NMDA f s -c • ρ NMDA b s d • ρ NMDA f s -d • ρ NMDA b s e • ρ NMDA f s -e • ρ NMDA b s f • ρ NMDA f s -f • ρ NMDA b
The different rates are given in Table M8.

NMDAr and age switch

The age-dependent expression ratio of the subtypes GluN2A and GluN2B (r age ) was obtained from hippocampal mice data 56 . We added noise to this ratio causing ∼1 NMDAr subunit to flip towards GluN2A or GluN2B (see Fig. M5E). The population of 15 NMDAr is divided in the two subtypes according to the ratio plotted in Fig. M5B as function of age. The ratio to define the number 57 compared to the decay time from the GluN2B channel estimated by our model (yellow) and data from Iacobussi's single receptor recording (purple) 55 . b, Comparison of our implementation of GluN2B:GluN2A ratio and the GluN2B:GluN2A ratio from the mice's CA1. c, Comparison of our implementation of NMDAr conductance change in response to the extracellular against data 58 . d, Forward and backwards temperature factors implemented to approximate NMDAr subtypes decay times at room temperature 55 and temperature changes observed in 57 . e, NMDAr subtypes number on our model given age. We add noise to have a smoother transition between different ages. f, Calcium concentration changes for causal and anticausal protocols in response to different aCSF calcium and magnesium compositions with fixed Ca/Mg ratio (1.5). Scale 50 ms and 5 µM.

NMDAr subtypes as function of age reads:

r age = 0.507 + 0.964 1 + e 0.099•(age-25.102 days) + N (0, 0.05)

N GluN2B = round N NMDA • r age r age + 1 N GluN2A = round N NMDA r age + 1 .
The round term in the two previous equations ensures that we have an integer value for the NMDAr subtypes, making the stair shaped curve seen in Fig. M5E.

NMDAr and temperature

We adjusted the GluN2A and GluN2B forward and backward rates to follow the temperature effects on NMDAr-mediated EPSP 57 (see Fig. M5A andD). Because GluN2B dominates the NMDAr-mediated EPSP, we fit its decay time on the NMDAr-mediated EPSP as function of temperature as reported by 57

NMDAr current and Ca 2+ -dependent conductance

NMDAr conductance is modulated by external calcium and is modelled according to the next equations using NMDAr subconductances A O1 and A O2 (GluN2A), and B O1 and B O2 (GluN2B).

γ NMDA = 33.949 + 58.388 1 + e 4•([Ca 2+ ] o -2.701 mM) pS B(V sp , [Mg] o ) = 1 1 + [Mg] o 3.57 • e -0.062•V sp NMDA = (B O1 + B O2 + A O1 + A O2 ) • B(V sp , [Mg] o ) • γ NMDA I NMDA = (Erev NMDA -V sp ) • NMDA
We now modify the conductance function γ NMDA reported by 58 . The reported NMDAr conductance at [Ca 2+ ] o = 1.8 mM is 53 ± 5pS. Here, we used the higher conductance 91.3 pS for NMDAr (for both subtypes) at [Ca 2+ ] o = 1.8 mM to compensate for the small number of NMDArs reported by 59 . Hence, we cover 58 data differently to account for this constraint: this caused a right-shift in the NMDA-conductance curve (Fig. M5C). The calcium influx Ca NMDA is modulated by the GHK factor, Equation ( 12), as function of the internal and external calcium concentrations and the spine voltage:

Ca NMDA = f Ca • Φ Ca • NMDA. ( 13
)
The combined effect of Magnesium 60 and extracellular Calcium concentration are displayed in Fig. M5F.

GABA(A) receptor

Markov chain

We used the GABA(A) receptor Markov chain (Fig. M6) presented in Bush and Sakmann 1990 63,64 and we estimated temperature adaptations using Otis and Mody 1992 65 measurements.

C 0 C 1 C 2 O 1 O 2 r b1 • [Gaba] r u1 r b2 • [Gaba] r u2 r o1 r c1 • ρ GABA b r o1 r c2 • ρ GABA b

GABA(A)r and temperature

Because the amplitude of GABA(A) current is altered by the GABAr shift 66 during development, we applied temperature changes only to the closing rates using a logistic the function ρ GABA b estimated by fitting 65 measurements (data comparison in the Fig. M7B andE). 191•(T -32.167) .

ρ GABA b = 1.470 - -1.279 1 + e 0.

GABA(A)r current and age switch

The GABA(A)r-driven current changes during development 3 passing from depolarizing (excitatory) to hyperpolarizing (inhibitory) 67 .

That is, the reversal potential of chloride ions permeating GABA(A)r shifts from above the membrane resting potential (inward driving force -excitatory) to below the membrane resting potential (outward driving force -inhibitory) 66 . Such effect mediated by chloride ions is associated with the KCC2 pump (K Cl co-transporter) which becomes efficient in extruding chloride ions during maturation 66 . To cover the GABA(A)r shift, we fit the chloride reversal potential (E Cl rev ) using the data published by 66 (Fig. M7C):

E Cl rev = -92.649 + 243.515 1 + e 0.091•(age-0.691 days)

I GABA = (O 1 + O 2 ) • (E Cl rev -V dend ) • γ GABA .
Table M9 presents the parameters to model the GABAr. 65 for GABA(A)r current at 25°C (same normalization as in panel B). f, Change in the polarization of GABA(A)r currents given the age driven by the E Cl rev .

VGCC -T, R and L type

Markov chain

A stochastic VGCC model was devised using the channel gating measurements from CA1 rat's (2-8 weeks) pyramidal neurons by Magee and Johnston 1995 at room temperature 70 . Our model has three different VGCC subtypes described by the Markov chains in Fig. M8: the T-type (low-voltage), the R-type (medium-to-high-voltage) and the L-type (high-voltage). The VGCC Markov chain derived from Magee and Johnston 1995 70 are composed of two gates (h,m) for T-(Fig. M9A andD) and R-types (Fig. M9B andE) and a single gate for L-type (Fig. M9C), as described in the equations below.

C 0 C 1 C 2 O R α R m (V sp ) • ρ V GCC f β R m (V sp ) • ρ V GCC b α R m (V sp ) • ρ V GCC f β R m (V sp ) • ρ V GCC b α R h (V sp ) • ρ V GCC f β R h (V sp ) • ρ V GCC b α R h (V sp ) • ρ V GCC f β R h (V sp ) • ρ V GCC b O L1 C 0 O L2 β L 1 (V sp ) • ρ V GCC b α L (V sp ) • ρ V GCC f α L (V sp ) • ρ V GCC f β L 2 (V sp ) • ρ V GCC b C 0 C 1 C 2 O T α T m (V sp ) • ρ V GCC f β T m (V sp ) • ρ V GCC b α T m (V sp ) • ρ V GCC f β T m (V sp ) • ρ V GCC b α T h (V sp ) • ρ V GCC f β T h (V sp ) • ρ V GCC b α T h (V sp ) • ρ V GCC f β T h (V sp ) • ρ V GCC b
R h = 100 h R in f (V sp ) = 1 1 + e Vsp+39 9.2 α R h (V sp ) = h R in f τ R h β R h (V sp ) = 1 -h R in f τ R h R-type m-gate rates β R m = 40 m R in f = 1 1 + e 3-10 8 α R m r = β R m • m R in f 1 -m R in f τ R m = 1 α R m + β R m m R in f = 1 1 + e 3-Vsp 8 α R m (V sp ) = m R in f τ R m β R m (V sp ) = 1 -m R in f τ R m T-type h-gate rates τ T h = 50 h T in f (V sp ) = 1 1 + e Vsp+70 6.5 α T h (V sp ) = h T in f τ T h β T h (V sp ) = 1 -h T in f τ T h
T-type m-gate rates

β T m = 1 m T in f = 1 1 + e -32+20 7 α T m r = β T m • m T * in f 1 -m T in f τ T m = 1 α T m + β T m m T in f = 1 1 + e -32-Vsp 7 α T m (V sp ) = m T in f τ T m β T m (V sp ) = 1 -m T in f τ T m L-type rates α L (V sp ) = 0.83
1 + e 13.7-Vsp 6.1

β L 1 (V sp ) = 0.53
1 + e Vsp-11.5 6.4

β L 2 (V sp ) = 1.86
1 + e Vsp-18.8 6.17

VGCC and temperature

We used the same temperature factor for every VGCC subtype, respectively ρ V GCC f and ρ V GCC b (see Fig. M9F), as follows: -36.279) .

ρ V GCC f = 2.503 - 0.304 1 + e 1.048•(T -30.668) ρ V GCC b = 0.729 + 3.225 1 + e -0.330•(T
The VGCC subtypes are differently sensitive to temperature, with temperature factors for decay times ranging from 2 71 to 50-fold 72 . It further complicates if T-type isoforms are considered. Indeed, they can have temperature factors that accelerate or slow down the kinetics. For instance, when passing from room to physiological temperatures, the isoform Ca v 3.3 has a closing time ∼50 % faster 71 and the isoform Ca v 3.1 becomes ∼15 % slower. To simplify, the same temperature factor is adopted to all VGCC subtypes.

VGCC currents

The VGCC currents are integrated to the dendritic spine and estimated using the GHK Equation ( 12), as follows:

I T = γ T • Φ Ca • O T (14)
I R = γ R • Φ Ca • O R (15)
I L = γ L • Φ Ca • (O L1 + O L2 ) ( 16
)
31/46 ) and both deactivation rates (β L 2 (V sp ) and β 1 2 (V sp )) for the L-type VGCC. d, Activation (α h (V sp )) and deactivation rates (β h (V sp )) for the T-type h-gate. e, Activation (α h (V sp )) and deactivation rates (β h (V sp )) for the R-type h-gate. f, Temperature factor applied to all the rates, forward change (ρ V GCC f ) for the α rates and backward change (ρ V GCC b ) for the β rates.

SK channel

The small potassium (SK) channel produces hyperpolarizing currents which are enhanced in the presence of intracellular calcium elevations. The SK channel current was based on the description of 62 as follows:

r(Ca) = Ca σ Ca σ + h σ SK ṁsk = r(Ca) • ρ SK f -m s τ SK • ρ SK b I SK = γ SK • (E SK rev -V sp ) • m sk • N SK .
We chose a temperature factor to decrease the decay time of hyperpolarizing currents by a factor of two when passing from physiological to room temperature. Despite that the ences for temperature effects on the SK channel are few, a report 74 suggests a left-ward shift in the SK half-activation when changing from 37°C (h SK = 0.38 ± 0.02 µM) to 25°C (h SK = 0.23 ± 0.01 µM) ; that is a 65% decrease.

Enzymes -CaM, CaN and CaMKII

Markov chain

To model the enzymes dynamics, we adapted a monomeric CaM-CaMKII Markov chain from 77 which builds over 78 . Our adaptation incorporates a simplified CaN reaction which only binds to fully saturated CaM. That is, CaM bound to four calcium ions on N and C terminals (see Markov chain in the Fig. M10). A consequence of the Pepke coarse-grained model is that calcium binds and unbinds simultaneously from the CaM terminals (N,C). We assumed no dephosphorylation reaction between CaMKII and CaN since 79 experimentally suggested that no known phosphatase affects CaMKII decay time which is probably caused only by CaM untrapping 79 . This was previously theorized by 80 's model, and it is reflected in Chang data 53,77 . The structure of the corresponding Markov chain is shown in Fig. M10.

Chang's data 77 provides a high-temporal resolution fluorescence measurements for CaMKII in the spines of rat's CA1 pyramidal neurons and advances the description of CaMKII self-phosphorylation (at room temperature). We modified Chang's model of CaMKII unbinding rates k 2 , k 3 , k 4 , k 5 to fit CaMKII dynamics at room/physiological temperature as shown by 53 supplemental files. Previous modelling of CaMKII 77,78 used a stereotyped waveform with no adaptation to model calcium. Our contribution to CaMKII modelling was to use calcium dynamics sensitive to the experimental conditions to reproduce CaMKII data, therefore, allowing us to capture physiological temperature measurements from 53 . Note that CaMKII dynamics has two time scales and we only capture the fastest one (after stimulation ceases, 60 s) and the relative amplitude of CaMKII between the different temperatures. The slowest one occurs at the end of the stimulus, close to the maximum (Fig. M11A); this can be caused by the transient volume increase in the spine as measured by 53 .

Table M12 shows the concentration of the enzymes and Table M13 shows the parameters to model enzymes reactions in the Table M12. Concentration of each enzyme.

Temperature effects on enzymatic-activity

We then included temperature factors in the coarse-grained model using Chang data 77 , as shown in Fig. M11. For CaMKII, we fitted the modified dissociation rates of the phosphorylation states k 2 , k 3 and k 5 to match the data on relative amplitude and decay time using the following logistic function: 511(T -45.475°C) .

ρ CaMKII b = 162.171 - 161.426 1 + e 0.
For CaN, we fit the 85 ' data at 25°C as seen in Fig. M12A. However, since CaN-CaM dissociation rates at physiological temperatures were not found, we set the temperature factor to CaN that fits the outcomes of the protocols we proposed to reproduce. A reference value from the CaN-AKAP79 complex 86 showed a Q 10 = 4.46 = (2.19 s -1 /9.78 s -1 ) which is nearly the temperature factor used in our model for CaM. Therefore, both the association and dissociation rates are modified using the following logistic functions: 78 The reaction description matches with the color: Releases 2Ca, consumes 2Ca, consumes mKCaM, releases mKCaM, releases CaM2C, CaM2N, CaM0, CaM2N, releases mCaN, consumes mCaN, phoshorylate K units to P units, phosphorylated states and dephosphorylation. , aplied to CaN association and dissociation rates. Our dissociation rate becomes 0.006 s -1 for 22°C (same temperature as in 87 ), 0.0012 s -1 ), at the same temperature the association rate becomes 9.45 •10 6 M -1 s -1 (46 • 10 6 M -1 s -1 in 87 ).

ρ CaN f = 2.
k 2N f k 2N b k 2N f k 2N b k 2C f k 2C b k 2C f k 2C b k K2N f k K2N b k K2N f k K2N b k K2C f k K2C b k K2C f k K2C b k CaM0 b k CaM0 f k CaM2N b k CaM2N f k CaM2C b k CaM2C f k CaM4 b k CaM4 f F • k 1 F • k 1 F • k 1 F • k 1 ρ CaMKII b • k 5 k 4 ρ CaMKII b • k 2 ρ CaMKII b • k 2 ρ CaMKII b • k 2 ρ CaMKII b • k 2 ρ CaN f • k CaN f ρ CaN b • k CaN b ρ CaMKII b • k 3

35/46

CaM0 + 2Ca ⇒ CaM2C CaM2N + 2Ca ⇒ CaM4 k 2C f = adapt(k 1C on , k 2C on , k 1C o f f , k 2C on ,Ca) 78 CaM0 + 2Ca ⇒ CaM2N CaM2C + 2Ca ⇒ CaM4 k 2N f = adapt(k 1N on , k 2N on , k 1N o f f , k 2N on ,Ca) 78 CaM2C ⇒ CaM0 + 2Ca CaM4 ⇒ CaM2N + 2Ca k 2C b = adapt(k 1C o f f , k 2C o f f , k 1C o f f , k 2C on ,Ca) 78 CaM2N ⇒ CaM0 + 2Ca CaM4 ⇒ CaM2C + 2Ca k 2N b = adapt(k 1N o f f , k 2N o f f , k 1N o f f , k 2N on ,Ca) 78 k 1C on = 5 • 10 6 M -1 s -1 1.2 to 9.6 • 10 6 M -1 s -178 k 2C on = 10 • 10 6 M -1 s -1 5 to 35 • 10 6 M -1 s -178 k 1N on = 100 • 10 6 M -1 s -1 25 to 260 • 10 6 M -1 s -178 k 2N on = 200 • 10 6 M -1 s -1 50 to 300 • 10 6 M -1 s -178 k 1C o f f = 50 s -1 10 to 70 s -178 k 2C o f f = 10 s -1 8.5 to 10 s -178 k 1N o f f = 2000 s -1 1 • 10 3 to 4 • 10 3 s -178 k 2N o f f = 500 s -1 0.5 • 10 3 to > 1 • 10 3 s -178
Coarse-grained model, KCaM-Ca reactions M13. Parameters for the coarse-grained model published in 78 and adapted by 77 and this work. 78 G andH) for temperature and age parameters obtained from an estimation of the body temperature regulation during development (or thermoregulation maturation, also called maturation of temperature homeostasis, estimated in Fig. S5G). Fig. S3 expands the presynaptic burst strategy hypothesized to recover the LTD in adult slices (Fig. C) for 900 pairing repetitions. Also, Fig. S3 tries to isolate the contribution of each age-dependent mechanism (NMDAr, GABAr, BaP efficiency switches) for 3 and 5 Hz predictions in 6 experiment. To this we fixed each of the three mechanisms coding for age in our model at P5 and P50, to observe how they shape the plasticity. Note the experiment in Fig. S4D-I

KCaM0 + 2Ca ⇒ KCaM2C KCaM2N + 2Ca ⇒ KCaM4 k K2C f = adapt(k K1C on , k K2C on , k K1C o f f , k K2C on ,Ca) 78 KCaM0 + 2Ca ⇒ KCaM2N KCaM2C + 2Ca ⇒ KCaM4 k K2N f = adapt(k K1N on , k K2N on , k K1N o f f , k K2N on ,Ca) 78 KCaM2C ⇒ KCaM0 + 2Ca KCaM4 ⇒ KCaM2N + 2Ca k K2C b = adapt(k K1C o f f , k K2C o f f , k K1C o f f , k K2C on ,Ca) 78 KCaM2N ⇒ KCaM0 + 2Ca KCaM4 ⇒ KCaM2C + 2Ca k K2N b = adapt(k K1N o f f , k K2N o f f , k K1N o f f , k K2N on ,Ca) 78 k K1C on = 44 • 10 6 M -1 s -1 78 k K2C on = 44 • 10 6 M -1 s -1 78 k K1N on = 76 • 10 6 M -1 s -1 78 k K2N on = 76 • 10 6 M -1 s -1 78 k K1C o f f = 33 s -1 78 k K2C o f f = 0.8 s -1 0.49 to 4.9 s -178 k K1N o f f = 300 s -1 78 k K2N o f f = 20 s -1 6

Fig. S3

. | Duplets, triplets and quadruplets for FDP, perturbing developmental-mechanisms for LFS and HFS in 6 , and age-related changes in STDP experiments 1,7,8 . Related to Fig. and . a, Mean synaptic weight change (%) for the duplet-FDP (2Pre50) experiment varying age. The panel shows showing that not only LTD is enhanced but also LTP. b, Mean synaptic weight change (%) for the triplet-FDP (3Pre50) experiment varying age. The panel shows that LTD magnitude further enhanced for adult rats and that a leftward shift of the LTD-LTP transition. c, Mean synaptic weight change (%) for the quadruplet-FDP (4Pre50) experiment varying age. The panel shows a further leftward shift on the LTD-LTP transition (compared to 3Pre50). d, Mean synaptic weight change (%) for the 1 Pre 900 at 30 and 3 Hz with 6 . The panel shows the fixed NMDAr at P5 (more GluN2B) causing an increase of LTD magnitude and a slight increase of LTP magnitude for adult rats compared to baseline (grey solid line). e, Same experiment as panel D but fixing BaP maturation at P5 (higher BaP attenuation). LTP is abolished, but LTD is not affected. This is because AP induced by EPSP attenuate too fast for 30 Hz not able to produce enough depolarization to activate NMDAr Mg-unblock. f, Same experiment as panel D but fixing GABAr maturation at P5 (excitatory GABAr) what causes only slightly enhances LTD (3 Hz) for adult rats. g, Same experiment as panel D but fixing NMDAr at P50 (more GluN2A). LTD appears with decreased magnitude for young rats compared to baseline (grey solid line). h, Same experiment as panel D but fixing BaP maturation at P50 (less BaP attenuation). LTP is enhanced for young rats because the BaP pairing with the slow closing GluN2B produces more calcium influx. i, Same experiment as panel D but fixing GABAr maturation at P50 (inhibitory GABAr) which does not affect the FDP experiment. j, Mean synaptic weight change (%) for 8 's single versus burst-STDP experiment for different ages. The data from Meredith (boxplots) were pooled by the age as shown in the x-axis. The solid line represents the mean, and the shaded ribbon the 2nd and 4th quantiles simulated by the model (same for panels A-F). k, Mean synaptic weight change (%) for 7 's STDP experiment in which the number of postsynaptic spikes increases. The x-axis marker from 14-21 indicates that only this interval was published without further specification. We use our model to estimate age related changes to 7 protocols. Note that the model does not cover the 1Pre2Post10 properly (model predicts only outcomes near the first data quantile). Notice that single and burst STDP leads to LTD, meanwhile 8 's to LTP or NC. l, Mean synaptic weight change (%) for 1 's STDP experiment which compares single versus burst STDP. The x-axis marker from 50-55 indicates that only a interval was published without further specification. We use our model to estimate age related changes to 1 protocols. It is noticeable that each STDP experiment has a different development. Fig. S5 shows multiple aspects related to temperature in STDP experiments and the temperature and age choices for the publications described in Table M1 compared to physiological conditions. We estimate how the rat's body temperature physiologically evolves in function of age using 9 and 10 's data. M1) showing a full LTD window. Our model also reproduces the fact that increasing temperature to 32-34°C the LTD is abolished (data not shown). b, Mean synaptic weight change (%) for 11 's STDP experiment for 1Pre2Post10, 70-100 at 5 Hz (see Table M1) showing a bidirectional window. c, Mean synaptic weight change (%) for 11 's STDP experiment for 1Pre2Post10, 20-30 at 5 Hz (see Table M1) showing a bidirectional window. We report that for 11 experiment done in room temperature the temperature sensitivity was higher than other experiments. d, Core temperature varying with age representing the thermoregulation maturation. This function (not shown) was fitted using rat 10 and mouse data 9 added by 1°C to compensate species differences 10 . The blue and white bars represent the circadian rhythm as shown in 9 . However, the "rest rhythm" for young rats (P5-14) may vary. e, Plot showing how far from being physiological are plasticity experiments done in physiological temperatures. Suggesting, there is scarcity of physiologically relevant data to model and understand plasticity. The dashed grey line is an approximation of the mean value from panel G. f, Initial conditions for CaN-CaMKII resting concentration for different [Ca 2+ ] o and temperature values. When [Ca 2+ ] o is changed temperature is fixed at 35°C, while temperature is changed [Ca 2+ ] o is fixed at 2 mM.

44/47

Chapter 6

Additional results and discussion

As shown in the article, the model can reproduce and predicts plasticity outcomes for experiments in different experimental conditions. Here, I will summarize other experiments that were not inserted in the main text thereby expanding the article scope and providing other preliminary results, such as: a synaptic version of CaN-CaMKII nonlinear filter by Fujii et al. 2013, the role of the temperature range in the frequency-dependent plasticity by O'Connor et al. 2005, the limits of nondynamical measures to predict plasticity as in Blackwell et al. 2019, a comparison between different types of presynaptic release and a deterministic version of the model, the possible interactions between laser-tissue warming and Ca 2+ dye, and two examples of in vivo-like plasticity.

CaN and CaMKII related experiments

Our model was based on the combined CaN and CaMKII activity in which the memory role was hypothesized and modelled by Lisman 1985;Lisman 1989. Lisman's core idea was that CaMKII would stay permanently activated through its self-phosphorylated activity and deactivated by CaN in a switch-like fashion. Also, the CaMKII activated states would surpass protein turnover by phosphorylating other CaMKII. Despite the innovation of such ideas [START_REF] Bear | In memoriam: John Lisman-commentaries on CaMKII as a memory molecule[END_REF], the memory maintenance role through self-phosphorylation (Chang et al. 2017;[START_REF] Michalski | The delicate bistability of CaMKII[END_REF] and phosphatase influence over CaMKII dephosphorylation (Otmakhov et al. 2015;Fujii et al. 2013) were not confirmed in ex vivo/in vivo experiments. However, a possible role of CaMKII in sustaining a transient "molecular memory" can be achieved by the CaMKII-GluN2B complex, which still needs to be investigated in vivo/ex vivo (Urakubo et al. 2014;Paul J Michalski 2014). Furthermore, fluorescence experiments by Fujii et al. 2013;Chang et al. 2017;Chang et al. 2019 show that CaMKII self-phosphorylation does not make it permanently activated as hypothesized by P. Michalski 2013 model. That is possibly due to the difference in the saturation profile of these enzymes in dendritic spines (model) and the soma (experiment). For instance, in Figure 6.1D, somatic CaMKII can accumulate activation for high frequency and pulses; in contrast, Figure 6.1A, the model shows that activity was slightly higher around 1 Hz, with a 20% difference between other frequencies. Such lack of frequency-dependence was noticed in Chang et al. 2017 CaMKII measurements. Thus, the synaptic version of the nonlinear in Figure 6.1C filter had an inverted outcome to what was measured by Fujii et al. 2013 in Figure 6.1F, with CaMKII less active for high frequency and pulse number than CaN for the synaptic version. The model suggests that the amplitude is a poor predictor for enzymatic activity since it can vary depending on where it is measured; CaN and CaMKII are more expressed in dendritic spines than the soma. Note that the model cannot evaluate the enzymatic dynamics in the soma measured by Fujii et al. 2013 since it is built just for dendritic spines. LTP was prevalent (Figure 6.2A under the red bottom line). Given the possible qualitative changes in the temperature range, precise temperature control would be more useful to identify frequency-dependent changes.

A synaptic version of the CaN and CaMKII nonlinear filter

O 'Connor et al. 2005 argument was that CaMKII controls LTP and CaN LTD. In Figure 6.2C, their data on inhibiting these enzymes showed how the FDP would be modified. However, given the qualitative temperature changes range that the model has identified, such experiments need to be better constrained. Our model does not reproduce data in Figure 6.2C since the CaMKII, or CaN inhibition effects could not be determined.

The limits of non-dynamical measures to predict plasticity (failed attempt)

Before developing the geometrical readout in the main publication, a strategy similar to the one in Blackwell et al. 2019 to predict plasticity was attempted. The simulated dynamics are collapsed into a single value (e.g. mean or median), and the prediction part would be based on the separability between outcome classes (LTP, LTD and no change) as shown in Figure 6.3. We tested the separability of different measures such as integral, max value, mean, median using Tigaret et al. 2016 protocols. We noticed that certain protocols were too overlapped (gray region) to be suitable for multiclass classifiers. In Blackwell et al. 2019 the plasticity prediction between pathological and healthy were done by averaging a protein in a given time interval, the shortcomings of this approach could be addressed by our geometrical threshold which also can be used to detect dynamical differences in protein dynamics without having the information collapsed in a single value.

How stochasticity affects the model 6.3.1 Averaged and discrete vesicle release

Part of the computational neuroscience models uses averages as input [START_REF] Tsodyks | Neural networks with dynamic synapses[END_REF]). However, the vesicle release in nature has an all-or-none behaviour (or discrete multivesicular release). Here we will test the differences between a discrete presynaptic release used in the article and an averaged presynaptic release. The left panel in Figure 6.4A shows the glutamate concentration between the discrete release that mimics a single vesicle and the average input estimated from our presynaptic model. Figure 6.4A illustrates the difference of the activation on these two approaches in the model, only changing the presynaptic input. Note in Figure 6.4A (Left) that since glutamate concentration is never zero for the averaged release, the enzymes elicited with the 1Pre1Post10 protocol from Tigaret et al. 2016 will be more activated as shown in Figure 6.4B (Right). Using an averaged stimulation as input may cause overstimulation since failures never occur. Also, the zero glutamate concentration is not reached, as shown in Figure 6.4A (Left). Thus, this experiment suggests that averaged inputs used in neuron networks may be outside the glutamate concentration a synapse usually experiences. The assumption that the neuron experiences an averaged concentration of glutamate may not be equivalent to stimulating a synapse with an averaged input.

Comparison with a deterministic version of the model

To show the relevance of stochasticity and to evaluate if it is possible to have a faster simulation, a deterministic version of the model was made and tested using the 1Pre1Post10 protocol from Tigaret et al. 2016. All Markov chains described in the methods 5.4 were implemented as a first-order chemical reaction. Figure 6.4B shows a comparison between the deterministic (purple lines) and stochastic version (yellow lines) using both discrete release and also a deterministic version using averaged release (green). Note that the deterministic and stochastic simulations using discrete release are slightly similar. It was identified that the source of this difference comes from the T-type flickering, which is continuously activated due to its low-voltage threshold (Magee and Johnston 1995). Meanwhile, other VGCCs, NMDAr and AMPAr, had little impact since their contribution are only at discrete times (pre and post spikes). Since we identified the T-type channel as the most relevant stochasticity source, we tested how the experimental conditions affect the versions of the model. Figure 6.4C shows how the CaN and CaMKII average baseline activity diverge for the model versions. When [Ca 2+ ] o increases, it amplifies the T-type long-term fluctuations on Ca 2+ and affects the Ca 2+ -binding enzymes resting concentrations.

Possible interactions between laser-tissue warming and

Ca 2+ dye

Fluorescent dyes revolutionized how we study neuronal functions (Knot et al. 2005); it allowed one to visualize Ca 2+ and proteins by exciting these dyes, for instance, with a laser. However, as with any other method that perturbs what is being observed, this method comes with drawbacks. Since the fluorescent dye binds to Ca 2+ , it modifies its diffusion properties (Maravall et al. 2000) and interactions with Ca 2+ -dependent pathways [START_REF] Zhang | Coupling of SK channels, L-type Ca 2+ channels, and ryanodine receptors in cardiomyocytes[END_REF]. For instance, as noted by Zhang et al. 2018 Ca 2+ dye disrupts interaction with the SK channel, which has a Ca 2+ and voltage-dependent activation. Also, plasticity protocols with and without the dye can display qualitative differences as shown by Mulkey and Malenka 1992. This experiment suggests that calcium dynamics observed with dye should carefully be considered (Maravall et al. 2000). For instance, in the Figure 6.5A-B, the protocols 1Pre2Post50 and 2Post1Pre50 from Tigaret et al. 2016 have similar Ca 2+ dye curves and produce different plasticity outcomes (without the dye). Our model better reproduces the dye curves when the temperature is increased and the SK channel contribution is set to zero [START_REF] Zhang | Coupling of SK channels, L-type Ca 2+ channels, and ryanodine receptors in cardiomyocytes[END_REF]. Also, in Figure 6.5C, the simulations predict that Tigaret et al. 2016 protocols have different shapes without the dye.

Furthermore, a concerning aspect of fluorescent dyes is the photobleaching or ablation caused by laser stimulation. Given the frequency, laser wavelength and duration, one can rapidly increase the temperature in neural tissue [START_REF] Picot | Temperature rise under two-photon optogenetic brain stimulation[END_REF]) and modify the kinetics and Ca 2+ source dynamics being observed (Podgorski and Ranganathan 2016;[START_REF] Ebbesen | Analysis of laser-induced heating in optical neuronal guidance[END_REF][START_REF] Oliver | Effects of temperature on calcium-sensitive fluorescent probes[END_REF]. Strategies to measure Ca 2+ fluorescence with high temporal resolution and manage heat transfer effects still need to be developed (Schmidt and Oheim 2018).

In vivo-like modelling

Since the experimental discovery of STDP (Markram et al. 2011) it has been hypothesized whether it is a relevant framework to understand in vivo spiking pat-terns. For instance, Froemke and Dan 2002 suggested that in vivo plasticity cannot be explained by only the delay time between pre and postsynaptic spikes since it disregards the history of previous spikes. Despite the possible lack of relevance for physiological learning, the STDP framework is widely used in neuromorphic computing sided by BCM. Graupner et al. 2016, using an irregular STDP with Poisson spikes simulated an in vivo-like firing, concluded that the plasticitydependence over the delay time is reduced. Another study using in vivo-like firing perturbed the delay time of the spike pairs (Cui et al. 2018). They show that perturbing the delay times reduces the LTP magnitude. Our model observed that depending on which perturbation is applied to the STDP firing structure, LTP could be either enhanced or abolished (see article Figure 7). Therefore, observing plasticity through one variable, such as the delay time, could be insufficient to explain in vivo firing. Also, the perturbing regular firing patterns to approximate in vivo-like firing is a limited approach enforced by data scarcity. However, new in vivo plasticity datasets (Ahnaou et al. 2020) offer a singular opportunity to test the robustness of current models for plasticity. For instance, "replay experiments" (Isaac et al. 2009;Bittner et al. 2017) or those only recorded with local field potential (Ahnaou et al. 2020). New sampling methods such as optogenetics STDP would allow the expansion of in vivo plasticity experiments (Anisimova et al. 2019).

Sleep stimulation patterns

Sleep is required for memory consolidation [START_REF] Whalley | Sleep to forget[END_REF] and the induction of long-term plasticity forms. However, the participation of characteristic neuronal oscillations is not fully understood. During sleep, the neural environment changes, as the Ca 2+ sources and its modulators (McCauley et al. 2020;[START_REF] Gordon | Effect of cage bedding on temperature regulation and metabolism of group-housed female mice[END_REF]. For instance, the astrocytes retract from the tripartite synapse modifying how glutamate is recycled during the sleep [START_REF] Haydon | Astrocytes and the modulation of sleep[END_REF]. Also, brain temperature oscillates, causing kinetic changes in the Ca 2+ sources [START_REF] Sela | Sub-minute prediction of brain temperature based on sleep-wake state in the mouse[END_REF]. A characteristic sleep oscillation, for instance, the non-rapid eyes movement phase (non-REM), displays activity spindles (or sharp-wave ripples, Sadowski et al. 2016) which induces LTP, as demonstrated by replay experiments in slices and in vivo [START_REF] Rosanova | Pattern-specific associative long-term potentiation induced by a sleep spindle-related spike train[END_REF]. It is hypothesised that sleep's primary function is to renormalise strengthened synapses by weakening (LTD) synaptic efficiencies that were strengthened during wakefulness [START_REF] Frank | Erasing synapses in sleep: is it time to be SHY?[END_REF][START_REF] Vivo | Evidence for sleep-dependent synaptic renormalization in mouse pups[END_REF]. Conversely, sleep deprivation shows reduced neuronal excitability, opposite to the strengthening hypothesis [START_REF] Borbély | Sleep deprivation decreases neuronal excitability and responsiveness in rats both in vivo and ex vivo[END_REF]. Our model could potentially be used to study replayed firing patterns related to sleep, possibly by including neuromodulators related to sleep. It can also be expanded by adding diffusion of enzymes, which could contribute to understand the interplay between neighbouring spines and investigate homeostatic plasticity [START_REF] Seibt | Primed to sleep: the dynamics of synaptic plasticity across brain states[END_REF]. Since sleep has a homeostatic role, the activation of a dendritic spine cluster during sleep can cause enzymes to diffuse out from the spine neck. For instance, Calcineurin can reach the dendritic shaft, but CaMKII has a shorter diffusion range, limiting the interaction with neighbouring spines (Fujii et al. 2013;Yasuda 2017). Yet, since sleep is differently important during development [START_REF] Vivo | Evidence for sleep-dependent synaptic renormalization in mouse pups[END_REF], coding for age should be a desirable model property when investigating plasticity related to sleep. An ideal way to study sleep patterns using the model would be with CA3-CA1 spikes recorded during behaviour, such as the data collected by Sadowski et al. 2016. Here, an empirical sleep-wave pattern was build based on the lowfrequency delta (1)(2)(3)(4) [START_REF] Siapas | Coordinated interactions between hippocampal ripples and cortical spindles during slow-wave sleep[END_REF]. Basically, a 100 pulses burst at 20 Hz presynaptic bursts separated by 1 Hz in 12 epochs as in Figure 6.6A. The experimental conditions used were physiologically plausible with 1.5 mM [Ca 2+ ] o and 1.2 mM [Mg 2+ ] o (Inglebert et al. 2020). To evaluate how physiological conditions constrain plasticity, age-dependent temperature regulation [START_REF] Wood | Treatment temperature and insult severity influence the neuroprotective effects of therapeutic hypothermia[END_REF]) was empirically estimated (Figure 6.6B). Figure 6.6C shows that a sleepwave pattern induces LTD during the first weeks of life, and later the same slowwave pattern induces less LTD (Figure 6.6C). The experiment in Figure 6.6C may reflect the synaptic pruning phase [START_REF] Semple | Brain development in rodents and humans: Identifying benchmarks of maturation and vulnerability to injury across species[END_REF]) in the rat development in which the synaptogenesis is up-regulated, and it is counterbalanced by synaptic elimination completed after P21-28 in rats [START_REF] Semple | Brain development in rodents and humans: Identifying benchmarks of maturation and vulnerability to injury across species[END_REF]. Note that further investigation, diversifying the firing patterns is necessary to tackle the development of sleep in rats since the oscillations (e.g. NREM, REM) are not the same through the first weeks of life [START_REF] Frank | Development of circadian sleep regulation in the rat: a longitudinal study under constant conditions[END_REF]).

Behavioral time scale synaptic plasticity (BTSP)

During behaviour, the hippocampus continuously maps the space. However, the plasticity rules forming spatial representation are not fully understood. An answer was proposed by Bittner et al. 2017. They discovered a synaptic rule for place field [START_REF] Hartley | Modeling place fields in terms of the cortical inputs to the hippocampus[END_REF] formation, which does not require immediate coincidence between pre and post-firing patterns, exceeding the canonical STDP delay time (Bi and Poo 1998) Bittner et al. 2017 BTSP curve. This preliminary result suggests that the same phenomenon may occur in rats. This result can motivate other behavioural time scale rules to be discovered. Bittner et al. 2017 mentions the possible signalling pathways evolved, and our model offers a dynamical explanation using CaN and CaMKII combined activity as shown in Figure 6.7A.

Chapter 7

Discussion and conclusion

The following discussion extends the one in the article and presents other arguments, limitations and questions the model could address.

Synaptic models can propose new questions and testable predictions, and in some cases, anticipate experimental advances on how learning and memory occur [START_REF] Sejnowski | The book of Hebb[END_REF]). Given the model scope, it is possible to accommodate the diversity of plasticity experimental outcomes or hypothesize over synaptic mechanisms interaction. In this thesis, a new model of synapse was proposed to unify usually neglected aspects affecting plasticity outcomes (as reviewed in chapter 3). Thus, a new synaptic rule is provided to link simulation and experiments using the combined activity of CaN and CaMKII interpreted with a geometrical readout.

It is not feasible for a model to include all differences between plasticity experiments since electrophysiology methods are highly heterogeneous (Tebaykin et al. 2018). However, an effort is necessary to model the differences which can contribute to understanding plasticity. For instance Letzkus et al. 2006 andSjöström et al. 2001 (Inglebert et al. 2020). With a model dedicated to identifying key experimental conditions, one can avoid possible contradictions and extend the number of covered experiments by the same model. The model in this thesis predicts that some stimulation protocols have sharper transitions (see plasticity maps in the article) depending on the experimental conditions (e.g. [Ca 2+ ] o inferior to 0.5 mM), thus, pointing out which parameters are more sensitive for ex vivo experiments and potentially improving the experimental design.

ex vivo experiments can be considered as extreme cases since the brain extracellular environment needs to be recreated. Usually, the firing patterns are different from in vivo and are highly correlated (i.e. regular). The model developed here can predict in vivo-like learning (e.g. Bittner et al. 2017) since it can be directly adjusted to follow physiological conditions. For instance, using Poisson spikes, i.e. simplified representation of in vivo firing ignoring bursts, the model found that an LTP window is maximal when the postsynaptic neuron fires at 10Hz for any presynaptic frequency (between θ and α rhythms) (Graupner et al. 2016).

Since the model codes for development aspects, it is possible to predict that such LTP window around the alpha-theta band emerges between the second and third week of life. This coincides with the beginning of rats' exploratory behaviour in which theta oscillations on EEG (4-10 Hz) are observed [START_REF] Vanderwolf | Hippocampal electrical activity and voluntary movement in the rat[END_REF][START_REF] Thiels | Weaning in rats: II. Pup behavior patterns[END_REF]). The model also predicts that two STDP protocols with firing patterns producing too weak or strong enzyme activation to generate LTP (Tigaret et al. 2016) can converge to the LTP region once they are jittered, suggesting that under irregular firing (with similar frequency), the specificity of spiking times or the delay are less relevant for in vivo-like plasticity.

Geometrical readout

The new synaptic rule presented in this thesis borrows the concept of orbits from dynamical systems to demonstrate the possibility of differentiating plasticity protocols using a stochastic model of postsynaptic Ca 2+ -cascade. Such synaptic induction rule states that the plasticity outcomes depends on the time spent by the joint CaMKII and CaN activity inside the plasticity regions and activate the LTP or LTD rates. This can be interpreted as an enzymatic version of the Ca 2+ duration hypothesis (Evans and Blackwell 2015) which claims that high and brief Ca 2+ transients induce plasticity. Meanwhile, long and moderate levels induce depression.

A Goldilocks hypothesis for plasticity (analogy with the Goldilocks principle) can be formulated from the geometrical threshold and can be tested through FRET probes: the enzymatic activity should not be excessively high or too low to induce plasticity, instead, it should reach and spend the right time inside the LTD or LTP plasticity region.

7.2

The role of stochasticity, is it a bug or a feature?

Neurons produce spikes that are not transmitted, resulting in synaptic failures (Rusakov et al. 2020). Hypothesis on the role of these failures are many: the failed spikes activate post to pre signalling, e.g. nitric oxide or retrograde endocannabi-128 CHAPTER 7. DISCUSSION AND CONCLUSION noid pathways [START_REF] Castillo | Endocannabinoid signaling and synaptic function[END_REF]; failed spikes form a self-communication within neuron's distant structures (dendrites, axon) (Short et al. 2017); that the absence of spikes is also part of the information between neurons. Examples of how failures are relevant or not is dubious, for instance, the hydra vulgaris thrive with low synaptic failure rates [START_REF] Dupre | Non-overlapping neural networks in Hydra vulgaris[END_REF]. Meanwhile, studies in rats suggest that spikes are the exception [START_REF] Borst | The low synaptic release probability in vivo[END_REF] due to high synaptic failure rate in vivo. However, brain regions have different synaptic failure rate given their function [START_REF] Joris | The calyx of Held: a hypothesis on the need for reliable timing in an intensity-difference encoder[END_REF], and present redundancy with multivesicular release [START_REF] Molnar | Human pyramidal to interneuron synapses are mediated by multi-vesicular release and multiple docked vesicles[END_REF]Rudolph et al. 2015). Also, various behaviours are altered when failures increase [START_REF] Prosser | The role of nervous systems in temperature adaptation of poikilotherms[END_REF]. The open question on why neurons produce spikes that will not result in a vesicle release requires more investigation. The model present here can capture the specificity of experimental conditions and firing adaptations (e.g. STD) in a stochastic manner and may be useful in raising new questions on the role of stochasticity. For instance, a hypothesis raised by the model, not discussed in the paper, is whether CaMKII autophosphorylation acts as a spike history reservoir against transmission failures. Once the failures start to be prevalent, the self-sustained activity prevents a total forgetting of the previously integrated spikes by the enzyme activation, allowing better integration of sparse spikes. A contrast with the mutant CaMKII without the autophosphorylation domain shows a CaMKII which decays faster and impairs learning. This hypothetical feature would allow that stimuli to be integrated even after the short-term depression predominate (more failures). Also, it could be useful during bursts in which vesicles needs time to be refiled every high-frequency event-cycle, for instance, as in the behavioural time scale plasticity or sleep-wave patterns (Bittner et al. 2017) shown in Figure 6.7.

Limitations

Considering the models in the mini-review of chapter 3, the model developed covers firing patterns structures, with extensive validation, and innovates in the experimental conditions. However, it falls short in the localization aspects, mainly for morphology, non-neuronal cells (e.g. microglia and astrocytes), presynaptic long-term plasticity and neuromodulators. The model is limited to only CA3-CA1 synapses from ex vivo studies. Meanwhile, abstract models, less constrained by data, can be more easily adapted to other brain regions.

The model uses a single set of parameters to predict plasticity. However, this may not be a valid assumption since plasticity dependent pathways are likely to be affected by experimental conditions. That was done for the sake of simplicity 7.4. FUTURE DIRECTIONS 129 and readability. A contrasting example is shown by Chindemi et al. 2020 which attributes to each synapse in a cortical neuron its own threshold system. Using our framework, one could investigate,using glutamate uncaging and CaN-CaMKII FRET probes, if the same synaptic rule is different or different synapses.

Another assumption is that enzymatic concentration is the same throughout the lifespan. However, it is suggested that for pup rats, the CaMKII concentration is reduced. An unchanging concentration is assumed in our model due to the lack of precise measurements. Also, despite the model robustness to experimental heterogeneities, a sensitivity analysis was not carried out to observe how robust the model is for parameter perturbations, such as variations in the numbers of receptors and ion channel (Marder and Taylor 2011). This is an important aspect to consider since stochastic effects on plasticity outcomes can have multiple sources: from the inherent flickering of receptors, vesicle release failures, different neurons, oscillations in the recording temperature.

Similarly to the fact that different parameters could equally fit neural data (Marder and Taylor 2011), another modelling approach suggests that the number of different subtypes of ion channels allows a richer repertoire of neuron dynamics [START_REF] Schneider | Biological complexity facilitates tuning of the neuronal parameter space[END_REF]. The variability of ion channels in the model developed is reduced with phenomenological modelling of BaP generation, and only three types VGCCs subtypes (see Jedrzejewska-Szmek et al. 2017). Also, such variability of expression is not fixed during the lifetime, Cizeron et al. 2020 and it was evaluated to be little differentiated during the infant and elder rats and more diversified during adulthood. Also, the model lacks validation material for older adult rats, which present different plasticity patterns from adults (Pinar et al. 2017).

Future directions

The model developed in this thesis can unfold different projects. However, a challenge of complex models is that they are difficult to replicate and scale-up. A deterministic version of the model was in order to gain performance and to evaluate the necessity of stochasticity, which showed to be most relevant for the presynaptic release and the resting concentrations of enzymes (see Figure 6.3). It is important to reduce the model complexity, without losing its capabilities, to make it more accessible. Another ongoing strategy is to make the model available online by a web tool in which experimentalist could use it to plan their experiments.

Another possibility is the application of the model to neuropathologies through 130 CHAPTER 7. DISCUSSION AND CONCLUSION Ca 2+ sources abnormalities. For instance, tau protein binds to the neurotransmitter vesicles and decrease release probability at high frequencies (Zhou et al. 2017). Also, the model can study APP intracellular interactions to evaluate the plasticity effects with an upregulated Ca 2+ influx. The study of disruption of a healthy synapse would make more sense if done with associated behaviours or in vivo-like firing patterns, such as those related to sleep, exploration and other behaviours. Also, a future direction for the model could be to include pathologies by having the contribution of different neuromodulators (e.g. acetylcholine Isaac et al. 2009) and glia cells (Nadkarni and Jung 2007). Also, having astrocyte contribution would be useful to understand how Ca 2+ -based communication in astrocytes modulates synaptic rules in neurons [START_REF] Verveyko | Modeling of astrocyte networks: towards realistic topology and dynamics[END_REF]. The limits between a pathological and healthy Ca 2+ influx should more precisely defined. Finally, a future direction is to formalize an algorithm to detect dynamical patterns. This could be relevant to explore multiple dynamical biophysical interactions automatically. This could be achieved by using the concept of time spent in a region (or volume) to understand how orbits produce different outcomes in a similar way to the supervised learning algorithms (Friedman et al. 2001). Also, an unsupervised learning algorithm could be potentially implemented using a geometrical threshold by attributing a criterion for class formation given the orbits distance or how they share the same space.

Figure 2 . 2 :

 22 Figure 2.2: Long term plasticity. A) High Ca 2+ entry caused by high-frequency stimulation. B) Induction of LTP using 50 Hz (adapted from (Dudek and Bear 1992)) showing the change in the excitatory postsynaptic potential (EPSP) slope (Byrne 2017). C) Receptor exocytosis and synaptic button enlargement. D) Low Ca 2+ entry caused by low-frequency stimulation. E) Induction of LTD using 3 Hz (adapted from (Dudek and Bear 1992)), same plasticity measure as in B. F) Receptor endocytosis and synaptic button shrinkage.

  Figure 2.2 depicts Dudek and Bear 1992's experiments inducing LTP and LTD, highlighting receptor endocytosis and exocytosis as the outcome for HFS and LFS stimulations, respectively. Dudek and Bear 1992 also blocked NMDA receptor to show how the long-term plasticity depends on it.

Figure 2 . 3 :

 23 Figure 2.3: Variations of the Ca 2+ hypothesis for plasticity. Left: Ca 2+ amplitude hypothesis. Middle: Ca 2+ duration hypothesis. Right: Ca 2+ location hypothesis.

Figure 2 . 4 :

 24 Figure 2.4: Ca 2+ -binding enzymes. The figure represents the CaMKII interaction with CaM, CaMKII phosphorylation of glutamatergic receptors and its role inducing AMPAr insertion (O'Connor et al. 2005). Also, it shows CaN interaction with CaM, CaN dephosphorylation of glutamatergic receptors and its role removing AMPAr (O'Connor et al. 2005).

Figure 2 . 5 :

 25 Figure 2.5: Histograms of experimental conditions in plasticity experiments using brain slices in the hippocampus and cortex.The data was extracted from hallmark plasticity studies reproduced by theoretical models(Graupner and Brunel 2012, Ebner et al. 2019, Jedrzejewska-Szmek et al. 2017, Inglebert et al. 2020). For articles which uses a interval, such as P14-21, the whole interval was counted. This histogram is based on the ones described inTebaykin et al. 2018.

Figure 3 . 1 :

 31 Figure 3.1: Examples of plasticity rules modulated by experimental conditions. On the Left, plasticity outcomes from LFS and TBS affected by age (see Dudek and Bear 1993). In the Middle, plasticity outcomes from STDP being affected by [Ca 2+ ]o and [Mg 2+ ]o (seeInglebert et al. 2020). On the Right, plasticity outcomes from TBS affected by temperature (seeKrelstein et al. 1990).

Figure 3 . 2 :

 32 Figure 3.2: Idealized experimental conditions implementation. On the Left, all implemented mechanisms have temperature control. In the Middle, All implemented ions have their reversal potential considering internal and external concentrations. On the Right, all implemented mechanisms are effects by age.

   show different experimental paradigms classified inTable 3.4.

Figure 3 . 3 :

 33 Figure 3.3: Idealized experimental firing patterns which could be covered by a single model. On the left, illustrations depict BCM related goals, data-driven comparisons, BCM-like modifications (sliding threshold and amplification), and BaP induced by EPSP. In the middle, data-driven comparison and different STDP curves (as in Wittenberg and Wang 2006). On the right, it illustrates a replay experiment in which spikes are recorded in CA3-CA1 in vivo and the firing patterns used as input in ex vivo.

  concentration for all ions external/internal concentration for one ions arbitrary concentration unit no control FIRING STRUCTURES FDP reproduces frequency-dependent data sliding threshold (if shown for the

4. 1 . 59 Figure 4 . 1 :

 15941 Figure 4.1: Increase of orbit separability with different dimensionalities. A) 20 pulses of two variables with different rise and decay constant representing different plasticity protocols. B) The same variables in panel B added combined with two other variables with different decay and rise times. C) Same as A and B combined to other two variables. D) Distance between curves a and b in panel A. E) Distance between curves a and b in panel B. F) Distance between curves a and b in panel C.

Figure

  Figure 4.2: Unidimensional and geometrical threshold.A) Variable with a slow decay (τ = 1.5s) integrating 30 pulses at 10 Hz with threshold at 2 a.u.. The color gradient (white to green) depicts the increment size in variable at each pulse B) Variable with a faster decay time (τ = 0.2s) integrating 30 pulses at 10 Hz with threshold at 2 a.u. . The color gradient (white to blue) depicts the increment size. C) Combined dynamics of variables in panels A and B showing a geometrical threshold (defined by the points[7.5,1.],[7.5,2.],[2.,2.] and[2.,1.]). D) Time above the threshold (t.a.t) for different increment sizes (white to green) for the variable in the panel A. E) Time above the threshold (t.a.t) for different increment sizes (white to blue) for the variable in the panel B. F) Time inside the polygon (t.i.p) for different increment sizes (white to red) for the variable in the panel B. G) Total time over the threshold for the variable in panel A changing the increment size. H) Same as G for panel B. G) Total time inside the polygon for orbits in the panel C.

64CHAPTER 5 .

 5 A STOCHASTIC MODEL OF HIPPOCAMPAL SYNAPTIC PLASTICITY WITH GEOMETRICAL READOUT OF ENZYMES -ARTICLE NMDAr conductance (Maki and Popescu 2014) which were phenomenologically implemented.

Fig. 1 .

 1 Fig. 1. | The synapse model, its timescales and mechanisms. a, Model diagram with the synaptic components, the pre and postsynaptic compartments. Inhibitory receptor bottom left. b, Stochastic dynamics of the different receptors and channels. Plots show the total number of open channels/receptors as a function of time. AMPArs and NMDArs are activated by glutamate, VGCC are activated by membrane potential, and GABArs are activated by GABA. The timescale of variable response increases from top to bottom panels. c, Dendritic spine membrane potential (left) and calcium concentration (right) as function of time for a single causal (1Pre1Post10) stimulus. d, Left: depletion of vesicle pools (reserve and docked) induced by 30 pairing repetitions delivered at 5 Hz 45(Methods). The same vesicle depletion rule is applied to both glutamate and GABA. Right: BaP efficiency as function of time. BaP efficiency controls the axial resistance between soma and dendrite in order to phenomenologically capture the distance-dependent BaP attenuation46,47 (Methods). e, Activated enzyme concentration for CaM, CaN and CaMKII, as function of time for the stimulus 1Pre1Post10, 30 pairing repetitions delivered at 5 Hz. Note that the vertical grey bar is the duration of the stimuli, 6 s.

Fig. 2 .

 2 Fig. 2. | The duration and amplitude of the joint CaN-CaMKII activity differentiates plasticity protocols. a, Activity of CaMKII (solid line) and CaN (dashed line) (µM) for two protocols. Experimentally the 1Pre2Post10 produces LTP, and 1Pre1Post10 produces no change (NC). Both are composed of 300 pairing repetitions delivered at 5 Hz. b, Joint enzymatic activity (CaN-CaMKII) for the protocols in panel a. The black dot indicates the initial resting activity and the arrows the trajectory direction as function of time. The grey points mark the time position (x-axis in panel a) for both protocols at 2, 10 and 60 s (when the stimulation stops). The black square is the zoomed region in panel c. c, The mean time spent (colorbar) for each protocol in panel b (100 samples for each protocol for panelsc, f and i). d, Same as in panel a, but for the LTP protocol, 1Pre2Post50 and, the NC protocol, 2Post1Pre50. Both are composed of 300 pairing repetitions at 5 Hz. e, Same as in panel b for protocols in panel d. f, The mean time spent (colorbar) for each protocol in panel e. g, Same as in panel a and d, but for two protocols with different frequencies and pulse repetitions. The LTD protocol, 2Pre50 900 at 3 Hz and, the NC protocol 2Pre10 300 at 5 Hz. h, Same as in panel b and e for protocols in panel g. i, The mean time spent (colorbar) for each protocol in panel h.

Fig. 3 .

 3 Fig. 3. | Read-out strategy an Tigaret et al. 2016 24 experiment. a, Illustration of the joint CaMKII and CaN activities crossing the plasticity regions. Arrows indicate the flow of time, starting at the black dot. Note that here time is hidden and one can only see the changes in enzyme concentrations. b, Region indicator showing when CaN and CaMKII crosses the LTD or LTP regions in panel a. Leaving the region activates a leaking mechanism that keeps track of the accumulated time inside the region. Such leaking mechanism drives the transition rates used to predict plasticity (Methods). c, Plasticity Markov chain with three states: LTD, LTP and NC. There are only two transition rates which are functions of the plasticity region indicator (Methods). The LTP transition is fast whereas the transition LTD is slow, meaning that LTD change requires longer time inside the LTD region (panel a). The NC state starts with 100 processes. d, Joint CaMKII and CaN activity for all Tigaret protocols (labelled in F). The stimulus ends when the trajectory becomes smooth. Corresponds to Fig. 2b,e and h, at 60 s. e, Region indicator for the protocols labelled in F.The upper square bumps are caused by the protocol crossing the LTP region, the lower square bumps when the protocol crosses the LTD region (as in panel d). f, Synaptic weight (%) as function of time for each protocol. The weight change is defined as the number (out of 100) of states in the LTP state minus the number of states in the LTD state (panel c). The trajectories correspond to the median of the simulations in panel g. g, Synaptic weight change (%) predicted by the model compared to data (EPSC amplitudes) from24 (100 samples for each protocol, also for panel h and i). The data (grey dots) was provided by24 (note an 230% outlier as the red asterisk), red bands indicate data means. h, Predicted mean synaptic weight change (%) as function of delay (ms) and number of pairing repetitions (pulses) for the protocol 1Pre2Post(delay), for delay between -100 and 100. LTD is induced by 2Post1Pre50 after at least 500 pulses. The mean weight change along each dashed line is reported in the STDP curves in panel i. i, Synaptic weight change (%) as function of pre-post delay. Each plot corresponds to a different pairing repetition number (legend). The solid line shows the mean, and the ribbons are the 2nd and 4th quantiles. The red dots are the data means estimated in24 , also shown in panel g.

Fig. 4 .

 4 Fig. 4. | Frequency dependent plasticity, Dudek and Bear 1992 9 dataset. a, Example traces of joint CaMKII-CaN activity for each of 9 's protocols. b, Region indicator showing when the joint CaMKII-CaN activity crosses the LTD or LTP regions for each protocol in panel a. c, Synaptic weight change (%) as function of time for each protocol, analogous to Fig. 3c. Trace colours correspond to panel a. The trajectories displayed were chosen to match the medians in panel e. d, Mean (100 samples) time spent (s) for protocols 1Pre for 900 pairing repetitions at 3, 10 and 50 Hz. e, Comparison between data from 9 and our model (1Pre 900p, 300 samples per frequency, TableM1). Data are represented as normal distributions with the mean and variance extracted from9 . Note that data from Dudek and Bear are given in field EPSP slope change. f, Prediction for the mean weight change (%) varying the stimulation frequency and pulse number (24x38x100 data points, respectively pulse x frequency x samples). The red dots show the 9 protocol parameters, the corresponding results are shown in panel e.

46 74CHAPTER 5 .

 465 Fig.7. | Jitter and spike dropping effects on STDP; and Poisson spike trains. a, Mean weight (%) for the jittered STDP protocols. The solid line is the mean, and the ribbons are the 2nd and 4th quantiles predicted by our model estimated using 100 samples (same for all panels). b, Mean weight (%) for the same24 protocols used at panel a subjected to random spike removal (sparsity %). c, Mean time spent (s) varying jittering. Contour plot shows 2Post1Pre50 and 1Pre1Post10 (300 at 5 Hz) with (grey contour plot) and with jittering (coloured contour plot). The circles and squares correspond to the marks in panel a. d, Mean time spent (s) varying sparsity. Contour plot in grey showing 0% sparsity for 2Post1Pre50 300 at 5Hz (see Fig.2f). The contour plots show the protocol with spike removal sparsities at 0% (NC) 30% (LTP), and 80% (NC). The triangles correspond to the same marks in panel a. e, Distribution of the 50 ms jittering applied to the causal protocol 1Pre1Post10, 300 at 5 Hz in which nearly half of the pairs turned into anticausal. The mean frequency is 5 ± 13.5 Hz. The protocol 2Post1Pre50 will have nearly half of the pairings turning into causal, making them have a similar firing structure and position inside the LTP region. f, Mean weight change (%) combining both sparsity (panel b) and jittering (panel a) for 2Post1Pre50, 300 at 5 Hz. g, Mean weight change (%) of pre and postsynaptic Poisson spike train delivered simultaneously for 10 s. The plot shows the plasticity outcome for different presynaptic firing rate (1000/frequency) for a fixed postsynaptic baseline at 10Hz. The upper raster plot depicts the released vesicles at 40 Hz and the postsynaptic baseline at 10Hz (including the AP evoked by EPSP). h), Mean weight change (%) varying the rate of pre and postsynaptic Poisson spike train delivered simultaneously for 10 s. The dashed line depicts panel g. 10/46

  (age, temperature, [Ca 2+ ] o , [Mg 2+ ] o ) with our NMDAr model, the main postsynaptic calcium source. -activation curve h adjusted to different [Ca 2+ ] o release probability sigmoid slope s = 2 adjusted to different [Ca 2+ ] o Vesicle release model (deterministic part) Ca pre attenuation recovery τ pre = 20 ms 50 -500 ms for with dye 17 therefore < 50 to 500 ms undyed (unbufered) deterministic jump attenuation recovery τ rec = 20 s ∼ 20 s 9 deterministic jump attenuation fraction δ ca = .0004 inactivation of pre calcium 18

Fig

  Fig.M1. | Presynaptic release. a, Presynaptic calcium in response to the protocol 1Pre, 300 at 5 Hz displaying adaptation. b, Release probability for the same protocol as panel A but subjected to the docked vesicles availability. c, Number of vesicles in the docked and reserve pools under depletion caused by the stimulation from panel A. d, Plot of the mean (300 samples) release probability (%) for different frequencies for the protocol 1Pre 300 pulses at [Ca 2+ ] o = 2.5 mM. Note that most of the frequencies are dominated by short-term depression, and the model also displays short-term facilitation (black curve, at 50 Hz). e, Release probability (%) for a single presynaptic spike given the [Ca 2+ ] o . Note that 14 model was multiplied by the experimentally measured release probability at[Ca 2+ ] o = 2 mM since their model has this calcium concentration as the baseline. Our model also does not cover the abolishing of release probability at [Ca 2+ ] o = 0.5 mM which can also be difficult to experimentally measure given the rarity of events13 .

84CHAPTER 5 .

 5 A STOCHASTIC MODEL OF HIPPOCAMPAL SYNAPTIC PLASTICITY WITH GEOMETRICAL READOUT OF ENZYMES -ARTICLE

Fig

  Fig.M2. | AP Evoked by EPSP. a, Model and data comparison for the distance-dependent BaP amplitude attenuation measured in the dendrite and varying the distance from the soma. The stimulation in panel A is set to reproduce the same stimulation as24 . Golding measurements have neurons that are strongly attenuated and weakly attenuated (dichotomy mark represented by the dashed line). However, in this work, we consider only strongly attenuated neurons. b, Attenuation of somatic action potential from37 and model in response to five postsynaptic spikes delivered at 100 Hz. The value showed for the model is the spine voltage with distance from the soma set to zero (scale 25 ms, 20 mV).c, Top panel shows the β soma used in Equation (5) to modify the axial conductance between the soma and dendrite. Bottom panel shows the age-dependent changes in the step of the resource-use equation, in Equation (6) that accelerates the BaP attenuation and decreases the sodium currents in the Equation (4). d, Probability of evoking an AP multiplied by the successfully evoked AP (p AP (V evoked ) • 1(evoked) for the protocol 1Pre, 300 at 5 Hz (2.5 mM Ca). d, Two-pool dynamics with the stimulation than panel D showing the vesicle release, the reserve and docked pools, and the evoked AP. e, Probability of evoking an AP for the protocol 1Pre 300 pulses at different frequencies (3 and 5 Hz have the same probability).

23 delay

 23 AP evoked by EPSP δ delay-AP = 15 ms 41 Table

Fig. M3 .

 M3 Fig. M3. AMPAr Markov chain with three sub-conductance states and two desensitisation levels. It includes parameters ρ AMPA f

  is made by changing the forward ρ AMPA f and backward ρ AMPA b rates in a temperaturedependent manner matching the decay time reported by 46 : -0.367•(T -28.976°C) .

Fig

  Fig. M4. | Effect of temperature in the AMPAr. a, Probability of AMPAr opening ( O2+O3+O4 N AMPA

Fig

  Fig.M5. | NMDAr changes caused by age, temperature and extracellular and magnesium concentrations in the aCSF.A Decay time of the NMDAr-mediated EPSP recorded from neocortical layer II/III pyramidal neurons (grey)57 compared to the decay time from the GluN2B channel estimated by our model (yellow) and data from Iacobussi's single receptor recording (purple)55 . b, Comparison of our implementation of GluN2B:GluN2A ratio and the GluN2B:GluN2A ratio from the mice's CA1. c, Comparison of our implementation of NMDAr conductance change in response to the extracellular against data58 . d, Forward and backwards temperature factors implemented to approximate NMDAr subtypes decay times at room temperature55 and temperature changes observed in57 . e, NMDAr subtypes number on our model given age. We add noise to have a smoother transition between different ages. f, Calcium concentration changes for causal and anticausal protocols in response to different aCSF calcium and magnesium compositions with fixed Ca/Mg ratio(1.5). Scale 50 ms and 5 µM.

  using a logistic functions ρ NMDA f and ρ NMDA b . The decay time comparison is shown in Fig. M5A. Then, we applied the same temperature factor ρ NMDA f and ρ NMDA b for GluN2A. The decay times of GluN2A and GluN2B are similar to the ones reported by Iacobucci 55 . The forward and backward factors are described as follows: ρ NMDA f = -1230.680 + 1239.067 1 + e -0.099•(T +37.631°C) , ρ NMDA b = 3.036 + 1621.616 1 + e -0.106•(T -98.999°C) .

Fig

  Fig. M6. | GABAr Markov chain model. Closed states (C 0 , C 1 and C 2 ) open in response to GABAr and can go either close again or open (O 1 and O 2 )

Fig. M8 .

 M8 Fig. M8. From left to right, R-, L-, and T-type VGCCs Markov chain adapted from Magee and Johnston 1995 70 . The R-(left scheme) and T-type (right scheme) have a single open state (red colour), respectively, O r and O T . The L-type VGCC (middle) has two open states, O L1 and O L2 .

Fig

  Fig. M9. | VGCC rates and temperature factors. a, Activation (α m (V sp )) and deactivation rates (β m (V sp )) for the T-type m-gate. b, Activation (α m (V sp )) and deactivation rates (β m ) for the R-type m-gate. c, Activation (α m (V sp )) and both deactivation rates (β L2 (V sp ) and β 1 2 (V sp )) for the L-type VGCC. d, Activation (α h (V sp )) and deactivation rates (β h (V sp )) for the T-type h-gate. e, Activation (α h (V sp )) and deactivation rates (β h (V sp )) for the R-type h-gate. f, Temperature factor applied to all the rates, forward change (ρ V GCC

  ρ SK b = 149.37 -147.61 1 + e 0.093•(T -98.85C) , ρ SK f = 0.005 + 2.205 1 + e -0.334•(T +25.59C)

  503 -0.304 1 + e1.048(T -30.668°C)

Fig

  Fig. M10. | Coarse-grained model of CaM, CaMKII and CaN adapted from77 and78 The reaction description matches with the color: Releases 2Ca, consumes 2Ca, consumes mKCaM, releases mKCaM, releases CaM2C, CaM2N, CaM0, CaM2N, releases mCaN, consumes mCaN, phoshorylate K units to P units, phosphorylated states and dephosphorylation.

Fig

  Fig. M11. | CaMKII temperature changes in our model caused by 1Pre, 30 at 0.49 Hz with Glutamate uncaging (no failures allowed), 1Mm Ca, 2mM Mg, P4-7 organotypic slices from mice hippocampus. a, CaMKII fluorescent probe lifetime change measured by 53 for 25°(blue) and 35°C (red) 53 . The decay time (τ) was estimated by fitting the decay after the stimulation (30 pulses at 0.49Hz) using a single exponential decay, y = a • e -t•b ; τ = 1\b. b, Simulation of the CaMKII concentration change (with respect to the baseline) at 25°in response to same protocol applied in the panel A. The simulations on the panels B, C, E and F show the mean of 20 samples. c, Same as in panel B but for 35°C. d, Estimated temperature change factor for the dissociation rates k 2 , k 3 and k 5 in the Markov chain at the Fig. M10. e, Change in the concentration of the CaMKII states (25°C) which are summed to compose CaMKII change in the panel B. f, Same as in panel E for 35°C with reference to the panel C.

Fig

  Fig. M12. | CaN temperature changes in our model caused by 1Pre, 100 at 20 Hz with Glutamate uncaging (no failures allowed), 2Mm Ca, Free Mg, 11-13 days in vitro. a, Simulated caN change (blue solid line) in response to the same stimuli of the CaN measurement from 85 RY-CaN fluorescent probe (green solid line). The decay time (τ) estimated from data (y = a • e -t•b ) is 94.83 s (dashed purple line) and for our model (solid purple line) is 82.66 s. b, Simulated CaN change for physiological temperature with decay time 54.44 s. Due to the lack of data, CaN kinetic change was set to fit plasticity on the protocols use in this work. c, Temperature change, ρ CaN f and ρ CaN b

  rate adaptation for the coarse-grained model adapt(a, b, c, d,Ca) = a•b c+d•Ca .

Fig

  Fig.M14. | Plasticity readout for the protocol 1Pre2Post10, 300 at 5Hz, from1 . a, CaMKII and CaN activity in response to protocol 1Pre2Post10. b, Enzymatic joint activity in the 2D plane showing LTP and LTD's plasticity regions. The black point marks the beginning of the stimulation, and the white point shows the end of the stimulation after 60 s. c, Region indicator illustrating how the joint activity crosses the LTP and the LTD regions. d, The leaky activation functions are used respectively as input to the LTP and LTD rates. The activation function has a constant rise when the joint-activity is inside the region, and exponential decay when it is out. e, The LTD rate in response to the leaky activation function, act D , in panel D. Note that this rate profile occurs after the stimulation is finished (60 s). The joint-activity is returning to the resting concentration in panel A. f, The LTP rate in response to the leaky activation function, act P , in panel D. g, Outcome of the plasticity Markov chain in response to the LTD and LTP rates. The EPSP change (%) is estimated by the difference between the number of processes in the states LTP and LTD, LT P -LT D. h, Normalized LTP and LTD rates (multiplied to their respective time constant, t D , t P ) sigmoids. The dashed line represents the half-activation curve for the LTP and LTD rates. Note in panel D that the leaky activation function reaches the half-activation K p = 1.3e4.

Fig

  Fig.S2shows variations of4 parameters for[Ca 2+ ] o , [Mg 2+ ] o , temperature and dendritic spine distance from the soma. Also, it shows the Poisson spike train protocol (as in Fig.G and H) for temperature and age parameters obtained from an estimation of the body temperature regulation during development (or thermoregulation maturation, also called maturation of temperature homeostasis, estimated in Fig.S5G).

Fig

  Fig. S2. | Varying experimental parameters in 4 and Poisson spike train during development. Related to Fig. and . a, Mean synaptic weight change for the FDP experiment varying the [Mg 2+ ] o , original [Mg 2+ ] o = 1.5 mM (dashed grey line). b, Mean synaptic weight change for the FDP experiment varying the [Ca 2+ ] o , original [Ca 2+ ] o = 2.5 mM (dashed grey line). c, Mean synaptic weight change for the FDP experiment varying the distant from the soma, original 200 µm (dashed grey line). Changing the distance from the soma modifies how fast BaPs evoked by EPSP will attenuate. Note that LTD is prevalent for a spine situated far from the soma. This could justify why spines distant from the soma are smaller in size since distance correlates with synaptic weight. d, Mean synaptic weight change for the FDP experiment varying the temperature, original temperature 35°C (dashed grey line). e, Mean synaptic weight change for the FDP experiment varying the pairing repetitions at 33°C showing how LTD is enhanced. f, Mean synaptic weight change for the FDP experiment varying the pairing repetitions at 37°C showing how LTD is abolished. g, Mean synaptic weight change for pre and postsynaptic Poisson spike train during 30 s for P5 and 34°C. The panel shows that there is weak and diffused LTP. h, Mean synaptic weight change for pre and postsynaptic Poisson spike train during 30 s for P15 and 35°C. The panel shows that there is a start of LTP window forming for slow postsynaptic rates (<1 Hz). i, Mean synaptic weight change for pre and postsynaptic Poisson spike train during 30 s for P20 and 35°C. The panel shows that a window forms around 10 Hz postsynaptic rate similar to what is shown by 5 and in Fig. H.

108CHAPTER 5 .

 5 Fig.S4presents modifications of 7 's STDP experiment and the reproduction of 3 data.

Fig

  Fig. S4. | [Ca 2+ ] o and [Mg 2+ ] o related modifications for 7 's experiment. Related to Fig. . a, Mean time spent for anticausal pairing, 1Post1Pre10, at different Ca/Mg concentrations. The contour plots are associated with the Fig. A, B and C. b, STDP and extracellular Ca/Mg. Synaptic weight change (%) for causal (1Pre1Post10, 100 at 0.3 Hz) and anticausal (1Post1Pre10, 150 at 0.3 Hz) pairings varying [Ca 2+ ] o from 1.0 to 3 mM (Ca/Mg ratio = 1.5). c, Varying frequency and extracellular Ca/Mg for the causal pairing 1Pre1Post10, 100 at 0.3 Hz. Synaptic weight change (%) for a single causal pairing protocol varying frequency from 0.1 to 10 Hz. [Ca 2+ ] o was fixed at 1.8 mM (Ca/Mg ratio = 1.5). d, Mean synaptic weight change (%) for 7 's STDP experiment showing how temperature qualitatively modifies plasticity. The dashed lines are ploted in panel B. e, Mean synaptic weight change (%) showing effects 0.5°C from panel A. Black and grey solid lines represent the same color dashed lines in panel A (30 and 30.5°C). The bidirectional curves, black and grey lines in panel A (dashed) and panel B (solid), becoming full-LTD when temperature increases to 34.5 and 35°C, respectively yellow and purple lines in panel A (dashed) and panel B (solid). Further increase abolishes plasticity. f, Mean synaptic weight change (%) for 3 's experiment in Free-Mg ([Mg 2+ ] o = 10 -3 mM for best fit) showing the different time requirements to induce LTP and LTD. For LTD, to simulate the NMDAr antagonist D-AP5 which causes a NMDAr partial blocking we reduced the NMDAr conductance by 97%. Note the similarity with Fig. S1F. g, Mean synaptic weight change (%) of 7 's STDP experiment changing [Ca 2+ ] o and Ca/Mg ratio. h, Mean synaptic weight change (%) of 7 's STDP experiment changing pre-post delay time and frequency. Note the similarity with Fig. S1C. i, Mean synaptic weight change (%) of 7 's STDP experiment changing pre-post delay time and age. Age has a weak effect on this experiment done at [Ca 2+ ] o = 2.5 mM.

Fig

  Fig. S5. | Age and temperature effects. Related to Fig. and . a, Mean synaptic weight change (%) for 11 's STDP experiment for 1Pre1Post10, 70-100 at 5 Hz (see TableM1) showing a full LTD window. Our model also reproduces the fact that increasing temperature to 32-34°C the LTD is abolished (data not shown). b, Mean synaptic weight change (%) for11 's STDP experiment for 1Pre2Post10, 70-100 at 5 Hz (see TableM1) showing a bidirectional window. c, Mean synaptic weight change (%) for 11 's STDP experiment for 1Pre2Post10, 20-30 at 5 Hz (see TableM1) showing a bidirectional window. We report that for11 experiment done in room temperature the temperature sensitivity was higher than other experiments. d, Core temperature varying with age representing the thermoregulation maturation. This function (not shown) was fitted using rat10 and mouse data9 added by 1°C to compensate species differences10 . The blue and white bars represent the circadian rhythm as shown in9 . However, the "rest rhythm" for young rats (P5-14) may vary. e, Plot showing how far from being physiological are plasticity experiments done in physiological temperatures. Suggesting, there is scarcity of physiologically relevant data to model and understand plasticity. The dashed grey line is an approximation of the mean value from panel G. f, Initial conditions for CaN-CaMKII resting concentration for different[Ca 2+ ] o and temperature values. When [Ca 2+ ] o is changed temperature is fixed at 35°C, while temperature is changed [Ca 2+ ] o is fixed at 2 mM.

Fujii

  et al. 2013 claims that these enzymes behave as a nonlinear decoder of frequency and pulse number (Figure 6.1D-F). They show that somatic CaMKII is sensitive to frequency, pulse number, and CaN only to the pulse by measuring enzymatic integral and amplitude. Using the model, we can test if such nonlinear filter observed in the soma of cultured neurons by Fujii et al. 2013 could be found in dendritic spines. For this, the same experimental conditions and stimulation protocol mentioned in the paper were used: glutamate uncaging (no failures), 25 • C, 2mM [Ca 2+ ] o free-Mg 2+ , P13 chosen arbitrarily since Fujii et al. 2013 uses cultured neurons, 200 µm from the soma, allowing BaP elicited by EPSPs. Figure 6.1A-C shows that the dendritic spine version of the nonlinear filter does not follow Fujii et al. 2013 findings.

6. 1 . 2

 12 The effects of the temperature range in the FDP Focusing on CaN and CaMKII contributions to plasticity, O'Connor et al. 2005 druginhibited CaMKII and CaN. They show that inhibiting CaMKII abolishes LTP, and inhibiting CaN abolishes LTD (Figure 6.2C). However, these molecules are not exclusively controlling LTP and LTD as mentioned by Coultrap and Bayer 2012; Kipanyula et al. 2016. Here, we would like to test how the wide temperature range affects qualitative changes in O'Connor et al. 2005 FDP ex vivo experiment. Since

Figure 6 . 1 :Figure 6 . 2 :

 6162 Figure 6.1: Synaptic version CaMKII and CaN nonlinear filter using the same experimental conditions estimated by the model (panels A-C) and Fujii et al. 2013 measurements (panels D-F). A) CaMKII maximal amplitude (normalized concentration) using different frequencies and pulse numbers. Low frequencies had maximal activation when compared to the data measures in the soma (panel D). B) CaN maximal amplitude (normalized values) using different frequencies and pulse numbers. The CaN profile was similar to what was measured in the soma (panel E). C) Difference between normalized CaN and CaMKII showing an inverted dominance to what was observed in the soma (panel F). D) Fujii et al. 2013 CaMKII amplitude measurement using the fluorescence change obtained in a FRET probe (also D-E). E) Fujii et al. 2013 CaN amplitude measurement using a fluorescent probe. F) Nonlinear filter found by Fujii et al. 2013 in the soma in which CaMKII is more active for high frequencies and pulses, and CaN for low frequencies.

Figure 6 . 3 :

 63 Figure 6.3: Comparison between measures derived from the simulated dynamics for the Tigaret et al. 2016 data. The grey region indicates an overlap in between different classes of plasticity outcomes (LTP, LTD and NC). A) Scatter plot between CaN and CaMKII integral. The axis was divided by 1e6 only to see axis ticks better (also in panel C). B) Scatter plot between CaMKII and CaN amplitude. C) Scatter plot between CaM amplitude and Ca integral.

Figure 6 . 4 :

 64 Figure 6.4: Comparison of stochastic and deterministic versions of the model. And a comparison between discrete and averaged release.. A) Left, discrete release (green) and averaged release. Right, enzymes activation using the different approaches. B) Comparison between a fully stochastic model (with discrete release), a deterministic version (with discrete release), and a deterministic version (with averaged release). C) Enzymatic baseline differences between stochastic and deterministic model versions given the external Ca 2+ .

Figure 6 . 5 :

 65 Figure 6.5: Dye implementation comparison. A) Tigaret Ca 2+ measurements of according to the protocols on top th row. B) Ca 2+ with dye, temperature increased and SK conductance set to zero. C) Ca 2+ without the dye.

Figure 6 . 6 :

 66 Figure 6.6: Model hypothesis for sleep under different age and temperature conditions. A) Stimulation pattern representing delta waves during sleep. B) Approximation of temperature oscillation during development based on Wood et al. 2016; McCauley et al. 2020. C) Mean weight change using the pattern in panel for different ages and temperatures.

Figure 6 . 7 :

 67 Figure 6.7: Behavioral time scale synaptic plasticity (BTSP). Plasticity states was changed set to 300. Bittner's parameters used to tune the model were 35 • C, 2 mM [Ca 2+ ]o, 1 mM [Mg 2+ ]o (0.7 mM gave better fit), P50, 200 µm from the soma. A) The geometrical threshold using the protocol from Bittner et al. 2017. B) A comparison between the model and data from the BTSP curve discovered by Bittner et al. 2017.

Table 3 .1: Plasticity experiments in the cortex (slices) reproduced by the models in Table 3.5.

 3 such as the astrocyte ion buffering properties[START_REF] Serrano-Gotarredona | STDP and STDP variations with memristors for spiking neuromorphic learning systems[END_REF]).

	L5 or L2/3 cortex -brain slices
	EXPERIMENT	PAPER	REP.	FREQ (Hz) AGE (DAYS) TEMP. ( • C) Ca (mM) Mg (mM)
	STDP + freq	Sjöström et al. 2001	1 min	0.1 -50	12 -21	32 -34	2.5	1
	STDP + freq	Sjöström et al. 2001	1 min	0.1 -50	12 -21	32 -34	2.5	1
	sub. t. depolarization
	1Pre1Post10 -15x
	L5-L5, L2/3-L5	Sjöström and Häusser 2006* -	50	14 -21	32 -35	2	1
	soma distances
	1Pre1Post10
	L5-L5, L2/3-L5	Sjöström and Häusser 2006* -	50	14 -21	32 -35	2	1
	soma distances
	sub. t depolarization
	STDP	Froemke et al. 2006	60 -100	0.2	10 -35	room	2.5	1.5
	STDP	Froemke et al. 2006	60 -100	0.2	10 -35	room	4	4
	bicuculine meth.
	5x 1Pre1Post10
	(STDP, at a given freq.	Froemke et al. 2006	30 -40	0.2	10 -35	room	4	4
	repeated 30-40 at 0.2 Hz)
	nx 1Pre1Post10
	(STDP, at a given freq.	Froemke et al. 2006	30 -40	0.2	10 -35	room	4	4
	repeated 30-40 at 0.2 Hz)
	xPre+yPost
	(STDP, at a given freq.	Froemke et al. 2006	30 -40	0.2	10 -35	room	4	4
	repeated 30-40 at 0.2 Hz)
	STDP (burst)
	(1Pre3Post ∆ t)	Nevian and Sakmann 2006	60	0.1	13 -15	32 -35	2	1
	Post at 50 Hz
	STDP (burst frequency)
	(1Pre3Post10 and 1Pre3Post-10)	Nevian and Sakmann 2006	60	0.1	13 -15	32 -35	2	1
	Post at X Hz
	STDP (burst frequency)
	(1PreXPost10 and 1PreXPost-10)	Nevian and Sakmann 2006	60	0.1	13 -15	32 -35	2	1
	Post at 50 Hz
	STDP (1Pre1Post)
	L2/3-L5	Letzkus et al. 2006	100 -200	1	21 -42	34 -35	2	1
	soma distances
	STDP (1Pre3Post)
	L2/3-L5	Letzkus et al. 2006	100 -200	1	21 -42	34 -35	2	1
	soma distances
	STDP	Shen et al. 2008*	10-15 (x5) 0.1	19 -26	room	2	1
	corticostriatal
	STDP	Pawlak and Kerr 2008	60	0.1	19 -22	31 -33	2.5	2.0
	corticostriatal
	STDP	Fino et al. 2010	100	1	15 -21	34	2	1
	corticostriatal

Table 3 .2: Plasticity experiments in cultured hippocampal neurons reproduced by the models in Table 3.5.

 3

	Hippocampus -brain slices
	EXPERIMENT	PAPER	REP.	FREQ (Hz)	AGE (DAYS)	TEMP. ( • C) Ca (mM) Mg (mM)
	STDP (warm) 1Pre1Post ∆ t	Wittenberg and Wang 2006 70 -100	5	14 -21	-34	2	1
	STDP (cold) 1Pre1Post ∆ t	Wittenberg and Wang 2006 70 -100	5	14 -21	-30	2	1
	STDP (cold) 1Pre2Post ∆ t	Wittenberg and Wang 2006 70 -100	5	14 -21	-30	2	1
	STDP (cold) 1Pre2Post ∆ t	Wittenberg and Wang 2006 20 -30	5	14 -21	-30	2	1
	STDP (burst)
	1Pre1Post ∆ t	Tigaret et al. 2016	300	5	50 -55		2.5	1.3
	picrotoxin
	STDP (single)
	1Pre2Post ∆ t	Tigaret et al. 2016	300	5	50 -55		2.5	1.3
	picrotoxin
	Doublets (duration)
	2Pre50, 2Pre10	Tigaret et al. 2016	300 and 900	5 and 3	50-55		2.5	1.3
	picrotoxin
	STDP (single, Ca) 1Pre1Post ∆ t	Inglebert et al. 2020	150, < delays 100, > delays	0.3	14 -20		1.3 -3	Ca/1.5
	STDP (single)+freq+Ca 2+	Inglebert et al. 2020	100, > delays	0.3 -10	14 -20		1.3, 1.8	Ca/1.5
	1Pre1Post10, X Hz		150, < delays
	STDP (burst)	Inglebert et al. 2020	100	0.3	14 -20		1.3, 1.8	Ca/1.5
	1PreXPost10, X=2..4
	STDP (single+age) 1Pre1Post ∆ t	Meredith et al. 2003	20	0.2	9 -45	-28	2	2
	STDP (burst+age) 1Pre2Post ∆ t	Meredith et al. 2003	20	0.2	9 -45	-28	2	2
	FDP	Dudek and Bear 1992	900	1 -50	∼35		2.5	1.5
	FDP+age	Dudek and Bear 1993	900	1	7 -35		2.5	1.5
	TBS+age	Dudek and Bear 1993	3 -4 ( 5) epochs	4Pre at 100 Hz	6, 14 and 17		2.5	1.5
				(10x at 5Hz)
	LFS (duration)	Mizuno et al. 2001	1 -600	1	12 -28		2.4	0
	LFS (duration)	Mizuno et al. 2001	1 -600	1	12 -28		2.4	0
	NMDA blocker
	dendritic cooperativity	Weber et al. 2016	50	3	49 -77	-35	1.25	1.3
	4 spines
	dendritic cooperativity
	4 spines	Weber et al. 2016	50	3	49 -77	-35	1.25	0.1
	low-Mg
	STDP
	astrocyte activity	Bonansco et al. 2011	60	1	25 -35	-22	2.5	2
	+picrotoxin

Table 3 .3: Plasticity experiments in the hippocampus (slices) reproduced by the models in Table 3.5.

 3

Table 3 .

 3

	reproduces frequency-dependent (quali) no FDP			single delay STDP no STDP			Poisson firing absent prediction without data				averaged release probability no presynaptic control			point neuron no morphology any morphological aspect		phen. excitability control by glia No glia		a single brain region Not specified one method (e.g. CA1 accute slices)		integrate-and-fire type no voltage			ionotropic receptors/ion channels only not specified			single Ca 2+ source no Ca 2+ (without specified source)			phenomenological	implementation of molecules absent	(e.g. unknown molecules)
	reproduces frequency-dependent data			single paired STDP	(multiple delay times)			Poisson firing, replay data	prediction without data				STP, presynaptic voltage or Ca 2+	discrete release			1D compartments	ball-and-stick	phen. glia control over aCSF	gliotr.	single region	more than one method		passive dendrite (if present)	at least K+ and Na+ channels			ionotropic	dose-dependent control			generic NMDA	generic VGCC				data-driven and phenomenological	implementation of molecules
	brain region)	rate stimulation (e.g. Poisson)	other direct stimulation experiment	burst STDP xPre(x1)-yPost(y1)-delay,	single paired STDP	STDP variations (e.g. ext. Ca 2+ , distance)	STDP jittered	Poisson firing, non-Poisson firing	replay data, in-vivo like	(e.g. sleep wave ripples)	regular with background noise	with data	presynaptic signaling (e.g. endocanabinoid, arc)	preVGCCs, preNMDA, preBuffers	presynaptic LTD/LTP	vesicles pool dynamics	STP (e.g. augmentation, facilitation)	discrete release	3D reconstructions	spines, dendrites, soma, axon	ion-buffering, aCSF filtrate, gliotransmission	recycling neurotransmission	Ca 2+ waves	slices, cultured, In-vivo	more than one brain region	active dendrite	differential channels expression	inhibitory/excitatory input integration	gap junctions, attenuation of voltage	ionotropic	metabotropic	dose-dependence	neuromodulators	extracelullar Ca 2+ , Pumps, dye inclusion	buffering	permeable AMPA	permeable Na+, NMDA subtypes, VGCC subtypes	ER	diffusion	interactions with receptors, ion channels	transporters (e.g. GLT-1)	turnover	messengers (e.g. g-protein)	buffering
					STDP				irregular firing structures		LOCALIZATION			presynaptic			morphology		glia		sampling methods	FAST AND SLOW VAR.	membrane potential			neurotransmitters & receptors				Ca 2+ sources						molecular dynamics

4: Table classifying different attributes of plasticity models.

Table 3 .

 3

5: Table classifying different attributes of plasticity models using the guide in the Table 3.5.

  Several models explored the combined role of CaMKII and CaN (or phosphatases and kinases) and their role in plasticity. For instance, Antunes et al. 2016b models AMPAr exocytosis and endocytosis triggered the enzymatic dynamics of CaN and CaMKII. CaMKII implementations are seen with different levels of detail, such as Pharris et al. 2019 which tries to model the combinatorial explosion of CaMKII states recently unachievable. Also, Ebner et al. 2019 implements a phenomenological version of the enzyme dynamics for pre and postsynaptic long-term plasticity using Clopath and Gerstner 2010 low-pass strategy to predict plasticity.

. For instance, Lisman 1989's model provided testable predictions for two problems simultaneously, memory storage and plasticity induction. His hypothesis that CaMKII would control memory maintenance was based on the dynamics of CaN and CaMKII. The former surpassing protein turnover through its self-sustained activity (illustrated by models

Miller et al. 2005; Graupner and Brunel 2007)

and could be reversed by CaN. Instead of the maintenance role,

Chang et al. 2017; Chang et al. 2019

suggests CaMKII takes part only in the induction of plasticity. Later, Otmakhov with Lisman and others showed experimental evidence suggesting that any known phosphatase does not dephosphorylate CaMKII

(Otmakhov et al. 2015)

. This empirical evidence was previously predicted by Michalski 2013's model. Regardless of physiological relevance, the concepts developed still influence models seeking different pathways to illustrate the memory storage hypothesis

(Urakubo et al. 2014

).

  . Variable-thresholding in neuron models date back to McCulloch and Pitts 1943. Experimentally, W. Singer in 1977 (as refereed by

  . Later, fast CaMKII measurements in dendritic spines done by Chang et al. 2017 were considered, marking a change of research direction since this new data described temperature effects on CaMKII (see article methods). Also, Chang et al. 2019 used their faster sensor in dendritic spines that differed from Fujii et al. 2013 in terms of CaMKII saturation. That raised the doubt whether the nonlinear filter found by Fujii et al. 2013 is valid for dendritic spines. The comparison of our model based on Chang et al. 2019 using Fujii et al. 2013 readout strategy is shown in Figure 6.1. The closest similarity between these papers was that dendrite measurements had a slower CaMKII saturation (Chang et al. 2017), which could explain Fujii et al. 2013 discrepancy. Chang et al. 2017 (supp. files) measurements displayed a slow decay in Ca 2+ fluorescence. This characteristic slow adaptation also affected CaMKII activation during glutamate uncaging (see article methods Figure

  attenuation was the transient spine volume which slightly decreased Ca 2+ concentration(Holst et al. 2021) and CaMKII activation. However, no data was available to model the transient volume since structural proteins dynamics are usually imaged in a time resolution of minutes or hours[START_REF] Bosch | Structural and molecular remodeling of dendritic spine substructures during long-term potentiation[END_REF].After reproducingTigaret et al. 2016, the challenge was to increase the range5.1. RESEARCH CONTEXT 63 of experiments the model can cover. However, data with similar experimental conditions were scarce. Given that the temperature could change CaMKII activation as shown inChang et al. 2017, the focus was to find the most similar plasticity experiments done at 35 • C.Dudek and Bear 1992 was the best candidate since it had a BCM-like curve and increased the number of LTD experiments to validate the model. After reproducing Dudek and Bear 1992's data, the initial strategy of averages and integrals was discarded since a long plasticity protocol (e.g. LTD) would scaled up the integral value due to noise. That was the motivation to find a readout in which the time would be relevant only when the enzymes reached a specific activation level, like a Goldilocks zone for plasticity.Parallel to adopting another prediction strategy, another experiment from Dudek and Bear 1993 had the same temperature and aCSF, but for different ages, showing that LTD and LTP had different requirements throughout developmental ages. Initially,Dudek and Bear 1993 motivated the inclusion of age control in the model by diversifying the NMDAr subtypes. To this, NMDAr shift was modelled using Sinclair et al. 2016's NMDAr subtypes expression profiles. In this way, we innovate on how to model age in plasticity models by adding a parameter controlling an age continuum, different from other computational neuroscience models (based on electrophysiology), which treated age as a class (e.g. young vs adult). Agerelated NMDAr changes were essential to calibrate the model since Tigaret et al. 2016 used older rats than

  by varying the age of STDP experiments done in different conditions.Effects of extracellular calcium and magnesium concentration on plasticity. The canonical STDP rule7 measured in cultured neurons, high [Ca 2+ ] o , and at room temperature, was recently found not to be reproducible at physiological[Ca 2+ ] o in CA1 brain slices38 . Instead, by varying the [Ca 2+ ] o and [Mg 2+ ] o they found a spectrum of STDP rules with either no plasticity

	or full-LTD for physiological [Ca 2+ ] o conditions ([Ca 2+ ] o < 1.8 mM) and a bidirectional rule for high [Ca 2+ ] o ([Ca 2+ ] o > 2.5 mM), shown in Fig. 6a-c (data).
	We attempted to reproduce 38 's findings by varying [Ca 2+ ] o and [Mg 2+ ] o with the following consequences for the model mechanisms (Methods). On the presynaptic side, [Ca 2+ ] o modu-lates vesicle release probability. On the postsynaptic side, high [Ca 2+ ] o reduces NMDAr conductance 79 , whereas [Mg 2+ ] o affects the NMDAr Mg 2+ block 80 . Furthermore, spine [Ca 2+ ] o influx activates SK channels, which hyperpolarize the membrane and indirectly modulate NMDAr activity

  • 10 -6 nS/µm 2 * see table legend 21 membrane leak conductance (only soma) g soma = 5.31 • 10 -3 nS/µm 2 3 • 10 -4 to 1.3 • 10 -3 nS/µm 222 47 to 2.1 • 10 3 nS (NeuroElectro:CA1) = 8.79 • 10 3 µm 2 π • D dend • L dend dendrite volume Vol dend = 4.4 • 10 3 µm 3 π • (D dend /2) 2 • L dend dendritic membrane capacitance C dend = 52.77 pF C m • A dend dendrite leak reversal potential g leakdend = 3.51 • 10 -2 nS g leak • A dend dendrite axial conductance g di f f = 50 nS R a • A dend Soma MΩ CS neck /(L neck • R a ) 50 to 550 MΩ (275 ± 27 MΩ) 32

		VALUE	REFERENCE
	Passive cable
	leak reversal potential membrane leak conductance (for spine and passive dendrite) g leak = 4 membrane capacitance E leak = -70 mV C m = 6 • 10 -3 pF/µm 2 axial resistivity of cytoplasm R a = 1 • 10 -2 GΩµm Dendrite	69mV from 20 1 • 10 -2 pF/µm 223 17 to 177 pF (NeuroElectro:CA1) 2 • 10 -3 GΩµm 24
	dendrite diameter dendrite length dendrite surface area A soma diameter D dend = 2 µm L dend = 1400 µm D soma = 30 µm soma area (sphere) A soma = 2.82 • 10 3 µm 2 soma membrane capacitance C soma = 16.96 pF soma leaking conductance g leaksoma = 15 nS Dendritic spine	same as 25 apical dendrites, 1200 to 1600 µm 26 21 µm 27 page 3 (4π/3) • (D soma /2) 3 ; 2.12 • 10 3 µm 228 C m • A soma 22 g soma • A soma
	spine head volume spine head surface spine membrane capacitance spine head leak conductance Dendritic spine neck	Vol sp = 0.03 µm 3 A sp = 4.66 • 10 -1 µm 2 C sp = 2.8 • 10 -3 pF g leaksp = 1.86 • 10 -6 nS	same as 29 4π • (3Vol sp /4π) 2/3 C m • A sp g leak • A sp
	spine neck diameter neck length neck cross sectional area neck resistance	D neck = 0.1 µm L neck = 0.2 µm CS neck = 7.85 • 10 -3 µm 2 g	0.05 to 0.6 µm 30 0.7 ± 0.6 µm 31 π • (D neck /2) 2

dend neck = 3.92 nS ≈ 255.1

  bottom), as follows.

	βage =	1 -I age τ age rec	-δ age • β age • I BaP (t),	δ age rec =	1.391 • 10 -4 1 + e 0.135(age-16.482 days) .

  86CHAPTER 5. A STOCHASTIC MODEL OF HIPPOCAMPAL SYNAPTIC PLASTICITY WITH GEOMETRICAL READOUT OF ENZYMES -ARTICLE

	NAME	VALUE	REFERENCE
	Soma parameters for Na+ and K+ channel
	sodium conductance potassium conductance reversal potential sodium reversal potential potassium BaP attenuation parameters	γ Na = 8 • 10 2 nS γ K = 40 nS Erev Na = 50 mV Erev K = -90 mV	0.32 nS/µm 233 see legend commentary 0.48 nS/µm 233 see legend commentary 33 33
	attenuation step factor (age)	δ age	see Equation (6) and Fig. M2B and C bottom 24,37
	attenuation step factor auxiliary attenuation step factor recovery time for the attenuation factor	δ decay = 1.727 • 10 -5 δ aux = 2.304 • 10 -5 τ rec = 2 s	adjusted to fit 24,37 adjusted to fit 24,37 adjusted to fit 24,37
	recovery time for the auxiliary attenuation factor τ aux rec = 2 s	adjusted to fit 24,37
	recovery time for the age attenuation factor	τ age rec = 0.5 s	adjusted to fit 24,37
	AP evoked by EPSP
	decay time for V evoke	τ

V = 40 ms

  88CHAPTER 5. A STOCHASTIC MODEL OF HIPPOCAMPAL SYNAPTIC PLASTICITY WITH GEOMETRICAL READOUT OF ENZYMES -ARTICLE

	NAME	VALUE	REFERENCE
	Glutamate parameters
	duration of glutamate in the cleft concentration of glutamate in the cleft glu amp = 1 mM glu width = 1 ms glutamate variability (gamma distribution Γ) glu cv = Γ(1/0.5 2 , 0.5 2 )	47 47 48
	glutamate signal	Glu	glu cv • glu amp for AMPAr, NMDAr and copied to GABA neurotransmitter
	AMPAr parameters
	number of AMPArs reversal potential subconductance O2 subconductance O3 subconductance O4 glu binding	N AMPA = 120 Erev AMPA = 0 mV γ A2 = 15.5 pS γ A3 = 26 pS γ A4 = 36.5 pS k 1	29 29 16.3 pS 45 28.7 pS 45 37.8 pS 45

Table M7 .

 M7 and 16. Calcium baseline concentration Ca ∞ = 50 nM 37 ± 5 to 54 ± 5 nM 17 Calcium decay time τ Ca = 10 ms 50 to 500 ms for with dye 17 therefore < 50 to 500 ms undyed (unbufered) Calcium diffusion D Ca = 0.3338 µm 2 ms -1 0.22 to 0.4 µm 2 ms -129,50 Calcium diffusion time constant τ CaDi f f = Postsynaptic calcium dynamics parameters.

	NAME	VALUE		REFERENCE
	Buffer and dye
	association buffer constant dissociation buffer constant buffer concentration	k Bu f f on k Bu f f o f f = 0.524 ms -1 = 0.247 µM -1 ms -1 Bu f f con = 62 µM	29 29 76.7 µM 29
	Calcium dynamics
	GHK	Vol sp 2D Ca •D neck	+ L neck 2D Ca = 0.5 ms	8 ms for a V sp = 0.7 µm 350
	temperature faraday constant gas constant Calcium permeability Calcium ion valence	T = 35°C F = 96.485 C mol -1 R = 8.314 J K -1 mol -1 P Ca = 0.045 µm ms -1 z Ca = 2	converted to Kelvin in the Equation (12) given the protocol 52 52 adjusted to produce 3 µM Calcium in response to a Glu release supplementary files from 53 52

  92CHAPTER 5. A STOCHASTIC MODEL OF HIPPOCAMPAL SYNAPTIC PLASTICITY WITH GEOMETRICAL READOUT OF ENZYMES -ARTICLE Fig. M7. | GABA(A)r current, kinetics and chloride reversal potential. a, States of GABA(A)r Markov chain at 25°C in response to a presynaptic stimulation.Opened = O 1 + O 2 , closed = C 0 +C 1 +C 2 .b, Model and data comparison65 for GABA(A)r current at 25°C. Even though data were recorded from P70 at 25°C and P15 at 35°C, we normalize the amplitude to invert the polarity and compare the decay time. This is done since the noise around P15 can either make GABAr excitatory or inhibitory as shown by E cl data in panel C. c, Chloride reversal potential (E Cl rev ) fitted to 66 data. Note we used both profiles from axon and dendrite age-depended E Cl rev changes since exclusive dendrite data is scarce. d, States of simulated from GABA(A)r Markov chain at 35°C in response to a presynaptic stimulation. e, Model and data comparison

	28/46

Table M10 .

 M10 Table M10 presents the parameters to model the VGCC channels. CaT = 12 pS same as 70 VGCC R-type conductance γ CaR = 17 pS same as 70 VGCC L-type conductance γ CaL = 27 pS VGCC parameters

	NAME	VALUE	ENCE
	VGCC
	VGCC T-type conductance γ same as 70 number of VGCCs 3 for each subtype 1 to 20 73

Table M11 .

 M11 TableM11presents the parameters to model the SK channel. SK channel parameters.

	96CHAPTER 5. A STOCHASTIC MODEL OF HIPPOCAMPAL SYNAPTIC PLASTICITY WITH GEOMETRICAL READOUT OF ENZYMES -ARTICLE
	NAME	VALUE	ENCE
	SK channel
	number of SK channels N SK = 15 SK conductance γ SK = 10 pS SK reversal potential E SK rev = -90 mV SK half-activation h SK = 0.333 µM SK half-activation slope σ = 6 SK time constant τ SK = 6.3 ms	10-200 75 76 62 62 4 in 62 62
			32/46

  • ms -1 fitted to cover all protocols in Table M1 rise constant inside the LTP region a P = 0.2 a.u. • ms -1 fitted to cover all protocols in Table M1 decay constant outside the LTD region b D = 2 • 10 -5 a.u. • ms -1 fitted to cover all protocols in Table M1 decay constant outside the LTP region b P = 1 • 10 -4 a.u. • ms -1 fitted to cover all protocols in Table M1 Plasticity Markov chain all protocols in Table M1 LTD region -upper diagonal (line equation in the 2D map) CaMKII = +1.64 •CaN + 20.20 fitted to cover all protocols in Table M1 LTD region -lower diagonal (line equation in the 2D map) CaMKII = -5.18 •CaN + 20.91 fitted to cover all protocols in Table M1

	LTD rate time constant LTP rate time constant hill coefficient LTP hill coefficient LTD half occupation LTP half occupation LTD Plasticity region (edges of the polygons)	t D = 1.8 • 10 4 ms t P = 1.3 • 10 4 ms h P = 2 h D = 2 K P = 1.3 • 10 4 a.u. K D = 8 • 10 4 a.u.	fitted to cover all protocols in Table M1 fitted to cover all protocols in Table M1 fitted to cover all protocols in Table M1 fitted to cover all protocols in Table M1 fitted to cover all protocols in Table M1 fitted to cover all protocols in Table M1
	LTP region (CaMKII) -top border	29.5	fitted to cover all protocols in Table M1
	LTP region (CaMKII) -bottom border	1.4	fitted to cover all protocols in Table M1
	LTP region (CaN) -right border	10.	fitted to cover all protocols in Table M1
	LTP region (CaN) -left border	6.35	fitted to cover all protocols in Table M1
	LTD region (CaMKII) -top border	29.5	fitted to cover all protocols in Table M1
	LTD region (CaMKII) -bottom border	1.4	fitted to cover all protocols in Table M1
	LTD region (CaN) -right border	1.85	fitted to cover all protocols in Table M1
	LTD region (CaN) -left border	6.35	fitted to cover

Table M14 .

 M14 Parameters to define the plasticity readout.

  is only to theoretically show how each age mechanism contributes to plasticity in Fig. . Also we compare predictions of between different STDP experiments across age.

	5.5. ARTICLE SUPPLEMENTARY RESULTS
	43/47

  (see Figure 6.7). Since our model has integrative properties granted by slow decaying enzymes, we would like to test if it can reproduce the Bittner et al. 2017 replay experiment. Although our model has not been designed for mice experiments, it can still reproduce

  differ by 0.5 mM [Ca 2+ ] o and are reproduced in Ebner et al. 2019's cortical model without any model adjustment for that. [Ca 2+ ] o is critical since it modifies the release probability and Ca 2+ influx, and a 0.5 mM [Ca 2+ ] o difference can cause qualitative changes in the hippocampus

68CHAPTER 5. A STOCHASTIC MODEL OF HIPPOCAMPAL SYNAPTIC PLASTICITY WITH GEOMETRICAL READOUT OF ENZYMES -ARTICLE

8/4672CHAPTER 5. A STOCHASTIC MODEL OF HIPPOCAMPAL SYNAPTIC PLASTICITY WITH GEOMETRICAL READOUT OF ENZYMES -ARTICLE

11/465.3. ARTICLE

78CHAPTER 5. A STOCHASTIC MODEL OF HIPPOCAMPAL SYNAPTIC PLASTICITY WITH GEOMETRICAL READOUT OF ENZYMES -ARTICLE

16/46 80CHAPTER 5. A STOCHASTIC MODEL OF HIPPOCAMPAL SYNAPTIC PLASTICITY WITH GEOMETRICAL READOUT OF ENZYMES -ARTICLE

5.4. ARTICLE METHODS

82CHAPTER 5. A STOCHASTIC MODEL OF HIPPOCAMPAL SYNAPTIC PLASTICITY WITH GEOMETRICAL READOUT OF ENZYMES -ARTICLE

23/465.4. ARTICLE METHODS

26/4690CHAPTER 5. A STOCHASTIC MODEL OF HIPPOCAMPAL SYNAPTIC PLASTICITY WITH GEOMETRICAL READOUT OF ENZYMES -ARTICLE

98CHAPTER 5. A STOCHASTIC MODEL OF HIPPOCAMPAL SYNAPTIC PLASTICITY WITH GEOMETRICAL READOUT OF ENZYMES -ARTICLE

100CHAPTER 5. A STOCHASTIC MODEL OF HIPPOCAMPAL SYNAPTIC PLASTICITY WITH GEOMETRICAL READOUT OF ENZYMES -ARTICLE

102CHAPTER 5. A STOCHASTIC MODEL OF HIPPOCAMPAL SYNAPTIC PLASTICITY WITH GEOMETRICAL READOUT OF ENZYMES -ARTICLE

104CHAPTER 5. A STOCHASTIC MODEL OF HIPPOCAMPAL SYNAPTIC PLASTICITY WITH GEOMETRICAL READOUT OF ENZYMES -ARTICLE

106CHAPTER 5. A STOCHASTIC MODEL OF HIPPOCAMPAL SYNAPTIC PLASTICITY WITH GEOMETRICAL READOUT OF ENZYMES -ARTICLE

110CHAPTER 5. A STOCHASTIC MODEL OF HIPPOCAMPAL SYNAPTIC PLASTICITY WITH GEOMETRICAL READOUT OF ENZYMES -ARTICLE

5.5. ARTICLE SUPPLEMENTARY RESULTS

The 0.2 mM difference makes the plasticity outcomes at high frequency closer to the data. The model has a strong attenuation to high frequencies due to fast vesicle depletion. The 0.2 mM difference increases Ca

2+ influx and release probability slightly. Other aspects in the age range were not investigated and could account for high frequency outcomes. Another possibility is that the high frequency outcomes LTP are washed out after few minutes given the in a age-dependent manner shown inCao and Harris 2012.

Acknowledgments

Acknowledgments

This work has been partially supported by the French government, through the UCAJEDI Investments in the Future project managed by the National Research Agency (ANR-15-IDEX-01), see ComputaBrain project. Funding was also provided by the UK Medical Research Council (MR/S026630/1) and from the European Union's Horizon 2020 Framework Programme for Research and Innovation under the Specific Grant Agreement No. 945539 (Human Brain Project SGA3). We thank Kim A. Blackwell, Serena M. Dudek, Matthew G. Gold, Jack R. Mellor, L. Mironov, Matthew C. Pharris, Gabriela K. Popescu and Gayle M. Wittenberg for advice and discussions.

Data and Software Availability

The code will be available on Github after peer-review process.

Competing interests

The authors declare no competing interests.

Additional Information

Supplementary information is available for this paper at http Correspondence and requests should be addressed to R.M. and Y.E.R. based on 63,65 Table M9. GABAr parameters.

Readout

We describe the readout mechanism which provides the plasticity event which takes place in the synapse. First, we define the following variables which are representative of "active CaMKII" and "active CaN":

Active CaN

CaN = CaN4

Active CaMKII

It is known that the calcium entry initiates a cascade of events that ultimately leads to short and long term plasticity changes. Specific concentrations of CaMKII and CaN trigger activation functions act D and act P when they belong to one of the two polygonal regions (P and D), termed plasticity regions:

To Specify the LTP/LTD rates, termed D rate and P rate , we use the activation functions, act D and act P , as follows:

The Markov plasticity chain (see Fig. M13) starts with initial conditions NC=100, LTD=0 and LTP=0. Fig. M14 shows how the readout works to predict plasticity for a single orbit. Table M14 shows the parameters to define the polygons of the plasticity regions (see Fig. M14

Supplementary files

Supplemental files present some experiments and predictions extending the notion of parameter sensitivity. Also, they show the effect of modifications in the experimental parameters in Table M1. For instance, Fig. S1 show variations on 1 's experiment.

Simulation of jump processes

This annexe provides information on the stochastic processes used in this thesis.

The model has both stochastic and deterministic parts. That is, the model is deterministic in between stochastic transitions (e.g. opening of AMPAr) satisfying the following flow:

in which x c are the continuous variables evolving deterministically in between discrete events and µ is a set of parameters. For example, continuous variables are x c = (V sp , V dend , V soma , Ca, • • • ) ∈ R c . At discrete times, transitions occurs with total rate R tot (x c , x d , µ) between states of the Markov chains described in the article methods 5.4. All the state variables are collected in x d , termed the discrete variables. For example, the states associated to the receptors and plasticity Markov chain x d = (AM P A, N M DA, V GCC, LT P, LT D • • • ) ∈ R d . This type of Markov process belongs to the class of Piecewise deterministic Markov processes (PDMP) and are extensively described in Davis 2018.

To simulate the synapse model as a stochastic process X, it is relevant to describe when the next jump will occur. The simulation of such process seems tricky at first due to the interdependence between the variables, that is, the x d affects the deterministic part, F (x c , x d , µ), and x c affects the total rate R tot (x c , x d , µ). For one to know when the process X jumps, the flow 1 must be determined. However, to know the flow 1, it is necessary to know when the process jumps, which is equivalent to determine the length of the continuous domain. Fortunately, the simulations of PDMP is well-reviewed, for example, in Graham and Talay 2013. There are two classes of simulation methods. The first one is the most precise and is called the rejection method. However, to work well, one needs to have 133 134 ANNEXE 1. ANNEXE: SIMULATION OF JUMP PROCESSES a tight bound λ on the total rate: ∀x d , x c , R tot (x c , x d , µ) ≤ λ. Given the intricacies of your model, we have not been able to find a good bound λ, meaning that we would have to reject a lot of "simulations", which increases the simulation time. The other alternative for simulation is the true jump method which takes an integral equation to find the jump time.

Hence, for the simulation of our model, we appeal to the tricky twist Veltz 2015 to the "true jump method", which circumvents the need to solve the integral equation by solving an ODE related to the flow 1 but of dimension c + 1. Using this new formulation to find the true jump opens the possibility to test different ODE solvers from package DiffferentialEquations.jl to chose the best one, which was in our case LSODA, that can handle both stiff and non-stiff equations [START_REF] Hindmarsh | LSODA, ordinary differential equation solver for stiff or non-stiff system[END_REF].

The simulation method is thus the following. Given a jump time