The first 100 nm domain, between the plasma membrane and the vesicles in the pre-synaptic terminal remains difficult to study, yet it seems that a displacement as small as a 10 of nanometers in the molecular organization can affect vesicular release. There are many examples, where a 10 nm precision has to be achieved in order to guarantee normal physiology function. This is the case for the apposition of pre- and post-synaptic terminal of neuronal synapses: this apposition is obtained by a set of redundant adhesion molecules, such as laminins that self-organize to maintain the synapse structure and stability in the central nervous system. A lack of the laminin β₂ subunit leads to a disruption of the hippocampal synapse structure, to a misalignment of the pre- and post-synaptic partners and to an increased post-synaptic density (PSD) size (Egles et al., 2007). In addition, mutations in PSD proteins are associated with neurological and psychiatric diseases (Sheng and Kim, 2011). Another example is autism spectrum disorders which have been associated with the mutations in genes encoding Shank2 and Shank3, PSD-93, and a mis-regulation of adhesion molecules neuroligin 3, neuroligin 4, and neurexin 1 (Durand et al., 2011; Sheng and Kim, 2011), affecting a precise geometrical apposition.

In the pre-synaptic terminal, some vesicles are concentrated in a region called the Active Zone (AZ), which is well-aligned with the PSD of the post-synaptic terminal. This apposition creates geometrical columns: a single column alignment was originally hypothesized and numerical simulations showed that it maximizes the synaptic current (Xie et al., 1997) while minimizing its fluctuations (Taflia and Holcman, 2011). Multiple nanocolumns were predicted in Freche et al. (2011) to sustain synaptic response and a more reliable transmission compared to several synapses containing a single column. Finally, these columns have recently been confirmed experimentally and observed at super-resolution (Tang et al., 2016; Biederer et al., 2017). This nanocolumn example shows that synaptic transmission uses a tens of nanometer precision for its organization and for example a misalignment of synaptic terminals is at the basis of several pathological disorders (Südhof, 2008).

Another example is the PSD, that cannot permanently retain glutamatergic receptors that are moving by random motion (Nicoll, 2017). After a long enough time, these receptors spread out, modifying the synaptic current (Freche et al., 2011; Fresche et al., 2013). In the absence of direct experimental approaches, studying the functional consequences of the nanometer precision in the domain between the membrane and the vesicles at the pre-synaptic terminal has recently benefited from three-dimensional modeling and numerical simulations. We focus this review on this paradigm shift of analyzing the diffusion of calcium ions and in particular about the nano-metric relation between the organization of calcium channels and vesicles and how it shapes the release probability, synaptic transmission, asynchronous release, and short-term plasticity.

Despite the fast advances of super-resolution microscopy, that allowed to reconstruct structural in vivo cell properties, or to follow calcium and voltage using genetically encoded indicators (Holcman and Yuste, 2015; Popovic et al., 2015; Beaulieu-Laroche and Harnett, 2018; Cartailler et al., 2018), it remains difficult to study the detailed molecular dynamics in nanometer domains at a time scale < 100 ms. Indeed, to understand how molecules interact in nanometer domains, the notion of concentration has to be abandoned because it does not make much sense due to the large fluctuations in the small number of molecules. However, molecular interactions can still be transformed into a cellular activation at the micrometer level, but the exact biophysical mechanisms remain in most cases unclear or controversial. Modeling and numerical simulations based on biophysical principles have emerged as orthogonal tools compared to experiments to describe molecular dynamics at this spatio-temporal scales (Holcman and Schuss, 2015).

At this intermediate level between the molecular and the cellular scale, physical modeling of diffusion is based on Brownian motion, which requires to specify an inherent time scale of simulations. Indeed the motion of molecules follows a random walk approximation expressed by the Euler's scheme for a trajectory X(t) at time t:

where η is a Gaussian variable, D(X) is the spatial dependent diffusion coefficient and Δt is the time scale to be chosen. It is usually a difficult choice. It should not be too small to avoid wasting simulation times and should not be too large compared to the small spatial scales involved in the microdomain, such as molecular binding sites.

In particular, taking into account in numerical simulations the region between vesicles and the plasma membrane has been particularly difficult to model because of its cusp-like geometry. It requires a specific mathematical treatment to estimate the mean time for a calcium ion after entering through a Voltage-Gated-Calcium-Channel (VGCC) to find a key calcium binding sensor, involved in triggering vesicular release (Guerrier and Holcman, 2014) (Figure 1). Such sites are Ca²⁺-binding proteins, located on synaptotagmins, that are involved in triggering directly or not vesicular fusion (Lee and Littleton, 2015). They are located precisely in this nanometric region below vesicles. Interestingly, spontaneous excitatory and inhibitory transmission are differently regulated by Ca2+ sensors (synapotagmin-1 and Doc2α/β a high-affinity Ca2+ sensors) (Courtney et al., 2018).

Calcium diffusion in the pre-synaptic terminal has traditionally been modeled as two- or three- dimensional diffusion (Modchang et al., 2010; Nadkarni et al., 2010; Stanley, 2016), but ignoring the three-dimensional complications of the vesicle shape. However, for auditory hair cells, Monte-Carlo simulations revealed (Graydon et al., 2011) that the spherical shape of the ribbon where vesicles are tethered, can generate a local Ca²⁺ microdomain that enhances vesicular fusion by trapping calcium ions (Graydon et al., 2011, Figure 3B). This spherical ribbon that aggregates vesicles is likely to create an intermediate microdomain for calcium dynamics between the pre-synaptic bulk and the boundary layer near the membrane, which should be further investigated. For other types of synapses, while estimating the time scale of calcium binding, and computing the vesicular release probability, one cannot ignore the specific three-dimensional organization of the first 100 nm, the region underneath vesicles and the position of calcium sensors (Rothman et al., 2017). The release probability not only depends on the binding of calcium ions to sensor proteins, located underneath the vesicle (Guerrier and Holcman, 2015), but also on the AZ organization: a sparse vesicular distribution vs. vesicular crowding, and channels clustered vs. uniformly distributed (Figure 1B) (Delvendahl et al., 2015; Nakamura et al., 2015). The major components in these dependencies being the distance between VGCC and vesicles coupled to the particular cusp-like geometry (Guerrier and Holcman, 2015).

The synaptotagmin family of molecules are Ca²⁺ sensors for vesicle fusion: following Ca²⁺ binding, activation of the SNARE-complex mediates membrane fusion. Although there are two specific sensors on synaptotagmin, the number of ions necessary for fusion was estimated around 4 or 5 (Kochubey et al., 2011; Schneggenburger et al., 2012). Indeed, for example synaptotagmin1 has two Ca2+ binding C2 domains (C2A and C2B) and three Ca2+ bind in C2A and two Ca2+ bind in C2B (Gruget et al., 2018).

The Hodgkin-Huxley model (Hodgkin and Huxley, 1952) can be used Guerrier and Holcman (2016) to generate a calcium influx current inside the pre-synaptic terminal. Following the opening of VGCC, this current corresponds to an entry per channel during a mean time of ≈1 ms for approximately 80 calcium ions, compatible with a previous estimation of 200 reviewed in Stanley (2016) or with the 45 ions per channels described in Keller et al. (2015): with a total of 12 channels, this would represent 540 ions entering during ≈0.1 ms. After entry, the calcium flux can be split into a ionic component that reaches the small region below the vesicle, and another one reaching the pre-synaptic bulk. The probability for a calcium ion to reach the binding region, defined as a small ribbon joining the vesicle and the plasma membrane, before the bulk, has been computed in models where vesicles are organized in a square lattice (Figure 1B) with length 2H, where H is of the order of the diameter of a vesicle, from 40 to 100 nm (Figure 1C). For a dense set of vesicles distributed on a square lattice, the splitting probability for a calcium ion (modeled as Brownian) to reach the ribbon before the bulk is (Guerrier and Holcman, 2015):

where r_v is the size of a vesicle, H is half the distance between two vesicles, A = 9.8, r(x) is the distance between the point source and the closest vesicle, and ε is the height of the small cylindrical ribbon (Figure 1C), where calcium sensors are located. This probability accounts for the particular geometry of the target and depends on the relative distance between the targets and the source points (Guerrier and Holcman, 2015).

The probability that a finite number T of calcium ions (T = 3, 4, and 5 ions) are bound at specific binding sites located between a vesicle and the synaptic membrane (section 3.1), when N calcium ions have entered through a cluster or a single VGCC located at a distance r from the center of the closest vesicle is defined by:

When the calcium unbinding events are too slow to be taken into account, the probability of vesicular release p_act(r, N) is thus the one to find at least T ions inside the cylindrical ribbon. The probability to find exactly k ions out of N follows the Binomial distribution B(N, p_s(r)), and the steady-state probability is:

where p_s(r) was computed in Equation (3). The maximal VGCC distance r_{max, pact}(N) to activate a vesicle with a probability p_act ≥ 0.8 is shown in Figure 1D. The probability p_act critically depends on the channels-vesicles distance, which can vary from a few to hundreds of nanometers. This fast decrease of the probability with the distance explains the large variability in the release probability as VGCC position can vary over time (Schneider et al., 2015).

The organization of vesicles in the AZ also influences the release probability: when vesicles are sparsely distributed (Figures 1B–D, H = 100 nm) and 100 ions entered through VGCC, then a 80% release probability p_act = 0.8 is reached when the distance between the vesicles and the channels is smaller than 24 nm, compared to the 20 nm radius of the vesicle. This result shows that the co-localization of VGCC with a vesicle is a key feature determining a high synchronous release probability. However, for a high vesicular crowding (Figures 1B–D, described by choosing the distance H = 35 nm) and when 100 ions are released instantaneously at a VGCC, then the probability p_act is higher than 0.9, regardless of the initial position of the channel, suggesting that vesicles are certainly released and leading to a synchronous release.

To conclude, a high crowding of vesicles should be associated with a high-release probability sustaining a synchronous release, while a sparse vesicle density might be associated with asynchronous release. Channels can be organized in clusters or uniformly distributed and this is also a major determinant governing release probability (Figures 1B–D). Indeed, the effect of channels clustering can be modeled by simply increasing the number of entering calcium ions. When vesicles are sparsely distributed, the 24 nm distance required to obtain a release probability p_act = 0.8 when 100 ions are entering through one channel, is increased to 61 nm for 500 ions. This effect results from the local geometry of the ribbon underneath the vesicle. When the number of ions is low, this maximal distance to guarantee p_act = 0.8 does not vary much when the activation threshold T increases from 3 to 5; however, for 500 ions, this distance changes significantly over 15 nm.

The maximal distance r_max,pact(N) between channels and vesicles to obtain a given release probability p_act depends on the number N of entering ions. For a fixed probability p_act, we plotted r_max,pact(N) in Figure 1D. For a sparse distribution of vesicles, characterized by a bulk located a distance H = 100 nm from the membrane, a vesicle is activated with a probability p_act = 0.8 (resp. p_act = 0.2), when 1,200 ions are entering at a distance 100 nm (resp. 450 ions), and 340 ions at a distance 50 nm (resp. 125 ions). This result has to be compared to the 200–500 nm diameter of the AZ (Südhof, 2012). Consequently, a sparse distribution of vesicles at the AZ requires a high number of entering calcium ions in order to trigger fusion, which can be achieved when channels are clustered. However, when channels are co-localized with vesicles, the activation probability p_act is significantly increased: indeed 450 ions are necessary for activation for p_act = 0.2 at a distance 100 nm. When the probability increases to 0.8, the distance reduces to 58 nm.

To conclude, a synapse with high release probability requires a nanometer precision of the channel location. However, this high requirement can be compensated by increasing the number of initial ions entering through VGCC clustering: with 2,000 ions, the maximum distance is relaxed to 140 nm. On the contrary, in a pre-synaptic terminal crowded at its surface with vesicles (characterized by H = 35 nm), very few initial ions are needed for an efficient release. Indeed, 50 ions are enough to activate a vesicle with probability 0.8, wherever the channels are located in AZ (Figure 1D).

To compute the time distribution of the release probability and to account for the calcium ions at the AZ and in the bulk, a full model of the pre-synaptic terminal is needed. The main challenge for such derivation is to account both for the stochastic regime governed by rare events of individual calcium ions arriving to a sensor binding site, and the continuous description of the calcium concentration in the bulk of the synapse (Guerrier and Holcman, 2017). The classical approach consists in using partial differential equations that often cannot take into account easily the specific AZ organization and in particular the geometry near vesicles. An important assumption of these approaches is the use of the forward rate for the calcium to sensors. This rate is often assumed and not derived (contrary to the approach described in section 3.1).

To compute the sensor activation, a different approach is to use Monte-Carlo or Brownian simulations to follow each ionic trajectory. But this approach is often computationally greedy to detect the rare events of calcium hitting a small target (Holcman and Schuss, 2015). Recently, a hybrid Markov-mass action model has been developed (Guerrier and Holcman, 2016), that combines a Markov chain to represent the stochastic events occurring at the AZ, with a mass-action laws model that represents calcium dynamics in the large bulk. The Markov chain and the mass-action model are coupled by the calcium ions coming from the bulk and binding to the sensor. The arrival time of such ions is Poissonian, with the rate computed taking into account the geometry of the vesicle, as discussed in section 3.1. This model is used to compute the time distribution of vesicular release (Guerrier and Holcman, 2016) and it shows that vesicular release is triggered by the binding of calcium ions that can originate either from the bulk or from VGCC.

The distribution of release time is bimodal although it is triggered by a single fast action potential (Figure 2). This simulation is initiated by three channels and each of them let a flow of 80 ions inside the cell during a time scale that was simulated from a Hudgkin-Huxley model (Guerrier and Holcman, 2016).

An example of specific simulation is as follows: the pre-synaptic terminal is a bulbous head of volume ≈1μm³. At the AZ, we positioned eight vesicles, distributed on a square lattice of surface ≈0.13μm² (Figure 1 and Guerrier and Holcman, 2016), so that the distance between two neighboring vesicles is 130nm, and each vesicle has a diameter of 40nm. The three calcium channels are uniformly distributed over the AZ, but remained from a distance of 6−10nm (we chose around 6nm here) from every vesicle. This distance corresponds to the radius of the red ribbon (Figure 3A) of height 1nm the calcium ions need to reach to simulate the binding to a sensor. For each simulation, the terminal undergoes three spikes at a fixed time interval Δt = 20−−150ms. Once a calcium ion enters the terminal, it can either reach a vesicle with probability p_s (Equation 3), or enter inside the bulk with probability 1−p_s according to the scheme of Figure 3C. We already discussed in subsection 2.2 how to chose the calcium diffusion coefficient: at the AZ it is 200μm²/s, and 20μm²/s in the bulk to account for crowding.

Inside the bulk, ions bind to buffers with a rate constant k0=5.6s-1 and unbind with rate k-1=500s-1. Calcium can be extruded by pumps with a rate kpump=0.88s-1 or they can leave the terminal with rate kes=6.1s-1. Finally, calcium can bind to the sensors located underneath a vesicle with a rate kS=0.3s-1. The number of buffer molecules in the bulk varies from 0 to 1,000. A calcium ion bound at the calcium sensor can unbind with a rate kU=2000s-1 (fast unbinding) or 5s⁻¹ (slow unbinding). Once 5 calcium ions are bound to a calcium sensor, then the vesicle fuses with the synaptic membrane and the vesicle spot becomes free for a new vesicle coming from the recycling pool to bind with a rate 1τDock (Figure 3D). An immediate refilling of vesicles at the AZ is obtained with a time τ_Dock = 0 ms, but other delay are possible such as τ_Dock = 50ms, or not refilling of vesicles with τ_Dock = ∞ as shown in Figure 4.

In that simulation framework, the first peak in the time distribution of vesicular release follows a single stimulation and the second one in Figure 2 (that is smaller in amplitude and wider) corresponds to the random arrival, over a much longer time period, of ions located in the synaptic terminal to small binding vesicular targets. To conclude, multiscale stochastic modeling approaches allow studying cellular events based on integrating discrete molecular events over long time scales from one milliseconds to seconds.

The lesson from modeling and numerical simulations of diffusion in the first hundred nanometers between docked vesicles and the plasma membrane is that this boundary layer is crucial for computing the vesicular release probability due to the critical position of vesicular sensors. This space is difficult to access experimentally and its role has been underestimated in short-term plasticity. However, experimental approaches using calcium chelator such as BAPTA or EGTA confirm the role of calcium ions traveling from far away compared to the ones entering directly through VGCC located underneath a vesicle to trigger release.

Another key feature relevant for short-term plasticity is the structural organization and the spatial correlation between the distribution of VGCCs and vesicles. Do vesicles contain the same amount of close VGCCs? What defines the exact location of VGCC underneath vesicles? Can this number fluctuate? What happens after vesicular fusion? How are VGCC redistributed? It is possible that VGCCs are constantly moving to find local optimal sites location (Schneider et al., 2015).

Short-term facilitation is classically thought as the accumulation of calcium in the synaptic bulk (calcium hypothesis of Katz and Miledi, 1968), due to various possibilities such as slow and fast buffers. But computational evidences (Guerrier and Holcman, 2017) reviewed here suggest that calcium accumulation at the sensor binding sites and not in the bulk is actually the determinant effect to pre-activate vesicular release (by binding a certain fraction of the sensor sites). A similar conclusion was reached in Bornschein et al. (2013) using fluctuation analysis, calcium imaging and numerical simulation analysis indicating that the residual calcium bound to the release sensors (see also Schneggenburger and Neher, 2000), after the first AP could cause Paired Pulse Facilitation at Purkinje neuron synapses.

A byproduct of facilitation should be asynchronous release, because the random calcium accumulation at sensor increases the time window when a vesicular can be released, due to ions arriving at random time from the bulk, thus leading to a high variability in the calcium arrival and the vesicular release times. However, synaptic facilitation requires a low concentration of calcium buffer, suggesting that for facilitating synapses, calcium buffers should be maintained at a low level. In general, the relation between facilitation and asynchronous release triggered by residual calcium can be mediated by the specific molecular composition of the vesicle: indeed synaptotagmin 7 (Syt7) (Turecek and Regehr, 2018) present at cerebellar granule cell synapses onto stellate cells and Purkinje cells (mice) plays such a role: in Syt7 KO, facilitation, and asynchronous release are smaller and shorter lived than in WT, although the residual calcium was unchanged. This is in contrast with synchronous released mediated by synaptotagmin I. It would be certainly interesting to model the distribution of various synaptotagmin on vesicles to determine when a vesicle will be more likely to be released asynchronously vs synchronously. In particular, is the distribution of these two SytI and 7 different on vesicles? Do they have very different calcium buffer affinity? These properties could be explored using numerical simulations.

In terms of calcium residual regulation, a low concentration level can be achieved by preventing ER or mitochondria to come in too close proximity of the AZ. In contrast, synchronous release is associated with a high calcium buffer concentration, preventing calcium ions to travel from far away. More specifically, two-dimensional numerical simulations (Delvendahl et al., 2015), modeling essentially the synaptic bulk, revealed that fixed endogenous buffers with low affinity, characterized by a low calcium-binding ratio, mixed with mobile buffers with high affinity, results in a fast AZ calcium clearance. This results in synchronous high-frequency transmission (at 200 Hz). But it remains unclear how calcium fluctuations is maintained low (Modchang et al., 2010; Nadkarni et al., 2010; Dittrich et al., 2013; Weinberg, 2015) especially near the vesicular calcium sensor.

Finally, how the rate of vesicular release can vary over 6 orders of magnitude for the same synapse (Kochubey et al., 2011) also remains enigmatic. The cusp geometry and rare binding events may hold the key to the solution of this spectacular modulation of the vesicular release rate (Guerrier and Holcman, 2015). Modeling vesicular trafficking and recycling at various synapses including ribbon synapses (Thoreson et al., 2016) should clarify the organization of the pre-synaptic terminal (Graydon et al., 2011, 2014) and the effect of VGCC trafficking (Schneider et al., 2015).

All authors listed have made a substantial, direct and intellectual contribution to the work, and approved it for publication.