Computational modeling of Adelta-fiber-mediated nociceptive detection of electrocutaneous stimulation

Abstract

Sensitization is an example of malfunctioning of the nociceptive pathway in either the peripheral or central nervous system. Using quantitative sensory testing, one can only infer sensitization, but not determine the defective subsystem. The states of the subsystems may be characterized using computational modeling together with experimental data. Here, we develop a neurophysiologically plausible model replicating experimental observations from a psychophysical human subject study. We study the effects of single temporal stimulus parameters on detection thresholds corresponding to a 0.5 detection probability. To model peripheral activation and central processing, we adapt a stochastic drift-diffusion model and a probabilistic hazard model to our experimental setting without reaction times. We retain six lumped parameters in both models characterizing peripheral and central mechanisms. Both models have similar psychophysical functions, but the hazard model is computationally more efficient. The model-based effects of temporal stimulus parameters on detection thresholds are consistent with those from human subject data.

Keywords: Computational models; Detection threshold; Nociceptive pathway; Stimulus detection; Stimulus parameters.

Fig. 1

Detection thresholds for each subject, study day, and combination of temporal stimulus parameters (Table 1). The larger solid circles and crosses present the mean detection thresholds for day 1 and day 2, respectively

Fig. 2

An electrode is attached to the skin of a subject to deliver pulse train stimulation. The dot-dashed concentric half circles represent the electric potential. Charging the nerve endings leads to traveling action potentials in the A$\mathit{\delta }$-fiber. The arrival of spikes at the presynaptic terminal triggers the release of the glutamate from the synapse resulting in an EPSC. The secondary neuron is charged, and the activity will converge upto the supraspinal part and lead to a binary response. Note that the number of signal channels is the number of secondary neurons, i.e., four in this diagram

Fig. 3

Illustration of the geometry of nerve endings under skin, a minimal depth is denoted by h. The endings with solid tips are recruited and those with empty tips are not recruited. The gray surface represents the recruited space, i.e., within critical distance $r_c$

Fig. 4

Activation of afferent fibers. a activation has a threshold nonlinear relation with amplitude A; b peripheral activity increases by increasing PW, eventually saturating

Fig. 5

Stochastic dynamics of the DDM using a pulse train current input with ($A=1$ mA, $\text{NoP}=2$, $\text{IPI}=50$ ms, $\text{PW}=0.42$ ms). a Noise-free dynamics (solid) with $\mathit{\sigma }=0$ and a stochastic realization (dashed). b $N=200$ Realizations of the stochastic dynamics (thin solid). Statistics of the potential are also shown, mean (solid), and mean plus or minus the standard deviation (dashed). See text for system parameter values

Fig. 6

Activities of secondary neurons using three different stimuli with the same $\text{PW}=0.42$ ms: $\text{NoP}=1$ (thick dashed); $\text{NoP}=2$ and $\text{IPI}=50$ ms (solid); $\text{NoP}=2$ and $\text{IPI}=150$ ms (dotted-dashed). a Lumped PSP stimulated by an electrical train of two pulses with amplitude $A=1$ mA, $\text{IPI}=50$ ms and $\text{PW}=0.42$ ms; b instantaneous firing rate; c the expected value of the number of spikes within a trial [0, T]; d psychophysical function value ${\mathit{\Psi }}_H$ depending on T

Fig. 7

Using 8 different combinations of temporal stimulus parameters. Temporal parameters used in panels (a–h) correspond to combinations A–H in Table 2, psychophysical function values ${\mathit{\Psi }}_D(A)$ using the DDM (solid lines) and best fitted ${\mathit{\Psi }}_H$ by the HM (dashed lines) with fitting error $E=0.0029$. The common parameters are set to ${\mathit{\alpha }}_1=0.5$ mA, ${\mathit{\tau }}_1=0.1$ ms, ${\mathit{\tau }}_2=50$ ms, and ${\mathit{\tau }}_s=1.5$ ms. The parameters corresponding to neuronal variability are different in these two models. In the diffusion model, we set the values of parameter as ${\mathit{\alpha }}_2=0.02$ A/s, $\mathit{\sigma }=0.05$ A/s, and $l=1$, while the fitting results in ${\mathit{\alpha }}_L=0.0220$ A/s, ${\mathit{\sigma }}_L=0.0021$ A/s, and ${\mathit{\lambda }}_L=0.4020$ kHz. The asymptotic behavior of the detection threshold with two independent pulses and its relation to the psychophysical curve with $\text{NoP}=1$ are illustrated by the thin dashed lines in panels (b) and (h)

Fig. 8

Distribution of the fitting error of the DDM by the hazard model with randomly chosen values of system parameters ${\mathit{\tau }}_1$, ${\mathit{\alpha }}_1$, ${\mathit{\tau }}_2$, ${\mathit{\alpha }}_2$, $\mathit{\sigma }$ and l in the DDM

Fig. 9

Dependence of detection thresholds on temporal parameters. a Experimental detection thresholds (solid circles and crosses) as in Fig. 1 and model-based thresholds (DDM: open circles, HM: open squares). The index refers to the combinations of Table 1 and see text for parameter values. The horizontal lines for $\text{NoP}=2$ indicate the asymptotic value $A_{2,50}$ for two independent pulses (DDM: solid line, HM: dashed line). b Non-monotone dependence of two-pulse threshold $A_{2,50}$ on IPI

Fig. 10

Model-simulated non-monotone effects of the IPI on the probability to detect when in the DDM (circles) and the HM (squares) for $\text{NoP}=2$. The amplitude is fixed as $A=0.35$ mA, and values of system parameters are identical to those used in Fig. 9