## Abstract

If, as is widely believed, information is stored in the brain as distributed modifications of synaptic efficacy, it can be argued that the storage capacity of the brain will be maximized if the number of synapses that operate independently is as large as possible. The majority of synapses in the brain are glutamatergic; their independence will be compromised if glutamate released at one synapse can significantly activate receptors at neighboring synapses. There is currently no agreement on whether “spillover” after the liberation of a vesicle will significantly activate receptors at neighboring synapses. To evaluate the independence of central synapses, it is necessary to compare synaptic responses with those generated at neighboring synapses by glutamate spillover. Here, synaptic activation and spillover responses are simulated in a model, based on data for hippocampal synapses, that includes an approximate representation of the extrasynaptic space. Recently-published data on glutamate transporter distribution and properties are incorporated. Factors likely to influence synaptic or spillover responses are investigated. For release of one vesicle, it is estimated that the mean response at the nearest neighboring synapse will be <5% of the synaptic response. It is concluded that synapses can operate independently.

- glutamate
- synapse
- AMPA
- NMDA
- transporter
- synaptic transmission
- diffusion
- tortuosity
- volume fraction
- porous medium
- synaptic cleft
- vesicle
- neurotransmitter

Spines of excitatory (glutamatergic) synapses have long been thought to be compartments capable of individual behavior, either electrical (Miller et al., 1985) or chemical (Zador and Koch, 1993). However, such compartmentation would be irrelevant if synaptic signaling at each spine was not independent of that at neighboring spines. Because information is stored in the brain as changes of efficacy at such synapses, their independent operation would have an important benefit: maximizing the information storage capacity of the brain. One prerequisite of independent signaling at spine synapses is that neurotransmitter (glutamate) released at one synapse should not significantly influence neighboring synapses, by activating (or desensitizing) their receptors.

Notwithstanding the evidence in favor of spillover or consistent with its occurrence in a number of experimental situations (Isaacson et al., 1993; Trussell et al., 1993; Barbour et al., 1994; Mennerick and Zorumski, 1994; Kullmann et al., 1996; Silver et al., 1996; Bergles et al., 1997; Clark and Barbour, 1997; Scanziani et al., 1997; Dzubay and Jahr, 1999; Carter and Regehr, 2000), little progress has been made in deciding whether glutamate release at a single synapse can activate receptors at a neighboring synapse. It is difficult to conceive of feasible experiments that could test this phenomenon, so modeling is still likely to give us a useful insight into the importance of this mechanism.

To evaluate the significance of spillover, it is necessary to compare receptor responses caused by spillover with those at the releasing synapse. Calculation of the spillover response obviously requires treatment of diffusion in the complex extrasynaptic structure of brain tissue. However, the manner of transmitter diffusion away from the synapse will also influence receptor activation within the synaptic cleft. Almost all previous models of synaptic activation have made no allowance for this, using an infinite-disk model that will dilute released transmitter much more slowly than the three-dimensional extrasynaptic tissue (Barbour and Häusser, 1997). Synaptic and spillover responses have been compared in one previous model that takes into account the extrasynaptic geometry (Rusakov and Kullmann, 1998b). The results of those calculations suggested that the ratio of spillover to synaptic responses could be as high as 75%.

Calculations are presented here that assess synaptic independence by comparing synaptic receptor activation and that caused by spillover, in a model in which the influence of the extrasynaptic space is taken into account. Analysis of the literature on diffusion measurements in brain tissue allows constraint of the glutamate diffusion coefficient. The influence of several factors that could alter synaptic and/or spillover responses is evaluated: glutamate receptors, glutamate binding sites, glutamate uptake, rate and amount of glutamate release, as well as the properties of the synaptic cleft. In contradiction to previous work, the present simulations predict that the effects of spillover are sufficiently minor to enable synapses to operate independently.

## MATERIALS AND METHODS

*Diffusion and receptor models.* A representation of the synaptic and extrasynaptic space was developed that required only one distance dimension, i.e., which had radial symmetry. It is illustrated in the diagram of Figure1*A*. The synaptic cleft—the volume between the presynaptic and postsynaptic elements—was represented by a flat cylinder (disk) within which transmitter diffused freely. The concentration of diffusate was assumed to be uniform in the direction perpendicular to the cleft. In the model synaptic cleft, two regions were distinguished for the purpose of placing receptors: a central region corresponding to the postsynaptic density (PSD) and a surrounding annulus. The extrasynaptic space was modeled as a porous medium (Nicholson and Phillips, 1981; Nicholson and Sykova, 1998) in which the standard diffusion equations apply, after slight modification. An apparent diffusion coefficient is used that is related to the free diffusion coefficient by a parameter called the tortuosity, represented by λ. In fact,*D*_{pm} (for porous medium) =*D*_{free}/λ^{2}. The tortuosity represents the requirement that diffusion must circumnavigate obstacles, following longer paths between two points than would have been the case in the absence of obstacles. It is also necessary to account for the fact that the diffusion in which we are interested is restricted to the extracellular space, which occupies only a fraction of the total brain volume. In a porous medium, this restriction is described by the volume fraction α. Typical values of volume fraction and tortuosity for brain tissue are α = 0.2 and λ = 1.6 (Nicholson and Sykova, 1998). The use here of these particular values will be justified carefully below (Results). How a smooth transition was effected between the disk and porous medium regions will be described shortly.

Diffusion in the synaptic cleft, porous medium, and transition region between the two was described by one general equation. Following an approach similar to that of Crank (1975), a diffusion equation with one distance variable was derived (although no claim to originality is made). Because the apparent diffusion coefficient in the porous medium will in general differ from that in the synaptic cleft, the form of the equation required for a non-constant diffusion coefficient was obtained:
Equation 1where *V*(*r*) is the volume available for diffusion within radius *r*. Putting *D* constant and substituting *V* = π*r*^{2}*h* or*V* = (4/3)π*r*^{3}(or even *V* = *kr, k* a constant), it is easy to verify that the above equation reduces to the standard forms for cylindrical or spherical (or plane) diffusion with a constant diffusion coefficient.

Little is known about the transition from disk-like synaptic cleft to porous medium extrasynaptic space, except perhaps the diameter of the synaptic cleft. The representation of the transition region is therefore somewhat arbitrary. A smooth transition with adjustable end points was constructed using an interpolating polynomial whose end points corresponded to the radii of the beginning and end of the transition region. In favor of a smooth transition, it can be argued that an abrupt transition would certainly be unrealistic. A smooth transition also helps reduce numerical errors. The interpolating function was chosen to have a value of 0 at its inner end point, of 1 at its outer end point, and zero first and second derivatives at both end points. These conditions ensured that no discontinuities occurred in Equation 1. The interpolating polynomial chosen was thus a quintic whose six coefficients were determined using the above boundary conditions. The volume function was thus:
Equation 2where *f*(*r*) is the interpolating function, and *a* and *b* are the radii of the inner and outer end points of the transition region. Also,
Equation 3(here, *h* is the cleft height: 20 nm) and
Equation 4(α is the volume fraction of brain tissue: 0.2). A similar transition of (apparent) diffusion coefficient between the cleft and the extrasynaptic space was also used:
Equation 5where, unless otherwise stated,*D*_{cleft} =*D*_{free}, the diffusion coefficient for unhindered diffusion (7.6 × 10^{−8}dm^{2}/sec = 7.6 × 10^{−6}cm^{2}/sec; in fact the value for glutamine; Longsworth, 1953), and *D*_{pm}=*D*_{free}/λ^{2}; λ = 1.6. The choices of diffusion coefficient and tortuosity will be justified in Results. Figure 1*B* allows a comparison of a typical choice of *V* with*V*_{cleft} and*V*_{pm}. A similar plot, of *D*with *D*_{cleft} and*D*_{pm} is shown in Figure 1*C*.

The model was populated with AMPA and NMDA-type glutamate receptors, glutamate transporters, and simple binding sites at concentrations described in Results. For each of the receptor types and the binding site, two concentrations were set independently: over the PSD and beyond the PSD. For transporters, two independent concentrations could also be set, but these were for the synaptic cleft (disk-like region) and for the transition region–porous medium. This was because most transporters are expressed in glial cell membranes, which are encountered only outside the region of apposition of the presynaptic and postsynaptic cells. The receptors were modeled according to published kinetic schemes. For AMPA receptors, the “set 1” rate constants of Jonas et al. (1993) were used. Although there is evidence that this model may be mechanistically inexact (Rosenmund et al., 1998), it well describes hippocampal AMPA receptor responses in patches. For NMDA receptors, the model of Lester and Jahr (1992) was used. The kinetics of the simple binding site were characterized by an association rate constant and a dissociation rate constant, the values of which will be given in Results.

*A glutamate transporter model with a rapid glutamate-trapping step.* In previous simulations of synaptic activation including transporters, only a very simple model was used (Barbour and Häusser, 1997; Rusakov and Kullmann, 1998b), in which transporters were represented as a binding site from which a slow reaction (transporter cycling) could eliminate bound glutamate. Recent work has suggested that an extension to the simple model is necessary. Recordings of transporter currents in the absence of anions that permeate the associated conductance have shown that a very brief transient is induced by rapid applications of glutamate (Bergles and Jahr, 1997; Auger and Attwell, 2000; Otis and Kavanaugh, 2000; Grewer et al., 2000). It is thought that this transient charge movement reflects a conformational change that traps the bound glutamate molecule, probably by translocation across the membrane. The rate of this step (of the order of 1000/sec) far exceeds the much slower overall cycling rate (14–50/sec). It seems likely that transporters with such glutamate-trapping behavior might be particularly effective (Auger and Attwell, 2000), so it is of interest to investigate the effect of such behavior in the present simulations.

A rapid glutamate-trapping step was accordingly incorporated into the simple transporter model. The new reaction scheme is shown in Figure2*A*. It is described by a second order differential equation. The two characteristic time constants of its solution assume simple forms for saturating and zero glutamate concentrations, allowing the four rate constants to be determined from as many simple experiments. The first piece of information is the apparent steady-state affinity for glutamate,
Equation 6The kinds of experiments that could be used to obtain the remaining kinetic information are illustrated in Figure 2, *B*and *C*, using simulated responses of the model to glutamate applications; the formulas relating the relevant time constants to the model rate constants are given in the figure. The kinetic determinations can to some extent be described by analogy with similar experiments on AMPA receptors. The time constants of “deactivation”, “desensitization”, and “recovery from desensitization” need to be determined (the first two in the absence of permeant anions). Of the four measurements required, only two have been determined accurately for hippocampal glial transporters (EAAT2/GLT1). These are the affinity (13 μm;Bergles and Jahr, 1997) and the recovery time constant (20 msec;Bergles and Jahr, 1998; Otis and Kavanaugh, 2000). The recovery time constant can be used to set *k*_{3} to 50/sec forthwith. Based in part on the few published experiments of the types shown in Figure 2, *B* and *C*, a range of plausible values was chosen for *k*_{+1}(Wadiche and Kavanaugh, 1998; Mennerick et al., 1999; Auger and Attwell, 2000) and *k*_{2} (Bergles and Jahr, 1997; Otis and Kavanaugh, 2000; Auger and Attwell, 2000), allowing *k*_{−1} to be obtained from Equation 6 for the steady-state affinity.

*Method of solution.* The coupled diffusion, receptor, and transporter equations were solved numerically using a space discretization by finite difference methods and integration by the first order explicit (Euler) method. Two spatial grids were used: the first extended from the point of release (*r* = 0) to a radius of 1 μm at 5 nm intervals; the second extended from there to 16 μm with 50 nm intervals. The first and second spatial derivatives of the concentration, volume, and diffusion coefficient were represented using the standard central difference approximations:
Equation 7
Equation 8where *X* is the function to be approximated (e.g., concentration). The derivatives of the volume and diffusion coefficient were calculated with a high accuracy (small δ*r*) at the beginning of each simulation, and the values were stored in a look-up table. At the outermost grid point, the concentration was clamped at zero, but this was sufficiently distant to influence only slightly the diffusion of transmitter in the region and time of interest. This can be shown by the close similarity between the receptor responses to numerically calculated glutamate transients and when the receptors were driven by analytical solutions for geometries where these exist (± <1% at their peaks; data not shown).

The apparent singularity in Equation 1 at *r* = 0 was dealt with according to Smith (1985). Briefly, noting that around*r* = 0, *V* = π*r*^{2}*h* and the diffusion coefficient was always constant, Equation 1 reduces to:
Equation 9The last term assumes the indeterminate form 0/0 at*r* = 0, but application of l'Hôpital's rule shows it to have a limit of*D*∂^{2}*C*/∂*r*^{2}as *r* = 0 is approached. Thus, Equation 9 is simplified to:
Equation 10
Equation 11When release was directly into a porous medium, the factor of 4 in Equation 11 was replaced by a factor of 6 as required for the case of spherical diffusion.

In the model the glutamate released was placed at *r* = 0 and at the next grid point. To predict the appropriate concentration, it would be necessary to know the effective volume represented by these two grid points. This problem was circumvented by an iterative adjustment of the initial concentrations to fit the model's output to theoretical predictions. When non-instantaneous release was modeled, a steady flux into the first two grid points was maintained for the appropriate time.

It was found convenient during the calculations to convert the distance units of all quantities, in particular the diffusion coefficient, to decimeters (10^{−1} m; 1 dm^{3} = 1 liter), since concentrations were then obtained directly.

For synaptic receptor responses, the average activation (*P*_{open}) of receptors over the PSD is shown in all figures. This is the average calculated over the area within the radius of the PSD (120 nm), obtained by integrating numerically:
Equation 12Most of the parameter values were obtained from experiments performed at room temperature, which can therefore be taken as the temperature in the model.

*Verification of the model.* A time step of 10 nsec was used for the calculations. Halving the distance between space grid points and decreasing the time step by a factor of slightly more than four (which is necessary for stability of the numerical method) had only a very small effect on the solution, as assessed by the peak receptor activation under typical conditions, which was altered by <0.1%. It was concluded that the approximate numerical solution was converging to the true solution.

The full model was solved numerically because no analytical solution is available. It thus poses a problem of verification. It was, however, possible to test the model partially in various ways, against previous models or analytical solutions.

First, the diffusion part of the model was tested. It was verified that the calculated solution for diffusion within an infinite disk reproduced to a reasonable accuracy the solution obtained by analytical methods (Crank, 1975; Barbour & Häusser, 1997). It may be seen from Figure 3*A*, which shows both the analytical solutions and the calculations, that the calculations give satisfactory results. Very brief discrepancies between the numerical and analytical solutions exist after release; this is presumably the result of the different initial conditions of release. Quite large deviations from the analytical solution are caused by the loss of glutamate at the edge of the model. This becomes significant after 10–20 msec (depending on the geometry) in the periphery of the model. It was shown above that these discrepancies do not greatly influence the receptor responses in the region and time of interest.

A general property of diffusion is that as time advances the initial conditions are “forgotten.” This provides one way of testing the calculations for composite geometries, in which a cylinder joins onto a porous medium (via the smooth transition), because at long times the solution should converge to that of the porous medium. This check is performed in Figure 3*B*. Although the concentration at the point of release initially follows the analytical solution for the infinite disk, the central concentration and the concentration at 500 nm display the expected behavior by approaching closely the analytical solution for the porous medium toward the end of the calculation.

It is depressingly easy in such calculations to generate or destroy diffusate. The total amount of glutamate present in each simulation was estimated by using Simpson's rule to evaluate the integral:
Equation 13where *C* is the concentration of glutamate in all forms. When the transition between disk and porous medium was too abrupt (<100 nm), significant amounts of glutamate were “lost.” For all the simulations shown, the value of the above integral (Eq. 13) did not deviate from the theoretical value by >0.5% between 1 μsec and 10–20 msec (depending on the geometry). The deviation at short times is attributable to the inaccuracy of the integration method when glutamate is present in only a few elements of the model. The deviation at longer times is attributable to the loss of glutamate at the edge of the model.

The implementation of the receptor models was checked by overriding the diffusion calculations of the simulation program and fixing the glutamate concentration as required. For the AMPA receptor model, numerous checks were performed using the values reported by Jonas et al. (1993) for the output of the model. In particular, the deactivation and desensitization time constants were verified, as well as the dose–response curve (EC_{50}, Hill coefficient, and *P*_{open,max}). The developers of the NMDA receptor model did not publish such values. The present model does, however, verify their stated*P*_{open,max} of ∼0.3. The output of the model conforms to classical NMDA receptor kinetics. In response to a 1 msec application of 1 mm glutamate, the peak*P*_{open} was 0.257. The 10–90% rise-time of the response was 9.9 msec, and the decay time constant was 85.4 msec. The responses of the transporter model to various test pulses were verified against the analytical solutions.

## RESULTS

### Synaptic receptors are not saturated and spillover responses are small

Most models of central synaptic transmission represent the synaptic cleft as a flat cylinder or disk. However, instead of limiting this representation to the true extent of the cleft, it is extended effectively to infinity. This introduces a potential error, because diffusion in the extracellular space around the cleft is probably quite different from that in the cleft. Most obviously, the volume of the extracellular space increases as the cube of the radius, whereas for a disk the volume only increases as the square. This means that eventually the extracellular space is expected to dilute released transmitter much more rapidly than the infinite disk. This behavior is likely to have an important influence both on synaptic responses and on responses caused by spillover at neighboring synapses. Both will be investigated here, because careful assessment of the significance of spillover responses requires comparison with synaptic responses under the same conditions.

A model with a more realistic geometry was developed, in which a disk and a porous medium are combined (see Materials and Methods, Fig. 1). The composite structure is intended to represent a “standard” central synapse—a hippocampal synapse. The model synapse (representing “macular synapses” in the terminology of Ventura and Harris, 1999) had a PSD of radius 120 nm (Takumi et al., 1999), and the disk-like region extended beyond the PSD to a radius of 180 nm (Ventura and Harris, 1999). The disk was joined via a transition region of 200 nm to a porous-medium extracellular space with a volume fraction of 0.2 and a tortuosity of 1.6. The length of the transition region was (arbitrarily) chosen to correspond approximately to the radius of the synaptic cleft. The influence of this choice will be examined below.

Transmitter diffusion and receptor activation were first investigated in the presence of very low concentrations of AMPA and NMDA receptors, to enable their responses to be monitored without altering transmitter diffusion. The diffusion coefficient within the synaptic cleft was set to 7.6 × 10^{−8}dm^{2}/sec, and the effective diffusion coefficient for the porous medium was therefore 3.0 × 10^{−8}dm^{2}/sec. Responses to instantaneous release of 5000 molecules of glutamate (the likely content of a small synaptic vesicle; Bruns and Jahn, 1995) at the center of the synaptic cleft were calculated. Figure 4 shows the mean synaptic *P*_{open} (averaged over the PSD) for AMPA and NMDA receptors, as well as their responses at a radius of 500 nm, a typical distance for a neighboring synapse, to give an idea of the likely spillover response.

Overall, levels of synaptic receptor activation are quite low (compared with maximum possible *P*_{open} values for the AMPA and NMDA receptor models of ∼0.8 and 0.3, respectively). AMPA receptors reach only 15% of the maximum*P*_{open} and the figure for NMDA receptors is 20%. The responses caused by spillover are considerably weaker than the synaptic responses, with AMPA receptor*P*_{open} reaching 0.8% of the synaptic response, and NMDA receptors attaining a somewhat higher level of activity (3.9% of the synaptic response). We can see immediately that it is difficult to generate significant receptor activation by spillover.

As suggested by the similar time courses of glutamate concentration at radii of the order of 500 nm (and greater) in the composite geometry and for spherical diffusion from a point source in a porous medium (Fig. 3*B*, *solid* and *dashed lines,*respectively), the receptor activities at 500 nm generated by these two glutamate transients are very similar (Fig. 4*C*). This suggests that the presence of the synaptic cleft per se has little influence on receptor activity at neighboring synapses. This is not the case, however, for synaptic responses, which, as will be shown later, are quite sensitive to the properties of the synaptic cleft.

The difference made by using the present synaptic geometry is illustrated by simulations run for the infinite disk morphology, using two different values of diffusion coefficient (Fig.4*A,B*). Even with unhindered diffusion, the slower dilution of glutamate in the infinite disk model results in a 1.8-fold increase of AMPA receptor activation and a 3.0-fold increase for NMDA receptors. An even greater increase is observed when a lower diffusion coefficient is used with the infinite disk (a combination found in many previous models). The geometry and diffusion coefficient are therefore critical in determining receptor activation.

In the composite geometry used, the length of the transition region between the disk and porous medium was chosen arbitrarily, so it is important to assess the possible error introduced by this choice. Extreme upper and lower limits on synaptic activation are given by the infinite disk and “open volume” geometries simulated in Figure 4,*A* and *B,* (the latter may approximate exposed synapses in culture). Plausible values of the length of the transition region are a few hundred nanometers. Very short transition regions pose technical problems (see Materials and Methods), whereas longer transition regions are constrained by the necessity for the volume fraction of the tissue to be 0.2. If too long a transition region is used, the transition regions from neighboring synapses will overlap frequently, and the tissue volume fraction will fall below the required value. The receptor activation occurring for transition regions ranging from 100 to 300 nm was calculated. The effects on the peak synaptic responses were less than ±5%, respectively for AMPA receptors and less than ±3% for NMDA receptors. The spillover responses were equally insensitive to the length of the transition region in this range (less than ±3%). These tests suggest that large errors are not introduced by the arbitrary choice of the length of the transition region.

It may seem surprising that NMDA and AMPA receptors are activated to a similar degree (∼15–20% of maximum*P*_{open}), because NMDA receptors are generally considered to have a much higher affinity for glutamate than AMPA receptors (note however that the definition of affinity for AMPA receptors is problematic because of their desensitization). The explanation of this paradoxical result is that the synaptic glutamate transient is too brief to allow equilibrium binding to be approached. This can be illustrated with the simplest possible reaction of an agonist binding to a single site. If the agonist concentration is initially zero and the sites are unoccupied, the fraction of sites bound following a step change of the agonist concentration (at*t* = 0) follows a time course described by:
Equation 14where *C* is the agonist concentration,*K*_{d} the dissociation constant (the ratio of the rate constants of dissociation,*k*_{off}, and association,*k*_{on}), and τ = 1/(*Ck*_{on} +*k*_{off}). The left-hand factor of Equation 14 (in curly braces) reflects the equilibrium binding position, whereas the right-hand term governs the exponential approach to equilibrium. If we set *Ck*_{on} ≫*k*_{off} (which is the same as*C* ≫ *K*_{d}), the above equation reduces to:
Equation 15which is independent of *k*_{off} and therefore of equilibrium affinity. It is interesting to point out that under these conditions (and assuming the agonist is removed before the sites are saturated) it is perfectly possible for a “low-affinity” site with a fast *k*_{on} to attain a higher occupancy than a “high-affinity” site with a slow*k*_{on}. Synaptic receptor activation appears to operate under conditions resembling those for which Equation15 applies, and the similar *k*_{on} values for AMPA and NMDA receptors lead to the approximately equivalent activation of the two receptor types.

The time-dependent binding behavior illustrated by Equations 14 and 15also explains why NMDA receptors at neighboring synapses are not saturated and are, in fact, hardly activated. Although glutamate concentrations that would saturate NMDA receptors can occur at neighboring synapses (e.g., 28 μm at 500 nm), these concentrations are not maintained for a sufficient length of time for equilibrium to be approached. In fact, the glutamate transient is so brief that receptor activation is unlikely.

### The glutamate diffusion coefficient is not a free parameter

The value of the glutamate diffusion coefficient is a critical parameter that ultimately determines the degree of synaptic and spillover receptor activation. If diffusion is slower, receptor binding has more time to approach equilibrium before the glutamate is diluted. Unfortunately, this issue is very confused in the literature, and the present choice will require careful justification. The problem can be divided into two: determining the diffusion coefficient within the synaptic cleft and within the extrasynaptic space. The latter will be treated first.

It is an experimental fact that many tracer molecules restricted to the extracellular space diffuse through brain tissue with effective diffusion coefficients some 2.5-fold smaller than in dilute aqueous solution (Nicholson and Sykova, 1998; Nicholson et al., 2000). This is equivalent to saying that their diffusion can be described using the diffusion equations modified for use with a porous medium with a tortuosity value of 1.6. If such measurements existed for glutamate, the appropriate effective diffusion coefficient could be used without further ado. Unfortunately, the presence in brain tissue of powerful glutamate uptake systems and various glutamate binding activities, coupled with the difficulty of detecting glutamate on a microscopic scale (glutamate-sensitive microelectrodes do not exist and glutamate cannot be oxidized during amperommetry), makes it very unlikely that accurate direct measurements of glutamate diffusion in brain tissue will be obtained. Thus, all modellers are faced with the problem of deducing the properties of glutamate diffusion from measurements made with other tracer molecules.

Tracer molecules of either positive or negative charge and of molecular weights bracketing that of glutamate all diffuse with tortuosity values of ∼1.6. It might therefore seem reasonable to assume that glutamate behaves similarly, but this course of action has rarely been followed, and a number of different arguments have been advanced in justification, invoking viscosity, nonspecific binding, or specific binding. These will be considered in turn.

The measurements of diffusion in brain naturally include the effects on the tracer molecule of any viscosity, so a different tortuosity value for glutamate need only be used if it is affected differently by the viscosity of the extracellular fluid. The effects of solution viscosity can depend on the size of the diffusate molecules (smaller molecules being less sensitive to viscosity), and it has been suggested that the diffusion of glutamate [molecular weight (MW) 146] would be more affected than that of typical tracer molecules such as tetramethylammonium (MW 60) (Rusakov and Kullmann, 1998b). In fact, because tracer molecules both larger and smaller than glutamate diffuse similarly (for review, see Nicholson and Sykova, 1998), it is unlikely that glutamate itself behaves differently. Amazingly, there exist molecules with molecular weights of up to 10^{6} that diffuse in brain tissue with tortuosities of 1.6 (Nicholson et al., 2000). However, some rigid macromolecules do exhibit somewhat higher tortuosities (up to 2.2). The values of tortuosity obtained for brain tissue could plausibly arise on geometrical grounds alone (Gardner-Medwin, 1980), and it is currently thought that the viscosity of the extracellular fluid is little different to that of normal saline (Nicholson and Sykova, 1998). It is concluded that a different tortuosity value for glutamate cannot be justified on the grounds of a different sensitivity to the viscosity of the extracellular fluid.

The principal argument against an important influence of nonspecific binding of glutamate and other molecules is that tracer molecules of different charges diffuse similarly (Nicholson and Sykova, 1998). The relatively low values of tortuosity determined for brain tissue also argue against a significant slowing of the diffusion of all tracer molecules. The absence of any such effect is not particularly surprising, because nonspecific binding sites would probably have low affinities for glutamate and would therefore need to be present in high concentrations to have a noticeable effect.

In contrast, it is certainly the case that glutamate is subject to specific interactions with numerous glutamate binding sites. These include the multitude of glutamate receptors and also several isoforms of glutamate transporters. The binding of glutamate to these binding sites will retard its diffusion. However, it is incorrect in a model to account for this by lowering the diffusion coefficient, because binding sites do not only retard diffusion, they also reduce the free concentration of diffusate. Although the slowed diffusion will tend to increase receptor activation, this effect will be opposed by the concomitant reduction of concentration. In fact, as will be shown in the next section, the introduction of glutamate “buffering” sites never increases receptor activation. Any binding sites considered potentially important should thus be modeled explicitly, without changing the diffusion coefficient of unbound diffusate. This also ensures correct modeling of the temporal and concentration-dependent behavior of the glutamate buffer.

Because empirical measurements of tortuosity are available, it is mostly of theoretical interest to estimate tortuosity by other means. Nevertheless, theoretical values have on occasion been preferred to the experimental determinations (Rusakov and Kullmann, 1998a,b). However, it should be noted that the *ab initio* calculations of geometrical tortuosity used to obtain those theoretical values (Rusakov and Kullmann, 1998a,b) are controversial and have attracted criticism from Chen and Nicholson (2000). The issue is amplified in a brief to this paper. For the simulations presented here, the empirically determined value of tortuosity has been used.

To summarize, if the extracellular space is to be modeled as a porous medium, which is the case here, the effective diffusion coefficient for glutamate should be obtained by setting the tortuosity equal to the value measured for small tracer molecules in the desired tissue. The influence of any important binding sites should be evaluated by explicit modeling, rather than by altering the diffusion coefficient for glutamate. The values of 1.6 for tortuosity and 0.2 for volume fraction that have been reported for many brain regions (Nicholson and Sykova, 1998) have been used here. Very similar values have been found for the rat hippocampus (Mazel et al., 1998). Significantly lower values of volume fraction for the hippocampus have been reported (α = 0.12; McBain et al., 1990) and used in previous modeling (Rusakov and Kullmann, 1998b). Mazel et al. (1998) suggest that the discrepancy arises because the influence of tissue anisotropy was not recognized in those studies. In any case, the disputed values concern the stratum pyramidale, which contains few excitatory spine synapses.

It is worthwhile to comment briefly on the problem of choosing the extracellular diffusion coefficient for a model in which the extracellular space is not treated as a porous medium but is instead represented explicitly, based on an electron microscopic reconstruction of a volume of neuropil. Diffusion in such models is conveniently investigated by Monte Carlo methods. As argued above, current thinking suggests that the diffusion coefficient for unhindered diffusion will correctly describe the behavior of glutamate within the extracellular space. However, it would in theory be necessary to confirm that the value chosen can reproduce the tortuosity measured in the tissue. This may not be trivial (and has not been reported), because it is not known how the volume and shape of the extracellular space are altered by fixation, although this compartment is suspected of shrinking.

Because synaptic clefts contribute a negligible fraction of the extracellular space, the properties of diffusion within them could be very different from those of the extracellular space without influencing diffusion measurements in bulk tissue. For this reason, little is known about diffusion in the synaptic cleft, and modellers still have some room for maneuver here. It is incorrect, however, to justify slowed diffusion in the synaptic cleft on the basis of the tortuosity of brain tissue, because the obstacles that make the brain tortuous are larger than the synaptic cleft and are outside it. It is the subject of speculation that the cleft may contain microscopic or submicroscopic obstacles to diffusion, and the possible effects of this will be investigated below. However, for most of the simulations here the choice of parsimony is made: it is assumed that glutamate diffuses freely in the synaptic cleft.

### Glutamate buffers can only reduce receptor activation

Given the large number of glutamate receptor subtypes and transporters and the paucity of information about the distribution and kinetics of some of these molecules, it will be useful to investigate in general terms the possible effects of binding sites, in particular to assess the consequences of omitting some binding sites from the model. Therefore, before populating the model with realistically modeled receptors and transporters (except NMDA and AMPA receptors at negligible concentrations, to monitor receptor activation in the various conditions), the effects of placing simple binding sites at the synapse and/or in the extracellular space were investigated. A series of simulations were run in which 100 μm (20 μm with respect to tissue volume) binding sites with a *k*_{on} of 10^{7}/(m sec) and*k*_{off} values of 10^{3}, 10^{2}, and 10/sec (giving *K*_{d} values of 100, 10, and 1 μm) were placed either in the synapse (over the PSD) and/or in the extracellular space (beyond the PSD). The results are illustrated in Figure 5.

In no case did the addition of binding sites lead to increased NMDA or AMPA receptor activation. It is sometimes argued that the slowed diffusion caused by glutamate binding will allow receptor activation to build up to higher levels. This is seen to be mistaken, for the simple reason that binding to the sites also reduces the free concentration of glutamate, which, of course, tends to reduce receptor activation. In fact, the time that a bound molecule spends not diffusing is exactly the time during which interaction with receptors (other than the bound site) is impossible. The reductions of free concentration and effective diffusion coefficient are therefore proportional. For receptor activation, the two effects exactly cancel under certain conditions (very brief glutamate transients). For binding sites with plausible kinetic properties, their effect will be time- and concentration-dependent and this, coupled with the receptor kinetic properties, causes a net reduction of receptor activation.

Several other important features emerged from these simulations. The first is that placing binding sites only in the cleft had remarkably little effect on receptor responses anywhere (data not shown). In addition, synaptic responses were in general rather difficult to alter, however binding sites were distributed. Importantly, in all the conditions modeled, the reduction of spillover responses after the addition of binding sites was greater than the corresponding decrease of the synaptic response, for both AMPA and NMDA receptors. This differential effect is plain in Figure 5.

Reasonable estimates are available for the synaptic and extrasynaptic concentrations of AMPA and NMDA receptors. It is thought that hippocampal synapses will contain no more than ∼200 AMPA receptors (Nusser et al., 1998), which is equivalent to 400 μmwithin the PSD radius of the model synapse. The synaptic NMDA receptor concentration is assumed to be some 15% (60 μm) of this on the basis of well clamped recordings of EPSCs in cerebellar granule cells (Silver et al., 1992). The concentrations of these receptors outside the synapse of interest can be estimated from binding experiments to be 1 μm for AMPA receptors and ∼150 nm for NMDA receptors (Barbour and Häusser, 1997). A simulation was run with these receptor concentrations, and the receptor responses were compared to the control case of negligible receptor concentrations. The differences were minor, for both synaptic and spillover responses. The largest effect was an 8% reduction of the spillover NMDA response (data not shown).

On this basis it is concluded that realistic AMPA and NMDA receptor concentrations have very modest effects on the glutamate transient and receptor activation in and around the synapse and were neglected in the following simulations (negligible concentrations were present to monitor their activation). Other receptors and binding sites (apart from glutamate transporters) will not be modeled. For most of them it is not expected that they will be present in sufficient concentrations to have an important effect on the responses under investigation here. Possible exceptions to this rule may include the kainate-binding protein present at high concentrations in avian and fish cerebella (Gregor et al., 1989; Barbour and Häusser, 1997). Metabotropic receptors may also be present in significant quantities, particularly in the cerebellum. However, the simulations above with the simple binding site allow these other receptors to be ignored in the knowledge that their most likely action would be to reduce spillover responses slightly, while leaving synaptic responses largely unaffected.

### Glutamate transporters selectively limit spillover

The simulations above further confirm the widely accepted notion that dilution by diffusion is sufficiently rapid to terminate the synaptic action of glutamate. However, transporter concentrations are so high (Lehre and Danbolt, 1998) that it is natural to ask whether they could nevertheless influence the synaptic glutamate transient. A potentially important role of transporters is to limit the spread of glutamate from its site of release, thereby reducing spillover responses. To evaluate these transporter functions, the model was populated with transporters. The reaction scheme developed in Materials and Methods was used, and a quantitative study of transporter expression in the hippocampus guided transporter placement and the choice of concentration (Lehre and Danbolt, 1998).

Lehre and Danbolt (1998) concluded that the EAAT2 (GLT-1; the predominant form in the hippocampus) and EAAT1 (GLAST) transporter isoforms are present at a total concentration of 25 μm(which corresponds to 125 μm when referred to the extracellular volume). Simulations were therefore performed with 25 μm of glutamate transporters outside the synaptic cleft. In contrast, none were placed in the synaptic cleft, because in the hippocampus they appear to be expressed essentially only in glial cells, which are absent from the cleft (Lehre and Danbolt, 1998). Three values of *k*_{+1} (association) were tested: 5 × 10^{6}, 10^{7}, and 5 × 10^{7}(m sec). These were combined with values for *k*_{2}(trapping) of 1000 and 2000/sec.

Including transporters, whatever their kinetic properties, led only to very modest reductions of synaptic responses (<2% for AMPA receptors, data not shown, and < 8% for NMDA receptors) (Fig.6*A*). This confirms that dilution by diffusion alone dominates the process of terminating the synaptic action of glutamate, at least at the kind of small isolated synapses modeled here.

Much more significant reductions of spillover responses were obtained after the addition of transporters to the model. For *k*_{+1} = 10^{7}and *k*_{2} = 1000, AMPA receptor responses at 500 nm were reduced by 32% (data not shown), whereas NMDA responses were reduced by 58% (Fig. 6*B*). The effect on spillover increased with distance: NMDA responses at 1 μm were reduced by 85%. Paradoxically, the results obtained (particularly for NMDA responses), showed almost no sensitivity to the rate constants chosen (over the combinations tested). The explanation for this counter-intuitive observation is that the dissociation rate constant was always determined by the choices of association and trapping rate constants, to preserve the apparent steady-state affinity. This dependence of the parameters introduces a compensatory action that ultimately leaves the effect of the transporters unchanged, suggesting that the present predictions may be relatively insensitive to uncertainties in the two poorly constrained rate constants.

These results show that transporters selectively limit spillover responses, especially over greater distances, but do not interfere with orthodox synaptic transmission. Although transporters have a significant effect in the above simulations, two uncertainties about the transport mechanism could reduce that effect. First, it is implicit in the model that in the absence of glutamate all the transporters soon accumulate in the state that can bind and rapidly trap glutamate. If the whole transport reaction is considered, it is not obvious that this should be the case, although the benefits of the system working in this way are plain. It is certainly possible, for instance, that at rest a fraction of transporters are in states requiring the binding of one or more sodium ions before glutamate can bind and be translocated, and this might be expected to reduce the effectiveness of the transporters on a synaptic timescale. In fact, the available evidence provides no support for this notion. There are two independent lines of evidence. In recordings of transporters (hEAAT2) expressed in oocytes, Wadiche et al. (1995) were able to resolve sodium-dependent charge movements after voltage steps. These currents were attributed to the binding and unbinding of one or more sodium ions. Inspection of their data show that near the resting potential of a hippocampal glial cell (−93 mV;Bergles and Jahr, 1997), the voltage-sensitive sodium binding would be nearly saturated. Raising the extracellular sodium concentration from the amphibian (∼100 mm) to the mammalian (∼145 mm) value will further decrease the small fraction of transporters without bound sodium.

An alternative point of view is suggested by an explicit model of the complete reaction cycle of a glutamate transporter (Otis and Jahr, 1998). According to their model, although 93% of transporters will have bound at least one sodium ion at rest, only a small fraction of transporters will have bound all three required for glutamate to bind. However, all the rates of sodium binding (and unbinding) are very rapid indeed (the sodium concentration is hardly limiting… ), so the states immediately preceding the glutamate-receptive state can be considered at equilibrium with it under all relevant conditions. The final effect will be to reduce the apparent affinity for glutamate. However, because the parameters of the model are determined by measurement of the apparent affinity, any such effect is automatically taken into account. The model of Otis and Jahr (1998) also leaves open the possibility that protons must bind before glutamate, with a rate that may not allow the approximation of equilibrium under all conditions. However, this issue has been reexamined by Auger and Attwell (2000), who show that proton binding probably occurs in a different part of the transport cycle.

The second problem is that glutamate transporters may form multimers, and it is currently unclear whether the functional unit of transport is a pentamer (Haugueto et al., 1996; Eskandari et al., 2000) or the monomer, which has been the assumption here so far. If a single functional transporter is a pentameric assembly, the true concentration of transporters in the model must be divided by five, which will reduce their effect. The results of a simulation with one-fifth of the previous transporter concentration have accordingly been added to Figure 6, *A* and *B*.

If nearby transporters are saturated by release of a vesicle of glutamate, they will become ineffective at limiting spillover for succeeding release events. Recordings of transporter currents in hippocampal glia have indicated that transporters are not saturated in this way (Diamond and Jahr, 2000). This behavior is reproduced in the model (Fig. 6*C*), which shows that the levels of transporter “occupancy” are in fact quite low, even at the exit from the synaptic cleft.

### Duration of glutamate release and timing of receptor activation

Thus far, release of transmitter has been modeled here as an instantaneous event. Obviously, release will in fact be spread out over time, although little precise information is available concerning the time course of release from a small synaptic vesicle. It is therefore worthwhile to evaluate by simulation the influence of the duration of release on receptor activation. The issue of release duration has been lent importance by the suggestion that an acceleration of initially slow release could explain the transition from “silent” synapses (Isaac et al., 1995; Liao et al., 1995), with only an NMDA receptor-mediated response detectable, to “normal” synapses, with an AMPA receptor-mediated component also observable (Choi et al., 2000).

It is instructive to calculate exactly when AMPA and NMDA receptors are activated after instantaneous transmitter release. Here, “activation” is defined as the binding of a second glutamate molecule, because in both receptor models this is the necessary and sufficient condition for the channel to open (with a certain probability). By following the rate at which binding of a second glutamate occurs, it is possible to determine the time course of receptor activation. The rate at which receptor activation occurs is therefore given by the probability of being in the state with one bound glutamate multiplied by the rate for glutamate binding, thus:
Equation 16In the notation of Jonas et al. (1993) the monoliganded state of the AMPA receptor model is *C1*, and *k*[glu] is given by *k*_{+2}*c* (where*c* is the glutamate concentration). In the AMPA receptor model, a second glutamate could, in theory, also bind to the monoliganded desensitized state (*C3*). This flux was not included when calculating receptor activation and always negligible under the conditions of the present simulations. In the notation ofLester and Jahr (1992), the monoliganded state of the NMDA receptor is “AR”, and the corresponding rate constant is*k*_{on}. If the activation rate (Eq. 16) is integrated with respect to time, it is possible to calculate the number of receptor activations occurring before a certain time. This analysis was performed for instantaneous release (in the presence of glutamate transporters), and the “activation integral” was then normalized to its maximum (Fig.7*A*). Receptor activation is terminated extremely quickly. Some 70 (NMDA) to 80% (AMPA) of activated receptors are activated within ∼35 μsec of release. For AMPA receptors, 95% of all receptor activations occur by 160 μsec, whereas for NMDA receptors the same point is reached in 760 μsec. This reinforces just how quickly dilution of transmitter by diffusion is able to terminate the synaptic action of glutamate. For instantaneous release, very little AMPA receptor activation therefore takes place after the peak of the EPSC (which occurs at 440 μsec). If release takes place over a time longer than a few tens of microseconds, it is clear that the time course of receptor activation will largely follow the time course of release.

To investigate the influence of the time course of transmitter release, the receptor responses to the same amount of glutamate released at various uniform rates were calculated (Fig. 7*B,C*). For instantaneous release, the ratio of peak NMDA response to peak AMPA response (in terms of fraction of maximal*P*_{open}) is 1.21. Releasing the glutamate over a relatively short time, 100 or 300 μsec instead of instantaneously, had little effect on the overall aspect of the AMPA (or NMDA) responses. Nevertheless, the time-to-peak of the AMPA response increased (from 440 μsec) to 500 and 640 μsec, respectively, whereas the amplitude decreased slightly, by 7 and 19%, with respect to the response to instantaneous release. However, as the glutamate release is prolonged, the NMDA/AMPA response ratio indeed increases dramatically, e.g., to 41.3 for release over 10 msec. The mechanism is simply that the slower dissociation of glutamate from the NMDA receptor permits the accumulation of the occupancy of those sites over much longer times than for the AMPA receptor.

These simulations confirm the theoretical validity of the mechanism proposed by Choi et al. (2000). Interestingly, the simulations generate an additional prediction with which to test their hypothesis. A slowing of release sufficient to render the AMPA response undetectable might also be expected to generate an observable delay of the rising phase of the NMDA response. Thus, if induction of LTP causes an acceleration of initially slow glutamate release (Choi et al., 2000), the rising phase of the NMDA component should be advanced in time.

Prolonging release duration can reduce the specificity of synaptic activation, because the spillover/synaptic ratio of AMPA responses falls somewhat (approximately twofold for release lasting 10 msec; data not shown). However, the effect is small, and synaptic activation is reduced to very low levels. The ratio for NMDA responses is hardly altered.

### Multivesicular release impairs synapse specificity

It has long been assumed that small central synaptic contacts are not capable of releasing more than one vesicle of transmitter per presynaptic action potential. The inaccessibility of most central synapses and the difficulty of isolating a single contact have made this a remarkably difficult hypothesis to test. Nevertheless, there is some evidence that multivesicular release can occur (Auger et al., 1998). The family of response curves in Figure8 shows the results of multivesicular release (and can also be used to assess the effects of variations in the vesicular content of glutamate). When five vesicles are released instead of one, synaptic AMPA *P*_{open}rises from 0.12 to 0.44, and NMDA*P*_{open} rises from 0.053 to 0.21. Note that, despite the increasing amplitudes, there are only minor kinetic changes, reflecting the rapidity of transmitter dilution compared with receptor gating kinetics. Spillover responses for both receptor types are increased even more than the synaptic responses. The AMPA*P*_{open} at 500 nm rises from 0.00069 to 0.016, whereas the NMDA response at the same distance is increased from 0.00093 to 0.019. Synaptic independence is therefore impaired.

The synaptic responses for multivesicular release are likely to be somewhat underestimated here, because all the vesicles are released at the same point, leading to some mutual occlusion of their responses near the point of release. Multivesicular release at real synapses will usually involve release at different points over the PSD. A limit on this error is provided by a saturated synaptic response. Thus, for five vesicles, the true responses could not be >80% (AMPA) and 38% (NMDA) greater than those shown. If the synaptic responses are underestimated, the relative significance of spillover will be overestimated in proportion.

### Obstructing diffusion in the synaptic cleft enhances synaptic independence

It was argued above that the value of the effective diffusion coefficient for the extracellular space, which is modeled as a porous medium, has been tightly constrained by experiment. In contrast, very little is known about the properties of diffusion in the synaptic cleft. It has been suggested that synaptic clefts might contain sufficient material to retard diffusion on a microscopic or submicroscopic scale, much in the same way that brain tissue retards diffusion on a larger scale. Such an effect (approximately twofold) has been deduced for diffusion of acetylcholine at the neuromuscular junction and attributed to the basal lamina (Land et al., 1984). There are also reports of material within central synaptic clefts (Ichimura and Hashimoto, 1988). However, because the arguments in favor of retarded diffusion are rather indirect, especially as concerns central synapses, this must be regarded as an open question. Nevertheless, the effect of retarding diffusion within the cleft is interesting and potentially beneficial to synapse operation. This can be judged from the simulations of Figure 9, where diffusion within the synaptic cleft is modeled as occurring in a porous medium with different values of volume fraction and tortuosity. The decreased volume fraction concentrates the released transmitter while the tortuosity causes retention of transmitter in the cleft. Both these changes lead to increased receptor activation at the synapse. What is particularly interesting is that the extrasynaptic responses are hardly affected. This means that introducing obstacles to diffusion within the cleft would produce more synaptic receptor activation for a given amount of transmitter, an economy of transmitter. It would also improve the spatial selectivity of synaptic transmission, because the ratio of the synaptic to spillover responses is increased.

The combinations of cleft volume fraction and tortuosity chosen are physically plausible combinations, in the sense that approximately equivalent combinations have been measured for anoxic brain tissue (Chen and Nicholson, 2000). Thus, unless the obstacles have a very special structure, very dense packing of obstacles is required to generate high tortuosities. It should be noted that the extreme values of cleft volume fraction and tortuosity modeled above represent an improbably dense packing of obstacles to diffusion, denser even than that, on a larger scale, of brain tissue.

These simulations make clear the sensitivity of synaptic responses to the properties of the synaptic cleft. Some possible variations in the details in the cleft geometry have not been modeled here, because quantitative data on cleft geometry have not been reported: a narrower cleft; a cleft that is wider centrally and narrows toward the PSD edge. Both of these situations would be expected to increase synaptic responses, without affecting spillover responses; synaptic independence would be enhanced.

### Receptor desensitization after synaptic activation

AMPA receptors undergo profound desensitization. It is possible that desensitization at the synapse could alter the ratio of spillover to synaptic responses, so the possible influence of this mechanism was investigated. A brief application of glutamate to a patch containing AMPA receptors can induce quite a large reduction of the amplitude of the response to a second application, despite little desensitization appearing to occur while channels are open during the first pulse (Colquhoun et al., 1992). In the AMPA receptor model used here this phenomenon is accounted for by a high probability of the receptor entering the desensitized state when only one glutamate is bound. Such observations have led to the expectation that desensitization should similarly provide a depressive mechanism during paired-pulse synaptic stimulation. Surprisingly, at least in some preparations, little or no desensitization could be observed (Stevens and Wang, 1995; Debanne et al., 1996; Hjelmstad et al., 1999). In contrast, Otis et al. (1996)have shown that at synapses in the chick cochlear nucleus, a compound EPSC will desensitize ∼40% of receptors.

Investigation of synaptic desensitization is often confounded by the profound presynaptic changes that occur during paired-pulse stimulation, so the present model was used to evaluate the expected contribution of receptor desensitization during paired pulse stimulation. The release of two vesicles separated by various intervals was simulated (Fig. 10), and the ratio of the AMPA receptor responses calculated. Desensitization has a modest effect (15% at a 10 msec interval). This is in fact expected, given the low levels of activation induced by one vesicle: receptors that do not bind glutamate cannot desensitize. Moreover, this represents a slight overestimate of the effect of desensitization, because normally a second vesicle would not be released exactly at the same point as the first and would therefore usually face even fewer desensitized receptors. These results provide a possible reconciliation between the desensitization behavior observed in patches and the apparent lack of desensitization at hippocampal synapses, without the need to postulate different kinetic properties for receptors in the two preparations. Nevertheless, even such a small effect might have been detectable in some of the preparations (Stevens and Wang, 1995; Hjelmstad et al., 1999), suggesting that synaptic AMPA receptors may desensitize more slowly than those in patches or that receptor occupancy is even lower than predicted here. In any case, the many different types of AMPA receptors probably have different desensitization behavior (Liu et al., 1999). Naturally, multivesicular release or the obstruction of diffusion would increase desensitization.

### Calculating the average spillover response

To summarize the various factors that can affect synaptic independence, the peak *P*_{open} values for AMPA and NMDA receptors are shown in Figure11, *A* and *B*, as a function of the distance from the point of release. Several simulations are plotted: a single vesicle in the absence of uptake, showing the effect of diffusion alone; a single vesicle in the presence of uptake; three vesicles in the presence of uptake; and a single vesicle with an obstructed synaptic cleft in the presence of uptake. AMPA*P*_{open} falls off very rapidly with distance in all the simulations, and it seems unlikely that activation of AMPA receptors by spillover glutamate will play an important role. The situation for NMDA receptors is less clear, however, and a more careful analysis is required to assess the effects of spillover.

To simplify comparison between different conditions, the curves for NMDA receptor *P*_{open} have been normalized to their respective average synaptic values (Fig.11*C*). To estimate average spillover responses, the calculated *P*_{open} values must be combined with knowledge of intersynaptic distances. Unfortunately, although a distribution of distances to nearest-neighbor synapses for hippocampus has been published (Rusakov et al., 1999), quite large discrepancies exist between the synaptic density reported (1.25/μm^{3}) and that observed in another study, in hippocampal slices (∼3.5/μm^{3}; Sorra and Harris, 1998). It will therefore not be possible to obtain a unique result. Estimated distributions of distances to nearest neighbor synapses were obtained as follows. It was first assumed that synapse positions are random (i.e., a Poisson process). However, in this distribution, synapses can be arbitrarily close, which is of course impossible. A region around the liberating synapse was therefore cleared of synapses within a certain radius (without altering the rest of the distribution). For consistency with the size of synapse modeled here, that radius should be 360 nm. However, synapses do vary in size, and shorter intersynaptic distances are found. To allow for this, a compromise boundary of 250 nm was chosen. The probability density function for the distance to the nearest neighbor (*nnd*) for this distribution is given by:
Equation 17where *r*_{0} = 250 nm as just described, and *N*_{V} is the density of synapses in the original Poisson distribution. This distribution was calculated for the two reported synapse densities (Fig.11*C*). Using these distributions, it is now possible to calculate the mean receptor activation for nearest neighbor synapses. In the order (1) one vesicle without uptake, (2) one vesicle with uptake, (3) three vesicles with uptake, and (4) one vesicle with uptake and an obstructed cleft, the mean spillover/synaptic NMDA response ratios were: (1) 0.07, (2) 0.04, (3) 0.10, and (4) 0.02 for the dense synapse packing; (1) 0.04, (2) 0.02, (3) 0.06, and (4) 0.01 for the lesser synapse density. The results show that NMDA receptor activation resulting from spillover can be considered unimportant, even in the absence of glutamate uptake, except perhaps if significant multivesicular release occurs.

## DISCUSSION

### Synapses can operate independently

The present synaptic model, which takes into account the three-dimensional extrasynaptic space and incorporates the correct diffusion coefficient for glutamate, predicts that synaptic receptors are far from saturated by a single vesicle of glutamate and that spillover responses are small in comparison. This suggests that typical small glutamatergic synapses can operate independently. However, whether such synapses do always operate independently, particularly as regards NMDA receptors, remains a finely balanced question. The above conclusions apply at room temperature (at which most of the parameters of the model were determined) for synaptic transmission that is isolated, both in space and time. The ability of NMDA receptors to integrate activation over tens or even hundreds of milliseconds means that it will be necessary to know more about the density and frequency of synaptic activations in the brain before deciding whether spillover has a role to play in activating NMDA receptors. It is therefore interesting that only a small fraction of the total number of synaptic inputs that a pyramidal cell receives need to be activated to excite it (Sayer et al., 1990; Larkman et al., 1991; Otmakhov et al., 1993). Moreover, *in vivo* recordings suggest that only a small fraction of hippocampal pyramidal cells are active at any one time (Henze et al., 2000). Together, these observations suggest that most synapses may be inactive most of the time under physiological conditions.

### Comparison with experiment

Because no parameters of the model were based on synaptic recordings, it is valid to compare the predictions with synaptic data. That release of one vesicle does not saturate AMPA receptors is now established for several synapses (Silver et al., 1996; Forti et al., 1997; Liu et al., 1999; McAllister and Stevens, 2000), although most studies have been performed in culture. The present predictions for AMPA receptor occupancy are somewhat lower than the 40% reported byLiu et al. (1999) or the upper limit of ∼50% of Silver et al. (1996). If this difference is confirmed, it seems likely that the geometry of the cleft or obstacles in it enhance synaptic responses somewhat. In the model, NMDA receptors are only activated to a slightly greater extent than AMPA receptors. Two-photon imaging of NMDA receptor-mediated spine calcium signals has provided direct evidence that the release of one vesicle is insufficient to saturate NMDA receptors (Mainen et al., 1999), a result also indicated by the variability and correlation of the AMPA and NMDA components of mEPSCs (McAllister and Stevens, 2000). Comparison of the NMDA receptor component of mEPSCs and responses to iontophoretic agonist applications at putative isolated synapses in culture (McAllister and Stevens, 2000) further suggested that an unexpectedly small fraction of NMDA receptors is activated by a single vesicle, a result that is consistent with the present simulations.

The present prediction of a small spillover effect may appear to contradict experimental results consistent with this phenomenon. However, the model is designed to simulate small central synaptic contacts operating in isolation, whereas the available experimental evidence concerns large synaptic contacts or widespread synaptic activity, for which the synergistic interaction of transmitter released from numerous sources is possible. Such effects may play interesting network roles physiologically (Scanziani, 2000). Recordings of glia have provided direct evidence for glutamate activating receptors outside the synaptic cleft (Mennerick and Zorumski, 1994; Bergles et al., 1997; Clark and Barbour; 1997; Dzubay and Jahr, 1999). However, receptor occupancy was not determined in those studies, and all of the observations were made when many vesicles were released simultaneously, making direct comparison with the present model difficult.

A powerful technique for determining the transmitter time course to which receptors are exposed is to measure the inhibition of a synaptic response by a rapidly unbinding (low-affinity) competitive antagonist. This technique was applied in a pioneering study using NMDA receptors (Clements et al., 1992), in which it was deduced that synaptically activated receptors were exposed to a glutamate transient that, when assumed to follow an exponential time course, was best approximated by a peak amplitude of 1.1 mm and a decay time constant of 1.2 msec. The same data were later fitted by a double exponential decay (Clements, 1996). An analogous experiment has since been performed for AMPA receptors (Diamond and Jahr, 1997). The calculated time courses are rather different from the theoretical predictions here (Fig.12*A*). This is not too surprising, because the expected glutamate concentration time course is not biexponential, but instead approximately proportional to*t*^{−3/2}. Because up to 99% of glutamate can leave the PSD region within 100 μsec (Fig.12*A*), before much antagonist has unbound (Clements et al., 1992), the method is also relatively insensitive to events occurring immediately after release. The slow decays that were reported are generally interpreted as showing that diffusion of glutamate out of the synaptic cleft is in some way slowed. However, Wahl et al. (1996)have already pointed out that such slowing would have to be quite significant, even in an infinite disk model. If release is assumed to be instantaneous in the present model, it is not possible to recover the slow time courses with plausible slowing of diffusion within the synaptic cleft (Fig. 12*A*). Although extrapolation of the present model to culture synapses may not be wholly accurate (synapse size could make a difference), another or an additional explanation is probably required to explain the difference. If release is prolonged, the slow time courses can be approximately recovered, but most of the glutamate released (>60%) must be released slowly (time constant, 300 μsec to 2 msec). These estimates were obtained by convolving the glutamate time course of the model with monoexponential or biexponential release time courses (data not shown). If diffusion is unhindered in the cleft, several vesicles are required for recovery of the published time courses. If diffusion is obstructed (as shown), less than one vesicle is required, but the kinetics of the release function are not faster than stated above. In the case of Clements et al. (1992), one possible alternative source of glutamate at late times could be neighboring synapses.

### Comparison with other models

Compared with the numerous infinite disk models of synaptic activation, the simulations here predict lower levels of receptor activation. This is expected on the basis of the faster dilution of transmitter in a three-dimensional space, coupled with the choice of diffusion coefficient.

In a recently published model, Rusakov and Kullmann (1998b) and Rusakov et al. (1999) make a very similar choice of geometry to the present simulations. Under certain conditions they predict quite significant extrasynaptic NMDA receptor activation via glutamate spillover after release of a single vesicle (e.g.,*P*_{open} of ∼0.03, 25 times the 0.0013 calculated here; Rusakov et al., 1999). A good deal of the discrepancy can be explained by those authors' use of a slower diffusion coefficient, and other different parameter choices also contribute. However, a discrepancy of about a factor of 2 for NMDA receptor activation at 0.5 μm remains even when their parameter values are adopted in the present simulations. Figure 12*B*replots a glutamate transient of Rusakov et al. (1999; their Fig. 3) and that produced when their published parameter values are used with the present model. Also plotted is the theoretical curve calculated for diffusion in a porous medium (i.e., without the synaptic cleft and transition regions). The theoretical curve is very close to the transient calculated in the present simulations. The time course ofRusakov et al. (1999) reaches almost double the concentration. A similar excess is observed for all the radii shown in their figure. It is difficult therefore to avoid the conclusion that their model contains >5000 glutamate molecules. It is not possible to identify precisely the source of this additional glutamate, but the description of the conditions imposed at the boundary of the synaptic cleft (Rusakov and Kullmann, 1998b) suggests that an imbalance of centripetal and centrifugal fluxes there might “create” glutamate. The larger glutamate transient will, of course, increase receptor activation and the difference probably accounts for the remaining discrepancy with the present model.

### Conclusion

It is extremely rare in biological modeling to be able to determine all the relevant parameters with the accuracy implied by the explicit calculations practiced on them. The author willingly acknowledges that this is still the case here. Critical issues remain: the properties of diffusion in the synaptic cleft; better information about the number of vesicles released, their glutamate content and how it is released; confirmation that synaptic receptors *in situ*behave as has been determined for extrasynaptic receptors in patches. Because models (such as this one) based on the porous medium approximations can only describe average diffusion behavior, simulations based on explicit tissue reconstructions will provide an important test of the present estimates and yield invaluable information regarding variations between synapses.

## Appendix

Rusakov and Kullmann (1998a,b) report a method for calculating the “geometrical” tortuosity of a porous medium in which diffusion is restricted to a network of thin lamellae. The method has attracted criticism from Chen and Nicholson (2000), who describe the difficulty of such attempts to calculate tortuosity *ab initio*. In addition, the derivation contains a number of problems that suggest it cannot be accurate. The following arguments will refer to the longer exposition in Rusakov and Kullmann (1998a).

Rusakov and Kullmann (1998a) define an “elementary geometrical tortuosity” (EGT; their Eq. 1) using a plane element of extracellular space oriented at an arbitrary angle to the imposed “direction” of diffusion. The tortuosity (following Chen and Nicholson, I define λ_{RK} as the tortuosity calculated by Rusakov and Kullmann) is then to be obtained by averaging the EGTs throughout the porous medium. An obvious difficulty with this first step of the calculation is that the EGT can assume very large or even infinite values for some angles (in their notation, β_{2} = 90°). Although this behavior does not prevent the subsequent average from remaining bounded, the large values are unrealistic, and the result therefore cannot be exact. A serious omission in the derivation is the absence of any treatment of the behavior of diffusion at branches of the network of lamellae, although the connection properties might be expected to have an important influence on diffusion through the structure.

For a porous medium that is by assumption uniform and isotropic, the unusual claim is made that tortuosity will be “different, depending on whether it is calculated with respect to unidirectional or radial diffusion”. Normally, in a uniform isotropic structure the value of the diffusion coefficient (and therefore the tortuosity) is expected to be constant and independent of both the form of the diffusion source and the coordinate system. The authors make their statement because they average the EGT over the structure in two different ways and obtain two different values. However, one method is erroneous; both methods should give the same result.

The EGT is defined on a plane element whose orientation can be defined by two angles. λ_{RK} is obtained by averaging many EGTs. The authors state that all orientations are equally likely in an isotropic medium (although this may be invalid for the orientations encountered along a path) and therefore in the average equal weight should be given to all orientations. This average can be calculated as an integral over the surface of a sphere (divided by its surface area). Accordingly, their Equation 1 is transformed into Equations 3 and 4 and integrated to yield λ_{RK} = π/2 (in fact, the integral is most simply evaluated using a surface of revolution). The authors then proceed to transform their Equation 3 to spherical coordinates (Eq. 5) and integrate it again, this time by a Monte Carlo method, to obtain their estimate of tortuosity for radial diffusion, 1.398, which differs from the value above. However, because the transformation to spherical coordinates amounts to nothing more than a change of variable in the integration, the result should be unaffected. In fact, in the Monte Carlo integration, all combinations of the two angles defining the plane element orientations, rather than the orientations themselves, are considered equally likely. This is a different integration that, not surprisingly, gives a different result. Put another way, the surface integral should use an elementary surface area of *dS* =*r*^{2}sin*θdθdφ* , but the Monte Carlo method implicitly uses *dS* =*r*^{2}*dθdφ*. Including the factor of sinθ enables recovery of the original result (λ_{RK} = π/2).

## Footnotes

This work was supported by the Centre National de la Recherche Scientifique (Unité Mixte de Recherche 8544). I thank Philippe Ascher and the members of the Laboratoire de Neurobiologie for their generous support and Steven Traynelis for helpful discussion. For their comments on this manuscript, I thank Etienne Audinat, Philippe Isope, Jacques Neyton, and Stéphane Supplisson (who also helped operate Mathematica).

Correspondence should be addressed to Boris Barbour, Laboratoire de Neurobiologie, Centre National de la Recherche Scientifique, Unité Mixte de Recherche 8544, Ecole Normale Supérieure, 46 rue d'Ulm, 75230 Paris Cedex 05, France. E-mail: barbour{at}ens.fr.