help button home button Biophys. J.
HOME HELP FEEDBACK SUBSCRIPTIONS ARCHIVE SEARCH TABLE OF CONTENTS
This Article
Right arrow Abstract Freely available
Right arrow Full Text (PDF)
Right arrow Alert me when this article is cited
Right arrow Alert me if a correction is posted
Services
Right arrow Similar articles in this journal
Right arrow Similar articles in PubMed
Right arrow Alert me to new issues of the journal
Right arrow Download to citation manager
Right arrow reprints & permissions
Citing Articles
Right arrow Citing Articles via HighWire
Right arrow Citing Articles via Google Scholar
Google Scholar
Right arrow Articles by Smart, J. L.
Right arrow Articles by McCammon, J. A.
Right arrow Search for Related Content
PubMed
Right arrow PubMed Citation
Right arrow Articles by Smart, J. L.
Right arrow Articles by McCammon, J. A.

Biophys J, October 1998, p. 1679-1688, Vol. 75, No. 4

Analysis of Synaptic Transmission in the Neuromuscular Junction Using a Continuum Finite Element Model

Jason L. Smart and J. Andrew McCammon

Department of Chemistry and Biochemistry, Department of Pharmacology, University of California at San Diego, La Jolla, California 92093-0365 USA

    ABSTRACT
Top
Abstract
Introduction
Methods
Results
Conclusions
References

There is a steadily growing body of experimental data describing the diffusion of acetylcholine in the neuromuscular junction and the subsequent miniature endplate currents produced at the postsynaptic membrane. To gain further insights into the structural features governing synaptic transmission, we have performed calculations using a simplified finite element model of the neuromuscular junction. The diffusing acetylcholine molecules are modeled as a continuum, whose spatial and temporal distribution is governed by the force-free diffusion equation. The finite element method was adopted because of its flexibility in modeling irregular geometries and complex boundary conditions. The resulting simulations are shown to be in accord with experiment and other simulations.

    INTRODUCTION
Top
Abstract
Introduction
Methods
Results
Conclusions
References

We are interested in the diffusional aspects of synaptic transmission. There have been many simulations studying the kinetic aspects of neurotransmitter/receptor binding and subsequent endplate currents (Wathey et al., 1979; Friboulet and Thomas, 1993; Khanin et al., 1994; Kleinle et al., 1996; Agmon and Edelstein, 1997; Edelstein and Agmon, 1997), but less work dealing with neurotransmitter diffusion in detailed models of the synapse (Bartol et al., 1991; Stiles et al., 1996; Bennett et al., 1997). Our simulations will address the issues of synaptic transmission from a continuum standpoint, employing finite element methods to solve the diffusion equation in a detailed model of the neuromuscular junction.

The advantages of continuum finite element methods include the ability to model complicated geometries, study long timescales with adaptive timesteps, include background concentrations of the neurotransmitter acetylcholine (ACh), include models for the individual clusters of acetylcholinesterase (AChE), and simulate the release of many synaptic vesicles (containing ACh) simultaneously. However, before taking full advantage of the method, it is necessary to establish its accuracy and reliability by comparing the results to other simulations, and to experiment.

The purpose of this work, therefore, is to validate the continuum approach for simulating neurotransmitter diffusion in the neuromuscular junction (NMJ). Additionally, this article examines how the different components of the vertebrate NMJ, such as the AChE, secondary folds, and release pore, affect the amplitude and timecourse of the transmitted signal.

Proper modeling of AChE and the acetylcholine receptors (AChR) is vital in reproducing the timescales and extent of ACh diffusion. Thus, after describing the computational methods used to simulate the diffusion of neurotransmitter, the AChE and AChR models will be discussed. Then, several test cases will be used to establish the reliability of the finite element model and the partial differential equation solver. Finally, the results of the simulations will be compared to other simulations and to experiment.

    METHODS
Top
Abstract
Introduction
Methods
Results
Conclusions
References

For all simulations, ACh diffusion will be governed by the diffusion equation
D∇<SUP>2</SUP>C=dC/dt (1)
where D is the diffusion constant of ACh in the NMJ, C is the concentration of ACh, and t is time.

From fits of kinetic and diffusional models to experimentally measured miniature endplate currents, the diffusion constant of ACh in the NMJ was estimated to be (Land et al., 1981, 1984) D = 4.0 × 10-6 cm2/s. Note that the diffusion constant accounts, in an averaged fashion, for the irregularities in the density of the synaptic medium.

Due to the irregular geometry of the NMJ, and the complexity of our AChE model (to be discussed below), it was useful to adopt the finite element method for the simulations (Becker et al., 1981; Martin and Carey, 1973; Brenner and Scott, 1994). The finite element method divides a simulation region into individual elements, and solves the given partial differential equation (PDE) at the nodes of those elements. Because the choice of element size and shape is highly flexible, the finite element method can handle complicated geometries.

In addition, the finite element mesh can be adaptively refined in regions where accurate solutions of the PDE are needed. Adaptive refinement increases the number and density of individual elements in areas where large local errors can occur, and thus reduces the errors associated with discretizing and solving the PDE. In the NMJ simulations, large gradients in the ACh concentration can lead to excessive discretization errors. Numerical techniques lacking adaptive refinement, such as standard finite difference techniques, therefore cannot be used for the NMJ simulations. The Test Cases section will further examine the reliability of the adaptive refinement algorithm.

All simulations discussed below use the finite element package Kaskade (Beck et al., 1995), with meshes consisting of tetrahedral elements. The actual geometries of the NMJ models vary from simulation to simulation, and will be described in detail later.

Neumann, or reflecting, boundary conditions were used on most boundaries of the NMJ model
∇C‖<SUB>0</SUB>=0 (2)
Surfaces with Neumann boundary conditions do not react with the diffusing neurotransmitter.

Acetylcholinesterase models

Background

The AChE populate vertebrate NMJ's at densities of 2000-3000 µm-2, both in the primary cleft and in the secondary folds (Salpeter et al., 1972, 1978; Salpeter, 1987; Anglister et al., 1995). They are found in clusters of three tetramers, attached by stalks to the postsynaptic membrane (Vincent and Wray, 1990; Zimmerman, 1993). In our NMJ model, the AChE clusters are represented by cubes, 16 nm on a side, spaced 64 nm apart. Within the continuum model, the proper choice of boundary conditions is vital in representing the true reactivity of these enzyme clusters with substrate. Although it is clear from our initial studies that mixed, or radiating, boundary conditions are the most appropriate, choosing the correct reactivity constant is essential. The mixed boundary conditions will be defined below.

The interactions between AChE and ACh will be modeled using the following kinetic scheme:
E+S → E+P (3)
where S represents our continuum ACh concentration. The corresponding rate equation is
<FR><NU>d[S]</NU><DE>dt</DE></FR>=<UP>−</UP>k[S][E]. (4)
The concentrations are in molar units, and the rate constant is in units of M-1 s-1.

To make contact between these standard kinetic relations and our continuum model, we must define the current of substrate into a single reactive "enzyme." For a single enzyme, modeled for simplicity as a sphere, the current of substrate flowing into the sphere is given by (Rice, 1985; Collins and Kimball, 1949)
I(R)=<UP>−</UP>4&pgr;R<SUP>2</SUP> D<FENCE><FR><NU>dC(r)</NU><DE>dr</DE></FR></FENCE><SUB><UP>R</UP></SUB> (5)
where R is the radius of the sphere, D is the diffusion constant, and C(r) is the substrate concentration at a distance r from the center of the sphere. The concentration is in molecular units.

The rate of substrate removal is then simply equal to the inward current (given by Eq. 5), times the enzyme concentration [E]
<FR><NU>d[S]</NU><DE>dt</DE></FR>=I(R)[E] (6)
which leads to the useful relation (compare Eq. 6 to Eq. 4)
I(R)=<UP>−</UP>k[S]<SUB>0</SUB> (7)
Here it is assumed that the substrate concentration is in excess of the enzyme concentration, and can be replaced by its initial value [S]0. [S]0 also corresponds to what would be measured at a large distance from the enzyme.

The simple relation defined by Eq. 7 provides the necessary connection between the experimental rate k and the ACh current into each enzyme in our continuum model. Although derived under the assumption of spherical symmetry, the relation is valid for an enzyme of any shape.

Experimental rates

There have been several experimental measurements of kcat and KM for the interaction of ACh with human AChE (Radic et al., 1992; Shafferman et al., 1995), yielding values of kcat/KM in the range of 2.5-3 × 109 M-1min-1. Values of kcat/KM for other species lie in a similar range. Assuming these experimental measurements obey our kinetic scheme, Eq. 3, the substrate concentration is governed by
d[S]/dt=<UP>−</UP><FR><NU>k<SUB><UP>cat</UP></SUB></NU><DE>K<SUB><UP>M</UP></SUB></DE></FR> [E][S]. (8)
kcat/KM is thus equal to the experimental rate constant k used below in determining the continuum reactivity constant, kact.

Motivation for the mixed boundary conditions

Mixed boundary conditions make the approximation that the reactivity of the enzyme is proportional to the probability that a reaction pair exists (Rice, 1985). If we equate the reactivity to the inward flux of substrate, we have an expression for our boundary condition:
4&pgr;R<SUP>2</SUP>D<FENCE><FR><NU>dC(r)</NU><DE>dr</DE></FR></FENCE><SUB><UP>R</UP></SUB>=k<SUB><UP>act</UP></SUB>C(R) (9)
where kact is the proportionality constant. The above statement of the boundary condition again assumes spherical symmetry. However, the 4pi R2 area term may be absorbed into the unknown kact to restate the boundary conditions in a form that is independent of the shape of the enzyme.

The above boundary condition will be used in the simulations to represent the reactive AChE clusters. It is therefore necessary to relate our unknown reactivity constant kact to the experimental forward-binding rate constant k.

Model of a single esterase

The direct approach to finding kact is to perform simulations of a single AChE molecule in a cluster, and determine what value of kact gives the experimental rate k. To do this, one AChE cluster was placed in a large simulation box, and the current of substrate into a portion of the cluster was measured. The steady-state current is given by the surface integral of the flux
I(r)=<UP>−</UP>D<LIM><OP>∫</OP></LIM><UP><B>dS</B></UP> · ∇C(r) (10)
I(r) may be equated to the experimental current calculated from Eq. 7
I(r)=<UP>−</UP>k[S]<SUB>0</SUB> (11)
The simulation thus involves varying k'act in
D<FENCE><FR><NU>dC(r)</NU><DE>dr</DE></FR></FENCE><SUB><UP>R</UP></SUB>=k′<SUB><UP>act</UP></SUB>C(R) (12)
until the measured current matches experimental results. Here the mixed boundary condition is written in its most general form, with the 4pi R2 term absorbed in the new constant k'act.

As in the full NMJ model, the AChE cluster is represented by a cube measuring 16 nm on a side. Since we are examining the current into a single AChE molecule, only 1/12 of the cube's surface will be reactive (using mixed boundary conditions). The remaining surface will have Neumann boundary conditions (dC/dr = 0). The full simulation region is 144 nm on a side, also with Neumann boundary conditions. The simulation begins with a constant substrate concentration of [S]0 = 1.66 M, though any initial concentration may be used. The simulation was performed with version 3.1 of Kaskade (Beck et al., 1995), using a 6000 tetrahedra mesh. Throughout the simulation, the substrate concentration at the outer edges of the box never dropped below 99.9% of [S]0, thus assuring us that the simulation region is sufficiently large.

Given [S]0 above, and the experimental rate constant k, Eq. 11 yields
I(r)=<UP>−</UP>k[S]<SUB>0</SUB>=<UP>−</UP>0.083 <UP>ns<SUP>−1</SUP></UP> (13)
By trial-and-error, the best value of k'act was
k′<SUB><UP>act</UP></SUB>=1.4×10<SUP><UP>−3</UP></SUP> <UP>nm ns<SUP>−1</SUP></UP> (14)
The predicted k'act will be used in the full NMJ simulations discussed below.

Acetylcholine receptor models

In the early stages of the NMJ simulations, the sole purpose of the receptor model is to define the postsynaptic region over which receptors are saturated by substrate. The area of saturation then gives some indication of the effects of junctional folds, AChE, and the release model, on the timecourse and extent of ACh diffusion. Since postsynaptic coverage is only meaningful if it is connected in some way with experimental observations, it is important to base our definition of coverage on the binding of ACh to receptors. Specifically, coverage will be defined as the area over which receptors are doubly ligated. In this way, postsynaptic coverage (i.e., receptor saturation) gives some indication of the amplitude of the miniature endplate current.

Model of a single receptor

The nicotinic AChR contains five domains, two of which (alpha  domains) bind ACh (Taylor et al., 1986). When ACh is bound to both alpha  domains, a central ion channel opens, yielding a current of ~2-5 pA (Dionne and Leibowitz, 1982; Colquhoun, 1990).

In our model, binding of ACh to a given receptor is determined by measuring the current into the receptor. The receptor is again modeled with mixed boundary conditions
D<FENCE><FR><NU>dC(r, t)</NU><DE>dr</DE></FR></FENCE><SUB><UP>R</UP></SUB>=k′<SUB><UP>act</UP></SUB>C(R, t) (15)
where the constant k'act determines the reactivity (and current). In order to determine k'act for a receptor, the same computational approach applied to AChE was used.

For simplicity, our AChR model contains only the two binding sites, modeled as reactive patches. The steady-state current into each AChR is given by
I(∞)=<LIM><OP>∫</OP></LIM> dx dy k′<SUB><UP>act</UP></SUB>C(x, y; ∞) (16)
integrated over both subunits (over both reactive patches), on the membrane surface.

The experimental current is based on the kinetic scheme (Dionne and Leibowitz, 1982)
2A+R <LIM><OP><ARROW>⇌</ARROW></OP><LL><SUB>k<SUB><UP>−</UP></SUB></SUB></LL><UL><SUB>2k<SUB><UP>+</UP></SUB></SUB></UL></LIM> A+AR <LIM><OP><ARROW>⇌</ARROW></OP><LL><SUB>2k<SUB><UP>−</UP></SUB></SUB></LL><UL><SUB>k<SUB><UP>+</UP></SUB></SUB></UL></LIM> A<SUB>2</SUB>R <LIM><OP><ARROW>⇌</ARROW></OP><LL><SUB>l<SUB><UP>−</UP></SUB></SUB></LL><UL><SUB>l<SUB><UP>+</UP></SUB></SUB></UL></LIM> A<SUB>2</SUB>R* (17)
where A represents ACh, R represents the nicotinic ACh receptor in the closed state, and R* represents the receptor in the open state. For lizard AChR, k+ = 3.3 × 107 M-1 s-1 (Land et al., 1984).

The experimental current into a single alpha  domain is
I<SUB><UP>expt</UP></SUB>(∞)=<UP>−</UP>k<SUB><UP>+</UP></SUB>[S]<SUB>0</SUB> (18)
where [S]0 is the initial substrate concentration.

To determine k'act for our model, a single receptor was placed at the base of a large simulation box. Only half of the receptor was active, to represent a single alpha  domain. The substrate current into the AChR subunit was then measured as a function of k'act. The correct value of k'act is the one for which the current is equal to the experimental result (Eq. 18). In our case, k'act = 9.0 × 10-4 nm/ns reproduced the correct experimental ACh current into the receptor. This constant was then applied to the AChR models implemented in the full simulations.

AChR in the full NMJ model

In the NMJ, the AChR density is highest at the tops of the folds, ~7,500-10,000/µm2 = 1 AChR/100-133 nm2. This gives a spacing of ~10-12 nm between each AChR molecule. The AChR density drops off sharply as one approaches the bottoms of the folds. For our simulations, the AChR are evenly spaced 11.3 nm apart at the tops of the junctional folds, and taper off down the sides of the secondary folds.

The current of substrate flowing into each AChR during the simulation is measured using
I(t)=<LIM><OP>∫</OP></LIM> dx dy k′<SUB><UP>act</UP></SUB>C(x, y; t) (19)
again integrated over both subunits. Here I(t) represents the time-dependent ACh current into a single receptor.

The precise definition of "postsynaptic coverage" is the area over which the total integrated current for each receptor is larger than two:
N(t′)=<LIM><OP>∫</OP><LL>0</LL><UL><UP>t′</UP></UL></LIM>dt I(t)>2 (20)
In other words, receptors with N(t') > 2 are considered to be in the doubly ligated state, and thus contribute to the coverage area.

Note that unbinding events occur on significantly longer timescales than the diffusion of neurotransmitter [~2 ms (Beeson and Barnard, 1990)], and are unlikely to affect substrate concentrations over the course of our NMJ simulations. After unbinding from receptors, it is believed that most ACh molecules are immediately hydrolyzed by AChE, and therefore make no further contributions to the observed endplate currents (Land et al., 1984; Salpeter, 1987). It is safe, then, to ignore unbinding events in our receptor model. Future simulations exploring longer timescales will include unbinding events.

An inside view of our standard NMJ model is shown in Fig. 1. The scale of the model is based on measurements of lizard NMJ's (Salpeter, 1987). The finite element mesh is comprised of 9248 vertices and 38772 tetrahedra.


View larger version (72K):
[in this window]
[in a new window]
  		FIGURE 1   Inside view of our standard NMJ model, containing a release cleft near the center of the image, and secondary folds at the bottom. The cubes represent the acetylcholinesterase clusters. The acetylcholine receptors (not shown) cover the postsynaptic membrane at the bottom of the image. For illustration, the finite element mesh is displayed on the surfaces of the model.

Test cases

Several test cases were used to verify that the diffusion equation was implemented correctly in Kaskade, and that the triangulation of the NMJ was sufficiently detailed to yield a robust solution.

Diffusion to a sphere

The diffusion equation was first solved for the case of a sphere with perfectly [C(r = R, t) = 0] or partially {D[dC(r, t)/dr]|R = k'actC(R, t)} absorbing boundary conditions, initially surrounded by a medium of constant concentration [C(r, 0) = C0, r > R]. These models are analogous to the interaction of ACh with a single AChE cluster. For both boundary conditions, the diffusion equation can be solved analytically. The numerical results compared well with analytic predictions (data not shown), demonstrating the method's accuracy in solving the diffusion equation and in handling complex boundary conditions.

Diffusion from a point source

In our model of the NMJ, the ACh is released from a small, concentrated volume near the presynaptic membrane. To establish whether our NMJ triangulation, and the adaptive mesh refinement capabilities of Kaskade, can sufficiently handle this type of release, diffusion from a point source was examined. Point source diffusion leads to large concentration gradients just after release, and therefore requires a robust finite element mesh and a finite element solver capable of adaptive mesh refinement. A similar test case was employed by Bartol et al. to test their Monte Carlo implementation (Bartol et al., 1991).

The test system was based on our NMJ model, with AChE, secondary folds, and the release cleft removed. ACh was released from a dense 103 nm3 region at the top-center of the model. This initial condition was chosen to replicate the true initial conditions of our full NMJ simulations. The radial profile of the concentration should be in rough accord with the analytic solution for diffusion from a point source in a hemisphere (Bartol et al., 1991; Rice, 1985)
C(r, t)=<FR><NU>N</NU><DE>[4&pgr;D(t+44)]<SUP>3/2</SUP></DE></FR><UP>exp</UP>[<UP>−</UP>r<SUP>2</SUP>/4D(t+44)] (21)
where N is the number of particles of ACh released (N = C0 · 103 nm3 = initial concentration · release volume). The 44-ns constant adjusts the analytic solution so that at time zero it has expanded roughly 7-8 nm radially, close to the initial volume of our simulated release.

The simulation and the analytic solution (Eq. 21) produce similar radial concentration profiles (data not shown), even though the initial conditions are slightly different. More importantly, the rate of expansion of the simulated neurotransmitter concentration is in good agreement with the analytic solution. The reasonable agreement establishes that the timescales for diffusion from a dense source are correctly reproduced by the finite element model.

    RESULTS
Top
Abstract
Introduction
Methods
Results
Conclusions
References

Transit time of acetylcholine across the synaptic cleft

After the release of ACh from vesicles at the presynaptic membrane, it takes ~5 µs for the bulk of ACh to reach the postsynaptic membrane (Edmonds et al., 1995). As an initial test of our NMJ model, the free diffusion time across the cleft was measured for a single quanta of ACh.

The area of the postsynaptic membrane over which the ACh concentration is >1% of the initial release concentration was measured as a function of time. The diffusion time across the cleft can be estimated as the time required for this postsynaptic coverage to reach a maximum. Note that in this case we are interested purely in the transit time, and therefore will not use the doubly ligated state of receptors as a measure of postsynaptic coverage. The 1% cutoff should give a more accurate indication of the transit time for the bulk of ACh, without including the kinetics of binding in the time estimates.

The NMJ model used for these simulations is within estimates for vertebrate NMJ's (Salpeter, 1987), and is similar to the model used by Bartol et al. (1991). The synaptic gap is 48 nm wide. Because we are only interested in ACh crossing times, secondary folds have been removed. Secondary folds will be included in the simulations discussed below. The ACh is released directly above the synaptic gap, localized in a 48 nm × 48 nm × 32 nm cleft. Integration of the initial concentration over the volume of the cleft yields ~13,000 ACh molecules, within estimates of the contents of skeletal muscle vesicles (Miledi et al., 1982).

The area of the postsynaptic membrane over which the ACh concentration is above 1% of its release concentration is plotted in Fig. 2. Simulations were performed with both active and inactive AChE to assess the effects of the reactive enzyme on the initial diffusion of ACh. From the plot it is clear that AChE has a strong effect on ACh concentration at the postsynaptic membrane. Inhibition of the esterases more than doubles the postsynaptic coverage, though the 1% criterion provides only an approximate measure of coverage. Additionally, the neurotransmitter concentration lingers for a significantly longer time at the postsynaptic membrane.


View larger version (20K):
[in this window]
[in a new window]
  		FIGURE 2   Postsynaptic area with an acetylcholine concentration >1% of the initial release concentration. Area is shown for both active esterases (solid line) and inactive esterases (dashed line).

The time required for the coverage to reach a maximum gives some indication of the diffusion time across the cleft. In Fig. 2, maximum coverage is reached in ~13 µs for inhibited AChE, and ~6 µs for active AChE. These measurements are in reasonable agreement with the aforementioned estimate. Note that a delay time of 6 µs corresponds to a mean distance of 120 nm for a diffusing particle in three dimensions (Chandrasekhar, 1943; Crank, 1975). Given the 48-nm synaptic gap and the 32-nm depth of the release cleft, the 6-µs delay time implies that the first neurotransmitter particles can travel anywhere from 50 nm to 110 nm along the postsynaptic membrane before the postsynaptic concentration reaches a maximum. It is therefore likely that some lateral diffusion occurs after quantal release, though most of the neurotransmitter is isolated in the region directly below the release cleft.

Although the transit time of ACh in the synaptic cleft is ~6 µs, miniature endplate current rise times are typically an order of magnitude larger. It is therefore evident that miniature endplate current rise times are dominated by receptor binding and isomerization, not by diffusion in the synapse (Matthews-Bellinger and Salpeter, 1978; Dwyer, 1981; Madsen et al., 1984). In our measurements of doubly ligated postsynaptic coverage discussed below, the coverage rise time is 85 µs, in good agreement with measured MEPC rise times.

Postsynaptic coverage

Estimates of the number of receptors activated by a single quanta of ACh range from 1000 to 2000 or more (Salpeter, 1987). Hartzell et al. (1975) have estimated that a quanta of neurotransmitter acts over 0.2 µm2 or less in snake NMJ's (Hartzell et al., 1975). Matthews-Bellinger and Salpeter have estimated that the area of receptor saturation due to a single quanta of ACh is ~0.3 µm2 in frog NMJ's (Matthews-Bellinger and Salpeter, 1978). In their case, saturation is defined as the area of the postsynaptic membrane containing a number of binding sites equal to the amount of ACh in a typical quanta. Land et al. (1981) have estimated that postsynaptic coverage extends 0.3 µm from the point of release for lizard intercostal muscle. In their case, postsynaptic coverage is defined as the region over which receptors are doubly bound. Since the measurements of Land et al. explicitly use the doubly ligated state of receptors, they provide a useful reference for our simulations.

The NMJ model previously described was used for these simulations along with a model containing secondary folds and a model with secondary folds, but inactive AChE. The release cleft was located in an active zone directly above a secondary fold (Kandel et al., 1991). Table 1 shows the number of doubly ligated receptors, total postsynaptic coverage, and coverage rise times for each NMJ model. As a reminder, postsynaptic coverage is defined as the maximum area over which receptors are doubly ligated. The coverage rise times are defined as the 20-80% time to peak for the postsynaptic coverage, and are indicative of the miniature endplate current (MEPC) rise times.

                              
View this table:
[in this window]
[in a new window]
  		TABLE 1   Comparison of the postsynaptic coverage, coverage rise times, and the number of doubly ligated AChR molecules, due to the release of a single quanta of acetylcholine

The complete NMJ model, containing secondary folds and active AChE, is in reasonable agreement with experimental measurements and estimates. 0.24 µm2 coverage indicates that receptors are saturated as far as 0.28 µm from the point of release, comparable to the estimates of Land et al. It is estimated that between 60% and 90% of the doubly ligated receptors will be in the open state (Salpeter, 1987). Thus, we can expect that anywhere from 1126 to 1689 receptors in our simulations will open. Again, this is within the expected range. Finally, the estimated rise time compares reasonably with reported MEPC rise times of 80-100 µs for lizard NMJ's (Land et al., 1981). Note that our rise time estimates do not account for open-receptor events, the last step in Eq. 17. Receptor isomerization will be included in future simulations, however.

When secondary folds are removed, postsynaptic coverage increases by 42%. In the simulations, the folds serve as buffer regions where excess ACh is hydrolyzed. The lower concentration of receptors in the secondary folds reduces coverage in these regions, while the constant presence of AChE leads to rapid hydrolysis. Therefore, the presence of secondary folds reduces the efficiency of binding to receptors, and thus the amplitude of the subsequent miniature endplate currents. Bartol et al. (1991) observed a similar phenomenon in their simulations, with a 34% increase in MEPC amplitude when folds were removed from their lizard NMJ model.

When AChE activity is blocked, postsynaptic coverage increases by 83%, with 3400 doubly ligated receptors. The 20-80% rise time extends to 290 µs. With the ACh now diffusing an average distance of 0.38 µm from the point of release, diffusion time in the cleft plays a larger role in determining rise times. Additionally, because it takes ~1 ms for postsynaptic coverage to reach a maximum, dissociation of ACh from receptors will likely play a larger role in determining MEPC rise times and amplitudes. Thus, although coverage increases by 83% when AChE is disabled, MEPC amplitudes should increase by significantly less, due to steady ACh dissociation from the receptor array. Indeed, Wathey et al. (1979) observed a 27% increase in MEPC amplitude in their one-dimensional continuum calculations when AChE was inhibited. Bartol et al. observed a 35% increase in MEPC amplitude in their Monte Carlo simulations (Bartol et al., 1991). Katz and Miledi (1972) observed a 55% increase in their observed miniature endplate potentials upon prostigmine inhibition of AChE.

Acetylcholine lifetime in the NMJ

Although miniature endplate currents decay over several milliseconds, the bulk of ACh is thought to be hydrolyzed, or to exit the NMJ via free diffusion, in as short as 200-500 µs (Colquhoun, 1990; Wray, 1990; Katz and Miledi, 1979). During our simulations, the total quantity of ACh in the junction was recorded as a function of time, and its mean lifetime in the junction was estimated. Fig. 3 plots the quantity of ACh in the junction relative to the initial quantity for active and inactive AChE. The lifetime of ACh is defined as the time required for the initial maximum concentration to fall to 10% of its original value.


View larger version (18K):
[in this window]
[in a new window]
  		FIGURE 3   Quantity of acetylcholine in the neuromuscular junction plotted as a function of time, relative to the initial quantity. The solid line represents active AChE, while the dashed line represents inactive AChE. In the inactive case, loss of ACh is due to diffusion out of the junction.

In the active AChE case, the lifetime is on the order of 215 µs, within the expected range. The effects of AChE are most pronounced in the early stages of diffusion, when the initial wave of ACh reaches the esterases. The curve representing active AChE shows a sharp decline at early times, while the curve representing inactive AChE is relatively constant over the first 5 µs. A comparison of the two curves suggests, then, that AChE has a strong impact on the initial burst of released ACh. Certainly this phenomenon would explain why substrate inhibition of AChE may be important biologically.

When AChE is deactivated, the ACh concentration is depleted by association with receptors, or by diffusion out of the junction. In order to simulate this phenomenon, the sides of the junction were extended to allow for ample room outside the proper "boundaries" of the original junction. The concentration of ACh within the original junction volume was still recorded over the course of the simulation. From Fig. 3 it is clear that the average ACh lifetime strongly depends on the presence of AChE, and increases to 900 µs when the esterases are inactive. Unlike the active AChE case, where the bulk of free ACh is hydrolyzed well before bound ACh begins to dissociate, the neurotransmitter in the inactive case remains in the synaptic cleft for over a millisecond. Prolongation of the synaptic event is therefore expected, due to repeated binding of the neurotransmitter to receptors.

The results of this section and the previous one suggest that AChE plays a role in both isolating neurotransmitter quanta in the immediate vicinity of release and in removing excess ACh before significant rebinding can occur.

ACh release models

In the initial stages of ACh release synaptic vesicles fuse to the presynaptic membrane, forming pores through which the neurotransmitter can diffuse. The fusion pores expand rapidly, reaching ~50 nm in size (Kandel et al., 1991). Earlier continuum models have suggested that diffusion through a vesicular pore is insufficient to account for the rapid rise time of miniature endplate currents (Khanin et al., 1994). More recent Monte Carlo simulations, however, have shown that diffusion through a rapidly opening pore can account for observed rise times (Stiles et al., 1996).

To further explore this issue, we have performed simulations using different release models. The different models are depicted in Fig. 4. The first model corresponds to rapid vesicular fusion on the timescale of diffusion of ACh. Here, the ACh is released in a shallow cleft at the top of the synaptic gap, and may diffuse directly into the synapse. The volume of the cleft is roughly equal to the volume of the synaptic vesicle used by Khanin et al. in their studies of neurotransmitter release (Khanin et al., 1994). The second model corresponds to diffusion through a rapidly expanding fusion pore. The pore size is fixed at 6 nm on a side, with a length of 9 nm, throughout the simulation; the pore length is identical to that of Stiles et al. (1996). The cross-sectional area of the pore (36 nm2) is equal to the average area over the first 100 µs of the pore used by Stiles et al. In their fast-opening pore model, the vesicle's contents emptied in ~100 µs, so that our 36 nm2 pore area provides a good estimate of their average pore area over the time that the bulk of ACh exited their vesicle model. The final release model corresponds to a slowly opening fusion pore. The pore size is fixed at 2 nm on a side, with a length of 9 nm. The cross-sectional area of the small pore is in good agreement with the average area of the slowly opening pore model of Stiles et al. In all three cases, the initial ACh concentration inside the vesicle is equal to 300 mM (Vincent and Wray, 1990).


View larger version (18K):
[in this window]
[in a new window]
  		FIGURE 4   Depiction of the three neurotransmitter release models. (A) Instantaneous release from a presynaptic cleft. (B) Release through a large pore, representing rapid pore expansion. (C) Release through a small pore, representing slow pore expansion. Figures not drawn to scale.

Fig. 5 shows the average ACh concentration in the vesicle as a function of time for the cleft and pore models. From the figure, it is clear that the release model, and hence the rate of pore expansion, has a strong impact on vesicle emptying times. In the cleft model, the average concentration falls by 90% in 9 µs, while in the large and small pore models, the concentration falls by 90% in 140 µs and 1 ms, respectively. In comparison, published estimates range from 100 µs to 500 µs for the release of ~90% of a vesicle's contents (Khanin et al., 1994; Stiles et al., 1996; Almers and Tse, 1990). Experimental measurements of the neurotransmitter serotonin in synaptic vesicles of cultured leech neurons found release times on the order of 100 µs (Bruns and Jahn, 1995). The large pore model, which corresponds to a radial expansion rate of 25 nm/ms (Stiles et al., 1996), agrees well with these estimates. The emptying times of the other two models lie outside the expected range.


View larger version (21K):
[in this window]
[in a new window]
  		FIGURE 5   Average vesicular concentration of acetylcholine for different release models.

Estimates of the average concentration of ACh at the postsynaptic membrane during a single quantal event ranges from 1 to 10 mM (Kleinle et al., 1996; Matthews-Bellinger and Salpeter, 1978; Land et al., 1980; Fertuck and Salpeter, 1976), much higher than the estimated Kd for AChR binding (Salpeter, 1987). For our three release models, ACh concentration on the postsynaptic membrane was measured as a function of time (Fig. 6). Here, the concentration has been averaged over a disk of radius 0.2 µm centered below the point of release. From the plot, it is clear that diffusion through a cleft or large pore leads to higher average ACh concentrations at early times. Both models lie within the 1 mM to 10 mM estimated range. Diffusion through the slowly opening pore model, however, leads to diffuse coverage. The ACh concentration lingers for several hundred microseconds on the postsynaptic membrane, but never rises above 0.25 mM. As will be discussed below, diffusion through the small pore significantly delays transmission of the synaptic signal.


View larger version (20K):
[in this window]
[in a new window]
  		FIGURE 6   Average acetylcholine concentration on the postsynaptic membrane for different release models.

The most important criterion for comparing the different release models is coverage of the postsynaptic membrane. Table 2 lists postsynaptic coverage, coverage rise time, the maximum number of open AChR channels, and the delay time for the three models. For the slowly opening pore model, postsynaptic coverage never rises above 0.12 µm2 in the first 500 µs after release. The rapidly opening pore model, however, allows for rapid diffusion of neurotransmitter out of the synaptic vesicle, and coverage close to that of the cleft model. Additionally, the coverage rise times for both the large pore and cleft models are in rough accord with MEPC rise time estimates for vertebrate NMJ's. Stiles et al. (1996) found that both instantaneous neurotransmitter release and passive diffusion through a rapidly opening pore gave similar MEPC rise times and amplitudes. Our results, along with theirs, lend support to the argument that diffusion out of a rapidly opening pore is sufficient to account for the observed amplitudes and rapid rise times typically measured in vertebrate NMJ's.

                              
View this table:
[in this window]
[in a new window]
  		TABLE 2   Comparison of the postsynaptic coverage, coverage rise time, number of doubly ligated AChR molecules, and initial delay time for different release models

The delay time in Table 2 is defined as the time required for the initial burst of ACh to reach the postsynaptic membrane. For the slowly opening pore model, the delay time alone is already longer than MEPC rise time estimates. Therefore, the slowly opening pore model, which corresponds to an expansion rate of 0.275 nm/ms (Stiles et al., 1996), can safely be ruled out as a possible mechanism for neurotransmitter release.

    CONCLUSIONS
Top
Abstract
Introduction
Methods
Results
Conclusions
References

In this work we have examined the effects of AChE, secondary folds, and different release models on the extent and timecourse of ACh diffusion. Our simulations have shown that coverage rise times (and hence MEPC rise times) are dominated by binding and isomerization events, rather than diffusion across the synaptic cleft. In our simulations, transit times were on the order of 6 µs, while coverage rise times were shown to be ~85 µs or greater.

From the simulations it is also clear that although AChE does have an impact on the initial burst of released ACh, its primary role appears to be during the later stages of synaptic transmission. We observed that the mean lifetime of ACh increases considerably when the esterases are inactive, and that postsynaptic coverage and coverage rise times increase modestly. Our simulations support previous observations that AChE clears ACh from the synaptic cleft before dissociation and subsequent rebinding of the neurotransmitter can occur (Land et al., 1984; Magleby and Terrar, 1975; Katz and Miledi, 1973). Additionally, the esterases keep the neurotransmitter from diffusing beyond the immediate vicinity of release.

However, our results predict that secondary folds play a smaller role in synaptic transmission. Due to the reduced surface density of receptors in the secondary folds, their presence has only a small effect on coverage, and therefore MEPC amplitude. Again the results are in accord with other simulations (Bartol et al., 1991), which suggest that secondary folds do not play a primary role in ACh binding or uptake, but serve some other purpose instead.

We have found that both release of ACh through a large pore, and instantaneous release from a presynaptic cleft, lead to similar degrees of postsynaptic coverage and similar rise times. In both cases, the rise times, coverage, and average concentration on the postsynaptic membrane are in accord with experimental measurements. Thus, our simulations suggest that passive diffusion of ACh through a sufficiently rapidly opening fusion pore can account for the measured timecourse and amplitude of endplate currents. However, release through a small (i.e., slowly opening) fusion pore leads to ACh concentrations on the postsynaptic membrane below published estimates. Additionally, the delay time of the transmitted signal is itself longer than experimental MEPC rise times. Therefore, a rapidly expanding fusion pore is vital for timely and efficient synaptic transmission.

By comparing our results to experiments, and to other simulations, we have demonstrated that continuum finite element methods are a viable alternative to Monte Carlo and analytic techniques for studying synaptic transmission. With the reliability of the method established, it is possible to move to more ambitious NMJ models, where experimental data, or data from other simulations, are sparse.

The primary advantages of the continuum method are its speed and expandability. Due to the discrete nature of ACh in Monte Carlo implementations, the size of individual timesteps in these simulations is often limited. In addition, the total time required for Monte Carlo simulations scales linearly with the number of particles in the system, and can be prohibitively large when many release sites are included. The continuum approach, however, utilizes adaptive timesteps to reduce computation time, and scales by the number of finite elements in the system, not by the number of ACh "particles." Therefore, with a sufficiently coarse finite element mesh, the continuum NMJ model may easily be expanded to study synaptic transmission on a much larger scale.

Future simulations will therefore examine much larger systems, including dozens, or hundreds, of interacting release sites. The role of AChE in isolating separate quanta will be further explored in these cases. In addition, potentiating effects between quanta will be examined. The eventual goal is the accurate reproduction of full endplate currents (EPC's), and possibly the study of desensitization effects due to rapid, repeated release. This will require the release of ~300 quanta, the quantal content of a typical EPC (Katz and Miledi, 1979). The continuum approach, with a sufficiently simplified NMJ model, can be used to study such large systems.

Finally, the full kinetic model for nicotinic receptors, shown in Eq. 17, will be implemented in future simulations. The addition of receptor kinetics is straightforward, and should involve little modification of the current code. Modeling of channel isomerization and unbinding events will allow for generation of full miniature endplate currents, and allow for direct comparisons to experimentally measured MEPC's.

    ACKNOWLEDGMENTS

We are very grateful for advice from Dr. Palmer Taylor, Dr. Randolph Bank, Dr. Adrian Elcock, Bodo Erdmann, and Erik Hom. We are also grateful to the referees for their valuable suggestions, which have led to an improved manuscript.

This work was supported in part by a grant from the NIH. J. S. was supported by a predoctoral fellowship from the UCSD NIH Molecular Biophysics Training Program.

    FOOTNOTES

Received for publication 30 December 1997 and in final form 8 June 1998.

Address reprint requests to Dr. Jason Smart, Dept. of Chemistry and Biochemistry, University of California at San Diego, La Jolla, CA 92093-0365. Tel.: 619-534-5843; Fax: 619-534-7042; E-mail: jsmart{at}ucsd.edu.

    REFERENCES
Top
Abstract
Introduction
Methods
Results
Conclusions
References

Biophys J, October 1998, p. 1679-1688, Vol. 75, No. 4
© 1998 by the Biophysical Society   0006-3495/98/10/1679/10  $2.00


This article has been cited by other articles:


Home page
 Proc. Natl. Acad. Sci. USAHome page
D. S. Dandy, P. Wu, and D. W. Grainger
Array feature size influences nucleic acid surface capture in DNA microarrays
PNAS, May 15, 2007; 104(20): 8223 - 8228.
[Abstract] [Full Text] [PDF]

Home page
 Biophys. JHome page
Y. Cheng, J. K. Suen, D. Zhang, S. D. Bond, Y. Zhang, Y. Song, N. A. Baker, C. L. Bajaj, M. J. Holst, and J. A. McCammon
Finite Element Analysis of the Time-Dependent Smoluchowski Equation for Acetylcholinesterase Reaction Rate Calculations
Biophys. J., May 15, 2007; 92(10): 3397 - 3406.
[Abstract] [Full Text] [PDF]

Home page
 Biophys. JHome page
Y. Song, Y. Zhang, C. L. Bajaj, and N. A. Baker
Continuum Diffusion Reaction Rate Calculations of Wild-Type and Mutant Mouse Acetylcholinesterase: Adaptive Finite Element Analysis
Biophys. J., September 1, 2004; 87(3): 1558 - 1566.
[Abstract] [Full Text] [PDF]

Home page
 Biophys. JHome page
Y. Song, Y. Zhang, T. Shen, C. L. Bajaj, J. A. McCammon, and N. A. Baker
Finite Element Solution of the Steady-State Smoluchowski Equation for Rate Constant Calculations
Biophys. J., April 1, 2004; 86(4): 2017 - 2029.
[Abstract] [Full Text] [PDF]

Home page
 Biophys. JHome page
K. Tai, S. D. Bond, H. R. MacMillan, N. A. Baker, M. J. Holst, and J. A. McCammon
Finite Element Simulations of Acetylcholine Diffusion in Neuromuscular Junctions
Biophys. J., April 1, 2003; 84(4): 2234 - 2241.
[Abstract] [Full Text] [PDF]

This Article
Right arrow Abstract Freely available
Right arrow Full Text (PDF)
Right arrow Alert me when this article is cited
Right arrow Alert me if a correction is posted
Services
Right arrow Similar articles in this journal
Right arrow Similar articles in PubMed
Right arrow Alert me to new issues of the journal
Right arrow Download to citation manager
Right arrow reprints & permissions
Citing Articles
Right arrow Citing Articles via HighWire
Right arrow Citing Articles via Google Scholar
Google Scholar