Computational analysis of a 9D model for a small DRG neuron

Abstract

Small dorsal root ganglion (DRG) neurons are primary nociceptors which are responsible for sensing pain. Elucidation of their dynamics is essential for understanding and controlling pain. To this end, we present a numerical bifurcation analysis of a small DRG neuron model in this paper. The model is of Hodgkin-Huxley type and has 9 state variables. It consists of a Na_v1.7 and a Na_v1.8 sodium channel, a leak channel, a delayed rectifier potassium, and an A-type transient potassium channel. The dynamics of this model strongly depend on the maximal conductances of the voltage-gated ion channels and the external current, which can be adjusted experimentally. We show that the neuron dynamics are most sensitive to the Na_v1.8 channel maximal conductance ($\overline {g}_{1.8}$). Numerical bifurcation analysis shows that depending on $\overline {g}_{1.8}$ and the external current, different parameter regions can be identified with stable steady states, periodic firing of action potentials, mixed-mode oscillations (MMOs), and bistability between stable steady states and stable periodic firing of action potentials. We illustrate and discuss the transitions between these different regimes. We further analyze the behavior of MMOs. As the external current is decreased, we find that MMOs appear after a cyclic limit point. Within this region, bifurcation analysis shows a sequence of isolated periodic solution branches with one large action potential and a number of small amplitude peaks per period. For decreasing external current, the number of small amplitude peaks is increasing and the distance between the large amplitude action potentials is growing, finally tending to infinity and thereby leading to a stable steady state. A closer inspection reveals more complex concatenated MMOs in between these periodic MMO branches, forming Farey sequences. Lastly, we also find small solution windows with aperiodic oscillations which seem to be chaotic. The dynamical patterns found here—as consequences of bifurcation points regulated by different parameters—have potential translational significance as repetitive firing of action potentials imply pain of some form and intensity; manipulating these patterns by regulating the different parameters could aid in investigating pain dynamics.

Appendix:

1.1 1.1 Na_v1.7 equations

$$ \alpha_{m_{1.7}} = \frac{15.5}{1+\exp\left( \frac{(V-5)}{-12.08}\right)} $$

(7)

$$ \beta_{m_{1.7}} = \frac{35.2}{1+\exp\left( \frac{V+72.7}{16.7}\right)} $$

(8)

$$ \alpha_{h_{1.7}} = \frac{0.38685}{1+\exp\left( \frac{V+122.35}{15.29}\right)} $$

(9)

$$ \beta_{h_{1.7}} = -0.00283 + \frac{2.00283}{1 + \exp\left( \frac{(V+5.5266)}{-12.70195}\right)} $$

(10)

$$ \alpha_{s_{1.7}} = 0.00003 + \frac{0.00092}{1+\exp\left( \frac{V+93.9}{16.6}\right)} $$

(11)

$$ \beta_{s_{1.7}} = 132.05 - \frac{132.05}{1+\exp\left( \frac{V-384.9}{28.5}\right)} $$

(12)

For x = m_1.7,h_1.7,s_1.7:

$$ x_{\infty} (V) = \frac{\alpha_{x} (V)}{\alpha_{x} (V) + \beta_{x} (V)}, $$

(13)

and

$$ \tau_{x} (V) = \frac{1}{\alpha_{x} (V) + \beta_{x} (V)} $$

(14)

The kinetics of Na_v1.7 were taken from Sheets et al. (2007) and Choi and Waxman (2011). m_1.7 corresponds to the activation gating variable, h_1.7 to the fast-inactivation gating variable, and s_1.7 to the slow-inactivation gating variable.

1.2 1.2 Na_v1.8 equations

$$ \alpha_{m_{1.8}} = 2.85 - \frac{2.839}{1+\exp\left( \frac{V-1.159}{13.95}\right)} $$

(15)

$$ \beta_{m_{1.8}} = \frac{7.6205}{1+\exp\left( \frac{V+46.463}{8.8289}\right)} $$

(16)

For x = m_1.8:

$$ x_{\infty} (V) = \frac{\alpha_{x} (V)}{\alpha_{x} (V) + \beta_{x} (V)}, $$

(17)

and

$$ \tau_{x} (V) = \frac{1}{\alpha_{x} (V) + \beta_{x} (V)} $$

(18)

$$ \tau_{h_{1.8}} = 1.218 + 42.043\times \exp\left( \frac{-(V+38.1)^{2}}{2\times15.19^{2}}\right) $$

(19)

$$ h_{{1.8}_{\infty}} = \frac{1}{1+\exp\left( \frac{V+32.2}{4}\right)} $$

(20)

The kinetics of Na_v1.8 were taken from Sheets et al. (2007) and Choi and Waxman (2011). m_1.8 and h_1.8 are similar fast activation and slow inactivation gating variables, respectively.

1.3 1.3 K equations

$$ \alpha_{n_{K}} = \frac{0.001265\times (V+14.273)}{1-\exp\left( \frac{V+14.273}{-10}\right)} $$

(21)

with $\alpha _{n_{K}} = 0.001265\times 10$ for V = − 14.273.

$$ \beta_{n_{K}} = 0.125\times \exp\left( \frac{V+55}{-2.5}\right) $$

(22)

$$ n_{K_{\infty}} = \frac{1}{1 + \exp\left( \frac{-(V + 14.62)}{18.38}\right)} $$

(23)

$$ \tau_{n_{K}} = \frac{1}{\alpha_{n_{K}} + \beta_{n_{K}}} + 1 $$

(24)

The kinetics of K channel were taken from Schild et al. (1994). It only includes an activating n_K gating variable.

1.4 1.4 KA equations

$$ n_{{KA}_{\infty}} = \left( \frac{1}{1 + \exp\left( \frac{-(V + 5.4)}{16.4}\right)}\right)^{4} $$

(25)

$$ \tau_{n_{KA}} = 0.25 + 10.04\times \exp\left( \frac{-(V + 24.67)^{2}}{2 \times 34.8^{2}}\right) $$

(26)

$$ h_{{KA}_{\infty}} = \frac{1}{1 + \exp\left( \frac{V + 49.9}{4.6}\right)} $$

(27)

$$ \tau_{h_{KA}} = 20 + 50\times \exp\left( \frac{-(V + 40)^{2}}{2 \times 40^{2}}\right) $$

(28)

The kinetics of KA channel were taken from Sheets et al. (2007). It consists of one fast activating and one slow inactivating gate.

It needs to be noted that “fast” and “slow” are relative terms for each individual channel. The time scales of each of these gates are demonstrated in Verma et al. (2020), where it was shown that (m_1.7,m_1.8), (n_KA), (h_1.7,h_1.8,h_KA), (n_K), and (s_1.7) belong to similar time scales (in increasing order).

The XPPAUT and MATCONT codes for this model can be found on ModelDB (http://modeldb.yale.edu/264591).