Abstract
Spreading depolarization (SD) is an important phenomenon in stroke and migraine. However, the processes underlying the propagation of SD are still poorly understood, and an elementary model that is both physiological and quantitative is lacking. We show that, during the onset and propagation of SD, the concentration time courses of excitatory substances such as potassium and glutamate can be described with a reaction–diffusion equation. This equation contains four physiological parameters: (1) a concentration threshold for excitation; (2) a release rate; (3) a removal rate; and (4) an effective diffusion constant. Solving this equation yields expressions for the propagation velocity, concentration time courses, and the minimum stimulus that can trigger SD. This framework allows for analyzing experimental results in terms of these four parameters. The derived time courses are validated with measurements of potassium in rat brain tissue.
Introduction
A spreading depression, or depolarization (SD), is a slow wave that progresses through neural tissue, depolarizing neurons and depressing EEG activity, as was discovered by Leao (1947). Spreading depressions are encountered in various pathologies in which metabolism and/or ionic homeostasis are compromised. SDs can be observed after stroke as peri-infarct depolarizations, enlarging the damaged tissue area (Nakamura et al., 2010). SDs can occur during global hypoxia and are the substrate of a migraineous aura as well (Dreier, 2011).
Typically, a slowly progressing (2–10 mm/min) wave is reported, in which neurons depolarize and ion concentrations in both the intracellular and extracellular space change drastically. After tens of seconds, the ion concentrations are gradually restored, neurons repolarize, and function returns. However, under unfavorable circumstances, no recovery of the ion concentrations takes place, hence the neurons no longer function and eventually die.
Experimentally, SD can be induced by any stimulus that excites or depolarizes the membrane potential. Notable examples are the application of glutamate or potassium, electrical stimulation, and hypoxia. Once triggered, an SD wave propagates independent of the amplitude or kind of stimulus (for review, see Somjen, 2001; Dreier, 2011; Grafstein, 2011; Lauritzen et al., 2011).
Even after 50 years of research, it is still not clear what the exact mechanisms of propagation of SD are. It is also not clear whether different forms of SD, for example, hypoxic versus normoxic, have different propagation mechanisms. We hope to contribute to solving the puzzle by providing analytical expressions that describe the dynamics of diffusible substances, such as potassium and glutamate, that lead to propagation and initiation of SD. To the best of our knowledge, such expressions have not been derived using physiological parameters.
Traditionally, extracellular diffusion of glutamate or potassium was hypothesized as the propagation mechanism for SD. However, various experiments suggest alternative mechanisms. In normoxic SD, for example, extracellular potassium does not increase until after the extracellular voltage changes (Herreras and Somjen, 1993). Therefore, diffusion through the opening of gap junctions is now thought to play a major role in SD propagation (Largo et al., 1997; Somjen, 2001; Herreras, 2005).
We show that, regardless of which substance or pathway of diffusion is most prominent, the contribution of diffusible substances to the propagation of SD can be described with the same general expressions.
We will discuss how the process of release and diffusion of an excitatory substance, such as glutamate or potassium, during the onset of SD can be described. This allows the construction of a reaction–diffusion equation that can solved analytically. This results in expressions for both the propagation speed and the stimulus strength needed to elicit SD, which indicates how susceptible tissue is to SD. The solutions that correspond to an SD wave are compared with measurements from literature and are shown to describe these very well. Furthermore, we give examples of how the results of experimental measurements and numerical simulations can be analyzed with these expressions.
Materials and Methods
In essence, the onset of an SD wave is a reaction–diffusion process that can be described by a reaction–diffusion equation. These equations express the changes in the concentration of a particular substance attributable to its generation, removal, and diffusion (Ghergu, 2012). The term diffusion is used here in the broadest sense. It can be either passive diffusion through Brownian motion or any other form of movement or transport through the tissue, as long as the substance is moved along its concentration gradient.
Reaction–diffusion models for SD with varying amounts of physiological detail have been proposed (Tuckwell and Miura, 1978; Tuckwell, 1981; Shapiro, 2001; Yao et al., 2011). However, the complexity of these models does not allow for analytical solutions. In a personal communication to Grafstein (1963), Hodgkin and Huxley already solved a reaction–diffusion equation for SD, calculating the concentration of extracellular potassium, with a cubic polynomial for the reaction rate. This approach produces an analytical solution for the propagation velocity of SD. We will use a similar approach here but use parameters and a reaction rate that are more closely related to the underlying physiology.
The case in which a cubic polynomial is used as reaction rate in SD has been discussed previously (Burěs et al., 1974; Wilson, 1999) and is more generally known as the Schlögl model (Schlögl, 1972). SD is mathematically analogous to other phenomena, such as a flame on a matchstick (Starmer, 2007). Therefore, the solutions that we will obtain for our choice of the reaction rate are a specific case of the general results for initiation (Idris and Biktashev, 2008) and propagation (Zeldovich and Frank-Kamenetskii, 1938) of waves in an excitable medium (Starmer, 2007).
Concentration dynamics of an excitatory substance.
We will discuss how the concentration of a single excitatory substance can be described, starting with the processes involved in the reaction. How multiple substances can interact during SD will be discussed briefly later (see Discussion, Multiple diffusing substances).
Because in experiments the application of extracellular potassium is often used to elicit SD and several concentration time courses are available from experiments in literature as well, we will use this as an example throughout. However, note that our approach is general and is not limited to this particular choice.
In physiological conditions, the concentration C of extracellular potassium is ∼3 mm. This is the resting concentration C0, at which there is an equilibrium between the efflux of potassium into the extracellular space from the neurons and glia and the influx of potassium through Na/K-pump activity. Physiological homeostasis mechanisms, such as increase of Na/K-pump activity and diffusion to the bloodstream, work to return the concentration to C0 after changes from its equilibrium value, for example, resulting from increased neural activity.
For simplicity, we assume that the concentration relaxes exponentially to this equilibrium, so that the rate Rr(C) at which the concentration decreases is proportional to the concentration difference from equilibrium:
where G is the constant of proportionality. Because we will mainly deal with excess concentrations, we will refer to G as the removal constant.
However, for pathologically high concentrations, another effect takes place as a result of the excitational effect of extracellular potassium. As an illustration, the time-averaged potassium release from a single, unstimulated, neuron is calculated with a model adapted from Zandt et al. (2011). Figure 1 shows this potassium efflux as a function of extracellular potassium. This release rate is negligible for low potassium concentrations. However, when [K+]e exceeds a critical concentration, the cell spontaneously generates action potentials (i.e., without synaptic input), which leads to a large efflux of more potassium (Kager et al., 2002; Zandt et al., 2011). Movement of this excess of potassium to neighboring tissue along its concentration gradient can generate a SD.
The average potassium efflux as a function of [K+]e from a single Hodgkin–Huxley-type (HH) neuron without input. Ion concentrations were kept constant at resting levels in the simulations, and only the extracellular potassium concentration was changed. The time-averaged leak and gated membrane currents were then summed. A critical concentration Ct can be observed at ∼20 mm, after which the neuron starts spiking and releases large amounts of potassium. The dotted line shows a step function approximation.
From the simulation shown in Figure 1, we further observe that the release rate of potassium, R(C), can be approximated with a step function of the extracellular potassium concentration C:
Here H denotes the Heaviside step function (H(x) = 1 when x > 0 and 0 otherwise), Ct denotes the concentration threshold above which neurons expel K+ and R0 (in millimolar per seconds) is the rate at which C increases as a result of this expulsion, which we will refer to as release rate. We will show below (see Using a sigmoid instead of a step function) that approximating the release with a simple step function is not an oversimplification but yields almost the same results as a sigmoidal function. This is important because neural tissue contains inhomogeneous populations of neurons, for which a sigmoidal reaction rate is a more realistic approximation than a step function.
Summing the expressions for the homeostasis process and the release process yields the total reaction rate:
The derivation of this reaction rate is generally valid. A pathological high concentration of an excitatory substance, when exceeding a certain threshold, can activate mechanisms that release more of this substance. As an example, Figure 2 shows the total reaction rate for extracellular potassium. Parameters and their typical values are presented in Table 1. We remark that the function has a shape similar to the polynomials used in other studies on SD (Grafstein, 1963; Reggia and Montgomery, 1994), with zeros at resting concentration, at a threshold Ct and at some large value above threshold. We can use Equation 3 to calculate waves resulting from the “generation” and diffusion of an excitatory substance.
The reaction rate, Rt = R(C) − G × (C − C0). Parameters are given in Table 1. The black dot denotes the physiological resting state. Far above threshold, the derived expressions are expected to be inaccurate, depicted by a dashed curve.
Parameters and their typical values
Note that Equation 3 can describe the depolarization of the neuronal membrane voltage too, because this is caused by an increase in concentration of excitatory substance, but it does not describe the recovery. However, the recovery process does not need to be modeled, because its contribution to the reaction rate is significant only at concentrations far above threshold. It can be shown that its influence on the initiation and the propagation of SD can hence be neglected, except for relatively small propagation speeds, which we will discuss below (see Discussion, Multiple diffusing substances).
In short, we argue that SD can be locally mediated through release of an excitatory substance, resulting from exceeding a concentration threshold. Positive feedback leads to net release of more excitatory substance, resulting in depolarization. Equation 3 describes the reaction rate of a single excitatory substance using the rate of concentration increase attributable to release R0, a concentration threshold Ct, and a removal constant G. Any effects from other substances or processes can be expressed in terms of their effect on R0, Ct, and G. As an example, extracellular potassium was used, but the derived expression is valid for any excitatory substance. In the next section, we will use the expression for the reaction rate to construct a reaction–diffusion equation that describes the propagation of SD.
A reaction–diffusion equation for SD.
To calculate the propagation velocity and the onset of an SD wave, a reaction–diffusion equation is constructed that describes the release and subsequent movement of excitatory substance through the tissue. Again, extracellular potassium is used as an example.
In general, a reaction–diffusion equation is written as follows:
C denotes the concentration of a diffusible substance, Rt its net production rate, and k its effective diffusion constant. The geometry of the extracellular space slows diffusion (Syková and Nicholson, 2008), whereas transporters may enhance movement. Furthermore, ions do not simply follow Fick's law of diffusion, because they are charged and therefore their diffusion depends on the co- and contra-diffusing ions, as excellently described by Shapiro (2001). Therefore, the effective diffusion constant needs to be determined experimentally.
Because the cortex can be considered as a two-dimensional sheet and the wave front is usually wide compared with the length of the wave, the propagation can be described in one dimension. When considering the initiation of spreading depression from a single point in space, two dimensions are needed and it is convenient to use polar coordinates.
In one dimension, Equation 4 reads
In the previous section, the total reaction rate was derived, using the expulsion of K+ by the neurons, R(C), and its removal from the extracellular space by neurons, glia cells, and diffusion to the blood, G × (C − C0). Substituting the reaction rate from Equation 3 into Equation 5 results in our main equation:
The equation for the two-dimensional case in polar coordinates is similar. Solutions of Equation 6 now describe the onset and propagation of SD.
Calculating a traveling wave solution.
Equation 6 is a non-ordinary partial differential equation that is typically difficult to solve analytically. However, exact expressions can be found for traveling wave solutions, the solutions that describe an SD. Numerical solutions show that these solutions exist, that they maintain their shape while traveling, and are stable against small perturbations. For more elaborate, mathematically rigorous, analyses of equations of this type for SD, we refer the reader to the work of Chapuisat and Joly (2011) and Dahlem et al. (2010).
Because the solution of Equation 6 is a traveling wave, it can be described as C(z) = C(x − vt), where v is the velocity, with the wave traveling in the positive x direction, assumed to be on the right. The solution can be translated arbitrarily, and we conveniently put
Note that
The solution can be split in two parts: (1) Cr(z), the part below threshold on the right; and (2) Cl(z), the part above threshold on the left, described by the following:
respectively. These are solved with the following constraints: continuity of the solution and its derivative,
a resting concentration far away from the wave on the right,
and the wave keeping a constant shape,
The constancy of the wave shape can be used to substitute 
When G = 0, it is easy to see that solutions to these equations are
with
Continuity of the derivative at z = 0 can then be used to determine v:
from which it follows that
When G ≠ 0, the physically relevant solutions to Equations 13 and 14 can be calculated by using an exponential Ansatz and solving the characteristic polynomial:
with
Again, continuity of the derivative is used to determine v:
with 
The last equality can be solved for λ− by substituting λ+ from Equation 28, resulting in
After elementary algebra, the velocity v is obtained:
The last approximation underestimates the solution no more than
Scaling the equation.
For numerical calculations and dimension analysis, it is convenient to scale C, x, and t. This yields an equation that depends on only one parameter, Ĝ:
In polar coordinates, the equation is scaled similarly. The scaled quantities, denoted with a hat, are listed in Table 2.
Quantities and their scaling
Calculating the critical stimulus strength.
In various pathological conditions, a large release of potassium or glutamate into extracellular space can elicit a spreading depression. Such a release can, for example, be caused by excessive activity of neurons or cell membrane damage. Again, we use extracellular potassium as an example. We assess the susceptibility of the modeled tissue to spreading depression by calculating whether or not a local release of potassium elicits an SD. This depends on the total amount of potassium released and on the duration of this stimulus. This is similar to determining the excitability of a neuron with an electrical stimulus. In that case, the current amplitude versus minimum pulse duration can be determined, which is characterized by a rheobase, a minimum charge, and chronaxy. The same framework can be applied to SD (Idris and Biktashev, 2008). We will limit our considerations to either a very short pulse stimulus or continuous stimulation by potassium injection, similar to determining the minimum charge and the rheobase of an electrical stimulus. We did not investigate the effect of the (spatial) width of the stimulus.
In the first scenario, a short release or injection of a certain amount of potassium is modeled with a delta pulse at t = 0 at the origin. The stimulus results in a very high concentration of substance around the origin, which diffuses away into a larger area. Around the origin, the concentration is above the threshold concentration and therefore potassium is released by the neurons. An SD is generated only when more extracellular potassium is generated than diffuses away. Provided Ĝ < ½, this is the case when the amount of injected potassium is large enough.
In the second scenario, a stimulus of infinite duration is modeled as a point source of potassium at the origin. The steady flux of potassium leads to a buildup of concentration around the source. If the stimulus is weak, an equilibrium is reached between the generation of potassium and its diffusion and removal. If the source is stronger, such an equilibrium does not exist and a traveling wave of high potassium is induced.
Calculating the initiation of SD in one dimension underestimates the effectiveness of diffusion and therefore polar coordinates need to be used. We use the scaled equation (Eq. 32), in polar coordinates. This allows us to calculate numerical results for different Ĝ and scale the results depending on k, ΔC, and R0, rather than simulating the full four-parameter space. The equation reads
Here,
is the amount of excitatory substance injected at t = 0 at the origin (in millimolar per square meters) and
is the flux of the point source of substance at the origin (in millimolar per square meters per seconds). [To avoid confusion, we point out that mm denotes millimolar = mmol/L (= mol/m3) and hence the use of “mm/L,” as sometimes encountered in literature on SD, is incorrect.] The scaling of these parameters can be inferred from their units. To obtain the amount of substance that would need to be injected in an experiment, this needs to be multiplied with the thickness of the cortex in which the concentration is increased. We further remark that, for a stimulus to be considered short, its duration should be much shorter than both ΔC/R0 and 1/G, typically ≪1 s. Furthermore, for a stimulus to be considered a point source, it should have a width smaller than k/v, typically <30 μm. This is not the case for most experimental stimuli. However, wider stimuli are expected to behave qualitatively similar to point sources, with the diffusion constant progressively losing influence with increasing stimulus width. Additional calculations are needed to investigate the effect of stimuli with finite areas, which was not further explored in this study.
No analytical solution was found for Equation 33, for neither P ≠ 0 nor S ≠ 0. Therefore, a numerical solution was obtained by using a simple finite elements calculation.
The elements were rings with inner and outer radii of rn and rn + Δr, with rn = n × Δr and n ranging from 0 to L/Δr, where L is the radius of the simulation area. The areas An and outer circumferences On were calculated, and the flux between neighboring elements was calculated as (C[n] − C[n + 1]))On/Δx. The outermost element was kept at C = C0, to create an absorbing boundary condition, although the simulation area was large enough to not let the solution reach the boundary. The changes in concentration were then calculated with an explicit time-stepping scheme.
Δx was decreased until no significant changes in simulation results were observed. The time step was chosen as 0.05/Δx2 to ensure convergence. Typical values used were L = 30, Δx = 0.1, and simulation duration T = 10. For various values of Ĝ, the smallest values for both S and P that were able to induce SD were then obtained empirically.
Using a sigmoid instead of a step function.
For tissue containing many neurons with inhomogeneous properties and varying inputs, a sigmoid is a better approximation to the release rate of potassium than a step function. To test whether this yields different results, the numerical calculations are done with a sigmoid instead of a step function. A sigmoid as wide as possible without having a significant slope at the origin was used. Figure 3 (left) shows this sigmoid together with the original step function. Figure 3 (right) compares the two waveforms. Even for this width of the sigmoid, the differences in shape and velocity are very small compared with the usually encountered variations in experimental measurements. Furthermore, numerical solutions show that the critical steady stimulus fS is smaller when R(C) is a sigmoid instead of a step function and the critical short stimulus fP is slightly larger, but the behavior of the functions is similar.
Comparison of SD waveforms for different expulsion functions. A sigmoid, R(Ĉ) = 1/(1 + exp(−(Ĉ − 1)/0.15)) − 1/(1 + exp(1/0.15)), is compared with the step function used throughout (left graph). The right graph shows how similar the two waveforms are. The propagation velocity when using a sigmoid is only 2% faster. The waveforms were calculated numerically from Equation 32, with Ĝ = 0.07.
From this, we conclude that a step function yields practically the same results as a sigmoidal function when calculating the propagation and initiation of SD.
Comparison with experimental measurements from literature.
As an initial validation of the model, its solution will be compared with the time courses of the extracellular potassium concentration measured during experimentally induced SDs. Two concentration time courses were obtained from literature. These measurements were chosen for their good time resolution, allowing us to obtain several data points during the onset of SD. The graphs were digitized manually from the figures using MATLAB (MathWorks).
A first set of data was taken from Adámek and Vyskočil (2011, their Fig. 1A), who elicited SD in rats of both sexes in the frontoparietal cortex. A filter paper soaked with 1 m KCl, applied from an opening in the skull, induced the SD. The potassium waveform was measured 5 mm away from the stimulus location. A second experimental result was obtained from Herreras and Somjen (1993, their Fig. 2), who elicited an SD in a female rat in the dorsal hippocampus by pressure injection of potassium. A propagation speed was obtained in both studies by measuring the arrival times of the SD wave at two electrodes. The values at the center of the observed velocity ranges were used. These were 60 and 100 μm/s, respectively.
For each dataset, the interval was selected in which the concentration first increases exponentially and then rises linearly. The solution to Equation 6 was fitted to the data points in this interval by a least-squares method with Ĝ, C0, ΔC, and R0, independent of k, which was set to yield the correct propagation speed.
Results
In the previous section, we derived a reaction–diffusion equation and calculated the concentration time courses and propagation velocity of SD as caused by an excitatory diffusible substance, such as potassium or glutamate. Furthermore, we assessed tissue susceptibility to SD by calculating the minimum amplitude of a stimulus that triggers SD.
Concentration time courses and propagation speed
We calculated a traveling wave solution to Equation 6, which describes an SD. Two examples of this solution with numerical values for extracellular potassium are shown in Figure 4. The concentration time course during SD is characterized by an exponential increase that lasts until the concentration reaches the threshold Ct, followed by an almost linear rise.
Traveling wave solutions of Equation 6. Shown are results for both no removal, G = 0 (dashed line), and a fairly large removal rate of G = 0.2/s (solid line). Parameters are given in Table 1. The propagation velocity of the waveforms was calculated as 35 and 24 μm/s, respectively.
When glial buffering is negligible, i.e., G = 0, the corresponding speed of the wave is
This equation shows that, when the removal rate of potassium is negligible, for example, after complete ischemia, the propagation speed increases as the square root of the effective diffusion constant k [as holds for any chemical wave (Showalter and JJ, 1987)]. The speed also increases as the rate of concentration rise R0 as a result of expulsion of potassium increases, for example, when the extracellular space is smaller. The speed of the SD wave also increases when the difference between resting and threshold concentration ΔC is decreased, for example, by hyperkalemia or a concurrent rise of glutamate, which lowers the threshold concentration.
Considering the scenario in which there is removal of extracellular K+, for example, by glia cells, bloodstream, and Na/K pumps, the conditions that result in the generation in an SD are slightly different. We find that a larger stimulus is needed to induce SD and SD will propagate more slowly. A useful quantity to calculate this is the normalized removal constant 
Note that this implies that, for fast normalized removal rates Ĝ > ½, the depolarization cannot propagate at all. Intuitively, one may suspect that the depolarization cannot propagate when diffusion is very slow either. This case, in which recovery cannot be neglected, is discussed below (see Discussion, Effects of recovery and limited total amount of released substance).
Critical stimulus for triggering spreading depression
To obtain a measure for the susceptibility of tissue to SD, we calculated both the minimum amount of excitatory substance in a short pulse P and the minimum flux delivered by continuous stimulation S that is sufficient to elicit SD. The critical amplitudes P and S are functions of Ĝ that scale with k, ΔC, and R0 and are given by
and
respectively. Here fp(Ĝ) and fS(Ĝ), shown in Figure 5, are functions that depend only on the normalized removal rate and were determined numerically. It can be observed for both stimulus types that SD is more difficult to induce in tissue that has a higher diffusion constant, larger concentration threshold ΔC, a smaller release rate R0, and a faster removal constant G.
The strengths of the critical stimuli, calculated as a function of Ĝ as described in Materials and Methods (Calculating the critical stimulus strength). The line drawn serves as a guide for the eye.
Comparison with experimental time courses
Figure 6 compares a modeled time course with the measurement obtained from Adámek and Vyskočil (2011). The obtained parameter values are C0 = 3.1 mm, R0 = 11 mm/s, ΔC = 9.4 mm, k = 3.4 × 10−9 m2/s, and G = 0.02/s. These values are in the physiological range, and the model correctly describes the shape of the measured time trace: an exponential rise followed by a steep almost linear rise.
Comparison of the solution of the model with experimental data, taken from Adámek and Vyskočil, 2011 (their Fig. 1A). The calculated and observed velocities of the waves are indicated in the legend. Dots are points digitally read out from the original graph.
Figure 7 compares a modeled time course with a second measurement (Herreras and Somjen, 1993). Here, we find C0 = 4.4 mm, R0 = 48 mm/s, ΔC = 9.0 mm, k = 1.9 × 10−9 m2/s, and G = 0.00/s. These values are also in the physiological range.
Comparison of the solution of the model with data taken from Herreras and Somjen, 1993 (their Fig. 2). The calculated and observed velocities of the waves are indicated in the legend. Dots are points digitally read out from the original graph.
Note that only the velocity ranges for the individual measured time courses are known with a typical spread of 20%, which leads to large errors, ±40%, in the determined diffusion constants.
Predictions of the model
Square root dependence of propagation speed on diffusion constant
Our model predicts a square root dependence of the propagation velocity on the effective diffusion constant k (Eq. 37).
This is in accordance with the findings of Shapiro (2001), who calculates with a detailed numerical model the propagation velocity of SD as a function of k for potassium. In his model, potassium diffuses mainly intercellularly through gap junctions rather than extracellularly. He finds that the propagation speed of SD increases as
SD cannot be induced when the concentration increase attributable to release is small
According to Equation 37, no SD can be elicited when G > R0/(2ΔC). Small rates of concentration rises, R0, therefore prevent SD. It is indeed observed that SD cannot be elicited in brains of young animals, in which the extracellular space is relatively large. Expelled potassium and neurotransmitters are diluted more and cause a smaller rise in extracellular concentrations, thereby protecting against SD (Somjen, 2004).
Tissue susceptibility to SD and its propagation speed are not necessarily correlated
Equations 37–39 predict that either decreasing R0 or increasing ΔC or G lowers the propagation speed and increases the critical stimulus amplitude for triggering SD. Decreasing the effective diffusion constant k also lowers the propagation speed but decreases the critical stimulus amplitude. Hence, the intuitive idea that SD is more difficult to induce in tissue in which it propagates more slowly is incorrect in cases in which the effective diffusion constant changes. This lack of correlation was indeed observed by Weimer and Hanke (2005). Their experiments in chick retina show that, during the recovery period that follows an SD, the threshold for electrical stimulation is increased, whereas cooling the tissue decreases this threshold. Both manipulations slow the propagation of SD.
Our model predicts from this measurement that cooling decreases the effective diffusion constant. It also quantifies this. Weimer and Hanke (2005) cooled the tissue from 30 to 20°C. This resulted in a reduction of propagation speed of 30%, in combination with a decreased threshold for electrical stimulation of 50%. These values both correspond with a decrease in diffusion constant of 50% (Eqs. 37, 38). In comparison, this cooling slows (passive) diffusion of potassium in water with only 20% (Bastug and Kuyucak, 2005). Hence, our model predicts that passive diffusion of extracellular potassium cannot have a primary role in the propagation of SD in these experiments: either potassium is actively transported, for example by the glia, or SD propagates through another mechanism, such as axonal conduction of action potentials or diffusion of glutamate.
We have to be careful in drawing conclusions from the electrical stimulation threshold alone because our model describes the stimulus amplitude in terms of substance release, which may not be proportional to the injected current and also may change during cooling.
Threshold concentration for SD is lower in female rat, whereas propagation speed is the same
Our model predicts that a lower concentration threshold increases the propagation speed (Eq. 37). This is in disagreement with observations from Adámek and Vyskočil (2011), who elicited SD by applying a cotton ball with 1 m KCl to the frontoparietal cortex. They measured the time course of the extracellular potassium concentration during the onset of SD close (1.5 mm) to the stimulus site. A measure for the concentration threshold was obtained by plotting this time course on a semilogarithmic scale, which shows a deflection point. This point was found at much lower concentrations in female rat cortex (5.0 ± 0.6 mm) than in male rat cortex (11.4 ± 0.4 mm) above the resting concentration (3.0 mm). The corresponding propagation speed for SD was determined to be 3.9 ± 0.3 mm/min in females and 3.6 ± 0.3 mm/min in males. The difference is not statistically significant.
Our simulations show that this deflection point in the semilogarithmic plot is not a reliable measure of the concentration threshold because it is influenced by the SD wave that is generated at the stimulus site. However, their measurements definitely suggest a concentration threshold that is lower in female rat but may be less pronounced than the deflection points suggest.
Because the lower concentration threshold for triggering does not have a large effect on the propagation speed, we predict that diffusing potassium has only a secondary role in the propagation of normoxic SD. The hypothesis that a lower concentration threshold for potassium leads to a higher susceptibility to SD is in agreement with the predictions of our model (Eqs. 38, 39).
Discussion
The main results of this work are the analytical expressions for the dynamics of diffusible excitatory substances during SD and how these contribute to the triggering and propagation. Expressions were derived for the propagation speed (Eqs. 36, 37) and the minimum amplitude of a stimulus that is sufficient to induce SD (Eqs. 38, 39). The theoretical concentration time courses during the onset of SD were compared with experimental recordings of extracellular potassium (Figs. 6, 7), and several predictions were made with the model by analyzing observations from literature.
Role of neural excitability
The most important assumption that we make to arrive at the expressions for the propagation and triggering (Eqs. 36–39) is that the release rate of a substance is a function of its own concentration that can be approximated with a step function or sigmoid. Although our results do not depend on the exact mechanism, we hypothesized that neuronal excitation mediates such a release rate of potassium and/or glutamate. However, it has been found that SD can still occur when TTX blocks sodium channels (Tobiasz and Nicholson, 1982; Müller and Somjen, 2000), showing that the generation of action potentials is not necessary for SD. In this case, depolarization of the neuronal membrane can still mediate SD. Also with TTX applied, a critical value for the neuronal resting membrane voltage is observed above which SD is triggered. This is a few millivolts higher than when TTX is not applied (Müller and Somjen, 2000). Hence, depolarization of the neuronal cell membrane by excitatory effects could mediate a pathological release of excitatory substances. This hypothesis can be tested by slowly increasing the glutamate or the potassium concentration in a large area and verifying that the threshold resting membrane voltages above which the depolarization is started are similar for both substances.
Validity of the derived expressions
The solution to Equation 6 was compared with two experimental time courses of extracellular potassium (Figs. 6, 7). For both measurements, the theoretical solution describes the measured data very well, until the concentration has risen to approximately two-thirds of the maximum concentration reached during the SD. Hereafter, recovery processes, which were not modeled, start to influence the concentration time course. In both cases, the obtained parameter values are in the physiological range.
Furthermore, the solution also deviates during the earliest phase of the SD wave from the dataset of Herreras and Somjen (1993) (Fig. 7). Rather than an exponential increase in concentration, the concentration initially rises approximately linearly. The start of this rise (t = −3.5 s) occurs simultaneously with the onset of prodromal spiking activity that Herreras and Somjen observed in extracellular recordings. This activity can cause the initially linear rise in [K+]e. Note that we only described the contribution of diffusing substances on the propagation of SD and not the possible contributions of other mechanisms, notably the propagation of neural activity through synaptic connections. This has been investigated, although not applied to SD, for constant ion concentrations (Amari, 1977), as well as numerically, for dynamic ion concentrations (Ullah et al., 2009).
Also, Herreras and Somjen observed a large change in extracellular potential, associated with neuronal depolarization, before the fast rise of [K+]e. Therefore, extracellular potassium seems to be a byproduct rather than the primary contribution to the propagation of SD in this experiment. Still, a critical point can be observed, after which a fast linear rise in potassium concentration occurs, preceded by an exponential increase, as caused by diffusion.
Although the model describes the experimental time courses of the concentration of extracellular potassium very well, this does not mean that diffusing extracellular potassium is the primary mechanism for SD propagation. For the release of a secondary substance, the same reaction–diffusion equation is valid. However, in that case, a different value for the concentration threshold is observed. We will discuss this in the next section.
Multiple diffusing substances
The release rate of one excitatory substance, such as potassium, can be influenced by the concentration of a second substance, such as glutamate. A rise in extracellular glutamate concentration lowers the concentration of potassium at which a certain membrane voltage is reached. If the release mechanism is mediated by the membrane voltage, this effectively lowers the potassium concentration threshold for potassium release. This results in a higher propagation speed (Eq. 37).
Two threshold concentrations can then be distinguished for a single substance: (1) the threshold for (experimentally) inducing SD, when no other substance is released yet; and (2) the observed threshold during propagation. The latter is lower because of the influence of a second excitatory substance. How close the threshold during propagation is to the threshold for inducing SD is a measure for the relative importance of a substance in the propagation of SD.
The model can be extended by adding an equation for a second substance to it and describing the concentration threshold of both substances as a function of concentration of the other.
Effects of recovery and limited total amount of released substance
In deriving Equation 6, it is assumed that cells in the tissue, once the concentration rises above the threshold, expel excitatory substance indefinitely. This assumption is justified in most cases, because only the front of the wave is important for an accurate description of propagation and initiation. Even for a diffusion constant close to zero, the model predicts that propagation of SD is possible, albeit slow. In that case, the concentration at a certain position in the tissue keeps rising, such that the concentration gradient eventually is steep enough to cause sufficient flux to the interstitial space of the neighboring tissue. This cannot occur in reality, because there is a natural limit on the concentration and the total amount of substance that can be released from the cells.
Therefore, when Equation 37 predicts relatively slow propagation velocities, these are either overestimated or propagation is not possible. The propagation velocity is considered slow, when the front of the waveform is affected by recovery effects.
To investigate this, recovery effects can be incorporated in the numerical model. This would also allow for investigating how consecutive SDs are induced by a prolonged stimulus. These effects can be incorporated, for example, by introducing a recovery variable in the model (Reggia and Montgomery, 1994). Such a variable slowly increases when the concentration is high and limits the release rate and enhances removal.
Conclusion
Concentration dynamics of diffusible excitatory substances, such as potassium and glutamate, were calculated during triggering and propagation of SD. For this purpose, a reaction–diffusion equation was constructed from four physiological parameters. Solutions to the equation describe the onset of SD and yield an analytical solution of the shape of the front of the wave, for example, concentration as a function of position and time, and the propagation speed (Eqs. 36, 37). A critical stimulus strength was calculated for both prolonged and short stimuli (Eqs. 38, 39).
The solution to our equations was validated with experimental measurements of extracellular potassium concentration and was found to describe the onset of the concentration time course very well. The estimates for the parameters are all in the physiological range. Predictions of the model were presented, and it was shown how the obtained expressions can be used to analyze several experimental findings. The model can be elaborated to describe multiple substances or recovery effects.
The presented expressions allow an analysis of the influence of an excitatory substance on SD propagation and triggering, in terms of concentration threshold, release rate, removal rate, and effective diffusion constant. The effects of other influences can be described by how they change these quantities. By comparing them with experiments or more elaborate numerical models, the derived expressions may be helpful in answering the question of which role diffusible substances play in various forms of SD.
Notes
Supplemental material for this article is available at https://www.researchgate.net/publication/234081422_simulation_and_fitting_of_SD_potassium_wave. MATLAB code was used for the simulations and the fitting. This material has not been peer reviewed.
Footnotes
- Correspondence should be addressed to Bas-Jan Zandt, Faculty of Science and Technology, University of Twente, Carré 3619, P.O. Box 217, 7500 AE Enschede, The Netherlands. b.zandt{at}utwente.nl
