Thermal block of action potentials is primarily due to voltage-dependent potassium currents: a modeling study

The squid giant axon varies in size, but we chose a standard length of 100 mm and diameter of 0.5 mm (500 μm) based on values reported in the literature (Rosenthal and Bezanilla 2000). We implemented this geometry in the NEURON simulation environment (Carnevale and Hines 2006). In general, the differential equations for the model (see below) were numerically integrated using a time step of 0.01 ms, and the voltages and ionic currents were measured in each segment of the model neuron. For studies in which compensatory current was injected to maintain the resting potential (see below), the time step was reduced considerably (to 0.2 μs) to ensure high accuracy.

To ensure stability of all state variables, the system was simulated for a total of 500 ms, and stimulating current was injected into the model after 250 ms. Depolarizing current (1 ms duration) was injected within the first segment of the first section of the model axons to initiate an action potential. For model axons of diameter greater than 10 μm, 2000 nA of depolarizing current was used. For model axons of diameter less than 10 μm, 100 nA of depolarizing current was used.

To explore how axon diameter affects its response to temperature, we varied axon diameters from 0.5 μm to 500 μm. The axon was divided into either two or three sections, each having different lengths and properties, i.e. temperature, resistivity, and ionic conductances (as described below). The number of segments within each section was chosen to allow for precise measurements of threshold block lengths over the entire length of the model axon. In the modified Hodgkin/Huxley model, the axial resistances varied with temperature. Thus, to construct regions of different temperature, we used the built-in NEURON function connect to ensure that each region had the appropriate axial resistance for its temperature.

Segmental currents over the j  different segments (where j  ranged from 1 to the total number of segments) were computed as

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{I}_{segment}}\left(\,j \right)=\frac{{{V}_{j-1}}-~{{V}_{j}}}{{{r}_{a,j-1,j}}}-\frac{{{V}_{j}}-{{V}_{j+1}}}{{{r}_{a,j,j+1}}},\nonumber \end{align} \tag{ 1 }$$

where V_j −1 corresponded to the voltage of the previous segment (if there was a previous segment), V_j +1 corresponded to the voltage of the next segment (if there was a next segment), r_a,j −1,j corresponded to the axial resistance between compartment j   −  1 and j , and r_{a,j ,j +1} corresponded to the axial resistance between compartment j  and j   +  1. Note that r_a  =  $\frac{{{R}_{a}}}{\pi {{d}^{2}}/4}*L$, where r_a is the axial resistance, R_a is the axial resistivity (in Ω cm), L is the length of the segment in μm, d is the axon diameter (in μm), and r_a's units of Ω cm μm μm⁻² yield 10⁴ Ω. Multiplying this by 10⁻² makes it possible to report the resistance in MΩ.

Temperature changes were applied to the model neurons in two ways. First, to determine the changes in currents that occur as an action potential encounters a change in temperature, we maintained one half of a model axon at the base temperature of 6.3 °C while varying the temperature of the second half of the axon. For these studies, the model axon consisted of two sections. By increasing the temperature of half of the axon (the second section), we ensured that transient currents induced by the initiation of an action potential, which occur at the beginning of the first section, were minimized. In the NEURON simulation environment, temperature is a global variable that is uniformly applied to every segment. We changed temperature to be a local variable so that we could access and modify the temperature of each segment independently. Ordinarily, NEURON uses the TABLE command in the .mod file to access a look-up table for all segments; by commenting out this command in the .mod file, the NEURON program evaluated the differential equations at every segment, rather than using the default parameters in the lookup table to decrease run time. This usually increased the simulation duration by a factor of about five.

To approximate a continuously changing temperature profile that might be induced by a thermal source such as a laser, we calculated a temperature profile, T(x), that interpolated smoothly between an initial temperature T₁ (6.3 °C) and a final temperature T₂ (25 °C) at all points x along the axon, such that the temperature began to change from its initial value T₁ at location a along the axon (10 mm before the midpoint, or 40 mm from the beginning of the axon) and reached temperature value T₂ at location b along the axon (10 mm after the midpoint, or 60 mm after the beginning of the axon). The temperature reached an intermediate temperature value of (T₂  +  T₁)/2 at the midpoint between points a and b, i.e. at (a  +  b)/2:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle T\left( x \right)={{T}_{1}}+\left( {{T}_{2}}-{{T}_{1}} \right)F\left( x,~a,~b \right),\nonumber \end{align} \tag{ 2 }$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle F\left( x,~a,~b \right)=\left\{\begin{array}{@{}cccccccccccccccccccc@{}} 0,~ & x~\leqslant a, \nonumber \\ 2{{\left( \frac{x-a}{b-a} \right)}^{2}}, & a~<x~\leqslant \frac{\left( a+b \right)}{2}, \nonumber \\ 1-~2{{\left( \frac{x-b}{b-a} \right)}^{2}}, & \frac{\left( a+b \right)}{2}~<x~\leqslant b \nonumber \\ 1, & x\geqslant b. \nonumber \\ \end{array} \right.\nonumber \end{align} \tag{ 3 }$$

Second, temperature change was applied along the central region of a model axon to determine the minimum length of elevated temperature that could block action potential conduction. The axon was modeled using three sections, and the length of the central section that was subjected to elevated temperatures was systematically varied until the action potential no longer propagated beyond the region of block. We classified an action potential as blocked if the potential at the end of the axon never exceeded  −60 mV, which is far below the threshold for initiating an action potential. Note that the resting potential of the model is  −65 mV, so that this is a fairly conservative criterion.

To capture the response of the neuron to changing temperatures, several features of the original Hodgkin/Huxley model became functions of temperature: (1) a temperature-varying Q₁₀ factor for the gating variables of the Hodgkin/Huxley model based on the experimental measurements of Rosenthal and Bezanilla (2000), (2) the peak potassium and sodium conductances varied with temperature, (3) a temperature-dependent sodium/potassium pump, and (4) the axial resistance as a function of temperature (a schematic comparing the original and modified Hodgkin/Huxley model is shown in figure 1(A)). After describing the detailed modifications that were made to the model, the entire set of equations is shown (below, model equations).

Zoom In Zoom Out Reset image size

The Q₁₀ factor.

In the standard Hodgkin and Huxley model (Hodgkin and Huxley 1952), temperature dependence is captured by altering the rate constants for the three gates controlling the voltage gated ion channels. The m gate (sodium activation), the h gate (sodium inactivation) and the n gate (potassium activation) are modified to account for changes in temperature. The corresponding modified rate constants are α_m and β_m for the m gate, α_h and β_h for the h gate, and α_n and β_n for the n gate. The values are changed by using the Q₁₀ relationship, which is

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{Q}_{10}}~={{\left( \frac{{{R}_{2}}}{{{R}_{1}}} \right)}^{\frac{10}{\left( {{T}_{2}}-{{T}_{1}} \right)}}}\nonumber \end{align} \tag{ 4 }$$

where temperature T₂ is greater than temperature T₁, rate R₂ is measured at T₂, and rate R₁ is measured at T₁. If both sides of the equation are raised to the power $\frac{{{T}_{2}}-{{T}_{1}}}{10}$, it is possible to solve directly for R₂ at temperature T₂ as

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{R}_{2}}={{R}_{1}}Q_{10}^{\frac{{{T}_{2}}-{{T}_{1}}}{10}}.\nonumber \end{align} \tag{ 5 }$$

Thus, for example, to modify the differential equation for the n gate (see below, Model Equations), since rate constants α_n and β_n are multiplied by the same factor to adjust them for temperature, the right-hand side of the entire differential equation for n is multiplied by the Q₁₀ term of equation (5). The same factor is used to modify the right-hand sides of the differential equations for the m and h gates. In the original Hodgkin and Huxley model, a fixed value of 3 was used for Q₁₀ for all temperature ranges, and the base temperature (T₁) was set at 6.3 °C (Hodgkin and Huxley 1952). Thus, the factor by which all the differential equations for gates m, n and h were multiplied simplified to ${{3}^{\frac{{{T}_{2}}-6.3}{10}}}$.

An earlier study of the change in the shape of the action potential (Hodgkin and Katz (1949), their table 2) showed that the Q₁₀ for the rates of fall (illustrated in their figure 1(B)), an approximate measure of the activity of the voltage-dependent potassium gates, changed in the temperature ranges 5 °C to 10 °C, 10 °C to 20 °C, and 20 °C to 30 °C. They also observed changes in the Q₁₀ for the rate of rise with temperature, though these changes were smaller. These data suggested that the Q₁₀ factors multiplying the differential equations for m, n and h were different and varied with temperature. We therefore used a systematic grid search to simultaneously find sets of Q₁₀ values for all three gates, varying them from a value of 1.0 (temperature independent) to 3.0 (the original Hodgkin/Huxley model value) in increments of 0.1. We excluded values that did not allow the model to generate an action potential. For those values that did generate an action potential, we determined the total amplitude of the action potential (measured from the resting potential to the peak of the action potential), by comparing these values to the amplitudes actually determined experimentally at that temperature (Hodgkin and Katz (1949), their table 1(c), p 243). Among those sets of Q₁₀ values producing total amplitudes within 5% of those measured experimentally, we then evaluated errors for maximum rates of rise and fall of the resulting action potentials (see below, equations (11)–(13)), choosing the Q₁₀ values that led to the smallest errors (supplementary figure 1(A) (stacks.iop.org/JNE/16/036020/mmedia)).

As a result of the systematic and exhaustive search, the multiplicative factor was found to be slightly different for each gate, and will be designated as ϕ_m(T), ϕ_h(T) and ϕ_n(T) for the m, h and n gates (respectively). Each of these are piecewise continuous functions:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{\phi }_{m}}\left( T \right)=\left\{\begin{array}{@{}lllllllccccccccccccc@{}} {{3}^{\frac{T-6.3}{10}}},{\hspace{104pt}} & 5\ {}^\circ {\rm C}~\leqslant T~\leqslant 10\ {}^\circ {\rm C}, \nonumber \\ {{3}^{\frac{10-6.3}{10}}}~{{3}^{\frac{T-~10}{10}}},{\hspace{71pt}} & 10\ {}^\circ {\rm C}<T\leqslant 15\ {}^\circ {\rm C}, \nonumber \\ {{3}^{\frac{10-6.3}{10}}}~{{3}^{\frac{15-10}{10}}}~{{2.8}^{\frac{T-15}{10}}},{\hspace{35pt}} & 15\ {}^\circ {\rm C}<T\leqslant 20\ {}^\circ {\rm C} \nonumber \\ {{3}^{\frac{10-6.3}{10}}}~{{3}^{\frac{15-10}{10}}}~{{2.8}^{\frac{20-15}{10}}}~{{2.7}^{\frac{T-20}{10}}}, & 20\ {}^\circ {\rm C}<T\leqslant 25\ {}^\circ {\rm C} \nonumber \\ \end{array} \right.\nonumber \end{align} \tag{ 6a }$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{\phi }_{h}}\left( T \right)=\left\{\begin{array}{@{}llllllcccccccccccccc@{}} {{3}^{\frac{T-6.3}{10}}}, & 5\ {}^\circ {\rm C}\leqslant T\leqslant 10\ {}^\circ {\rm C}, \nonumber \\ {{3}^{\frac{10-6.3}{10}}}~{{2.9}^{\frac{T-10}{10}}}, & 10\ {}^\circ {\rm C}<T\leqslant 15\ {}^\circ {\rm C}, \nonumber \\ {{3}^{\frac{10-6.3}{10}}}~{{2.9}^{\frac{15-10}{10}}}~{{3}^{\frac{T-15}{10}}}, & 15\ {}^\circ {\rm C}<T\leqslant 20\ {}^\circ {\rm C} \nonumber \\ {{3}^{\frac{10-6.3}{10}}}~{{2.9}^{\frac{15-10}{10}}}~{{3}^{\frac{20-15}{10}}}~{{3}^{\frac{T-20}{10}}},~ & 20\ {}^\circ {\rm C}<T\leqslant 25\ {}^\circ {\rm C} \nonumber \\ \end{array} \right.\nonumber \end{align} \tag{ 6b }$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{\phi }_{n}}\left( T \right)=\left\{\,\,\begin{array} {@{}llllllcccccccccccccc@{}} {{3}^{\frac{T-6.3}{10}}}, & 5\ {}^\circ {\rm C}\leqslant T \leqslant 10\ {}^\circ {\rm C}, \nonumber \\ {{3}^{\frac{10-6.3}{10}}}~{{2.8} ^{\frac{T-10}{10}}}, & 10\ {}^\circ {\rm C}<T\leqslant 15\ {}^\circ {\rm C}, \nonumber \\ {{3}^{\frac{10-6.3}{10}}}~{{2.8}^{\frac{15-10}{10}}}{{2.4}^{\frac{T-15}{10}}}, & 15\ {}^\circ {\rm C}<T\leqslant 20\ {}^\circ {\rm C} \nonumber \\ {{3}^{\frac{10-6.3} {10}}}~{{2.8}^{\frac{15-10}{10}}}{{2.4}^{\frac{20-15}{10}}}{{2.3}^{\frac{T-20}{10}}},~ & 20\ {}^\circ {\rm C} <T\leqslant 25\ {}^\circ {\rm C}{\rm .} \nonumber \\ \end{array} \right.~\nonumber \end{align} \tag{ 6c }$$

Sodium–potassium pump.

In the original Hodgkin and Huxley model, the sodium/potassium pump was not incorporated, because their initial studies suggested that the resting membrane potential was unchanged by temperature (Hodgkin et al 1952). In the studies of Rosenthal and Bezanilla (2000) in Loglio pealei, however, the resting potential varied with temperature and season, suggesting that components of the resting potential were temperature dependent. Studies in other molluscan species, such as Aplysia californica (Carpenter and Alving 1968) and Anisodoris nobilis (Gorman and Marmor 1970), demonstrated that a temperature-dependent sodium–potassium pump significantly contributed to changes in the resting potential with temperature. Rakowski et al (1989) characterized an electrogenic sodium–potassium pump in the squid giant axon that exchanged 3 sodium ions for 2 potassium ions. We used the equations of van Egeraat and Wikswo (1993) for the pump currents, based on the study of Rakowski et al (1989), and incorporated temperature-dependence based on the studies in Anisodoris (Gorman and Marmor 1970). We estimated the Q₁₀ for the pump conductance from its change in permeability ratio of P_Na/P_K from 0.028 at 4 °C to 0.0068 at 18 °C (Gorman and Marmor 1970), yielding an estimated Q₁₀ of 1.88 (using equation (1)), which was then used to scale the peak pump conductance, g_NaKinitial. Given that the pump exchanges three sodium ions for two potassium ions, the pump equations are:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{g}_{{\rm NaK}}}(T)={{g}_{{\rm NaK}initial}}~{{1.88}^{\frac{\left( T-6.3 \right)}{10}}},\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{I}_{{\rm Na}\_pump}}=3~{{g}_{{\rm NaK}}}\left( T \right)\left( {{V}_{m}}-{{E}_{pump}} \right),\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{I}_{{\rm K}\_pump}}=-2{{g}_{{\rm NaK}}}\left( T \right)\left( {{V}_{m}}-{{E}_{pump}} \right),\nonumber \end{align} \tag{ 7 }$$

where I_{Na_pump} and I_{K_pump} are the sodium and potassium currents due to the pump, respectively, V_m is the membrane potential, and E_pump is the reversal potential of the pump.

Peak sodium and potassium conductances.

In the original Hodgkin and Huxley model, the peak sodium and potassium conductances, ${{\bar{g}}_{{\rm Na}}}$ and ${{\bar{g}}_{{\rm K}}}$, were assumed to be constant and independent of temperature. As reviewed above, more recent work suggests that these values may vary with temperature.

To estimate the peak conductance values, rates of rise and rates of fall were extracted from figure 3 of Rosenthal and Bezanilla (2000) using the program DataThief (www.datathief.org). We extracted the center of the triangle indicating the mean data value as well as the values of the top and bottom error bars to provide a measure of the errors associated with each data point. In the model, the peak of the derivative of the rising phase was used as an estimate of the peak rate of rise, and the peak of the derivative of the falling phase was used as an estimate of the peak rate of fall (supplementary figure 1(B)).

We used a systematic grid search to determine the peak sodium and potassium conductance values that minimized the difference between the rates of rise and fall observed experimentally (Rosenthal and Bezanilla 2000) and those observed in the model (see below for details of how errors were minimized; supplementary figure 1(C)). The data at each of the temperatures reported by Rosenthal and Bezanilla were then fit by exponentials:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{g}_{{\rm K}max}}\left( T \right)=~1.60~{{e}^{-{{\left( \frac{T-27.88}{12.85} \right)}^{2}}}},\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{g}_{{\rm Na}max}}\left( T \right)=~0.42~{{e}^{-{{\left( \frac{T-31.83}{31.62} \right)}^{2}}}},\nonumber \end{align} \tag{ 8 }$$

where the maximum potassium and sodium conductances, g_Kmax and g_Namax, have units of S cm⁻² (supplementary figure 1(D) shows plots of these functions).

Axial resistance.

Earlier work showed that the product of the axial resistance and the membrane capacitance could be estimated by the time constant of the rise of the foot of the action potential and its conduction velocity (Jack et al 1975). In our study, the 'foot' was defined as the region of the simulated action potential between the time points when the action potential begins (i.e. the time that the potential increased by 1% above the resting membrane potential; another 0.1 ms of data prior to this (100 points) was added to improve the exponential fit) and halfway between this point of onset and the point at which the derivative of the action potential reached its maximum value (supplementary figure 1(B), region indicated by the blue line). The resulting data was fit to an exponential, whose time constant, τ_foot, was then used for estimating the axial resistance (see below). The time constant of the model action potential shown in supplementary figure 1(B) was 179.9 μs.

We calculated conduction velocity as a function of temperature by measuring the time it took the action potential to travel between two points along the axon separated by 8 mm on either side of the midpoint of the axon. For all conduction velocity studies, the model axon was brought to the same uniform temperature along its entire length. The time of arrival of the action potential was determined as the instant when the action potential reached half its maximum value at a specific point on the axon. Conduction velocity was then calculated by dividing the distance between the two points on the axon by the time it took the action potential to travel between the points.

Having measured the time constant, τ_foot, and the conduction velocity, θ, from the model, we then determined the product of the membrane capacitance per unit area (C_m, in μF cm⁻²) and the axial resistivity (R_a, kΩ cm) from the equation

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{R}_{a}}{{C}_{m}}=\frac{r}{2{{\tau }_{foot}}{{\theta }^{2}}}\nonumber \end{align} \tag{ 9 }$$

where r is the radius of the axon (Jack et al 1975, pp 116–8). To derive this equation, it is assumed that the conduction velocity of the action potential is large relative to the electrotonic velocity, i.e. the passive propagation of signals along the axon; in addition, it is assumed that the extracellular resistance is negligible relative to the axial resistance.

To illustrate the computation for a specific example: at 5 °C, the time constant of the foot of the simulated action potential, τ_foot, was 179.9 μs (supplementary figure 1(B)); the conduction velocity of that action potential, θ, was 11.3 m s⁻¹. The simulated axon diameter is 500 μm, so that its radius is 250 μm. The resultant value of R_aC_m is $\frac{250\ {{10}^{-6}}\ {\rm m}}{2\left( 179.9\ {{10}^{-6}}~{\rm s} \right){{\left( 11.3\ {\rm m}\ {{{\rm s}}^{-1}} \right)}^{2}}}~$, which equals 0.0054 s m⁻¹ or 0.0054 Ω Farads m⁻¹. This value is similar to the value shown in figure 5(B) of Rosenthal and Bezanilla (2000) for axons at 5 °C from squid acclimatized to the cooler temperature of May. (Please note that the units shown in their figure 5, ΩFarads cm, do not appear to be correct.) Rosenthal and Bezanilla also determined experimentally that the capacitance of the membrane did not vary significantly from 1 μF cm⁻² in both warm and cold seasons. Thus, the equation for the axial resistance, R_a(T), in Ω cm, as a function of temperature (in °C; the coefficient of determination R² for the fit is 0.94) is

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{R}_{a}}\left( T \right)=56.84~{{e}^{-0.03\ {\rm T}}}.\nonumber \end{align} \tag{ 10 }$$

To fit parameters to the actual data, we minimized the squared deviations between model measurements for a given set of parameter values and the original data for the rates of rise and fall at each temperature. To capture changes both in the rate of rise and the rate of fall, we used a systematic grid search, varying the peak potassium and sodium conductances from their initial values in the Hodgkin/Huxley model in increments of 0.02 S cm⁻², and testing for two failure conditions: (1) spontaneous activity in the absence of any excitatory input, which indicated that the potassium conductance was too small to maintain a stable resting potential, and (2) failure to generate an action potential in response to excitatory current, indicating that the potassium conductance was large enough to overwhelm the inward sodium current.

The peak rate of fall was determined in the model action potential by finding the maximum of the derivative of the falling phase (supplementary figure 1(B)). The error in the rate of fall was computed using equation (11):

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{E}_{fall~}}\left( {{T}_{i}} \right)=\left| \frac{{{F}_{exp}}\left( {{T}_{i}} \right)-{{F}_{est}}\left( {{T}_{i}} \right)}{{{F}_{exp}}\left( {{T}_{i}} \right)} \right|,\nonumber \end{align} \tag{ 11 }$$

where ${{E}_{fall~}}\left( {{T}_{i}} \right)$ is the error at temperature T_i during the falling phase of the action potential, ${{F}_{exp}}\left( {{T}_{i}} \right)$ is the rate of fall measured experimentally at temperature T_i, and ${{F}_{est}}\left( {{T}_{i}} \right)$ is the estimated rate of fall from the model at temperature T_i.

We also computed the error in the experimentally measured rate of rise and that produced by the model. The peak rate of rise was determined from the derivative of the rising phase of the model action potential (supplementary figure 1(B)). Errors in the rate of rise were computed using equation (12):

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{E}_{rise~}}\left( {{T}_{i}} \right)=\left| \frac{{{R}_{exp}}\left( {{T}_{i}} \right)-{{R}_{est}}\left( {{T}_{i}} \right)}{{{R}_{exp}}\left( {{T}_{i}} \right)} \right|,\nonumber \end{align} \tag{ 12 }$$

where ${{E}_{rise~}}\left( {{T}_{i}} \right)$ is the error at temperature T_i in the rising phase, ${{R}_{exp}}\left( {{T}_{i}} \right)$ is the rate of rise measured experimentally at temperature T_i, and ${{R}_{est}}\left( {{T}_{i}} \right)$ is the estimated rate of rise from the model at temperature T_i.

Total error was then computed by summing the squares of the errors at each temperature over the entire set of temperatures that were matched, which were 5 °C, 7.5 °C, 10 °C, 12.5 °C, 15 °C, 17.5 °C, 20 °C, and 25 °C (based on the measurements of Rosenthal and Bezanilla (2000)). Thus, the final computed error was

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{E}_{total}}=&~\sum\nolimits_{{\rm all}~{{T}_{i}}}{{{\left( \frac{{{R}_{exp}}\left( {{T}_{i}} \right)-{{R}_{est}}\left( {{T}_{i}} \right)}{{{R}_{exp}}\left( {{T}_{i}} \right)} \right)}^{2}}} \nonumber \\ & +~\sum\nolimits_{{\rm all }~ {{T}_{i}}}{{{\left( \frac{{{F}_{exp}}\left( {{T}_{i}} \right)-{{F}_{est}}\left( {{T}_{i}} \right)}{{{F}_{exp}}\left( {{T}_{i}} \right)} \right)}^{2}}}.\nonumber \end{align} \tag{ 13 }$$

The total error is plotted on a logarithmic scale (supplementary figure 1(E)) for several variations of the model. The results demonstrate that if only the sodium/potassium pump is temperature dependent (supplementary figure 1(E), second set of bars), or if only the peak sodium and potassium conductances are varied as a function of temperature (supplementary figure 1(E), third set of bars), error is not reduced as much as if both features of the model vary with temperature (supplementary figure 1(E), fourth set of bars). Since Rosenthal and Bezanilla (2002) observed that at 29.5 °C the action potential in Loligo pealei was blocked, we used this information to determine the value for g_Kmax at this temperature.

Similar optimization was performed for the temperature dependence of the axial resistivity. The simulations were initialized with the standard Hodgkin/Huxley model value of membrane resistivity of 35.4 Ω cm and the conduction velocity at a given temperature T_i was evaluated as described above. The resultant conduction velocity was compared with the published results, and membrane resistivity was changed until the modeled conduction velocity was within 1% of the experimentally published values. No other parameters were varied to adjust the conduction velocity.

In the original experimental data (Rosenthal and Bezanilla 2000), measurements of the rates of rise and fall had accompanying error bars. We therefore measured not only the mean values shown in the data plots, but the values of the error bars, and propagated those values through the model. This allowed us to provide some estimate of the errors of the parameters, which are the maximum sodium and potassium conductances, the axial resistivity, and the product of R_aC_m. Supplementary figure 1(F) illustrates this for one parameter at one temperature (17.5 °C): the range of measured values for the rates of fall are shown along the x axis; the function relating the rate of fall to the peak potassium conductance is plotted, and the corresponding range in the peak potassium conductance parameter is indicated along the y  axis. The arrow points to the measured mean value of the parameter.

These are the model equations that were used in each segment of the NEURON model:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \begin{array}{@{}lllllccccccccccccccc@{}} {{C}_{m}}~\frac{d{{V}_{m}}}{dt}=\!\!\!\!&-{{g}_{L}}\left( {{V}_{m}}-~{{E}_{L}} \right)-{{g}_{{\rm K}max}}\left( T \right){{n}^{4}}\left( {{V}_{m}}-{{E}_{{\rm K}}} \right) \nonumber \\ &-~{{g}_{{\rm Na}max}}\left( T \right){{m}^{3}}h\left( {{V}_{m}} \right. \left. -{{E}_{{\rm Na}}} \right)-{{I}_{{\rm N}{{{\rm a}}_{pump}}}}-{{I}_{{{{\rm K}}_{pump}}}}+{{I}_{segment}}, \end{array}\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \frac{dm}{dt}={{\phi }_{m}}\left( T \right)\left( {{\alpha }_{m}}\left( 1-m \right)-{{\beta }_{m}}m \right),\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \frac{dh}{dt}={{\phi }_{h}}\left( T \right)\left( {{\alpha }_{h}}\left( 1-h \right)-{{\beta }_{h}}h \right),\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle \frac{dn}{dt}={{\phi }_{n}}\left( T \right)\left( {{\alpha }_{n}}\left( 1-n \right)-{{\beta }_{n}}n \right),\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{\alpha }_{n}}\left( V \right)=\frac{-0.01\left( V+55 \right)}{\left( {{e}^{\frac{-\left( V+55 \right)}{10}}}-1 \right)}\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{\beta }_{n}}\left( V \right)=0.125\left( {{e}^{\frac{-\left( V+65 \right)}{80}}} \right)\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{\alpha }_{m}}\left( V \right)=\frac{-0.1~\left( V+40 \right)}{\left( {{e}^{\frac{-\left( V+40 \right)}{10}}}-1 \right)}\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{\beta }_{m}}\left( V \right)=4\left( {{e}^{\frac{-\left( V+65 \right)}{18}}} \right)\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{\alpha }_{h}}\left( V \right)=0.07\left( {{e}^{\frac{-\left( V+65 \right)}{20}}} \right)\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{\beta }_{m}}\left( V \right)=\frac{1}{\left( {{e}^{\frac{-\left( V+35 \right)}{10}}}+1 \right)}\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{g}_{{\rm K}max}}\left( T \right)=~1.60~{{e}^{-~{{\left( \frac{T-27.88}{12.85} \right)}^{2}}}},\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{g}_{{\rm Na}max}}\left( T \right)=~0.42~{{e}^{-~{{\left( \frac{T-31.83}{31.62} \right)}^{2}}}}\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{g}_{{\rm NaK}}}\left( T \right)={{g}_{{\rm NaK}initial}}~{{2}^{\frac{\left( T-6.3 \right)}{10}}}\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{I}_{{\rm Na}pump}}=3{{g}_{{\rm NaK}}}\left( T \right)\left( {{V}_{m}}-{{E}_{pump}} \right)\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{I}_{{\rm K}pump}}=~-2{{g}_{{\rm NaK}}}\left( T \right)\left( {{V}_{m}}-{{E}_{pump}} \right)\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{I}_{segment}}=~\frac{{{V}_{j-1}}-~{{V}_{j}}}{{{r}_{a,j-1,j}}}-~\frac{{{V}_{j}}-{{V}_{j+1}}}{{{r}_{a,j,j+1}}}\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{\phi }_{m}}\left( T \right)=~\left\{\begin{array}{@{}llllllcccccccccccccc@{}} {{3}^{\frac{T-6.3}{10}}}, & ~5\ {}^\circ {\rm C}\leqslant T\leqslant 10\ {}^\circ {\rm C}, \nonumber \\ {{3}^{\frac{10-6.3}{10}}}{{3}^{\frac{T-~10}{10}}}, & 10\ {}^\circ {\rm C}<T\leqslant 15\ {}^\circ {\rm C}, \nonumber \\ {{3}^{\frac{10-6.3}{10}}}{{3}^{\frac{15-10}{10}}}{{2.8}^{\frac{T-15}{10}}}, & 15\ {}^\circ {\rm C}<T\leqslant 20\ {}^\circ {\rm C} \nonumber \\ {{3}^{\frac{10-6.3}{10}}}{{3}^{\frac{15-10}{10}}}{{2.8}^{\frac{20-15}{10}}}{{2.7}^{\frac{T-20}{10}}}, & 20\ {}^\circ {\rm C}<T\leqslant 25\ {}^\circ {\rm C} \nonumber \\ \end{array} \right.\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{\phi }_{h}}\left( T \right)=\left\{\begin{array}{@{}llllllcccccccccccccc@{}} {{3}^{\frac{T-6.3}{10}}}, & ~5\ {}^\circ {\rm C}\leqslant T\leqslant 10\ {}^\circ {\rm C}, \nonumber \\ {{3}^{\frac{10-6.3}{10}}}{{2.9}^{\frac{T-10}{10}}}, & 10\ {}^\circ {\rm C}<T\leqslant 15\ {}^\circ {\rm C}, \nonumber \\ {{3}^{\frac{10-6.3}{10}}}{{2.9}^{\frac{15-10}{10}}}{{3}^{\frac{T-15}{10}}}, & 15\ {}^\circ {\rm C}<T\leqslant 20\ {}^\circ {\rm C} \nonumber \\ {{3}^{\frac{10-6.3}{10}}}{{2.9}^{\frac{15-10}{10}}}{{3}^{\frac{20-15}{10}}}{{3}^{\frac{T-20}{10}}}, & 20\ {}^\circ {\rm C}<T\leqslant 25\ {}^\circ {\rm C} \nonumber \\ \end{array} \right.\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{\phi }_{n}}\left( T \right)=~\left\{\begin{array}{@{}lllllccccccccccccccc@{}} {{3}^{\frac{T-6.3}{10}}}, & ~5\ {}^\circ {\rm C}\leqslant T\leqslant 10\ {}^\circ {\rm C}, \nonumber \\ {{3}^{\frac{10-6.3}{10}}}{{2.8}^{\frac{T-10}{10}}}, & 10\ {}^\circ {\rm C}<T\leqslant 15\ {}^\circ {\rm C}, \nonumber \\ {{3}^{\frac{10-6.3}{10}}}{{2.8}^{\frac{15-10}{10}}}{{2.4}^{\frac{T-15}{10}}}, & 15\ {}^\circ {\rm C}<T\leqslant 20\ {}^\circ {\rm C} \nonumber \\ {{3}^{\frac{10-6.3}{10}}}{{2.8}^{\frac{15-10}{10}}}{{2.4}^{\frac{20-15}{10}}}{{2.3}^{\frac{T-20}{10}}}, & 20\ {}^\circ {\rm C}<T\leqslant 25\ {}^\circ {\rm C} \nonumber \\ \end{array} \right.~\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{R}_{a}}\left( T \right)=56.84~{{e}^{-0.03\ {\rm T}}}\nonumber \end{align}$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{r}_{a}}\left( T \right)=\frac{{{R}_{a}}\left( T \right)}{\pi {{r}^{2}}/4}.\nonumber \end{align}$$

The simulation environment NEURON was scripted using Python to search for parameter values while minimizing errors, as described above. The Python package 'numpy' was used for numerical analysis, and the package 'matplotlib' was used for plotting data. Fits of data to curves were done using MATLAB's Curve Fitting toolbox. All code used in the modeling and optimization for this paper are available for download under the appropriate licenses from the GitHub repository: (https://github.com/mohitganguly/IRBlock_Vanderbilt).

We first determined whether adding temperature-dependence to components of the original Hodgkin/Huxley led to qualitatively new behavior, or primarily altered the quantitative response of the modified model to temperature. Since most temperature manipulations generate a smooth change from an initial temperature to a final temperature, we generated an action potential along a model axon whose first half was at the control temperature of 6.3 °C, and then smoothly changed the temperature of the second half of the model axon to 25 °C (see equations (2) and (3), Methods; figures 2(A) and (B)). By generating the temperature change near the middle of the axon, we ensured that initial currents at the point of action potential initiation were gone by the time the action potential encountered the change in temperature. The consequence of a gradual change in temperature was a clear increase in net hyperpolarizing (outward) current for both the original Hodgkin/Huxley model (figures 2(C)–(F), position 2) and for the modified model (figures 2(G)–(J), position 2) at the center of the region undergoing a change in temperature. In the region at the control temperature, the relative ratio of outward to inward current was more similar for both models (figures 2(C)–(F) and (G)–(J), position 1). In the region that was at the new, higher temperature in the modified model, the ratio of outward to inward current was also increased (figures 2(G)–(J), position 3), although the currents were much smaller.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

To better understand the development of the net hyperpolarizing current as a function of temperature, we confined the temperature change to a single step change (figure 3), since in figure 2, we observed that the changes in currents primarily occur where the axon shows significant changes in temperature. Both the original Hodgkin/Huxley model and the modified model qualitatively showed similar responses: as temperature increased from 5 °C to 25 °C, the sodium inward current became shorter (figures 3(B) and (F), for Hodgkin/Huxley model and for modified model, respectively), due to the increased speed of sodium ion channel inactivation. Over the same temperature range, the potassium outward current became faster and larger (figures 3(C) and (G), for Hodgkin/Huxley model and modified model, respectively), due to the increased speed of potassium ion channel activation. As a consequence, the net current (I_total) became increasingly outward (hyperpolarizing) with increasing temperature as the membrane was depolarized by the action potential (figures 3(D) and (H) for Hodgkin/Huxley model and modified model, respectively). For both models, much more outward charge was transferred relative to inward charge (figures 3(E) and (I)).

In both models, increased temperature alters both the sodium and the potassium ion channels, so both could contribute to thermal block. To determine the contributions of the voltage-dependent ion channels to thermal block, we changed the temperature of a model neuron in its center while selectively removing ion channels computationally from that central region, and then initiating an action potential at one end of the axon. We recorded the voltage changes at the beginning of the central region at which the temperature change occurred, and at the end of the model axon (figures 4(A) and (K), points marked 1 and 2). In both models, when all ion channels were present, we determined a minimum length and temperature at which action potential propagation was blocked (35 °C applied to a 5.6 mm length at the center of the axon for the Hodgkin/Huxley model, 29.5 °C applied to an 0.9 mm length at the center of the axon for the modified model).

Zoom In Zoom Out Reset image size

As the action potential reached the region of changed temperature, its amplitude was significantly reduced (solid red lines in figures 4(B) and (L) for the Hodgkin/Huxley and modified Hodgkin/Huxley models, respectively; dashed blue lines show the action potential when the central region is at the control temperature of 6.3 °C). No action potential was recorded at the end of the model axon (solid red lines in figures 4(C) and (M) for the Hodgkin/Huxley and modified Hodgkin/Huxley models, respectively; dashed blue lines show the action potential when the central region is at the control temperature—the action potential propagates through the central region unchanged to the end of the axon).

To demonstrate that all of the voltage-gated ion channels were necessary for thermal block, we set the conductances of both the sodium and potassium voltage-gated ion channels to zero within the central region. Removing ion channels can change the resting potential of the central region because the voltage-dependent ion channels contribute significantly to the resting potential. Since, at rest, the membrane is primarily permeable to potassium ions, removing all ion channels will tend to depolarize the membrane. In turn, a depolarized membrane will alter the states of voltage-dependent ion channels on either side of the central region. To eliminate this potentially confounding factor, we applied compensatory current to the central region to ensure that it remained at the control resting potential after ion channels were removed computationally (see methods, compensatory current, and supplementary figure 2).

When the action potential encountered the central region, its amplitude decreased and the waveform broadened (figures 4(D) and (N), dashed blue lines). By the time the action potential reached the end of the model axon, it was identical to the control action potential (figures 4(E) and (O), dashed blue lines). Removing all the voltage-dependent ion channels did not prevent passive currents from propagating through the region in which the ion channels were blocked, and these passive currents re-activated an action potential that propagated to the end of the model axon.

When the central region was subjected to an increased temperature, the Hodgkin/Huxley model showed no changes (figures 4(D) and (E); the solid red lines corresponding to the increased temperature condition are completely covered by the dashed blue lines). The lack of response is predictable, as all of the temperature-dependent aspects of the model have been removed. In the modified Hodgkin/Huxley model, both the axial resistance (figure 1(D)) and the sodium/potassium pump depend on temperature. When the central region was subjected to an increased temperature, the action potential showed a small increase in width relative to the control (figure 4(N), dashed red line), but the action potential at the end of the axon was essentially unchanged (figure 4(O); the dashed red line is covered by the control solid blue line). Thus, both models clearly require the voltage-dependent ion channels to generate thermal block in response to increased temperature.

To determine if voltage-gated sodium ion channels are necessary for thermal block, we set the conductance of the voltage-gated sodium ion channels to zero. We did not modify the conductance of the voltage-gated potassium ion channels. Again, compensatory current was used to ensure that the resting potential was unchanged by computationally removing the voltage-gated sodium ion channels. In both models, when the central region was kept at the control temperature of 6.3 °C, the action potential peak was reduced as it reached the central region (figures 4(F) and (P) for the Hodgkin/Huxley and modified models, respectively, solid blue lines). Passive currents were able to propagate through the region and re-activate an action potential, which arrived at the end of the axon essentially unchanged (figures 4(G) and (Q), respectively, solid blue lines). When the central region was subjected to increased temperature, the peak of the action potential was significantly reduced in both models as it reached the region of increased temperature (figures 4(F) and (P), respectively, dashed red lines). The action potential was completely blocked by the end of the axon (figures 4(G) and (Q), respectively, dashed red lines). These results demonstrate that voltage-gated sodium ion channels are not necessary, whereas voltage-gated potassium ion channels are alone sufficient to induce thermal block.

To determine if voltage-gated potassium channels are necessary for thermal block, we set the conductance of these channels to zero, but did not modify the conductance of the voltage-gated sodium ion channels. Again, we used compensatory current to ensure that the resting potential was unchanged despite the computational removal of the voltage-gated potassium ion channels. In both models, when the central region was kept at the control temperature of 6.3 °C, the action potential broadened slightly as it reached the central region (figures 4(H) and (R) for the Hodgkin/Huxley and modified models, respectively, solid blue lines), but propagated normally to the end of the axon (figures 4(I) and (S), respectively, solid blue lines). When the central region was subjected to increased temperature, the action potential narrowed slightly as it reached the central region (figures 4(H) and (R), respectively, dashed red lines), but propagated through and reached the end of the axon essentially unchanged (figures 4(I) and (S), respectively, dashed red lines). Under these conditions, voltage-gated potassium ion channels are necessary for thermal block, and voltage-gated sodium channels are not sufficient alone to induce thermal block.

To generalize these results, minimum block lengths for all four conditions (all voltage-gated ion channels intact, all voltage-gated ion channels blocked, voltage-gated sodium ion channels blocked, or voltage-gated potassium ion channels blocked) were measured in both models over a wide range of temperatures, with compensatory current applied to ensure no change from the resting potential after computational removal of ion channels. Since voltage-gated sodium ion channels contribute depolarizing current and voltage-gated potassium ion channels contribute hyperpolarizing current, the block lengths when all ion channels are intact should lie midway between the block lengths for either of the two voltage-gated ion channels alone. If one of the voltage-gated ion channels dominates, then the block lengths when all ion channels are intact should lie closer to the block lengths for that ion channel alone. The data suggest that, at the lowest temperature that can cause block in the Hodgkin/Huxley model, both voltage-gated sodium and potassium ion channels contribute to the block, but as the temperature is raised, the block length is dominated by the voltage-gated potassium ion channels (figure 4(J)). In contrast, for the modified Hodgkin/Huxley model, at all temperatures capable of inducing block, the block length associated with the potassium channels dominates (figure 4(T)).

As another way of testing the role of the different voltage-dependent ion channels in thermal inhibition, we changed the temperature-dependence of the gating variables m, h and n in the modified Hodgkin Huxley model. By setting the Q₁₀ of n to 1, we ensured that the voltage-dependent potassium ion channels would not open more rapidly as the axon temperature increased. After this modification, application of increased temperature did not induce block of the action potential (data not shown). In contrast, after setting the Q₁₀ of m and h to 1, eliminating the temperature-dependence of the voltage-dependent sodium ion channels, increased temperature could still induce block (data not shown). Because it is not currently feasible to experimentally eliminate the temperature-dependence of ion channels, we did not pursue this approach further.

Applying compensatory current to an axon may be experimentally challenging, so we determined whether qualitatively similar results were observed if the central region's resting potential was allowed to change after computationally removing voltage-dependent ion channels in either the original Hodgkin/Huxley model or the modified model (supplementary figures 3(A) and (F)). As expected, in both models, computational removal of all voltage-dependent ion channels from the central region caused a small depolarization of the central region and regions around it (supplementary figures 3(B) and (G)); computational removal of all voltage-dependent potassium ion channels caused a larger depolarization (supplementary figures 3(C) and (H)); and computational removal of all voltage-dependent sodium ion channels caused a small hyperpolarization (supplementary figures 3(D) and (I)). Despite these changes in potential, the qualitative results we obtained for length of block needed as a function of temperature were very similar to those obtained using compensatory current (supplementary figures 3(E) and (J); compare figures 4(J) and (T)): at higher block temperatures, the length of the block was clearly dominated by the voltage-dependent potassium ion channels. In the Hodgkin/Huxley model, at lower block temperatures, the voltage-dependent sodium ion channels clearly played an important role in determining block length (supplementary figure 3(E)). In the absence of compensatory current, the central region of the axon depolarizes once the voltage-gated potassium ion channels are blocked. In turn, this depolarization reduces the inward current through the voltage-dependent sodium ion channels because they inactivate in response to depolarization. The depolarizing current is also partially shunted by the increased conductance when the voltage-dependent sodium ion channels are activated. The net result is the same: the block lengths when all ion channels are intact lie between block lengths for either of the two voltage-gated ion channels alone.

Since neurons generally fire repetitively in response to stimuli, we tested whether a thermal block that was sufficient to stop an action potential could also block repetitive firing (figure 5(A)). Thermal block sufficient to block a single action potential was able to block repetitive firing, whether induced by a steady depolarizing current (figure 5(B)), or by individual depolarizing current pulses (figure 5(C)).

Zoom In Zoom Out Reset image size

Unmyelinated fibers range widely in size across phylogeny; some are as small as 0.2 μm (e.g. unmyelinated C fibers in a muscle nerve or a cutaneous nerve; Gardner and Johnson 2013). We sought to determine the minimal thermal block lengths as a function of fiber diameter. A mathematical analysis of the cable equation predicted that modalities acting on the surface of axons would scale with the square root of the diameter of the axon (Lothet et al 2017, supplemental material). In the mathematical model, a generalized form of the cable equation was derived in which the dependence of membrane capacitance and resistance and axial resistance on axon diameter was made explicit. The equation was then rewritten to re-scale length by the square root of the axon diameter, and in this new coordinate system, the dependence of all the terms on axon diameter could be factored out, making the new equation independent of axon diameter. In turn, this immediately implied that if one axon has diameter d₁, and could be blocked by a high-temperature region whose minimum length was L, then an axon of diameter d₂ would be blocked by a high temperature region of minimum length $L\sqrt{\frac{{{d}_{2}}}{{{d}_{1}}}}$. These mathematical predictions have been confirmed by the experimental studies in Lothet et al (2017), which demonstrated that thermal inhibition could selectively block small-diameter axons before large-diameter axons in both an invertebrate nerve (from Aplysia californica) and the vertebrate vagus (from the musk shrew Suncus murinus). The conduction velocity and length constant of an unmyelinated axon also scale with the square root of axon diameter (Hodgkin 1954, Debanne et al 2011). These observations suggest that the spread of current away from the region of increased temperature should also scale as the square root of the axon diameter. Thus, if voltage changes were plotted against the distance scaled by the length constant of the axon, which is proportional to the square root of the axon diameter, the resulting voltage changes should be independent of axon diameter.

The length constant of an unmyelinated axon is defined from the solution of the cable equation after time-dependent changes in voltage have died out while a steady current is injected into the cable, yielding

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{V}_{m}}\left( x \right)=~{{V}_{0}}{{e}^{-x/\lambda }},\nonumber \end{align} \tag{ 14 }$$

where V_m(x) is the membrane potential V_m at point x along the cable, V₀ is the voltage at the point of current injection, and λ is the length constant (Jack et al 1975). Note that when the original voltage V₀ has fallen by 1/e of its original value, x  =  λ. Thus, to determine the length constant of an axon with a region of increased temperature, we injected a steady hyperpolarizing current, and determined the distance at which the hyperpolarization had fallen to a value of 1/e of its value at the point of injection. For a 500 μm axon, this value was 3.96 mm. The length constant of the same axon at its original temperature was 3.93 mm (a difference of 1%). To determine how the current spread scaled with axon diameter, an action potential and its blocked version were also determined for a 10 μm diameter axon, and then rescaled by the length constant (the length constant was 0.56 mm; note that 3.96 is equal to 0.56 times the square root of the ratio of the diameters, i.e. 0.56  ×  $\sqrt{\frac{500}{10}~}$  =  3.96, as predicted by the square root scaling law).

To determine the spatial extent of an action potential in the absence of block, an action potential was initiated in response to a stimulating current (2000 nA for 1 ms) in an axon at its original temperature. At the instant the peak of the action potential reached the middle of the axon, we recorded the voltage along the entire axon (figure 6(A), dark blue line). Note the large spatial extent of the action potential. It is also worth noting the difference in shape from the action potentials shown in figure 4: in that figure, we record the change in voltage over time at the initial point of the block region and at the end of the axon, so we look at the time profile at a single point in space. In contrast, here we record the change in voltage over the entire length of the axon (i.e. over all space) at a single instant in time. To determine the effect of the spatial extent of the thermal block on the action potential, an action potential was initiated while the minimum length along the axon necessary for thermal block (1.12 mm) had its temperature increased to 29.5 °C, starting at the middle of the axon. At the identical time after the stimulating current initiated the action potential, we recorded the voltage along the entire axon. Note the large spatial extent of the block (figure 6(A), dashed green line). We rescaled these results by the axonal length constant so they could be compared to results obtained in axons of different diameters (figure 6(B)). We then repeated the simulation using a 10 μm diameter axon (figure 6(C)); when the results were rescaled by the length constant, they were essentially identical to those observed in the 500 μm diameter axon (figure 6(D); compare figure 6(B)).

Zoom In Zoom Out Reset image size

The scaling of thermal block suggested that it would act first on smaller-diameter axons, rather than on larger-diameter axons. Once axons are scaled by the square root of their diameters, the spatial spread of voltage at a fixed point of time becomes identical. In turn, this implies that the voltage changes experienced by the voltage-dependent ion channels will also be identical, and thus the currents induced by these depolarizations will be identical. As a consequence, this numerically confirms the mathematical model of Lothet et al (2017) (supplemental material), which demonstrated that the length of block would scale with the square root of axon diameter. In contrast, extracellular currents are known to act on larger-diameter axons before smaller-diameter axons, since these currents scale as the square of the axon diameter (Rattay 1986). To further test the predictions of the mathematical analysis (Lothet et al 2017), we determined the minimum lengths of the central region that needed to be increased in temperature to induce thermal block as a function of axon diameter. As predicted by the mathematical model, the block lengths scaled linearly with the square root of the axon diameter (figure 7). The scaling effects were observed for the Hodgkin/Huxley model (figure 7(A)), as well as for the modified Hodgkin/Huxley model (figure 7(B)) even at diameters that are similar to those of unmyelinated C fibers in vertebrates (figure 7(C)). Furthermore, as temperatures were increased, the minimum block length decreased for both models (figure 7).

Zoom In Zoom Out Reset image size

When an action potential is initiated in the initial segment of the axon, the currents generated are larger than those generated by the propagating action potential (Kole et al 2008). Since the thermal blocks described above were applied to a propagating action potential, thermal block of action potential initiation might be different (schematics of experiment shown figures 8(A) and (D), respectively). To determine how thermal block of initiation scales with axon diameter, we first determined the minimum current needed to initiate an action potential for a given diameter axon. Whenever an action potential is initiated using a very short current pulse, there is little time for activation of voltage-dependent ion channels, and thus the current needed to reach the voltage threshold is inversely proportional to the passive input resistance of the cable. In turn, this implies that the minimum current needed to initiate an action potential would scale as the square root of the axon diameter cubed (Jack et al 1975, p 32). Indeed, the minimum currents needed to initiate an action potential scaled in this way for both the Hodgkin Huxley and the modified model (figures 8(B) and (E), respectively). When thermal block was applied to either model after it was stimulated by the appropriate threshold current to initiate an action potential, the inhibition block temperature again scaled with the square root of axon diameter (figures 8(C) and (F), respectively).

Zoom In Zoom Out Reset image size

These results suggest that thermal block in unmyelinated axons, especially at higher temperatures, are primarily due to the activation of voltage-dependent potassium ion channels in response to depolarization, which induce a strong hyperpolarization and thus block action potential initiation and propagation, and can act preferentially to block smaller-diameter axons before larger-diameter axons.