## Abstract

Immobile and mobile calcium buffers shape the calcium signal close to a channel by reducing and localizing the transient calcium increase to physiological compartments. In this paper, we focus on the impact of mobile buffers in shaping steady-state calcium gradients in the vicinity of an open channel, i.e. within its “calcium microdomain.” We present a linear approximation of the combined reaction–diffusion problem, which can be solved explicitly and accounts for an arbitrary number of calcium buffers, either endogenous or added exogenously. It is valid for small saturation levels of the present buffers and shows that within a few hundred nanometers from the channel, standing calcium gradients develop in hundreds of microseconds after channel opening. It is shown that every buffer can be assigned a uniquely defined length-constant as a measure of its capability to buffer calcium close to the channel. The length-constant clarifies intuitively the significance of buffer binding and unbinding kinetics for understanding local calcium signals. Hence, we examine the parameters shaping these steady-state gradients. The model can be used to check the expected influence of single channel calcium microdomains on physiological processes such as excitation–secretion coupling or excitation–contraction coupling and to explore the differential effect of kinetic buffer parameters on the shape of these microdomains.

- Ca
^{2+}microdomains - Ca
^{2+}diffusion - Ca
^{2+}buffers - buffer kinetics
- diffusion modeling
- synaptic transmission

Calcium is involved in a multitude of fast intracellular signal transduction mechanisms ranging from excitation–contraction coupling to synaptic transmission (Augustine et al., 1985; Cheng et al., 1993; Bruns and Jahn, 1995; Clapham, 1995). To achieve a high bandwidth of signal transmission at specific sites, the cell needs to localize the calcium signals in time and space. Mobile calcium buffers are elegant tools to achieve this temporal and spatial functional compartmentalization (Roberts, 1994): mobile calcium buffers act as calcium shuttles or sinks to produce steep gradients in a close neighborhood of channels. In these “calcium microdomains,” [Ca^{2+}] readily reaches many tens of micromolar levels to activate low-affinity processes. Complementary to this, the microdomains dissipate very rapidly by virtue of the mobilities of the buffers.

Recent experimental studies (Eilers et al., 1995; Yuste and Denk, 1995) have used imaging techniques to observe the temporal and spatial dynamics of [Ca^{2+}] in different cell types. Unfortunately, currently available imaging technology does not simultaneously provide a sufficiently high resolution in time and space. Temporal resolution is sacrificed to get a decent spatial resolution, or vice versa. This dilemma, in concert with the insight that “Ca^{2+} signaling takes the local route” (Augustine and Neher, 1992b), has in turn triggered a huge body of simulation studies that attempt to calculate the expected time course of calcium concentration increases and their spatial extent (Simon and Llinas, 1985; Stern, 1992; Nowycky and Pinter, 1993; Roberts, 1994;Klingauf and Neher, 1997) in the presence of multiple buffers.

From these simulations, we have learned about the impact of cellular buffers on calcium signaling, but each simulation represents only one point in a multidimensional parameter space. Nevertheless, it would be helpful if we could make a compromise between the computational and the analytical complexity of the buffered diffusion problem. One approach is to linearize the corresponding differential equations. It was first done for one mobile buffer (Neher, 1986), which was assumed not to change its concentration by binding to calcium, a reasonable approximation if a high concentration of a mobile chelator is present. Later, Stern (1992) considered a channel as a point source in an isotropic medium and approximated the steady-state calcium concentration profiles surrounding the source in terms of nonlinear integral representations. More recently, Pape et al. (1995) calculated the concentration profiles of Ca^{2+} and a mobile buffer for “small saturation levels” of the buffer. Within the limit of no spatial gradients of the mobile buffer (by virtue of its high concentration), they achieved the results of Neher (1986) andStern (1992).

In the present paper, we extend this approach to an arbitrary number of mobile buffers. We derive analytical solutions for the steady-state profiles of calcium and calcium-bound buffers surrounding a channel as a point source. It is shown that in the range of tens of nanometers, the effectiveness of buffered solutions depends strongly on chemical kinetics and diffusional mobility of the buffers. The results are used to discuss relative effects of BAPTA and EGTA (Borst and Sakmann, 1996) in blocking synaptic transmission.

## MATHEMATICAL MODEL

We will be considering the concentration profiles of calcium and calcium-bound buffers around a single calcium channel. In particular, as we will see in the sequel we can focus on steady-state gradients because the concentration profiles stabilize rapidly in the vicinity of a channel after its opening. Thus, the first step is to show why we can only consider mobile buffers to be present in the course of our analysis.

**I. Immobile buffers do not affect the standing calcium gradients around an open calcium channel.** Let us consider the interaction of calcium ions with N different buffer species (which we shall denote by B_{1}, B_{2}, … , B_{N}) in a reaction cell according to the kinetic scheme:
Equation 1The parameters k_{i} and k_{−i} are the kinetic rate constants of this reaction in units of 1/(M.sec) and 1/sec, respectively. Furthermore, we denote by K_{i} the dissociation constant of the i-th buffer (under the specific pH and temperature conditions) given by K_{i} = k_{−i}/ k_{i}. If there are no sources or sinks for the buffers B_{i} (which we assume to be the case here), the notion of the total concentration of the i-th species in the reaction volume is a well defined conserved quantity, i.e.: [B_{i}]_{T} ≡ [B_{i}] + [CaB_{i}]. For the sake of simplicity, we shall further make use of the following abbreviations in the rest of this paper: y_{i} ≡ [CaB_{i}], y_{N+1} ≡ [Ca^{2+}], x_{i} ≡ [B_{i}], i = 1, … , N.

Taking into account the diffusive mobility of all molecules, say in an isotropic reaction medium, we need to introduce the diffusion coefficients of calcium ions, D_{N+1} ≡ D_{Ca}, as well as that of the free form of B_{i}, D_{Bi}, and of its calcium-bound form, D_{CaBi}. The spatiotemporal evolution of the concentrations is then incorporated into a system of partial differential equations given by:
Equation 2
where Δ is the Laplace operator. Assuming that D_{Bi} = D_{CaBi} ≡ D_{i}, we can add the first and last equation of the system (2) to get the following equation for the new variable z_{i} = x_{i} + y_{i}:
Equation 3which is the standard diffusion equation for the total i-th buffer. If the buffer is homogenously distributed in the cell at time zero, i.e. z_{i}(t = 0, r⃗) ≡ [B_{i}]_{T} for all space coordinatesr⃗, the unique solution to Equation 3 is just the constant z_{i}(t, r⃗) ≡ [B_{i}]_{T}. This means that the total buffer remains spatially constant in the sample if this was the case at time zero. Assuming this in the following, we can now simplify Equation2 to get the following system:
Equation 4
describing the dynamics of concentration changes attributable to diffusion and reaction in three dimensions. Several studies have investigated the temporal evolution of this system for specified geometries, initial and boundary conditions with or without simplifications. The objective often was to calculate the concentration of free calcium at locations that are not experimentally accessible because of resolution limitations of the currently available imaging technology, e.g. the calcium concentration within a few tens of nanometers from the cytosolic mouth of a calcium channel.

Let us now assume that the first M < N buffer species are mobile and the last N–M are immobile, i.e., have vanishing diffusion coefficients: D_{i} = 0 μm^{2}/sec for i = M + 1, … , N. In steady-state, all of the time derivatives in Equation 4 are zero, giving rise to 0 = k_{i} · y_{N+1} · ([B_{i}]_{T} − y_{i}) − k_{−i}· y_{i}, i = M + 1, … , N, which reduces Equation 4 to:
Equation 5
Equation 5, however, is exactly the system one has if there were only the first M buffers, i.e., the mobile buffers, present in the cell. The conclusion is in studying the steady-state concentrations, all the fixed buffers can be neglected because they do not influence these concentration profiles (Stern, 1992). The effect of the immobile buffers is to prolong the time needed to reach steady-state, but once this is achieved, the standing gradients are completely identical to those one would see if there were no fixed buffers present at all. The reason for this is simply the fact that the fixed buffers cannot be replenished by means of diffusion and, hence, get saturated. Because of this, without loss of generality, we may assume that only mobile buffers are present whenever we look at steady-state profiles.

**II. The reaction–diffusion system can be linearized if the “buffer saturations,” i.e., the increase in the concentrations of calcium-bound buffers on channel opening, are small.** Let us now come back to the general system (Eq. 4) to establish the conditions under which a linearization can be carried out. The reaction–diffusion system is nonlinear because of the products of the concentrations according to the law of mass action. If we denote the resting concentrations by ȳ = (ȳ_{1}, … ,ȳ_{N+1})^{t}, any deviations from the resting concentrations, say by virtue of calcium injection through an open channel, can be written as: y_{i} =ȳ_{i} + δy_{i}, i = 1, … , N, i.e., y = δy +ȳ, which suggests linearization of Equation 4 aroundȳ. We can write Equation 4 as a vector field equation:
Equation 6where the reaction vector field f(y) incorporates all the nonlinearity and is given by:
D is just the diagonal matrix of the diffusion coefficients:
The concentration increases (above resting values) of the calcium-bound buffers, so-called “buffer saturations,” are maximal in steady-state, at any point in space. Hence, if they are small in steady-state, they can only be smaller during the transient evolution of the gradients. For small saturations, we can approximate the nonlinear contribution of the reactions to the evolution of the gradients with linear expressions. This is exactly like approximating an exponential function of time, say with a time constant τ, with a linear function. Such an approximation will be quite good if we confine ourselves to a time window of length τ. The linear function must be chosen to have the same slope as the slope of the exponential at the point where the two curves overlap. Hence, we need to compute the derivative of the reaction term f(y) in the system (Eq. 6), which is accomplished by standard Taylor series expansion of f. As a result, the nonlinear term f(y) in (Eq. 6) is substituted by a linear matrix product A · δy giving:
Equation 7which is a simple linear matrix partial differential equation. A is here the matrix of partial derivatives of f evaluated at the resting point ȳ, given by:
Equation 8with τ_{i} = 1/(k_{−i} + k_{i}[Ca^{2+}]) being the reaction time constant for binding of calcium to the buffer B_{i} (Bernasconi, 1976) and κ_{i} the “binding ratio” of B_{i} defined as κ_{i} ≡ (∂[CaB_{i}]/∂[Ca^{2+}]) = ([B_{i}]_{T} · K_{i}/([Ca^{2+}] + K_{i})^{2}) (Zhou and Neher, 1993). But what conclusions can we draw from the linearized system (Eq.7)?

**III. The steady-state solution to the linearized reaction–diffusion problem is determined purely by three factors: the mean reaction times of the buffers with calcium (τ _{i}), the buffering power of each buffer species (κ_{i}), and the mobility of the buffers (D_{i}). The solution is given by a sum of exponentials divided by the distance from the channel. Every buffer can be assigned a characteristic length-constant, which is indicative of its kinetically limited buffering capability close to the channel.**The structure of the system (Eq. 7) and the matrix A already show the determinants of the linearized problem, the diffusive mobility of the buffers, and their binding kinetics. If a buffer is very mobile, i.e., has a very high diffusion coefficient, then only its kinetics should govern its influence on the calcium signal. The more buffer one has, the higher the chance of a calcium ion to be bound by the buffer and thus the shorter the mean time that elapses until a calcium ion is captured by a buffer molecule. This intuitive notion is expressed in the ratio κ

_{i}/τ

_{i}. The binding ratio κ

_{i}is a measure of the capacity of a buffer to bind calcium at a given concentration [Ca

^{2+}]. For instance, κ = 100 means that ∼1% of a given calcium load will appear as free calcium and 99% will be bound by the buffer. κ

_{i}/τ

_{i}can easily be shown to be equal to [B

_{i}] · k

_{i}, a quantity that we refer to as the “buffer product” and is the reciprocal mean time required for a calcium ion to be bound by B

_{i}.

Before proceeding with these lines of thought, however, let us outline how the system (Eq. 7) can be solved to get analytical solutions. Because we are interested in the concentration profiles within a very small area, i.e., within the microdomain of a channel, we can consider the channel mouth to be embedded in a hemisphere and in doing so reduce the 3-D problem to a 1-D radial diffusion problem. Transformation of the system (Eq. 7) to spherical coordinates gives us: Equation 9where r now represents the distance from the channel mouth.

Finally, bringing our calcium source, i.e., the open channel, into play, we put it at the origin given by r = 0 and assume that it has a constant flux Φ (in mol/sec) of calcium ions and represents no sources or sinks for other buffers. Appendix outlines the derivation of the transient solution to the system (Eq. 9) as a function of time after channel opening and distance from channel. There, we also demonstrate numerically that the transient solution rapidly approaches steady-state within a few hundred nanometers from the channel. For this reason, we will focus only on the steady-state gradients in the following sections.

Within the microdomain of the channel, calcium is carried either as calcium-bound form of one of the buffers or in its ionic free form, and the total flux of calcium at each distance r from the origin must be constant and equal to Φ at steady-state. Thus, we have the following calcium flux conservation constraint:
Equation 10where each term in the sum is just the contribution of each calcium-carrying species to the total calcium flux Φ. So, we need to solve (Eq. 9) in steady-state subject to the constraint given by Equation 10. This is outlined in Appendix . The main conclusion that can be drawn is that the gradients around the channel can be described as a sum of exponentials multiplied by 1/r according to:
Equation 11The vectors u_{i} as well as the length-constants of the exponentials, 1/
, are determined by the matrices A and D, hence by buffer kinetics τ_{i}, buffer binding power κ_{i}, and buffer diffusivity D_{i}. Our objective in the next sections will be to investigate the behavior of this solution under different conditions and for different distances from the channel.

**IV. For distances from the channel that are big compared with the length-constants, i.e., if the mean diffusion time for calcium is much bigger than the mean reaction times of the buffers, the steady-state concentrations are determined purely by equilibrium buffer properties, namely the binding ratios (κ _{i}) and the diffusion coefficients (D_{i}).** A calcium ion, in the absence of any buffers, will on average need a time t to diffuse a distance r from the channel, which is proportional to r

^{2}/D

_{Ca}, i.e., t ∝ r

^{2}/D

_{Ca}. Hence, as one moves away from the channel, the mean diffusion time eventually gets much bigger than the mean buffer reaction times τ

_{i}. For such distances, one would expect that the buffers would be in chemical equilibrium with local calcium. As a consequence, the binding kinetics of a buffer should not affect the concentration profiles but its calcium affinity and concentrations should.

To see that our solution (Eq. 11) confirms this prediction, we evaluate it for long distances, i.e., for r
≫ 1, i = 1, … , N. Then, the first N terms in Equation 11 vanish, and δy is approximately given by δy ≈ (1/r)a_{N+1}**u**_{N+1}. Equation 10 finally results in an expression for a_{N+1} as (see Appendix for the definition of **u**_{N+1}):
Equation 12
This confirms our intuition that the concentrations are determined merely by equilibrium properties if we are at distances that are bigger than the length-constants. In particular, δ[Ca^{2+}] ≈ [Φ/4πr (∑_{i=1}^{N}κ_{i}D_{i} + D_{Ca})]. Thus, the increase in the concentration of free calcium scales linearly with the single channel current Φ and falls according to a 1/r -law with an apparent calcium diffusion coefficient given by D_{app} = ∑_{i=1}^{N}κ_{i}D_{i} + D_{Ca}.

**V. For distances from the channel that are short compared with the length-constants, i.e., if the mean diffusion time for calcium is much shorter than the mean reaction times of the buffers, the calcium concentration behaves like in the unbuffered scenario.** We have just demonstrated that buffer kinetics is insignificant far away from the channel. We shall now proceed to show that very close to the channel, no buffer is capable of shaping the calcium concentration profile because of its finite reaction time. In Appendix , we calculate the increase in the concentration of free calcium, δ[Ca^{2+}], as the difference between the expected calcium increase in the absence of any buffers (Φ/4πrD_{N+1}) and the scaled sum of the concentrations of calcium-bound buffers above rest (∑_{j=1}^{N}D_{j}/ D_{N+1}δ[CaB_{j}]). Approaching the channel, the first term grows with 1/r, whereas the sum is finite (because δ[CaB_{j}] < [B_{j}]_{T}). Thus, sufficiently close to the channel, the first term will dominate such that calcium is approximated by Φ/4πrD_{N+1}. It can be shown that this happens for r ≪ r_{min} ≡ min 1/
. This describes the behavior of calcium close to the channel. The next issue is the near-channel behavior of the buffers. Clearly, the buffer saturations can only be finite because the total buffer concentration is finite. But what determines the saturation level of a buffer?

**VI. The maximal buffer saturations depend heavily on the kinetics of interaction of calcium with the buffers. The faster a buffer, the higher its chance to bind calcium close to the channel and the bigger its saturation. This saturation scales linearly with the single channel current and can be calculated explicitly to check the validity of the small saturation assumption.** We have mentioned that our approach is based on the assumption of small buffer saturations. Hence, one needs to be able to check the validity of this assumption before applying the theory. This in turn raises the question of whether it is possible to estimate the degree of buffer saturation, given the experimental conditions such as single channel current and buffer characteristics. Because the concentration deflections from rest are maximal in steady-state and the higher, the closer we are to the source, the limit for r → 0 gives us an upper bound for the maximal expected buffer saturations. For this sake, let x = (δy_{1}(0) δy_{2}(0) · δy_{N}(0))^{t} denote the vector of buffer saturations at the source. Then, using Equation EAI.10, one arrives at the following identity:
Equation 13
Thus, the buffer saturations can be calculated analytically as the solution to the above linear algebraic system of equations. It can be shown that the buffer saturations can also be calculated according to:
Equation 14
In other words, we compute the (positive) square root of B, perform the matrix product −
·**u**, and the j-th component of this vector gives the concentration deflection from rest of the calcium-bound form of the j-th buffer at the source in steady-state.

To develop an intuitive understanding of the parameters that determine the maximal buffer saturation, it is instructive to calculate the saturation for the case of one buffer. Then, the matrix B has the eigenvalues μ_{1} = 1/τ_{1}D_{1} + κ_{1}/τ_{1}D_{Ca}, μ_{2} = 0. Using Equation 14, we arrive at:
Equation 15Hence, the buffer saturation at the source is a product of three factors: one involves only equilibrium buffer properties, κ_{1}/4π(κ_{1}D_{1} + D_{Ca}); the other,
, is just the inverse length-constant of the buffer proportional to the kinetic term 1/
; and the last factor is just the single channel flux Φ. Consequently, the following conditions lead to large saturations: large current, short buffer length-constant (fast kinetics), and low diffusion coefficient. Comparison with Equation 12 also reveals that the buffer saturation at the source is identical with the expression for the buffer concentration deflections for long distances from channel, evaluated at r = 1/
.

Because our approach is based on the small saturation assumption and the buffer saturation is maximal at the source, Equation 13 can be used to check, as a sufficient condition for the validity of the method, whether the expected saturation is indeed going to be small.

**VII. For distances from the channel that are comparable with the length-constants, the action of the buffers is to produce a “relay race diffusion”: a buffer takes calcium from one with a shorter length-constant and hands it over to one with the longer length-constant.** So far, in sections IV–VI, we have investigated the behavior of the steady-state gradients for distances much longer or shorter than the characteristic buffer length-constants and developed a formalism that allows us to check for the buffer saturations, and hence the applicability of our approach. In this section, the objective is to look at distances from channel, which are of the same order of magnitude as the buffer length-constants, and explore the main properties of our linearized solution to reveal some insights into what the characteristic buffering length means. This is best demonstrated by studying the fluxes of different calcium-carrying species. We choose to illustrate the general principles for buffer conditions that match whole-cell recording situations in bovine chromaffin cells in the sequel.

We assume the presence of 0.5 mm of a fast, slowly mobile endogenous buffer (Zhou and Neher, 1993) as well as 2 mmATP in the cell. Figure 1A illustrates the differential effect of the successive addition of different buffers. We have plotted δ[Ca^{2+}] · r (computed according to Equation EAI.10) over the distance r to eliminate the “ 1/r -law” and show the pure buffer effect. Starting with 2 mm ATP as the only buffer, one sees that it is buffering even at distances within nanometers of the channel in accordance with its length-constant of 10 nm (*topmost curve* in Fig. 1A). The addition of 0.5 mm of a low-affinity, poorly mobile, fast endogenous calcium buffer shows up in the range between 20 and 200 nm and is of rather small amplitude. With this buffer background, we then add either 2 mm EGTA or 2 mm BAPTA (Table1). Note that both buffers have very similar binding ratios (4600 and 4300, respectively) and thus give rise to similar concentration deflections for large distances from the source according to Equation 12. Under the specified conditions, the length-constant for EGTA is 419 nm and for BAPTA 28 nm. In the case of BAPTA, the buffering is completely dominated by BAPTA above 30 nm, giving rise to steep gradients. In the case of EGTA, one clearly sees the multiphasic decay: under 100 nm, ATP and the endogenous buffer are doing the job, whereas EGTA produces only small, almost linear gradients between 100 and 400 nm. At very small distances, i.e., within nanometers of the channel which is below all the length-constants, none of the buffers can act significantly and all curves overlap and approach Φ/4πD_{Ca}. To demonstrate the relative contribution of the different buffers in the above system, we have plotted the individual fluxes of all calcium-carrying species in Figure 1B. They were calculated using the explicit representation of Equation EAI.10 together with Equation 10 and can easily be seen to be given by (Φ_{i}: flux of the i-th species, i = 1, … , N + 1):
Equation 16
The total flux corresponds to 3.1 × 10^{6}calcium ions/sec and is spatially constant. Initially, ATP carries calcium away from the source and has almost 42% of the calcium flux (1.3 × 10^{6} ions/sec) at 50 nm. As one moves away from the source, the calcium ions are progressively taken over by the endogenous buffer, which has its maximal flux between 200 and 300 nm. In accordance with its large length-constant, EGTA slowly takes over most of the calcium load, corresponding to its big binding ratio. This differential buffering within distances in the range of the length-constants is a kinetic effect. Equilibrium considerations here are inadequate. On the other hand, as seen in Equation 12, for distances larger than 2 μm, the relative fluxes are governed purely by the relative binding ratios and diffusion coefficients. To further illustrate the kinetic notion of buffer length-constant, we have calculated the spatial changes of the individual buffer fluxes, i.e., the spatial derivative of the fluxes, which is given by the simple expression:
Equation 17
Figure 1C plots these flux changes, normalized to the respective peak changes, as a function of distance from the channel. In the case of calcium, the absolute value of the flux changes is plotted to get positive values. As one would expect intuitively, the individual curves for the buffers peak at the length-constants of the buffers. This clarifies another interpretation of the notion of length-constant: it is the point in space where a buffer maximally changes its contribution to carrying calcium. The widths of the curves also correlate with the length-constants: the smaller the latter, the faster the changes and the smaller the half-maximal widths. The changes in the flux of free calcium, of course, are maximal in absolute value where the first buffer, i.e., the one with the shortest length-constant, maximally changes its flux. Thus, the free calcium curve peaks with the ATP curve. In a sense, these spatial flux changes constitute a “spectroscopy of the buffered diffusion problem”: as we move away from the source, by looking at the buffer flux changes we get a unique fingerprint of the individual buffers present (and their contribution to chelating calcium close to the channel) because they appear, one after the other, as individual peaks in the flux changes.

**VIII. Application: a single calcium channel cannot control release of neurotransmitters at the calycal synapse in the MNTB.** Finally, as an application, we now want to demonstrate how these theoretical considerations can be used to draw conclusions about the nature of transmitter release at a fast central synapse, namely the calyx of Held. The issue we examine is the question of whether the opening of a single calcium channel is sufficient to trigger phasic transmitter release. Because this spatial and temporal domain is not accessible to direct measurements, one needs to indirectly interfere with the transmission process to reveal some of its properties. For the calyx-type synapse in the rat medial nucleus of the trapezoid body,Borst and Sakmann (1996) have performed kinetic competition experiments by dialyzing the terminal with different concentrations of fast and slow exogenous calcium buffers of similar affinities, namely BAPTA and EGTA. They have observed the extent to which these buffers are able to compete with the endogenous calcium sensor that triggers rapid release and reduce transmission by reducing the free calcium concentration at the sensor. Their basic observation is that 10 mm EGTA is as potent in blocking synaptic transmission in these terminals as 1 mm BAPTA. We will conclude that no reasonable position for a release site can be found, which within the framework of our model would predict such a result.

Assuming a power law between the free calcium concentration and transmitter release, the result of Borst and Sakmann (1996) implies that the secretion apparatus is located at a mean distance from the channel in which both chelators can give rise to similar calcium concentrations, if a single calcium channel would be able to control release by virtue of its proximity to the calcium sensor in the microdomain. Thus, we have searched for a location, i.e., a distance from the calcium channel, in our model where similar calcium concentrations are present with 1 mm BAPTA or 10 mm EGTA. This position was found to be >1.1 μm away from the channel where the calcium concentration was <20 nmabove resting values, even with single channel currents as high as 10 pA. Studies of fast synapses, like the bipolar terminals of the retina (Heidelberger et al., 1994), demonstrate that the calcium sensor for triggering release has a rather low affinity in the micromolar range. In neuroendocrine chromaffin cells, the threshold for activating secretion is ∼500 nm (Augustine and Neher, 1992a). In addition, if 20 nm free calcium concentration above rest could trigger fast release, the cell would need to control calcium with the precision of a few, which seems to be impossible. Because these calcium concentrations cannot be responsible for fast release, the small buffer saturation regimen cannot provide a scenario for secretion in these terminals. Hence, the microdomain of a single channel cannot govern the release process unless we face significant buffer saturation on channel opening in the microdomain. The latter, in turn, is only possible with single channel currents much above 50 pA, which is a completely unreasonable proposition. Thus, a single calcium channel cannot by itself control release at a nearby release site (by virtue of close proximity between the channel mouth and the calcium sensor) in these terminals, arguing against microdomain-dominated release. Only clusters of channels can contribute to the buildup of high enough calcium concentration domains to trigger the transmission. Once clusters of channels produce sufficiently high calcium concentrations, BAPTA gets locally saturated because of its short length-constant (as we saw before), whereas EGTA is much less saturated. Eventually, this BAPTA saturation depletes it to an extent that EGTA and BAPTA can exert the same influence on secretion. With increasing source strength, the main advantage of BAPTA over EGTA in buffering calcium close to the source, namely its fast kinetics, loses importance, because buffer saturation leads to a point at which availability of free buffer is the limiting factor (for instance, in the limiting case of 100% saturation, BAPTA cannot buffer at all). Then, the binding ratio (κ) is the critical determinant of free buffer availability, and thus equilibrium considerations and buffer properties dominate the buffered diffusion (see also discussion of the “rapid buffer approximation” in the next section). Note that under the conditions of the above experiments, the exogenous binding ratio is so high (2100 and 22,000, respectively) that the buffer length-constant is determined purely by κ/(τ · D_{Ca}), which is equal to [B] · k_{on}/ D_{Ca}, the buffer product divided by D_{Ca} (see definition of μ_{1} in VI). Consequently, the effect of two buffers can be compared by multiplying the free buffer concentrations with the on-rates to obtain the buffer products as good estimates of the length-constants.

## DISCUSSION

There is an ongoing discussion regarding the spatial relation between presynaptic calcium channels and the secretion apparatus at synapses and its impact on the free calcium concentration at the calcium sensor responsible for triggering exocytosis. It is debated whether secretion is triggered by the opening of a single calcium channel, which is somehow associated with its prefusion complex, or whether the simultaneous action of many calcium channels leads to a sufficiently high calcium level to trigger exocytosis (Roberts et al., 1990; Augustine et al., 1991; Stanley, 1993; Llinas et al., 1995). The first scenario, the “single channel domain,” is supported by recent findings indicating molecular interactions of calcium channels with components of the famous SNARE-complex (Bennett et al., 1992; Rettig et al., 1996; Sheng et al., 1996). The other scenario, the “domain overlap regimen,” is suggested by experiments showing that even slow buffers of EGTA type are capable of reducing synaptic transmission at moderate concentrations (Borst and Sakmann, 1996) or the existence of synapses with low release probabilities but a high number of docked vesicles (Rosenmund et al., 1993). Unfortunately, in general, it is not possible to measure directly the free calcium concentration at the relevant release sites. The solution to this problem has been to use extensive numerical simulations of the buffered diffusion problem.

This paper outlines a linearization of the general reaction–diffusion problem that is aimed at determining the parameters shaping the calcium gradients close to a calcium channel in analytical and intuitive terms. Our work is an extension of the results of Pape et al. (1995). It is a submicroscopic theory that can only be usefully applied if we are considering distances up to a few hundred nanometers from the calcium source. The theory is specifically tailored for use in a temporal and spatial domain that is not accessible to imaging, with the objective to obviate time-consuming simulations of diffusion–reaction equations.

Because standing gradients rapidly evolve close to a channel (Appendix), we focus on steady-state concentration profiles of calcium and calcium-bound buffers. It is noted that fixed buffers do not affect the steady gradients. This is because they cannot be replenished by diffusion of free buffer. Other mechanisms affecting calcium signaling, such as active transport of calcium ions, can only have little effect within this range because it is impossible to achieve such a high density of pumps and exchangers within a small area to counteract the calcium influx within <1 msec. Hence, we can focus on the impact of mobile buffers. The main feature of our approach is the assumption that the saturation of the involved buffers, i.e., the incremental increase in the concentration of calcium-bound buffers on channel opening, is small enough to permit a linearization of the problem. This is mostly equivalent with sufficiently small single channel currents or sufficiently high buffer concentrations. Note that we do not require the increase in free calcium concentration to be small (which in general is not the case as one approaches the source). The reason for this freedom is the fact that we compute the free calcium concentration as the difference between the totally expected calcium increase (in the absence of any buffers) and the scaled sum of the calcium-bound buffer concentrations (see Eq. EAI.10). If this applies, one can deduce many properties of our solution to the reaction–diffusion problem.

(1) The increase in [Ca^{2+}] above resting values, δ[Ca^{2+}], as well as that of the calcium-bound buffers, is proportional to the source strength given by the single channel current.

(2) To each buffer present, one can assign a buffer length-constant that is determined by equilibrium (dissociation constant and total buffer concentration) and kinetic (on-rate for calcium complexation) properties of the buffer as well as its diffusional mobility. It mainly represents the average distance a calcium ion will diffuse before it is captured by the buffer and the average distance the buffer will diffuse until it gets into local equilibrium with calcium.

(3) The shape of the gradients surrounding a channel is determined by a multiexponential law (with rates given by the inverse length-constants) superimposed on a “ 1/r -law,” which one expects in the absence of any buffers.

(4) For distances large compared with the length-constants, the concentration deflections above resting values are determined purely by equilibrium buffer properties, i.e., the on-rates of the buffers have no influence on the concentrations, because there is enough time for the reactions to reach chemical equilibrium.

(5) On the other side, for distances short compared with the length-constants, the buffers are not able to bind calcium. This is because calcium ions will diffuse to regions of lower [Ca^{2+}] within the time the buffers need on average to bind calcium. As a consequence, δ[Ca^{2+}] will behave as it does in the “no buffer regime,” showing the pure 1/r dependence.

(6) The individual buffer saturations for a specified single channel current and defined buffer conditions depend heavily on the individual length-constants. The bigger the length-constants, the smaller the buffer saturations. Maximal buffer saturations can be computed on the basis of the buffer parameters to check in advance whether the theory will be applicable.

(7) The solution to the problem can be written analytically in terms of the buffer parameters. Numerical computations can be reduced to standard eigenvalue problems involving spectral decomposition of matrices. The latter can be done most effectively using existing program packages and thus eliminates the need for expensive discretizations of the involved differential equations.

(8) Finally, in case all diffusion coefficients are identical, the temporal evolution of the gradients in the vicinity of an open channel can be written as a convolution product of a purely reactive term and a purely diffusive term, which permits effective computation of the time courses.

If the main assumption is fulfilled, one can use the present approach to estimate [Ca^{2+}] at different distances from the channel and explore the differential effect of buffers on synaptic transmission. But when is the main assumption fulfilled? By tolerating, say not more than 20% increase in [CaB], the following rule of thumb can be stated: with 100 μm of an EGTA-type buffer, 4 pA single channel current is allowed, and with 100 μm of a BAPTA-type buffer, 0.3 pA calcium current is allowed. It should be noted, however, that these are really conservative estimates. From simulations of reaction kinetics and comparison of the linearized solution with the full numerical solution (not shown here), we actually expect the linear approximation to be valid over a much wider range of buffer saturations, and, of course, as the concentration of mobile buffers increases, higher single channel currents are tolerable.

An alternative approach to the problem of understanding the steady gradients surrounding a channel is the “rapid buffer approximation (RBA)” (Smith et al., 1996). Here, one assumes that the buffers are so fast that they are in chemical equilibrium at every point in time and space. In other words, if the calcium gradient is given by an exponential with a length-constant L and the reaction time constant of the buffer is τ, then one requires τ ≪ L^{2}/(2D_{Ca}). This simplifies the general reaction–diffusion system to one nonlinear partial differential equation. It has its own merits if one is not facing sharp gradients. For instance, if we have 1 mm BAPTA, then τ = 4 μ sec and L must be >150 nm for the validity of RBA; however, close to the channel, i.e., between 10 and 150 nm, [Ca^{2+}] falls exponentially with the length-constant of BAPTA, which is 30 nm. This implies that we cannot apply RBA this close to the source. Figure2 illustrates this point by comparing directly the RBA solution according to Smith (1996) with our solution. Clearly, RBA is underestimating [Ca^{2+}] (and overestimating [CaBAPTA]) because it assumes chemical equilibrium of the reactions where, as we have seen, no equilibrium is attainable. On the other hand, for long distances, our solution matches the RBA as is visible in Figure 2 for r > 150 nm. Indeed, it can be shown that Equation 10 of Smith (1996) is identical to our Equation 12. Thus, our approach is valid for small saturation levels (or sharp calcium gradients). On the other hand, if one is facing high single channel calcium currents, which almost completely saturate the mobile buffer over extended regions, RBA would be appropriate to explore the microdomain properties.

Finally, let us outline how one could handle channel clustering. The source in our considerations so far has been a point source of calcium. Nevertheless, there is evidence that in some systems, tens or hundreds of calcium channels can be clustered (Tucker and Fettiplace, 1995). We cannot treat them as a single channel and accordingly scale up the single channel current because we would then ignore the spatial arrangement of the channels. Fortunately, the linearized scheme has an intrinsic feature that still allows us to cope with these circumstances, namely the well known superposition principle. We can compute calcium and buffer concentrations in response to the opening of a single channel in the cluster and then consider a specific spatial arrangement of an array of channels and sum up the individual concentration deflections by virtue of one channel opening to get the profiles in response to the activation of the whole channel cluster. This would allow us to rescue our intuition within this model to complicated channel arrangements without any need to discretize differential equations.

## Appendix

In this section, we solve the system (Eq. 9) in steady-state, subject to the condition (Eq. 10) to derive an analytical solution for the standing gradients surrounding an open calcium channel in its microdomain. Equation 9 in steady-state, i.e., withδ˙y = 0, can be written as:
Equation AI.1This is a system of ordinary differential equations that can easily be simplified by the transformation δz ≡ r · δy to give:
Equation AI.2a linear second-order system with constant coefficients given by B = (b_{ij})_{i,j=1}^{N+1}. Furthermore, integrating Equation 10 from r to ∞ and realizing that at infinity, the concentration deflections above resting values are zero, i.e., δy_{i}(r = ∞) = 0, we can rewrite (Eq. 10) as:
Equation AI.3Equation EAI.3 is the transformed constraint subject to which Equation EAI.2 needs to be solved. We make the Ansatz δz = exp(− Lr)**u** for some matrix L ∈ R^{(N+1)×(N+1)} and a vector**u** ∈ R^{N+1}. Differentiating δz with respect to r twice according to Equation EAI.2, we find that the following identity must be valid:
Equation AI.4
L is called the (positive) square root of the matrix B, L ≡
, and is determined by the square root of the eigenvalues of B and the corresponding eigenvectors (Golub and Van Loan, 1989). The vector **u**incorporates the constraint (Eq. EAI.3) and must satisfy the identity:
Equation AI.5It is easy to see that the vector **u** is given by:**u** = Φ/(4πD_{N+1}) · (0 0 · 0 1)^{t}. The last equation effectively reduces the dimensionality of the homogenous problem (Eq. EAI.2) by projecting it from the N + 1 -dimensional Euclidean space to an N -dimensional submanifold. Mathematically, this is equivalent to a nonhomogenous problem in the N -dimensional space, the latter having the advantage of closed-form solutions. Incorporating Equation EAI.5 into the system (Eq. EAI.4), we can rewrite the latter in the N -dimensional space spanned by δz_{1}, … , δz_{N} in terms of a new matrix differential equation given by:
Equation AI.6
Equation AI.7
Explicitly, C and ϕ are given by:
The solution to Equation EAI.6, which is finite at r = ∞, is given by:
Equation AI.8with some **v** ∈ R^{N}. The vector **v** is specified by the condition that for small distances r from the channel, all of the flux is carried by calcium, i.e.:
Equation AI.9and is given by **v** = C^{−1}ϕ. Inserting this in Equation EAI.8 with δz ≡ r · δy finally gives us the desired closed-form solution to the linearized, steady-state buffered diffusion problem as:
Equation AI.10
We can also arrive at a more intuitive representation of this solution by using a spectral decomposition of the matrix B = −D^{−1} · A · B has the spectral representation as B = S · Λ_{B} · S^{−1}, where Λ_{B} is the diagonal matrix composed of the eigenvalues μ_{i} of B and S = (**u**_{1}**u**_{2} ·**u**_{N}**u**_{N+1}) with**u**_{i} the eigenvector of B corresponding to the eigenvalue μ_{i}. The general solution to the system (Eq. EAI.2) can then be written as:
Equation AI.11
where a_{i} are some real numbers that incorporate the boundary condition (Eq. EAI.5). It is important to note that μ_{N+1} = 0 is an eigenvalue of B with the corresponding eigenvector**u**_{N+1} = (κ_{1} κ_{2}· κ_{N} 1)^{t}, which is just the vector of equilibrium binding ratios. Futhermore, it is easy to see that the vector **a** = (a_{1} a_{2} · a_{N} a_{N+1})^{t}satisfies the equation S · **a** =**u**[**u** as in (EAI.5)] and is uniquely determined by the same. Thus, our solution is a sum of exponentials with unique coefficients determined by the calcium flux condition, divided by the distance r. For N = 1, we get:
Equation AI.12
which is the result of Pape et al. (1995), written in other terms.

## Appendix

The purpose of this section is to derive the transient solution of the linearized buffered diffusion problem for the case of an arbitrary number of mobile buffers with equal diffusion coefficients but different kinetics. The condition of equality of the diffusion coefficients has been used in many simulation studies (Nowycky and Pinter, 1993) because most exogenous mobile buffers (such as EGTA and BAPTA but also ATP) and calcium ions have estimated intracellular diffusion coefficients ∼200 μm^{2}/sec. Here, it is used merely to simplify the algebra of our analysis.

We start with Equation 9 describing the dynamics of the linearized reaction–diffusion problem and substitute δy = exp (At)δz, with A representing the chemical relaxation matrix according to Equation 8. This substitution eliminates the reaction term in Equation 9 and results in a pure diffusion equation given by:
Equation AII.1where d is now the diffusion coefficient of the buffers, i.e., d = D_{1} = … = D_{N+1}. Next, we perform a Laplace transform of EquationEAII.1 with respect to t. If δ̂z denotes the Laplace transform of δz, and realizing that δz(t = 0, r) = δy(t = 0, r) = 0 for all r, we get:
Equation AII.2the solution of which is easily seen to be:
Equation AII.3for some vector B(s) ∈ R^{N+1} yet to be determined. This is achieved by considering the calcium flux condition. Again, at r = 0 and for t > 0, we have a constant flux Φ of calcium without any sources or sinks for the buffers. Hence, we have the flux of all calcium-carrying species given by the vector:
Equation AII.4
Equation AII.5Equation EAII.3 gives us the Laplace transform of δz, hence we transform Equation EAII.5 and compute (∂/∂r)δ̂z to determine B(s):
Equation AII.6
The last equation gives together with Equation EAII.3 the Laplace transform of the desired result:
Equation AII.7The inverse transform of e^{−r
} is g(t, r) = dr · e^{−r2}^{/4td}/(2td ·
) and the inverse transform of (A + sId)^{−1} · (0 . 0 1)^{t} is just e^{−At}(0 . 0 1)^{t}. This together with the convolution theorem of Laplace transformation finally results in the following representation of δz(t, r) as the convolution product:
Equation AII.8
The computation of δy(t, r) can be done efficiently by performing a spectral decomposition of the chemical relaxation matrix A. Because A is uniquely determined by the kinetic properties of N different mobile calcium buffers, it has N + 1 distinct eigenvalues, including zero. Hence, it has the spectral representation as:
Equation AII.9where Λ_{A} is the diagonal matrix of eigenvalues of A, i.e.,
Equation AII.10and the i-th column of T is the eigenvector of A corresponding to the eigenvalue λ_{i}. EquationEAII.8 can then be rewritten as:
Equation AII.11This implies that the computation of δy(t, r) consists of the following steps. First, we compute A, its eigenvectors and eigenvalues, as well as T^{−1}(0 0 . 1)^{t}. Then, we perform N + 1 numerical integrations of the integrals ∫_{0}^{t} g( u , r)e^{λi}^{u}d u and subsequently the matrix products.

Let us consider the case of one mobile buffer, i.e., N = 1, to estimate the time needed to achieve steady-state. The eigenvalues are: λ_{1} = −(1 + κ_{1})/τ_{1}, λ_{2} = 0, and the corresponding eigenvectors are e_{1} = (−1 1)^{t}, e_{2} = (κ_{1}1)^{t}. Then:
Equation AII.12Finally, using Equation EAII.11, we get the following representation of the transient concentrations after channel opening:
Equation AII.13
Figure 3 depicts the expected time courses of relaxation of δ[Ca^{2+}] and δ[CaB] toward steady-state, according to Equation EAII.13, at different spatial positions for an EGTA-type buffer, present at 2 mm total concentration. We have computed the ratio of the instantaneous concentration of δ[Ca^{2+}] (A) and δ[CaB] (B) over the expected steady-state concentration. These fractions of the steady-state concentration values, which are independent of the single channel calcium flux, are plotted as pseudocolored images as a function of time and distance from channel. Superimposed on the images are the lines where a constant fraction (indicated by the numbers) of the steady-state concentrations are achieved. As is visible in A, even at 500 nm distance, 90% of the steady-state calcium concentration is reached within 0.5 msec. We can state that within a few hundred nanometers from the channel, calcium gradients reach steady-state within <1 msec.

## Footnotes

We thank Dr. Christian Rosenmund for helpful comments on this manuscript.

Correspondence should be addressed to Dr. Erwin Neher, Department of Membrane Biophysics, Max-Planck-Institute for Biophysical Chemistry, Am Fassberg 11, D-37070 Göttingen, Germany.