Abstract
The rate of exocytic events from both neurons and nonneuronal cells exhibits fluctuations consistent with fractal (selfsimilar) behavior in time, as evidenced by a number of statistical measures. We explicitly demonstrate this for neurotransmitter secretion atXenopus neuromuscular junctions and for rat hippocampal synapses in culture; the exocytosis of exogenously supplied neurotransmitter from cultured Xenopus myocytes and from rat fibroblasts behaves similarly. The magnitude of the fluctuations of the rate of exocytic events about the mean decreases slowly as the rate is computed over longer and longer time periods, the periodogram decreases in powerlaw manner with frequency, and the Allan factor (relative variance of the number of exocytic events) increases as a powerlaw function of the counting time. These features are hallmarks of selfsimilar behavior. Their description requires models that exhibit longrange correlation (memory) in event occurrences. We have developed a physiologically plausible model that accords with all of the statistical measures that we have examined. The appearance of fractal behavior at synapses, as well as in systems comprising collections of synapses, indicates that such behavior is ubiquitous in neural signaling.
 quantal secretion
 vesicular exocytosis
 miniature endplate currents
 fractal
 longterm correlation
 lognormal process
 Xenopus neuromuscular junction
 myocyte autoreception
 fibroblast
 hippocampal synapse
Communication in the nervous system is mediated by action potentialinitiated exocytosis of multiple vesicular packets (quanta) of neurotransmitter (Katz, 1966). Even in the absence of such action potentials, however, many neurons spontaneously release individual packets of neurotransmitter (Fatt and Katz, 1952). A packet may contain from 7000 to 10,000 molecules of acetylcholine (ACh), if we use the neuromuscular junction as an example (Kuffler and Yoshikami, 1975). On arrival at the postsynaptic membrane, the ACh molecules induce elementary endplate currents (EECs), which take the form of nonstationary two (or multi) state on–off sequences (Sakmann, 1992). Current flows when the ACh channel is open (i.e., when its two binding sites are occupied by agonist) and ceases when the channel is closed. A postsynaptic miniature endplate current (MEPC) comprises some 1000 EECs (Sakmann, 1992). It was shown by Del Castillo and Katz (1954) that superpositions of MEPClike events comprise the postsynaptic endplate currents elicited by nerve impulses.
It generally has been assumed that the sequence of MEPCs forms a memoryless stochastic process (Fatt and Katz, 1952). However,Rotshenker and Rahamimoff (1970) discovered that exocytosis in the frog neuromuscular junction can exhibit correlation (memory) over a period of seconds, provided that extracellular Ca^{2+} levels are elevated above their normal values. Cohen et al. (1974b)subsequently found such correlation (event clustering) even in the absence of elevated Ca^{2+} levels.
In this paper we study the statistical properties of exocytic events over a far larger range of time scales than previously examined. MEPCs were recorded from innervated myocytes in Xenopusnerve–muscle cocultures and from rat hippocampal neurons in cell culture. MEPCs from nonneuronal preparations also were examined: the quantal secretion of ACh from isolated myocytes (autoreception) and from rat fibroblasts, both exogenously loaded with ACh (Dan and Poo, 1992; Girod et al., 1995).
We direct particular attention toward those statistical measures that reveal the presence of memory. Our analysis reveals that the time sequences of the MEPCs, and therefore of the underlying exocytic events, exhibit memory that appears to decay away slowly in both neuronal and nonneuronal cells. This longduration correlation is present over the entire range of time scales investigated, which stretches to thousands of seconds. The occurrence of an MEPC, therefore, makes it more likely that another MEPC will occur at some time thereafter. The analysis of long MEPC data sets reveals that the rate of exocytic events behaves in a manner consistent with a fractal process, exhibiting fluctuations over multiple time scales. Fractals are objects that possess a form of selfsimilarity: parts of the whole can be made to fit to the whole by shifting and stretching. The hallmark of fractal behavior is powerlaw dependence in one or more statistical measures over a substantial range of the time (or frequency) scales at which the measurement is conducted (Lowen and Teich, 1995; Thurner et al., 1997). Because multiscale fluctuations are at the heart of this behavior, selecting short data segments that exhibit minimal fluctuations will dilute whatever fractal characteristics might be present in a given data set, as we illustrate. The classic work of Fatt and Katz (1952) is revisited in light of these findings.
MATERIALS AND METHODS
Xenopus nerve–muscle cocultures. Cultures were prepared by following previously reported methods (Spitzer and Lamborghini, 1976; Anderson et al., 1977; Tabti and Poo, 1991). In brief, the neural tube from 1dold embryos (stage 20–24; Nieuwkoop and Faber, 1967) was dissociated in Ca^{2+}/Mg^{2+}free Ringer’s solution supplemented with EDTA, plated on clean glass coverslips, and incubated at 20–22°C for 1 d before recording. Recording and culture medium consisted of 50% (v/v) Leibovitz’s medium (Life Technologies, Gaithersburg, MD), 1% (v/v) fetal calf serum (Life Technologies), and 49% (v/v) Ringer’s solution (115 mm NaCl, 2 mm CaCl_{2}, 2.5 mm KCl, and 10 mm HEPES, pH 7.3).
Hippocampal cultures. Hippocampal cultures were prepared by following protocols previously described (Goslin and Banker, 1991), but with the following modifications. Hippocampi from newborn Sprague Dawley rats (postnatal day 0) were dissociated and plated at low density, without glial support cells. Cells from 7 to 14dold cultures were used for experiments. The recording solution consisted of (in mm) 140 NaCl, 5 KCl, 1 CaCl_{2}, 1 MgCl_{2}, 10 HEPES, 24 dglucose, and 0.01 TTX (Sigma, St. Louis, Mo), pH 7.4.
Fibroblast cultures. Parenteral 3Y1 cells, a line of rat skin fibroblasts (Sternberg et al., 1993) were kindly provided by Paul Greengard (Rockefeller University, New York, NY). The cells were cultured in DMEM and 10% fetal calf serum at 37°C with 5% CO_{2} and were used 2–3 d after subculture. For detection of ACh secretion from fibroblasts, the latter were treated with trypsin–EDTA solution (Life Technologies); the suspension of the fibroblasts was added to the recording chamber that contained culturedXenopus myocytes.
Electrophysiology. Recordings of miniature endplate currents at Xenopus neuromuscular junctions were made with the gigaohm seal, nystatin perforatedpatch technique (Horn and Marty, 1988). Conventional heatpolished patch pipettes were filled with pipette saline containing (in mm) 150 KCl, 1 NaCl, 1 MgCl_{2}, 10 HEPES, pH 7.4, and nystatin (Sigma). Nystatin was added at a final concentration of 460 μM (final concentration of DMSO, 1%). Nystatin stock (46 mm; 1 mg/ml in DMSO) was prepared before each experiment, stored at room temperature in a lightproof container, and used for up to 6 hr after preparation. A conventional gigaohm seal was formed by pressing the pipette gently against the myocyte and providing light suction. Spherically shaped myocytes were selected.
For all other experiments, including hippocampal cells and nonneuronal cells, the conventional gigaohm seal wholecell recording method was used (Hamill et al., 1981). For experiments in which MEPClike events were recorded from an isolated myocyte loaded intracellularly with ACh (autoreception), the intrapipette solution containing (in mm) 150 KCl, 1 NaCl, 1 MgCl_{2}, and 10 HEPES, pH 7.4, was supplemented with 20 mm AChCl (Dan and Poo, 1992). The intrapipette solution used for recording from hippocampal neurons consisted of (in mm) 155 Kgluconate, 1 MgCl_{2}, 10 HEPES, 1 sodium gluconate, 5 MgATP, 0.5 NaGTP, and 0.1 leupeptin, pH 7.4. The internal solution in the patch pipette for fibroblasts contained (in mm) 105 Kgluconate, 1 Nagluconate, 10 HEPES, 5 MgATP, and 0.5 NaGTP, pH 7.4, supplemented with 50 mm AChCl. The detached fibroblast was patched before being manipulated into contact with an isolatedXenopus myocyte. The recording pipette for the myocyte contained (in mm) 150 KCl, 1 NaCl, 1 MgCl_{2}, and 10 HEPES, pH 7.4.
In all experiments, currents were recorded at room temperature under voltageclamped conditions (holding voltage V _{h} = −70 mV for myocytes, −80 mV for fibroblasts, and −65 mV for hippocampal neurons), with an Axopatch 1D amplifier (Axon Instruments, Foster City, CA). The currents were filtered at 1 kHz, digitized, and stored on videotape for later playback. Computer analysis was performed with the SCAN program kindly provided by Dr. J. Dempster (Strathclyde University, UK). This analysis resulted in a series of exocytosis event times, which were quantized to 55 msec: multiple events falling within a single 55 msec window were registered as a single event. To guarantee sufficient statistical accuracy for estimating parameters of the fractal behavior, we retained only data sets with N ≥ 400 events for analysis.
RESULTS
Spontaneous vesicular exocytosis from neuromuscular junctions, hippocampal synapses, and nonneuronal cells
Within 1 d of plating, Xenopus embryonic spinal neurons establish functional synaptic transmission with cocultured myocytes, exhibiting stable spontaneous and evoked release characteristics (Chow and Poo, 1985; Xie and Poo, 1986). Spontaneous pulsatile membrane currents resembling MEPCs (Fig. 1) developed in a wholecell voltageclamped myocyte within several minutes of forming a tight seal with a pipette containing normal intracellular solution supplemented with nystatin. The configuration was permitted to stabilize for 10 min before experiments were begun. The observed current pulses were virtually identical in their rate and amplitude distributions to events recorded by the conventional wholecell recording configuration.
These pulsatile current events represent spontaneous exocytosis of AChcontaining synaptic vesicles at the developing neuromuscular junction, because their rate and amplitude distribution are not affected by the addition of TTX (data not shown; but see Xie and Poo, 1986). The large amplitude variability presumably results from immature filling of the synaptic vesicles (Evers et al., 1989). (Vesicles containing an unusually small amount of ACh could result in a current event with a magnitude that lies below the threshold of detectability of the recording system, thereby escaping detection.) In the lowdensity cultures we use, each myocyte is innervated by a single neuron. Thus, the neurotransmitter release that is detected arises from a single synapse, which, in general, comprises a number of individual release sites.
Similar spontaneous pulsatile inward currents were observed from isolated myocytes exogenously loaded with ACh (Dan and Poo, 1992). MEPClike events appeared within several minutes after establishing the wholecell configuration. The average amplitude and rate of these events increased with time thereafter, ceasing to exhibit systematic changes after ∼10 min. Analysis of spontaneous release from the cells was begun after a 15 min loading period to ensure that stability had been established. In this preparation the inward currents result from the spontaneous, quantal release of ACh packets from the myocyte and the subsequent detection of this release by activation of its own surface ACh receptors (autoreception).
Spontaneous quantal ACh secretion from exogenously loaded fibroblasts also was observed (Girod et al., 1995) with the help of the ACh detection system inherent in Xenopus myocytes. Wholecell voltageclamped recordings from a myocyte in contact with the fibroblast displayed transient inward currents resembling MEPCs. The frequency and amplitude of these events also increased gradually with time during a period of ∼30 min. The events were recorded after this period to permit stability to be achieved. It appears that the secretion events result from the sequestration of ACh into vesicles and the subsequent exocytosis of those vesicles. It is likely that the vesicles are components of the constitutive endocytosis or membrane recycling pathways.
The analysis and comparison of the secretion in the three forms of vesicular release discussed above support the view that the basic machinery for secretion, both in neuronal and nonneuronal cells, is relatively similar (Girod et al., 1995).
Inward currents also were observed in 7 to 14dold hippocampal neurons cultured at low density. These currents are also the result of spontaneous neurotransmitter release, because they are apparent even in the presence of TTX. In this case, however, the neurotransmitter was most likely glutamate, because the addition of the AMPA receptor blocker CNQX and the NMDA receptor blocker AP5 abolished the currents. Unlike the spontaneous currents observed at neuromuscular junctions inXenopus cell culture, it is likely that these MEPCs arise from excitation by several presynaptic neurons. Even at low densities it is not possible to trace exact pre and postsynaptic neuronal cell pairs. Moreover, by 7 d in culture, an extensive network of neurites forms, making it likely that each neuron receives input from a multitude of neurons.
We have analyzed the statistical patterns of the spontaneous secretion events generated in all of these preparations and have formulated a suitable model for the exocytic behavior: the fractal lognormalnoisedriven doubly stochastic Poisson process (FLNDP). We begin by examining various statistical measures of the sequence of MEPCs observed in the Xenopus neuromuscular junction.
Intereventinterval histogram (IIH)
The solid curve in Figure 2 is a semilogarithmic plot of the MEPC IIH, arising from spontaneous vesicular release activity observed in a typical Xenopusneuromuscular junction. The IIH is a measure of the relative frequency of the times between successive events. The sequence of events from which this histogram was constructed has a duration L = 8164 sec and contains N = 2644 interevent intervals, thereby exhibiting a mean interevent interval E[t] = 3.09 sec and a mean rate λ = 1/E[t] = 0.324 events/sec.
Traditional mathematical descriptions of vesicular exocytosis generally assume that the MEPC sequence forms a renewal process. Renewal processes are memoryless; successive intervals are all independent and are drawn from a single distribution. They are, therefore, characterized completely by the IIH, which is an estimate of the intereventinterval probability density p(t).
The simplest renewal model is the homogeneous Poisson point process (HPP; Cox and Lewis, 1966). The HPP is characterized by a single constant quantity, its rate λ, which is the number of vesicular release events expected to occur in a unit time interval. The HPP intereventinterval probability density functionp(t) behaves as a decreasing exponential functionp(t) = λexp(−λt),t ≥ 0, where t is the interevent interval and E[t] = 1/λ is the mean interevent time. Incorporating the effects of dead time (absolute refractoriness) or sick time (relative refractoriness) in the process preserves the exponential tail of the intereventinterval distribution while suppressing p(t) for shorter times. Because the exponential function plotted on semilogarithmic coordinates is a straight line, the HPP model clearly does not provide a good fit to the experimental IIH presented in Figure 2. Cohen et al. (1974a) reached this same conclusion from their studies of the frog neuromuscular junction.
A reasonable fit can be obtained if a slightly more complex renewal process, the gamma renewal process (GRP; Cox and Lewis, 1966), is used. This is the approach taken by Hubbard and Jones (1973). The theoretical IIH is then the twoparameter gamma density functionp(t) =t ^{a−1}exp(−t/τ)/Γ(a)τ^{a}, where τ is a characteristic time, a is a parameter known as the order of the gamma process, and Γ(·) is the gamma function. The mean interevent interval for this distribution is E[t] = aτ; its variance is Var[t] = aτ^{2}. (The particular case a = 1 corresponds to the exponential distribution, illustrating that the HPP is a special case of the GRP.) The gamma IIH, with its mean and variance set equal to those of the data, is shown as the dashed curve in Figure 2. The fit is very good. Thus, were the MEPC sequence renewal in nature, it would be describable by a GRP, and nothing more need be said about it.
However, we demonstrate that the MEPC sequence is not renewal using three statistical measures that are sensitive to the presence of memory in a point process: vesicular releaserate measurements over different time scales, the Allan factor (AF), and the periodogram (PG). The dependencies among the interevent intervals evidenced by these measures reveal that a fractalrate stochastic point process (FRSPP) (Teich et al., 1996a; Thurner et al., 1997) represents the sequence of MEPCs. A fractalrate process is required because the IIH does not follow a powerlaw form and hence does not exhibit scaling. This indicates that the sequence of exocytic events itself does not form a fractal in time. However, the rate of event occurrencesis consistent with scaling behavior, leading to the FRSPP model for the MEPC sequence.
Indeed a particular FRSPP, the FLNDP, provides an excellent representation for all of the statistical measures of the exocytic events that we have investigated, as discussed at the very end of this section. The associated threeparameter simulated intereventinterval density function is displayed as the dotted curve in Figure 2; it clearly provides an excellent fit to the IIH, even a bit superior to that of the GRP. Because the MEPCs do not form a renewal process, however, it is clear that the IIH alone is inadequate for choosing among alternative models.
Selfsimilarity of vesicular release rates
Perhaps the simplest measure of neuronal activity is the estimate of the rate: the number of events registered per unit time (Teich, 1992). For vesicular release, even this straightforward measure is consistent with the presence of fractal properties; the magnitude of the fluctuations of the rate (e.g., its standard deviation) decreases more slowly, as the counting time used to compute it increases, than would be expected for independentevent counts.
In Figure 3 a we illustrate the vesicular release rate for the same Xenopus neuromuscular data as those analyzed in Figure 2. Two different counting times were used to compute the rate: T = 25 sec (solid curve) and T = 250 sec (dashed curve). The total duration of the solid curve is 775 sec (31 consecutive samples, each of 25 sec), whereas that of the dashed curve is 7750 sec (31 consecutive samples, each of 250 sec). Evidently, increasing the averaging time by a factor of 10 reduces the magnitude of the fluctuations only slightly. This indicates that longduration fluctuations are present in the train of vesicular release events, consistent with a fractal rate.
This behavior derives from dependencies among the interevent intervals, as confirmed by using a surrogate data set; the fractal properties of the rate estimate are destroyed by shuffling (randomly reordering) the intervals. This operation removes the dependencies among the intervals while exactly preserving the intereventinterval histogram. With all dependencies among the intervals eliminated by shuffling (aside from those inherent in retaining the same IIH), the resulting surrogate essentially behaves as a renewal point process. The vesicular release rate for these same data, after such shuffling, is illustrated in Figure 3 b. The T = 250 sec shuffled data (dashed curve) now exhibits noticeably smaller fluctuations than does the T = 25 sec shuffled data (solid curve). This more rapid reduction in the magnitude of the fluctuations with larger averaging time is typical for nonfractal rates.
The quantitative behavior of the magnitude of the rate fluctuations with counting time is considered more conveniently in terms of the AF, which is discussed next.
Powerlaw behavior of the AF
A highly useful measure that is sensitive to correlations in a point process is the AF, a relative variance based on a particular wavelet transform (Teich et al., 1996a). Its definition and properties are presented in Appendix . The AF calculated at a particular counting time T provides a quantitative measure of the variability exhibited by the rate estimates displayed in Figure 3 a. For general well behaved processes, the AFA(T) is a function of the counting time T; the unique exception is the HPP, for whichA(T) = 1 for all counting timesT. Useful values of the counting time T range from onehalf of the minimum interevent interval to approximately onetenth of the duration of the recording.
The solid curve in Figure 4 is the AF for the same sequence of Xenopus neuromuscular junction MEPCs, the IIH of which is shown in Figure 2 and the rate functions of which are shown in Figure 3. The AF is seen to increase steadily for counting times greater than ∼10 sec, exceeding a value of 100 at a counting timeT = 400 sec. For sufficiently large counting times, the AF is well approximated by a straight line on this doubly logarithmic plot, so that it is well fit by an increasing powerlaw function of the counting time, A(T) ∝T ^{αA}, with α_{A} ∼1.5 for this particular junction. A monotonic, powerlaw increase indicates the presence of fluctuations on many time scales. The quantity α_{A} is identified as an estimate of the fractal exponent of the point process (Lowen and Teich, 1995; Thurner et al., 1997).
A plot of the AF alone does not reveal whether its substantial magnitude arises from the distribution of the interevent intervals (the IIH) or from their ordering. This issue is addressed by plotting the AF for the shuffled intervals. AFs constructed from shuffled data retain information about the relative sizes of the intervals, but all correlations and dependencies among the intervals are destroyed by the shuffling process, as discussed earlier. The curves designated “shuffled data” in Figure 4 illustrate the AFs obtained by this method; the long dashes indicate the mean value obtained from 100 shufflings, and the short dashes delineate a range of ± 1 SD about this mean value. The lack of observed fractal behavior in the AFs of the shuffled data indicates that it is the ordering of the intervals that gives rise to the powerlaw growth of the AF for those records. This confirms that the original data cannot be modeled as a renewal process.
The longdashed curve behaves very much like the AF for a GRP, approaching a maximum value close to 1/a = 4.01 for large counting times, as expected from theory (Teich et al., 1997). Indeed, AF simulations for the GRP that best fits the IIH (longdashed curve in Fig. 2) closely resemble the longdashed curve in Figure4. It is clear, therefore, that the GRP provides a good fit to the AF only for the shuffled data and therefore cannot possibly describe the original (unshuffled) data. In contrast, the simulated FLNDP AF shown in Figure 4 (dotted curve) follows the original AF much more closely; moreover, a shuffled version of the FLNDP also provides a good fit to the shuffleddata AF (shuffled simulation not shown).
Powerlaw behavior of the PG
Fractal variability of the exocytic events also manifests itself in other statistical measures, perhaps the most familiar of which is the PG, which is an estimator of the power spectral density. Much as for continuoustime processes, the power spectral density computed for the exocytic events reveals how power is concentrated in various frequency bands. In Figure 5 we present a doubly logarithmic PG plot for the same Xenopus data sequence (solid curve) examined in Figures 2, 3, 4. For low frequenciesf (corresponding to long time scales T), the PG is well fit by a decreasing powerlaw function of the frequency,S(f) ∝f ^{−αS}. Thus, MEPC activity exhibits 1/ftype noise. The quantity α_{S} provides an alternative means of estimating the fractal exponent of the vesicular exocytosis process. For this particular neuron, the PG yields α_{S} ∼1.6, which is in close accord with the value α_{A}∼1.5 obtained from the AF, as expected (Teich et al., 1996a; Thurner et al., 1997).
The PG computed from a shuffled version of the data (dashed curve), in contrast, is quite flat at low frequencies, providing further evidence that it is the ordering of the intervals, rather than their relative magnitudes, that is responsible for the fractal aspects of the rate of vesicular activity. Again, the simulated FLNDP PG shown in Figure 5 (dotted curve) provides a good fit to the experimental PG, and shuffled versions of the FLNDP also lead to results that accord with the shuffleddata PG (shuffled simulations not shown).
Because the PG is the Fourier transform of the joint coincidence rate (a measure of correlation used for a process of events), the results presented here are not inconsistent with those obtained by Rotshenker and Rahamimoff (1970; their Fig. 1), who showed excess correlation to 5 sec in preparations subjected to extracellular Ca^{2+}levels above their normal values. Evidence for excess power at low frequencies in the absence of elevated Ca^{2+} also was provided by Cohen et al. (1974b). In our case, however, the powerlaw form for the PG shown in Figure 5 reaches down to ∼1.24 × 10^{−4} Hz, indicating that excess correlation extends to at least 8164 sec (the reciprocal of 1.24 × 10^{−4}), which is the full length of the data set. We see evidence of powerlaw behavior in the PG of data fromall of the preparations we have examined, including theXenopus neuromuscular junction (with and without added KCl).
Alteration of the exocytosis pattern induced by depolarization with KCl
The addition of KCl to the bath solution depolarizes the nerve terminal, thereby resulting in an elevation of the cytosolic Ca^{2+} concentration and a consequent increase in the rate of spontaneous vesicular exocytosis (Katz, 1962). Normalized IIHs for the Xenopus neuromuscular junction in the presence of added KCl are shown in Figure 6 a(dashed curve corresponds to 10 mm KCl;dotted curve corresponds to 20 mm). The IIH in the absence of added KCl, which is a normalized version of the solid curve in Figure 2, is presented for purposes of comparison (solid curve). Normalized IIHs were used to facilitate direct comparison of the curves. Each recording was obtained from a different neuromuscular junction because long data sets are required to obtain accurate statistics, and individual preparations typically do not remain viable long enough to permit more than a single recording to be obtained from a given preparation. The IIH curves all follow an exponential decay for long interevent intervals; however, the value near t = 0 is suppressed in the presence of KCl, and the coefficient of variation for the intervals is decreased.
The same data used to construct the IIHs in Figure 6 a were used to generate the AF curves in Figure 6 b. The three AFs show evidence of powerlaw behavior at long counting times, although the minimum value of the AF is reduced under stimulation. The increased rate of vesicular release serves to regularize the process and thereby leads to a reduction of the Allan variance over time scales where refractoriness is operative (Lowen and Teich, 1997). The presence of KCl also seems to reduce the strength of the longterm correlation present in the data (see also Cohen et al., 1974b), thus increasing the counting time at which powerlaw behavior becomes apparent.
Data selection: dilution of fractalrate behavior
Fatt and Katz (1952) were the first to investigate the statistical behavior of sequences of MEPCs, finding that the IIH was exponentially distributed. For a renewal process, this implies that the sequence can be described by a memoryless HPP. The IIHs measured in subsequent studies often have been variants of the exponential, which sometimes has fostered the (erroneous) notion that the event sequences are describable by processes akin to the Poisson, such as the gammarenewal (Hubbard and Jones, 1973) or deadtimemodified Poisson (VereJones, 1966) processes.
Fatt and Katz (1952) were very careful to note that the segment of data they selected for analysis (their Figs. 11–13) was sufficiently short (duration L = 176.8 sec comprising N = 800 MEPCs; E[t] ∼0.221 sec) so as to exclude, as they put it, the “occasional occurrences of short highrate bursts” of events, and to avoid “progressive changes of the mean,” present in their data. The observation of fractalrate behavior requires long data sets, and burstiness and apparent trends are at its very core, existing as natural components of exocytic behavior. We therefore would like to believe that the MEPCs observed by Fatt and Katz indeed did exhibit fractalrate fluctuations but that these researchers removed most traces of it by selecting relatively short segments of data for analysis and moreover by choosing precisely those segments that exhibited minimal fluctuations.
Indeed, an analysis by Cox and Lewis (1966, page 220) of even the special segment selected by Fatt and Katz reveals a departure from Poisson behavior that takes the form of a “relatively longterm effect.” Similar conclusions were reached by Cohen et al. (1973,1974a,b) and by Van der Kloot et al. (1975).
We proceed to explicitly demonstrate the consequences of selecting such segments of data with our own measurements. We choose the 20 mm KCl data set (16358 events; duration 4451 sec; E[t] ∼0.272 sec) because it has a large number of events and its rate is comparable with that of Fatt and Katz’s classic data set.
The PG for our full data set is displayed as the solid curve in Figure7. We now select a segment from the center of the full data set (N = 800 events, excising both the 7779 events preceding it and the 7779 events after it; L = 162 sec; E[t] ∼0.203 sec) that is comparable with the segment analyzed by Fatt and Katz in number of events, data duration, and mean interevent interval. The PG for this truncated MEPC segment is shown as the dashed curve in Figure 7. Because of its limited length, it is clear that the lowest frequency available to this PG isf _{min} = 1/L = 1/162 ∼6 × 10^{−3} Hz, the lefthand endpoint of the dashed curve. Although the agreement of the two curves is reasonable over the range where they coexist, it is plain that the fractal behavior in the full data set cannot be accessed in the truncated version. Similarly, the AF of the truncated MEPC data set cannot be estimated reliably forT > 16 sec, assuming that a minimum of 10 samples is required for statistical accuracy. The dotted curve in Figure6 b reveals that the AF begins to depart from simple renewal behavior only for counting times >16 sec. Moreover, selecting aparticular short segment of data on the basis of lack of variability (i.e., lack of burstiness or lack of progressive changes of the mean) serves to reduce further any manifestations of fractalrate behavior.
We conclude that fractalrate behavior, although it well may have been present in the original data set collected by Fatt and Katz, could not be discerned in the 176.8 sec segment that they analyzed.
Powerlaw behavior of the AF for rat hippocampal synapses,Xenopus myocytes, and rat fibroblasts
Although our attention thus far has been directed principally toward the Xenopus neuromuscular junction, we also have observed vesicular exocytosis consistent with fractalrate behavior from other neuronal and nonneuronal preparations, as mentioned earlier. Figure 8 displays the AFs for spontaneous vesicular release from two neuronal and two nonneuronal cells: theXenopus neuromuscular junction (solid curve; reproduced from Figs. 4, 6 b), the rat hippocampal synapse in the presence of TTX (longdashed curve),Xenopus–myocyte autoreception (shortdashed curve), and the exogenously loaded rat fibroblast brought into synaptic contact with a Xenopus myocyte (dotted curve). The curves presented in Figure 8 are representative of the data sets that were sufficiently long to merit analysis (N ≥ 400); these comprise four Xenopusneuromuscular junctions, four hippocampal synapses, oneXenopus myocyte, and two fibroblasts.
For large counting times, all of the eleven vesicularexocytosis data sets examined to date exhibit AFs that increase with the counting time, in a form consistent with the presence of fractal behavior. Fractal exponents estimated from the AF plots were in the range α_{A} = 0.1–2.7 (mean = 1.23), whereas those from the PG plots were in the range α_{S} = 0.2–4.0 (mean = 1.79). The fractal exponents calculated from the AF and PG were in general agreement (Teich et al., 1996a; Thurner et al., 1997); the overall correlation coefficient was +0.62, with substantially superior agreement for the longer data sets. Improved correlation would, no doubt, be obtained if the recordings were of yet longer duration so that asymptotic powerlaw behavior could be attained. Nevertheless, fractalrate activity seems to be ubiquitous in exocytic events.
Biophysical origin of the fractal behavior
A number of possible origins exist for the observed fractal behavior. One plausible scenario is that exocytosis is governed by fractal Ca^{2+}ion channel activity and that this activity ultimately derives from 1/ftype fluctuations of the membrane voltage. We proceed to provide a biophysical description of this process. The mathematical formulation, which leads to the FLNDP, is developed in Appendix .
It generally is accepted that voltagegated Ca^{2+}ionchannel openings are responsible for vesicular exocytosis (Zucker, 1993). For a fixed membrane voltageV near the resting potential, calcium flow is negligible. Occasionally, however, random thermally induced channel openings occur, which often lead to spontaneous exocytic events for nearby vesicles. Such spontaneous behavior is almost completely memoryless and is therefore well modeled by an HPP, with rate λ given by the Arrhenius equation (Berry et al., 1980) λ = 𝒜exp{−[E _{A} −zFV]/R𝒯}, as given in Equation B1. Here 𝒜 is a rate constant (often called the frequency factor),E _{A} is the constant activation energy associated with the ionchannel opening, z is the valence of the charge involved in the channel opening, F is the Faraday constant (coulombs/mol), R is the thermodynamic gas constant, and 𝒯 is the absolute temperature. Channelopening events that do not lead to exocytosis are accounted for in the values of 𝒜 and E _{A}. According to this picture, different fixed membrane voltages V lead to spontaneous exocytic patterns that differ only in their average rates; all are HPPs. These rates are exponential functions of the membrane voltage, as prescribed by the Arrhenius equation. Indeed, the lowvariability sections of data selected for analysis by Fatt and Katz (1952) likely would be associated with regions for which the rate is relatively constant so that the sequence of events could be well approximated by an HPP, as they found.
However, the membrane voltage is not fixed but, rather, varies randomly in time. Denoted V(t), it has a gaussian amplitude distribution and a 1/ftype spectrum (Verveen and Derksen, 1968). The rate λ(t) of the Poisson process, therefore, also varies in time, as described by Equation B2, which shows that the rate is the exponential transform of the voltage. Because the latter has a gaussian (normal) amplitude distribution, the rate is described by a lognormal amplitude distribution (Saleh, 1978) with a closely related spectrum. (A lognormal random variable is one for which the logarithm is normally distributed.) The rate process therefore is called fractal lognormal noise (FLN) and is described in Appendix . In short, the channel openings are characterized by a doubly stochastic Poisson process (Cox and Lewis, 1966) with a rate that is FLN: the FLNDP.
The net result is that the calciumflow events, and therefore the exocytic events, are described by the FLNDP, the properties of which are provided in Appendix . Unlike the HPP, the FLNDP has memory. Thus, the fluctuating membrane voltage imparts fractal correlations to the rate of exocytic events so that the observation of a short (long) interevent interval, for example, signifies a locally high (low) rate λ(t), which in turn indicates that the next interevent interval is also likely to be short (long). The FLNDP is clearly a nonrenewal process. However, the sequence of events generated by the FLNDP model does not form a proper fractal in time because it does not itself scale. Rather, the rate of event occurrences scales so that the FLNDP belongs to the family of fractalrate stochastic point processes.
Analytical predictions and computer simulations based on this model were compared with the exocyticevent data for a variety of statistical measures. We performed 100 simulations of the FLNDP, using parameter values obtained from the data set that is displayed in Figures 2, 3, 4, 5. IIH plots were computed for each simulation; these differed only in the seeds used for the random number generator. These 100 plots were averaged together to yield an aggregate IIH plot. This same process was used to construct AF and PG plots, using the same simulations (and therefore the same parameters and random seeds).
The resulting theoretical curves, denoted FLNDP, are presented as the dotted curves in Figures 2, 4, and 5. The theoretical results shown in all of these figures used a single set of parameters derived from the experimental data. Agreement with the data is excellent over all time scales. The slight deviations between the behavior of the FLNDP and the experimental data evidenced in the AF (Fig. 4) and in the PG (Fig. 5) would diminish no doubt were refractoriness included in the simulation. Indeed, it is known that refractoriness produces a dip in the AF for counting times in the vicinity of the refractory period (Lowen and Teich, 1997), the very region where the agreement is least satisfactory in Figure 4. Moreover, the AF and PG calculated using shuffled FLNDP simulations are in excellent accord with those of the shuffled data (longdashed curves) in Figures 4 and 5, respectively. We conclude that, aside from its physiological plausibility, the FLNDP provides an excellent mathematical model for characterizing sequences of MEPCs observed in our experiments.
It is, however, possible that the fractalrate behavior manifested in exocytic data derives from other mechanisms. Wide ranges of conformational states and time scales seem to be ubiquitous in large proteins (Liebovitch and Tóth, 1990) so that fractalrate exocytic behavior could originate from fractal behavior of the specialized proteins directly involved in vesicular exocytosis rather than via mediation by calcium. Or, fractal Ca^{2+}ionchannel openings could lead to average intracellular Ca^{2+}ion concentrations that behave in a fractal manner, in turn modulating global docking and transportprotein behavior (Zucker, 1993). Even in these cases, however, the FLNDP model would provide a useful mathematical description of the vesicular release process, albeit with a different biological interpretation.
DISCUSSION
As indicated above, activity consistent with fractalrate behavior is present in every in vitro spontaneous vesicularsecretion preparation of sufficient length (N≥ 400 events) that we have examined, neuronal and nonneuronal alike. Because vesicular exocytosis at the synapse shares many features and proteins in common with exocytosis and intracellular trafficking in all eukaryotic cells (Bennett and Scheller, 1993), it may be that such behavior is present in these systems as well.
Fractal and fractalrate behavior are also present in excitable tissue recordings for various biological systems in vivo, from the microscopic to the macroscopic (Bassingthwaighte et al., 1994; West and Deering, 1994). Examples include the openings and closings of ion channels (Läuger, 1988; Millhauser et al., 1988; Liebovitch and Tóth, 1990; Teich et al., 1991; Lowen and Teich, 1993a–c); patterns of actionpotential firings in the auditory system (Teich, 1989, 1992; Teich et al., 1990; Powers and Salvi, 1992; Kumar and Johnson, 1993; Kelly et al., 1996; Lowen and Teich, 1996, 1997), visual system (Turcott et al., 1995; Teich et al., 1996a,b, 1997), somatosensory cortex (Wise, 1981), and mesencephalic reticular formation (Grüneis et al., 1993); and even the sequence of human heartbeats (Kobayashi and Musha, 1982; Saul et al., 1988; Turcott and Teich, 1993, 1996). In almost all of these cases, the upper limit of the observed time over which fractal correlations exist is imposed by the duration of the recording. The appearance of fractalrate behavior at synapses, as well as in systems comprising collections of synapses, indicates that such behavior is ubiquitous in neural signaling.
The connection between fractalrate fluctuations and information encoding and transmission in neurons, if there is one, remains unclear. Fractal noise exhibits larger fluctuations at lower frequencies and thereby generally renders difficult the detection of the slowest, most gradual changes in a signal. Thus fractal exocytic activity could represent a fundamental source of noise ubiquitous in living cells, to which natural systems must adapt. However, many natural signals are themselves fractal (Voss and Clarke, 1978), and it may be that fractal activity in neurons provides some advantages in terms of matching the detection system to the expected signal (Teich, 1989, 1992).
Fractalrate activity also represents a form of memory, because the occurrence of an event at a particular time increases the likelihood of another event occurring later, with the strength of this memory persisting for some time. Fractalrate synaptic activity therefore provides a distributed network for memory and may provide a mechanism for potentiation.
Although it is difficult to ascribe definitively the observed longterm correlation in neurotransmitter release events to fractalrate behavior for a single data set of the limited size available in these experiments, it is gratifying that the FLNDP model, which relies on only three free parameters (intereventinterval mean and variance, and fractal exponent), provides good agreement with the exocytic data. That activity consistent with fractalrate behavior exists in all of the 11 data sets that we have examined indicates that the FLNDP is a useful model for describing our observations. Further confidence in the use of this model is engendered by the success of its close relative, the refractorinessmodified fractal binomialnoisedriven doubly stochastic Poisson process (FBNDP), for modeling actionpotential activity of primary afferent auditorynerve fibers in the cat (Lowen and Teich, 1995, 1996; Thurner et al., 1997). For voltage fluctuations small enough so that the exponential transform in Equation B2 can be approximated by a linear function, the FLNDP becomes nearly the same as the FBNDP. Moreover, the FLNDP and another related process, the fractal binomialnoisedriven doubly stochastic gamma process (FBNDG), have been used successfully to model retinalganglioncell and lateralgeniculatecell actionpotential activity in the visual system of the cat (Teich et al., 1997).
In conclusion, it is clear that traditional renewal models treating vesicular exocytosis as a memoryless stochastic process are wholly inadequate for representing many of its salient features. Rather, a new class of models that rely on fractalrate stochastic point processes is required.
DEFINITION AND PROPERTIES OF THE AF
The Allan factor (AF) is defined as the ratio of the Allan variance to twice the mean of the event count (Teich et al., 1996a). The Allan variance, in turn, is the average variation in the difference of adjacent counts (Allan, 1966). To compute the Allan variance at a specified counting time T, the data record of durationL first is divided into L/T contiguous counting windows, each of duration T. Much as in the procedure used to calculate the rate estimate, the number of eventsZ
_{k}(T) falling within the kth window is registered for all indicesk corresponding to windows lying entirely within the data record. The difference between the count numbers in a given window [i.e., Z
_{k}(T)] and the one after it [Z
_{k+1}(T)] is then computed for all k. The mean square of this quantity, E{[Z
_{k+1}(T) −Z
_{k}(T)]^{2}} , is the Allan variance. Dividing the Allan variance by twice the mean yields the AF:
An increase in the AF near a specific time scale occurs if event clusters of that particular scale are present in the data; the AF then will reach a plateau beyond the largest time scale present. Such behavior would be manifested, for example, by a Bartlett–Lewis cluster process such as that used by Cohen et al. (1974b) to model the MEPC sequence in the frog neuromuscular junction. A fractalrate stochastic point process, in contrast, generates a hierarchy of clusters of different durations, which leads to an AF plot that continues to rise as each cluster time scale is incorporated in turn. The net result is an AF that rises in powerlaw manner with increasing counting time T (straight line on a doubly logarithmic plot). For such processes the AF begins to rise above its asymptotic value of unity at a counting time that depends on the relative strength of the fractal component of its rate.
It is important to note that the random fluctuations inherent in any finite data set lead to AF plots that exhibit variability about the values predicted for exactly defined (nonrandom) point processes (Thurner et al., 1997). For data sets of sizes comparable with those used in this paper, such fluctuations can prove significant for any single plot, and conclusive proof of fractal behavior is not always possible. However, given a number of data sets of this size, uncertainty is greatly reduced. A similar argument applies for the PG.
The AF is preferred to the count index of dispersion (Fano factor), a similar measure constructed from the ordinary variance–time curve, for the analysis of fractalrate stochastic point processes because of its greater generality and freedom from bias (Lowen and Teich, 1996; Teich et al., 1996a; Thurner et al., 1997). A particular advantage of the use of the AF (or, equivalently, the Allan variance–time curve, which contains the same information) over the ordinary variance–time curve (Cox and Smith, 1953; Cohen et al., 1974b) lies in its insensitivity to linear trends, a result of the fact that it relies on a firstorder difference. This is a very important feature, which stems from the close relation of the AF to wavelet theory and, in particular, to the Haar wavelet. Generalizations of the AF, based on other wavelets, are insensitive to higher order trends (Teich et al., 1996a).
MATHEMATICAL FORMULATION OF THE FLNDP
For a membrane voltage near the resting potential, calciumion flux into the cell generally remains well below the value needed to trigger an exocytic event (Zucker, 1993), because voltagegated Ca^{2+}ion channels tend to remain closed under these conditions. Random thermal fluctuations, however, occasionally cause a few nearby channels to open in the same time frame, although the membrane voltage greatly favors the closed state. It turns out that only a few such channels need open to provide sufficient calcium to trigger vesicular exocytosis (Zucker, 1993). Were the membrane voltage fixed, this spontaneous behavior would be essentially memoryless; knowledge of all previous event occurrence times would yield no additional information about the future beyond that given by the average rate λ of the events. The mathematical model for discrete events in this case is the HPP, which is specified by a single parameter: the expected rate λ of nearby, nearly simultaneous Ca^{2+}ionchannel openings and, therefore, of the associated exocytic events.
Transitionstate theory (Berry et al., 1980) describes the dependence of this expected rate λ on various parameters of the ion channels (Hille, 1992) and predicts that it is given by the Arrhenius equation:
However, the resting voltage V of an excitabletissue membrane is not fixed but, rather, exhibits fractal (1/ftype) fluctuations with a gaussian amplitude distribution (Verveen and Derksen, 1968). We therefore replaceV by V(t) in Equation B1 to accommodate these voltage fluctuations:
Thus, three parameters completely describe the membrane voltage: μ_{V} ≡ E[V], c, and α. (The membrane voltage variance, ς_{V} ^{2} ≡ Var[V], is expressible in terms of these three parameters.) We could, alternatively, have used a formulation in terms of the autocorrelation function R _{V}(τ) ≡ E[V(t)V(t+τ)], which is the Fourier transform of the power spectral density.
Because the voltage is gaussian (normal), and the rate is the exponential transform of the voltage in accordance with Equation B2, the rate has a lognormal amplitude distribution [the exponential of a gaussian is defined to be lognormal (Saleh, 1978)]. The random process describing the rate λ(t) therefore is given the appellation fractal lognormal noise (FLN); its amplitude distribution is lognormal, and its spectrum is derived from 1/ftype noise. The mathematical model for the discrete events in this case is therefore a Poisson process with a rate described by FLN, i.e., the FLNDP.
It remains for us to determine the relevant statistics of λ(t) for the FLN described by Equation B2. These, in turn, will be used to determine the statistics of the exocytic events themselves, which are described by the FLNDP. We proceed by defining an auxiliary process X(t), which is a normalized version of the membrane voltage V(t):
We now turn to the autocorrelation function of the rate:
Finally, we consider the process of the exocytic events themselves, which we denote N(t), as in Appendix . These events are described in terms of a Poisson process driven by a fractal lognormal rate function λ(t). The sequence of the coincident calciumflow events, and therefore of the vesicular release events themselves, is then described by a FLNDP, as promised. The statistics of the FLNDP process are computed readily from those of the FLN rate. If the results for a general Poisson process (Saleh, 1978;Lowen, 1996) are used, the power spectral densityS
_{N}(f) of the events becomes:
If we assume further that the rate λ(t) [or equivalently the voltage V(t)] exhibits fluctuations that are slow in comparison with the average rate of channel openings E[λ], closedform expressions for the moments of the times tbetween channel openings also can be obtained, again with the help of general Poissonprocess theory (Saleh, 1978):
Footnotes

This work was supported by Grants from the Whitaker Foundation to S.B.L., from the National Institutes of Health (NS31923) to Mm.P., and from the Office of Naval Research (N0001492J1251) to M.C.T. We thank Conor Heneghan and Eric Schwartz for helpful suggestions.
Correspondence should be addressed to Professor Malvin C. Teich, Department of Electrical and Computer Engineering, Boston University, Boston, MA 02215.