Predictive analytics identifies key factors driving hyperalgesic priming of muscle sensory neurons

Hyperalgesic priming, a form of neuroplasticity induced by inflammatory mediators, in peripheral nociceptors enhances the magnitude and duration of action potential (AP) firing to future inflammatory events and can potentially lead to pain chronification. The mechanisms underlying the development of hyperalgesic priming are not well understood, limiting the identification of novel therapeutic strategies to combat chronic pain. In this study, we used a computational model to identify key proteins whose modifications caused priming of muscle nociceptors and made them hyperexcitable to a subsequent inflammatory event. First, we extended a previously validated model of mouse muscle nociceptor sensitization to incorporate Epac-mediated interaction between two G protein-coupled receptor signaling pathways commonly activated by inflammatory mediators. Next, we calibrated and validated the model simulations of the nociceptor’s AP response to both innocuous and noxious levels of mechanical force after two subsequent inflammatory events using literature data. Then, by performing global sensitivity analyses that simulated thousands of nociceptor-priming scenarios, we identified five ion channels and two molecular processes (from the 18 modeled transmembrane proteins and 29 intracellular signaling components) as potential regulators of the increase in AP firing in response to mechanical forces. Finally, when we simulated specific neuroplastic modifications in Kv1.1 and Nav1.7 alone as well as with simultaneous modifications in Nav1.7, Nav1.8, TRPA1, and Kv7.2, we observed a considerable increase in the fold change in the number of triggered APs in primed nociceptors. These results suggest that altering the expression of Kv1.1 and Nav1.7 might regulate the neuronal hyperexcitability in primed mechanosensitive muscle nociceptors.

Introduction

Acute pain arising from musculoskeletal injuries sometimes persists beyond the typical tissue-recovery time and transitions to chronic pain. The mechanisms underlying this transition from acute to chronic pain are not well understood, hampering the development of new treatment options (Heinricher, 2016; Pozek et al., 2016). However, we do know that primary afferent nociceptors in muscle tissue recognize noxious stimuli (e.g., heat and mechanical forces) and respond to them by firing action potentials (APs) via the action of numerous membrane proteins and intracellular molecules. These nociceptors can be sensitized by inflammatory mediators, such as prostaglandin E₂ (PGE₂), and develop long-lasting hyperexcitability (i.e., an increase in the magnitude and duration of AP firing) to a subsequent challenge, a phenomenon known as hyperalgesic priming (Gold and Flake, 2005; Reichling and Levine, 2009; Gold and Gebhart, 2010). Priming can occur at either the peripheral or central terminals (in the spinal cord) of nociceptors, and it is believed to initiate the transition from acute to chronic pain (Gold and Gebhart, 2010; Gangadharan and Kuner, 2013). Thus, the “primed” state represents a kind of “pain memory” that could potentially be reversed to treat chronic pain (Kandasamy and Price, 2015). However, the specific neuroplastic changes in the nociceptor that constitute this memory and maintain its priming are not well understood. We know that inflammatory mediators, such as PGE₂ and tumor necrosis factor-α, activate a number of different signaling pathways within the neuron, including G protein-coupled receptors (GPCR), nuclear factor-κB, and mitogen-activated protein kinases (Niederberger and Geisslinger, 2008; Wuertz et al., 2012). However, only a few molecular mechanisms, such as protein kinase C (PKC) activity or PKC_ε-dependent switching between the distinct GPCR pathways, are believed to lead to hyperalgesic priming (Reichling and Levine, 2009). Moreover, we still do not know the contribution of the different membrane proteins to the increased neuronal response in primed nociceptors. Because the AP response originates at the peripheral nociceptors before it is transmitted to other neurons, knowledge of the key neuroplastic modifications contributing to the priming of these nociceptors and of the AP responses to multiple inflammatory insults is necessary to improve our understanding of the pain transition process and to help identify potential therapeutic targets to prevent it. In addition, targeting the peripheral nervous system is advantageous because it provides an opportunity to avoid side effects from targeting neurons in the central nervous system, e.g., targeting opioid receptors present in spinal cord neurons to treat pain can paradoxically contribute to the development of priming in many different nociceptors (Khomula et al., 2019, 2021).

While numerous behavioral animal models of hyperalgesic priming have helped identify certain key molecules and pathways, including PKC activation by exchange protein directly activated by cAMP (Epac), among others (Aley and Levine, 1999; Hucho and Levine, 2007; Ferrari et al., 2013a; Araldi et al., 2016a,b), it is difficult to know whether targeting them can particularly regulate the responses of primary afferent nociceptors because pain behavior results from the aggregation of responses by many different neurons located in the peripheral and central nervous systems (Dubin and Patapoutian, 2010). On the other hand, in vivo models to assess primary afferent nociceptors, especially muscle afferents, are challenging because these neurons are heavily embedded in peripheral tissues (e.g., skin and muscle), making their nerve endings hard to access (Mense, 2010). In addition, quantifying the effects of priming through experiments involves blocking, knocking out, or enhancing individual ion channels or receptors in different inflammatory scenarios, some of which may not be feasible. Furthermore, these experiments are complicated by the presence of ongoing inflammation due to endogenous mediators, preventing the distinction between mechanisms of ongoing acute pain and the neuroplastic changes that specifically underlie the increase in the amount and frequency of AP firing and its prolongation during multiple inflammatory events.

Computational modeling can complement traditional experimentation in the search for key proteins or processes that could potentially regulate the response of primed nociceptors to future inflammatory events. Through computational modeling, we can compute the effects of various neuroplastic modifications, such as knocking out or overexpressing a given protein, or the effects of blocking or activating a molecular process on AP generation in response to a combination of different types of noxious stimuli and inflammatory mediators in a systematic and time-efficient manner. In fact, previous computational models of pain signaling in nociceptive neurons yielded insights into the roles of specific ion channels. For example, these models have characterized the contributions of different Na⁺ and K⁺ channels, such as TRPA1, Piezo2, Kv1.1, Nav1.6, Nav1.7, Nav1.8, and the Ca²⁺-activated K⁺ channel, to AP generation in neurons innervating the trigeminal nerve (Tanaka et al., 2016), gastrointestinal tract (Chambers et al., 2014), urinary bladder (Mandge and Manchanda, 2018), muscle tissue (Nagaraja et al., 2021), and other nonspecific afferent neuron dorsal root ganglions (Amir and Devor, 2003; Baker, 2005; Tigerholm et al., 2014). However, none of these analyses incorporated neuronal sensitization. The computational models that do account for inflammation-induced sensitization (Suleimanova et al., 2020), including the model recently developed by our group (Nagaraja et al., 2023), only represent one inflammatory event and do not account for the priming of muscle nociceptors and its effect on the AP response to subsequent inflammatory events. In addition, with one exception (Nagaraja et al., 2021, 2023), none of the previous models focused specifically on the signaling of muscle nociceptors. Given the high degree of variability exhibited by neurons depending on the physiological tissue they innervate, computational models must incorporate the specific membrane proteins known to be present on the pertinent neuron type. Therefore, to understand the mechanisms of hyperalgesic priming of muscle nociceptors and its effect on their responses to subsequent inflammatory events, we must use a computational model that represents membrane proteins and intracellular mechanisms specific to muscle nociceptors as well as the pathways activated within them by commonly encountered inflammatory mediators.

In this study, we primarily investigated the key transmembrane proteins and intracellular signaling processes that regulate the response of a mouse muscle nociceptor to mechanical forces during two subsequent inflammatory events. To this end, we extended a previously developed and validated mathematical model of inflammation-induced sensitization in a mechanosensitive muscle nociceptor (Nagaraja et al., 2023) to incorporate the kinetics of an Epac-mediated pathway activated by inflammatory mediators (i.e., PGE₂) and a mechanosensitive membrane protein, TRPV4. The enhanced model represents 15 ion channels, two ion pumps, one ion exchanger, four endoplasmic reticulum (ER) membrane proteins, 29 intracellular components (including Ca²⁺ buffering proteins, kinases, enzymes, and second messenger molecules), and 46 associated processes. Upon model validation, we performed a global sensitivity analysis (GSA) by simulating the responses of 50,000 neurons to mechanical forces before and after the administration of an inflammatory mediator at two subsequent time points to quantify the contribution of the different modeled proteins and signaling processes to the change in AP firing magnitude after each inflammatory event. From this analysis, we first identified five ion channels and two processes as key regulators of the increase in neuronal AP firing in a primed muscle nociceptor (i.e., in the nociceptor after the second inflammatory event), and then investigated the effects of modifying these key regulators by introducing specific neuroplasticity after the first inflammatory event and evaluating the increase in the magnitude of neuronal AP firing after the second inflammatory event. These analyses allowed us to generate experimentally testable hypotheses regarding the role of ion channels and molecular processes in the priming of mouse muscle nociceptors and their response to mechanical stimuli.

Materials and methods

Computational model

We extended a previously developed and validated model of inflammation-induced sensitization of muscle nociceptors to include intracellular interactions that contribute to hyperalgesic priming. The existing model simulates the effects of inflammatory mediators, such as PGE₂, on neuronal responses via the independent activation of two distinct GPCR-mediated pathways, resulting in the activation of PKC and protein kinase A (PKA) within the neuron (Nagaraja et al., 2023). One of the GPCRs activates adenylyl cyclase, which catalyzes the activation of adenosine 3′,5′-cyclic monophosphate (cAMP). Under normal conditions, cAMP primarily activates PKA. However, after inflammation, cAMP also activates Epac, which in turn can activate PKC, thus providing feedback between the two GPCR-mediated signaling pathways (Hucho and Levine, 2007). To account for this feedback, we extended the model by adding descriptions of Epac kinetics (Figure 1A) and a new ion channel, i.e., TRPV4, based on evidence of its contribution to mechanical sensitivity in nociceptors (Ho et al., 2012; Mickle et al., 2015). The equations used to describe the Epac signaling and TRPV4 kinetics are provided in the Supplementary material.

FIGURE 1

Figure 1. (A) Implementation of two parallel inflammation-induced G protein-coupled receptor (GPCR) signaling pathways in our nociceptive muscle neuron model. Shown are the four modeled neuronal transmembrane proteins, Nav1.8, Nav1.7, Kv1.1, and TRPA1, whose activation and inactivation kinetics are modified by protein kinases C (PKC) and A (PKA), which are activated by the two prostaglandin E₂ (PGE₂)-initiated intracellular signaling pathways in the model. The arrows of the ion channels indicate the direction of flow of the ions through the corresponding channel. In the first pathway, phosphorylation of the GPCR activates subunits G_αq, β, and γ of the receptor. The G_αq subunit activates membrane-bound phospholipase C (PLC) and phosphatidylinositol 4,5-bisphosphate (PIP₂) to produce diacylglycerol (DAG), which in turn activates PKC. In the second pathway, phosphorylation of the GPCR activates subunits G_αs, β, and γ of the receptor. The G_αs subunit activates membrane-bound adenyl cyclase (AC) that activates cAMP, which in turn activates PKA. We also implemented the cAMP activation of Epac, which increases PKC activation (blue dashed arrow), thus providing feedback between the two inflammation-induced pathways. In addition, we added a mechanosensitive ion channel, TRPV4, previously not present in the model. Finally, we modeled the phosphorylation of Nav1.8, Nav1.7, Kv1.1, and TRPA1 by PKC and PKA, which modified their activation and inactivation kinetics. (B) Schematic showing the simulation inputs and outputs. In each simulation (one performed using the nominal parameter set and another using 50,000 distinct parameter sets generated for the global sensitivity analysis), we applied a step input of 5 or 100 mN mechanical force for 60 s at the simulation time point of 1 h. Next, at the 47- and 95-h simulation time points, we applied a step input of an inflammatory mediator (IM) at 100 nM for 30 min immediately followed by a 60-s step input of either 5 or 100 mN mechanical force. For each simulation that ran successfully, we used the time course of the membrane potential output to calculate the number of action potentials (APs) fired after the application of the mechanical force input once before (AP1) and twice after (AP2 and AP3) exposure to IM. Finally, we calculated the fold change in the number of APs fired after the first (FC1) and second (FC2) exposure to the IM, by dividing AP2 and AP3, respectively, by AP1. We classified the simulations where FC1 and FC2 were > 1 as primed neurons and those where FC1 and FC2 were ≤ 1 as non-primed neurons.

Model simulations, inputs, and outputs

Our extended model of inflammatory pain consists of 59 ordinary differential equations (ODEs) and 141 parameters. Each equation represents one variable in the model, where a variable represents activation or inactivation factors of 18 transmembrane proteins; the intracellular concentrations of K⁺, Na⁺, Ca²⁺, and inositol trisphosphate (IP₃); the endoplasmic reticulum (ER) Ca²⁺ concentration; the active and inactive subunits of the two GPCRs; three membrane-associated enzymes; the concentrations of 18 intracellular proteins and second messenger molecules; and the membrane potential (V_m). We provide a list of the model variables, their descriptions, and initial values in Supplementary Table 1. Using a lumped Hodgkin–Huxley-type formalism (Hodgkin and Huxley, 1952), we calculated changes in V_m at a given time point based on the changes in the currents of all neuronal transmembrane proteins listed above as follows:

$\frac{\text{d}{V}_{\text{m}}}{\text{d}t}=(I_{\mathrm{Nav1}.8}+I_{\mathrm{Nav1}.9}+I_{\mathrm{Nav1}.7}+I_{\text{Piezo}}+I_{\mathrm{ASIC3}}+I_{\text{TRPA}1}\hspace{0.17em}$

$\hspace{0.17em}+I_{\text{TRPV}4}+I_{\text{TREK}}+I_{\mathrm{Kv7}.2}+I_{\mathrm{Kv1}.1}+I_{\text{BKCa}}+I_{\text{Ka}}+I_{\text{Kleak}}$

$\hspace{0.17em}+I_{\text{CaT}}+I_{\text{CaL}}+I_{\text{PMCA}}+I_{\text{NaK}}+I_{\text{NCX}})/C{}_{\text{m}}(1)$

where C_m denotes the membrane capacitance and I represents the current through the different transmembrane proteins (described by the subscripts). The model simulates the time course of changes in the variables over 96 h in response to a mechanical force applied before and after two subsequent exposures to an inflammatory mediator (Figure 1B). We used 141 parameters to describe all the modeled mechanisms, and provide a list of the model parameter numbers (used to keep track of the parameters in our simulations), names, values, descriptions, units, and sources of the computational or experimental study from which we adapted or derived their values in Supplementary Table 2. We modified a subset of the model parameters (designated as “modified” in Supplementary Table 2) to match the inflammation-induced changes in AP firing magnitude from literature data (see “Model calibration and validation” section below). In all simulations, we maintained the extracellular concentrations of Na⁺, K⁺, and Ca²⁺ as well as the volume of the nociceptor nerve ending and its C_m at constant values (Nagaraja et al., 2023). We provide the ODEs and other equations describing all the modeled mechanisms, as well as the Nernst potentials and ionic balances for the intracellular concentrations of Na⁺, K⁺, and Ca²⁺ in the Supplementary material.

To drive the model, we provided as input a rectangular pulse with a mechanical force of either 5 or 100 mN (Figure 1B). We selected these values to represent inputs from both the innocuous and pain-inducing ranges of mechanical forces for mouse nociceptors (Cavanaugh et al., 2009; Deuis et al., 2017). We applied three pulses for a period of 60 s each, starting with the first pulse at the 1-h simulation time point. Next, at both the 47-h and 95-h simulation time points, we provided as input a 30-min rectangular pulse of an inflammatory mediator at a concentration of 100 nM, immediately followed by another 60-s rectangular pulse of 5 or 100 mN mechanical force. We chose the 30-min period to evaluate the inflammation-induced sensitization of the neuron’s response to mechanical forces based on literature data, which showed that the threshold reduction due to inflammation peaked between 30 min and 1 h in rat and mouse neurons (Hendrich et al., 2013). At the end of each simulation, our model generated a 96-h time course for each of the 59 model variables.

In all our computational analyses, we focused on the V_m time course, from which we calculated the total number of APs generated following the application of each pulse of mechanical force before and after the addition of the inflammatory mediator. We defined an AP as a V_m spike that crossed over 0 mV and reached a minimum value of +10 mV. We used the MATLAB function FINDPEAKS to identify the APs and to record their heights, widths, as well as the simulation time points at which they occurred. Next, we calculated two fold-change values (FC1 and FC2) for AP firing by dividing the total number of APs generated in response to the application of the force after the first and second inflammatory events (47 h and 95 h; Figure 1B), respectively, by the corresponding value before inflammation. The fold changes in the magnitude of AP firing following the two inflammatory events calculated using the nominal parameter set represented the baseline inflammation-induced sensitization. In addition to input forces of 5 and 100 mN, we calculated FC1 and FC2 values for force inputs of 10, 20, and 50 mN. We performed all computations in the software suite MATLAB R2015b (MathWorks, Natick, MA, USA) and solved the model equations using the MATLAB solver ODE15s with default tolerance levels.

Model calibration and validation

Model calibration

To ensure that our model accurately captured the inflammation-induced increase in neuronal excitability (i.e., AP firing) as measured experimentally, we calibrated a subset of 10 out of the 141 model parameters using experimental data from the literature (Bar et al., 2004). To perform the calibration, we modified the values of 10 parameters associated with Epac and TRPV4 activation and inactivation, and Epac-mediated PKC activation (designated as “modified” in Supplementary Table 2), such that the simulated fold changes in the magnitude of AP firing to mechanical forces of 20 and 100 mN after two subsequent exposures to the inflammatory mediator were within ± 1.96 standard error of the corresponding experimental measurements in rat spinal neurons (Bar et al., 2004). We defined the model’s “nominal parameter set” as the final parameter values obtained after performing the calibration procedure.

Model validation

To validate our model, we compared its predictions (using the nominal parameter set) with literature data not used for model calibration. First, we compared the predictions of fold change in the number of APs fired in response to each of two mechanical forces of 100 mN after an exposure to 100 nM PGE₂ at two subsequent time points with the corresponding experimental data from rat neurons (Hendrich et al., 2012). Second, using experimental data from the same study, we compared the predictions of the reduction in mechanical threshold (i.e., the minimum amount of force needed to elicit an AP) following the addition of 100 nM PGE₂ at the two time points.

Sensitivity analysis

We performed a sensitivity analysis to identify proteins and intracellular molecules and their associated signaling processes whose alterations were key for the development of hyperalgesic priming in the neuron and regulation of subsequent inflammation-induced changes in neuronal AP firing. First, we assessed the model’s robustness by performing a local sensitivity analysis (LSA), where we varied the model parameters near their nominal values (± 1%) as previously described (Nagaraja et al., 2014). Second, we performed a global sensitivity analysis (GSA) to account for the known heterogeneity in the expression of the various proteins at different nerve endings of muscle nociceptors. This analysis allowed us to account for the known variability in the conductance, activation, and inactivation gating factors of the same membrane proteins under different stimuli (Gold and Gebhart, 2010).

To this end, we simulated 50,000 distinct nociceptors and stimulated them with two inflammatory events followed by a mechanical force. First, we generated 50,000 unique parameter sets by randomly selecting parameter values from a fourfold range (twofold in each direction) around the nominal parameter values. To generate the random parameter sets, we used Latin hypercube sampling (MATLAB function LHSDESIGN) (Nagaraja et al., 2014). Next, we performed two sets of simulations using the 50,000 parameter sets, where we drove each simulation using a mechanical force of either 5 or 100 mN applied once before and twice 30 min after the addition of an inflammatory mediator at the 47- and 95-h simulation time points (Figure 1B). We interrupted and eliminated the time course simulations of V_m that did not reach the 96-h time point within 15 min of computation time (wall-clock) or that required time steps smaller than 1 × 10^–12 s. We used this lack of convergence in the simulations to flag parameter sets that resulted in non-physiological kinetic behavior. Accordingly, we only used the simulations that ran to completion to calculate the fold changes FC1 and FC2 described above in the number of APs fired (in response to either 5 or 100 mN force) after each of the two inflammatory events. We used the calculations of APs, FC1, and FC2 to further eliminate the simulations in each force group where the number of APs fired in response to a mechanical force in the absence or presence of an inflammatory mediator was zero. Finally, we eliminated the simulations in which the AP fold-change values, FC1 and FC2, after each of the two inflammatory events were equal to one. We assumed that these simulations represent “non-nociceptive neurons” because their responses to noxious mechanical forces did not change after two consecutive inflammatory events, i.e., they did not have the ability to sensitize. Using the results from the remaining simulations, we performed two analyses, a partial rank correlation coefficient analysis and a parameter distribution analysis.

Partial rank correlation coefficient analysis

For this analysis, we calculated the Spearman’s partial rank correlation coefficient (PRCC) and the associated p-values between the primary outputs (i.e., the AP firing fold changes after each inflammatory event, FC1 and FC2) and each of the 141 model parameter values. The values of the PRCC varied between −1 and +1, with large absolute values reflecting a high impact of the model parameter on the model output (i.e., FC1 and FC2). The sign of the PRCC indicated the positive or negative directionality of the correlation between the model parameter and the output. A PRCC with a p-value < 0.01 indicated that it was significantly different from zero. Upon completion of this analysis, we obtained two sets of 141 PRCC values for each force input (i.e., 5 and 100 mN), along with their associated p-values.

Parameter distribution analysis

For this analysis, we first separated the simulations for each of the two input forces of 5 and 100 mN into two groups. We defined a group of simulations in which both AP fold changes FC1 and FC2 were > 1 as “primed” neurons and a group in which both FC1 and FC2 were ≤ 1 as “non-primed” neurons. Next, we generated histograms of the parameter value distributions for each of the model’s 141 parameters using the MATLAB function HIST, with 50 bins partitioning the interval between the minimal and maximal values for each model parameter in the two groups of simulations (i.e., primed and non-primed neurons). We calculated the percentage of the simulations for each distribution curve by dividing the number of simulations in which a given parameter’s value fell within the range of a bin by the total number of simulations in that group. Then, for each model parameter, we quantified the area of overlap between the sensitized and non-sensitized neuron group distributions by calculating the Bhattacharyya coefficient, which varied between 0 and 1, representing no and 100% overlap, respectively, as previously described (Mitrophanov et al., 2015). A small overlap area indicated that a parameter (and the protein it represents) was consistently over (or under) expressed in a primed neuron relative to a non-primed neuron and was therefore more likely to be associated with inflammation-induced neuronal hyperalgesic priming than parameters with larger distribution overlap areas.

We utilized the results from the PRCC and the parameter distribution analyses to identify key transmembrane proteins whose alterations could regulate the hyperalgesic priming in muscle nociceptors after inflammation. Using the results from the PRCC analysis, for each force input, we divided the set of 141 PRCCs into five clusters using a k-means clustering algorithm (MATLAB function KMEANS) (Nagaraja et al., 2017). We considered the model parameters in the cluster that had the highest absolute PRCC values and had p-values ≤ 0.01 as key contributors to hyperalgesic priming induced after inflammation. Using the results from the parameter distribution analysis, we first ranked the absolute values of the 141 Bhattacharyya coefficients in ascending order and designated the parameters within the top five lowest values as key for hyperalgesic priming of the neuron. Finally, we combined the model parameters identified as key in both analyses and labeled the proteins/molecules or the intracellular signaling processes represented by those parameters as key for AP-response regulation during multiple inflammatory events.

Ion channel modification simulations

For each transmembrane protein identified as key for FC1 and FC2 changes, we performed simulations to mimic specific modifications associated with that protein which are likely to occur following the first inflammatory event that “primes” its response to future inflammatory events. In these simulations, we either reduced or increased the expression of a transmembrane protein at 48 h into the simulation (i.e., after the first inflammatory event). To simulate an increase or decrease in a protein’s expression, we multiplied or divided, respectively, the current in Eq. (1) corresponding to that protein by a factor of 2. First, we simulated the protein modifications described above using a model with the nominal parameter set, which represented an average nociceptive muscle afferent neuron. Then, to verify that we could reproduce the effects of different modifications in a population of neurons, we repeated the simulations for every modification using parameter sets from simulations that were classified as primed and non-primed neurons. As in the GSA, we interrupted the simulations that did not reach the 96-h time point of the V_m time course within 15 min of computation time (wall-clock) to flag parameter sets where a modification resulted in non-physiological kinetic behavior. Of the simulations that converged successfully, we calculated the mean ± 1 standard error (SE) of the AP fold changes (FC1 and FC2) after each inflammatory event for every modification. Finally, we performed a Wilcoxon rank sum test to compare the mean values of each simulation with a modification to the corresponding simulation without any modification.

Results

The model captured inflammation-induced changes in the activation threshold and in the AP firing response to mechanical stimulation

To ensure that the model accurately captured the effect of multiple inflammatory events on the magnitude of the AP firing response to mechanical forces, we calibrated the model using data obtained from electrophysiological measurements in rat spinal neurons with peripheral terminals in the knee joint (Bar et al., 2004). The calibration procedure resulted in fold changes in the number of APs fired by the neuron in response to an innocuous (i.e., 20 mN; Figure 2A) and a noxious (100 mN; Figure 2B) mechanical force input after two subsequent applications of inflammatory mediators to fall within one SE of the experimental data. We designated the final set of model parameter values obtained after this calibration procedure as the nominal parameter set. Figure 2C shows the simulated APs for the case where we stimulated the neuron by a 20 mN force before (i.e., at the 1-h time point) and at 30 min after the addition of PGE₂ at the 47- and 95-h simulation time points. In addition, based on the results of the LSA (performed using the nominal parameter set), we found that V_m was not overly sensitive (sensitivity indices > 100) to any of the model’s 141 parameters, suggesting that the model was stable and robust to small perturbations (± 1%) of its nominal values.

FIGURE 2

Figure 2. Model calibration. We calibrated the model’s predictions of the increase in action potential (AP) firing magnitude by fitting them to experimental data from rat spinal cord neurons in response to mechanical forces applied after the addition of an inflammatory mediator. Open bars (N = 6) show mean experimental data ± 1 standard error fold change in the magnitude of AP firing in response to innocuous (A) and noxious forces (B) after two subsequent exposures to inflammatory mediators (Bar et al., 2004). Solid bars show the results of model fitting to the experimental data. We used innocuous and noxious forces of 20 and 100 mN, respectively, to stimulate the neuron in the model. (C) Simulated time course trajectory of the membrane potential output in response to a 20 mN force applied before (at the 1-h time point) and after two subsequent inflammatory events (at the 47- and 95-h time points). The inset shows the shape and number of APs generated in each response.

To validate the model, we simulated the AP response (based on the nominal parameter set) to a mechanical force of 100 mN applied before (i.e., at the 1-h time point) and at 30 min after the administration of 100 nM PGE₂ at the 47- and 95-h simulation time points and determined the number of APs fired after each stimulus (Figure 3A). We first compared the model predictions of the fold changes (FC1 and FC2) with the corresponding experimental data derived from rat gastrocnemius muscle neurons (Figure 3B). In addition, we compared the model predictions of the percentage reduction in the mechanical threshold (i.e., the minimum force required to elicit an AP from a neuron) 30 min after the addition of 100 nM PGE₂ at the two time points in the simulation with the corresponding data in rat gastrocnemius muscle neurons (Hendrich et al., 2013). In the model, before the addition of the inflammatory mediator, the minimum force needed to elicit an AP was 1 mN. In both comparisons, our predictions were within ± 1.96 SE of the validation experimental data for all the comparisons, indicating that we cannot differentiate between the average results of the experiment and the model predictions (Figure 3C, error bars show 1 SE of experimental data).

FIGURE 3

Figure 3. Model validation. We validated the model by comparing the simulations of inflammation-induced action potential (AP) firing increase and mechanical threshold reduction with the corresponding experimental data. (A) Simulated time course trajectory of the membrane potential output in response to a 100 mN force applied before (at the 1-h time point) and after two subsequent inflammatory events (i.e., at the 47- and 95-h time points). The inset shows the number of APs generated and their magnitude. (B) Mean experimental data ± 1 standard error (SE) fold change in the number of APs fired in response to a noxious force after exposure to two subsequent inflammatory events in rat gastrocnemius muscle neurons (open bars, N = 10) (Hendrich et al., 2013). (C) Mean experimental data of the percentage reduction ± 1 SE in the mechanical threshold induced by two subsequent exposures to 100 nM PGE₂ (open bars, N = 6). Solid bars in (B,C) show the corresponding model predictions. (D) Model predictions of the AP firing fold change in response to five different force stimulations after two inflammatory events. The fold changes after the first inflammatory event are referred to as FC1 and those after the second inflammatory event as FC2.

Upon validation, we used the model to establish the baseline inflammation-induced increase in AP firing after two inflammatory events. The AP fold-change values for the individual forces after the first inflammatory event (FC1) were 11 for 5 mN, 10 for 10 mN, 8 for 20 mN, 2.1 for 50 mN, and 1.8 for 100 mN (Figure 3D, solid bars). The AP fold-change values for the individual forces after the second inflammatory event (FC2) were 66 for 5 mN, 43 for 10 mN, 21 for 20 mN, 1.75 for 50 mN, and 3.5 for 100 mN (Figure 3D, open bars). The fold-change values after either inflammatory event decreased as the intensity of the mechanical force stimulus increased.

For the next set of simulations, we selected 5 and 100 mN to represent a mechanical force input from the innocuous and pain-inducing range, respectively, for mouse nociceptors. In addition to simulating two inflammatory events using the nominal parameter set, we repeated the simulations using 50,000 virtual neurons. Of the 50,000 simulations performed for both the 5 and the 100 mN mechanical force inputs, 48,588 and 48,651, respectively, ran successfully. After two additional filtering procedures (see section “Sensitivity analysis” in “Materials and methods”), we extracted 4,310 and 4,027 simulations from the 5 and 100 mN groups, respectively, to use for further analysis. Lastly, in both groups, we classified the simulations as either “primed” or “non-primed” neurons. For the 5 mN simulation set, we identified 3,191 and 1,119 simulations as primed and non-primed neurons, respectively, and for the 100 mN simulation set, we identified 2,101 and 1,926 simulations as primed and non-primed neurons, respectively. Figures 4A, C shows the number of APs fired before (pink circles) and after the first (green circles) and second (blue circles) inflammatory events for the primed neurons in both the 5 and 100 mN groups, where the insets show the average APs fired (horizontal black bars). For both the 5 and 100 mN force inputs, the average fold changes in the number of APs fired in response to the second inflammatory event compared to those fired in the absence of inflammation (FC2; Figures 4B, D, solid bars) were higher than those for the first inflammatory event (FC1; Figures 4B, D, open bars). Overall, these results show that our simulations were able to generate a population of neurons in which hyperalgesic priming occurs in response to both innocuous and noxious stimuli.

FIGURE 4

Figure 4. Baseline inflammation-induced sensitization and hyperalgesic priming of muscle nociceptors. Shown are the number of action potentials (APs) fired in response to a mechanical force of (A) 5 mN and (C) 100 mN applied before (pink circles) and after the first (green circles) and second (blue circles) exposures to an inflammatory mediator. The insets show the mean (horizontal black lines) number of APs fired. (B,D) show the corresponding mean ± 1 standard error (SE) of the AP fold-change values after the first (FC1) and second (FC2) inflammatory events in response to 5 and 100 mN mechanical forces, respectively.

Key proteins and processes of inflammation-induced AP response regulation

To identify the transmembrane proteins that strongly regulated the AP response (specifically the number of APs generated) following the addition of an inflammatory mediator to a pre-inflamed muscle nociceptor across many different nociceptor-signaling conditions, we used two distinct analyses (PRCC and parameter distribution; see section “Sensitivity analysis” in “Materials and methods”). We used the AP fold-change values (FC1 and FC2) and the respective parameter values for each simulation group to perform the PRCC and parameter distribution analyses. For the simulations in which we used a mechanical force of 5 mN as input, the PRCC analysis results showed that for the primed neurons, the model parameters associated with Kv7.2 yielded high and statistically significant correlations (p < 0.01) with the FC1 values (Figure 5A, solid black bar), and the model parameters associated with proteins NaK and TRPA1 as well as molecular processes of G_αq-coupled receptor phosphorylation and phosphorylation of Nav1.8 and Nav1.7 yielded high and statistically significant correlations (p < 0.01) with the FC2 values (Figure 5B, solid black bars). Moreover, for the non-primed neuron simulations, the PRCC analysis results showed that the model parameters associated with proteins NaK, Piezo2, and Nav1.7 yielded high and statistically significant correlations (p < 0.01) with the FC1 values (Supplementary Figure 1A, solid black bars), and the model parameters associated with proteins TRPA1 and Nav1.7 and G_αq-coupled receptor phosphorylation yielded high and statistically significant correlations (p < 0.01) with the FC2 values (Supplementary Figure 1B, solid black bars).

FIGURE 5

Figure 5. Partial rank correlation coefficient (PRCC) analysis identified key proteins and processes for action potential (AP) regulation. The bars show the PRCCs of the 141 model parameters with fold changes in the total number of APs generated after the first (FC1) and second (FC2) inflammatory events in response to a mechanical force of 5 mN (A,B) and 100 mN (C,D) in the primed neuron simulations. The PRCCs above their respective thresholds (dotted horizontal lines) that were statistically significant (i.e., p < 0.01; p-values depict the probability of seeing the observed correlation if no correlation exists) are indicated by solid black bars, and the labels of the bars show the ion channels/ion pumps or the rates of intracellular processes that these parameters describe in the model. We used the Spearman’s rank correlation method to compute the PRCC values, and obtained the p-values of the correlation using a Student’s t-test. Because the data were not normally distributed, we used a large-sample approximation while performing this test. The analyses are based on 2,438 and 1,605 simulations classified as primed neurons for applied mechanical forces of 5 and 100 mN, respectively.

For the simulations in which we used a mechanical force of 100 mN as input, the PRCC analysis results showed that the model parameters associated with Kv7.2 yielded high and statistically significant correlations (p < 0.01) with the FC1 values (Figure 5C, solid black bars), and the molecular processes of G_αq-coupled receptor phosphorylation and activation, phosphatidylinositol 4,5-bisphosphate (PIP₂) activation, and Nav1.8, and Nav1.7 phosphorylation by PKA and PKC yielded high and statistically significant correlations (p < 0.01) with both the FC1 and FC2 values (Figures 5C, D, solid black bars). For the non-primed neuron simulations, the PRCC analysis results showed that the model parameters associated with TRPA1 yielded high and statistically significant correlations (p < 0.01) with the FC1 values (Supplementary Figure 1C, solid black bar), and the model parameters associated with the activation of proteins Nav1.8 and Nav1.7 yielded high and statistically significant correlations (p < 0.01) with the FC2 values (Supplementary Figure 1D, solid black bars).

In the parameter distribution analyses, we calculated the 141 Bhattacharyya coefficients to determine the overlap between the distributions of the normalized values of the model’s 141 parameters in the primed and non-primed neuron groups for both the 5 and 100 mN simulations. While none of the parameters had a considerably small overlap area between the two groups of 141 parameters, the parameters that demonstrated the five lowest values in the 5 mN group were associated with activation or inactivation of ion channels Nav1.7, Kv1.1, and TRPA1 and the rates of Nav1.8 and Nav1.7 phosphorylation by PKC and PKA post-inflammation. In the 100 mN simulation group, the parameters that demonstrated the five lowest values were associated with activation or inactivation of ion channels TRPA1, TRPV4, and Nav1.7 and the rates of TRPA1 phosphorylation by PKC and PKA post-inflammation. Supplementary Tables 3, 4 provide a list of the 141 Bhattacharyya coefficients for all model parameters for the 5 and 100 mN simulation groups, respectively. Figure 6 shows representative examples of two such parameters from each simulation set. For the parameters representing Nav1.7 activation (in the 5 mN force simulation group) and TRPA1 activation (in the 100 mN force simulation group), a larger percentage of the simulations in the primed neuron group fell in the lower range of their normalized values compared to those of the non-primed group (Figures 6A, C, solid vs. dashed lines), indicating that Nav1.7 and TRPA1 activation might be downregulated in primed neurons. Conversely, for the parameters representing Kv1.1 activation (in the 5 mN force simulation group) and Nav1.7 inactivation (in the 100 mN force simulation group), a larger percentage of the simulations in the primed neuron group fell in the higher range of their normalized values compared to those in the non-primed neuron group (Figures 6B, D, solid vs. dashed lines), indicating that these processes might be upregulated in primed neurons. Finally, we combined the results of both analyses and identified five ion channels (Kv7.2, Kv1.1, Nav1.7, Nav1.8, and TRPA1) and two molecular processes (G_αq-coupled receptor phosphorylation and Nav1.7 and Nav1.8 phosphorylation) whose modifications could potentially regulate hyperalgesic priming of mouse muscle nociceptors.

FIGURE 6

Figure 6. Parameter distribution analysis identified key proteins and processes for action potential (AP) regulation. Shown are the distributions of parameter values across the simulations of primed (solid lines) and non-primed (dashed lines) neuron groups representing (A) Nav1.7 activation and (B) Kv1.1 activation with a 5 mN mechanical force used as input, and (C) TRPA1 activation and (D) Nav1.7 inactivation with a 100 mN mechanical force used as input. The x-axis indicates the normalized parameter values, and the y-axis represents the percentage of simulations in each neuron group in which the parameter values fell within a particular range (described in the “Materials and methods” section).

Analysis of modification of model-identified key proteins and molecular processes

To quantify the extent to which different inflammation-induced neuroplastic modifications in membrane proteins affect subsequent neuronal response, we performed simulations where we modified each model-identified key protein, one at a time, after the first inflammatory event in the primed and non-primed neuron groups for both the 5 and 100 mN force inputs. Specifically, for the three proteins Nav1.7, Nav1.8, and TRPA1, we simulated the effect of their overexpression, and for Kv1.1 and Kv7.2, we simulated the effect of reducing their expression on AP firing. We used these specific modifications based on our model results and literature evidence (Nicol et al., 1997; Blair and Bean, 2002; Kwan et al., 2009; Hameed, 2019). We also performed simulations where we combined all protein modifications to occur simultaneously. We induced these modifications after the first inflammatory event, i.e., at the 48-h simulation time point. We then compared the mean AP fold changes (FC1 and FC2) due to each of these modifications to the corresponding value in the simulations with no modifications (i.e., the baseline sensitization) (Figure 7).

FIGURE 7

Figure 7. In silico analysis identified the relative contributions of specific neuroplastic changes in model-identified key ion channels to action potential (AP) generation in primed and non-primed neurons. We simulated an individual twofold expression increase in Nav1.7, Nav1.8, and TRPA1; an individual twofold expression decrease in Kv7. 2 and Kv1.1; as well as simultaneous modification of Kv7.2, Nav1.7, Nav1.8, and TRPA1 (gray bars) after the first exposure to an inflammatory event in primed and non-primed neuron sets identified during global sensitivity analysis. We compared the means and one standard error (SE) of the AP fold changes (FC1 and FC2) in response to 5 and 100 mN forces, after two subsequent inflammatory mediator exposures derived from simulations implementing the six modifications involving the individual and combined increase or decrease in expression of proteins with corresponding simulations with no modifications (i.e., the solid bars) in the (A,B) primed and (C,D) non-primed neuron groups. Solid bars in all panels indicate the means and one SE of the magnitude of AP fold change in simulations with no modification. Because the data from all the groups were not normal (which we established by performing a Kolmogorov–Smirnov test), we used a Wilcoxon rank sum test, which is typically used to compare the means of two independent samples when we cannot assume normality, to determine if any of the protein modification significantly changed AP firing. An asterisk (*) indicates that the mean AP fold change due to a particular modification was significantly different from the simulations with no modification, with p ≤ 0.01.

Of the simulations with protein modifications, we identified 2,438 and 868 that ran successfully in the primed and non-primed neuron groups, respectively, when 5 mN force was used as input. When 100 mN was used as the force input, we identified 1,605 and 1,800 simulations that ran successfully in the primed and non-primed neuron groups, respectively. When either an innocuous (i.e., 5 mN) or a noxious (i.e., 100 mN) force was used as input, in the set of primed neurons (i.e., neurons in which subsequent inflammatory events increased AP firing without any modification), modifying TRPA1 (Figures 7A, B, vertically dashed bar for FC2), Kv1.1 (Figures 7A, B, horizontally dashed bar for FC2), and Kv7.2 (Figures 7A, B, white bar for FC2) expression as well as the combined modification of all four key proteins (Figures 7A, B, gray bar for FC2) significantly increased the average AP fold change after the second inflammatory event compared to the simulations with no modifications (Figures 7A, B, black bar for FC2). Of note, while the simulations captured the extent to which distinct neuroplastic modifications induced by hyperalgesic priming affected the neuron’s response to subsequent inflammatory events, the effect was quantitatively higher (for all modifications) when an innocuous force was used as input compared to a noxious force (Figures 7A, B, FC2 values).

In the set of non-primed neurons (i.e., neurons in which subsequent inflammatory events decreased AP firing), when either an innocuous force (i.e., 5 mN) or a noxious force (i.e., 100 mN) was used as input, all the implemented neuroplastic modifications significantly increased the average AP fold change after the second inflammatory event (Figures 7C, D) compared to the simulations with no modifications (Figures 7C, D, black bar for FC2). With the exception of Kv1.1 and Kv7.2 expression reduction (Figure 7D, open and horizontally striped bars for FC2), as observed in the primed neuron group, the AP fold changes induced by the different protein modifications after exposure to the second inflammatory event were higher for the innocuous force input (Figure 7C) than for the noxious force (Figure 7D). Overall, based on our simulation results, hyperalgesic priming involving Nav1.7 and Kv1.1 neuroplasticity (compared to Nav1.8, Kv7.2, and TRPA1) appeared to have the largest effect on the sensitivity of muscle nociceptors to mechanical forces during multiple inflammatory events.

Discussion

Hyperalgesic priming in peripheral nociceptors causes them to become hyperexcitable to subsequent injuries or inflammatory events when they have experienced a prior sensitization event, thus preserving a “memory” of the prior event (Woolf and Ma, 2007; Reichling and Levine, 2009). Hyperalgesic priming could potentially lead to the transition of acute to chronic pain, however, its underlying mechanisms are not well understood (Parada et al., 2005; Kandasamy and Price, 2015) although such knowledge might hold the key to developing therapies to prevent or reverse pain chronification. Inflammatory mediators present during injury sensitize nociceptors by decreasing their activation thresholds and increasing their responsiveness (i.e., AP firing) to stimuli. This increase reflects changes in the expression and function of the nociceptor’s membrane proteins. However, we still do not know which ion channels, pumps, and exchangers are key for the generation of AP responses after multiple inflammatory events. In this study, we used computational analysis to identify such key proteins and molecular processes. We first extended our validated, multi-compartment computational model of a mouse muscle nociceptor to incorporate Epac-mediated PKC activation (Parada et al., 2005) and the kinetics of a mechanosensitive ion channel, TRPV4. The extended model accounts for the activity of 18 membrane proteins, four ER membrane proteins, and 29 second-messenger molecules, including proteins, molecules, kinases, as well as Na⁺, K⁺, and Ca²⁺ ions, and describes 46 intracellular processes, including the activation and inactivation of the various membrane proteins as well as intracellular proteins and molecules by inflammatory mediators. We then calibrated and validated the model using experimental data capturing the effects of two subsequent inflammatory events on the AP responses in rat nociceptors. In agreement with experimental observations, our simulations of two subsequent inflammatory events predicted an increase in AP firing and a reduction in a nociceptor’s mechanical activation threshold with each inflammatory event (Figures 3B, C).

To identify key regulators of the AP response, we used the model to simulate pain signaling responses to mechanical forces and two inflammatory events in 50,000 unique muscle nociceptors intended to represent the heterogeneity in protein expression and activity and the numerous plausible protein interactions that can occur in vivo. We found that modifications to five ion channels (Kv7.2, Kv1.1, Nav1.7, Nav1.8, and TRPA1) and two molecular processes (G_αq-coupled receptor phosphorylation and Nav1.8 and Nav1.7 phosphorylation) strongly regulated the increase in the total number of APs fired by the neurons in response to both an innocuous (5 mN) and noxious (100 mN) mechanical force during inflammation. Moreover, by separately simulating the specific modifications in the expression of each of the five proteins, we showed that modifying the expression of Nav1.7 or Kv1.1 after the first inflammatory event significantly regulated the magnitude of AP firing. Therefore, Kv1.1 and Nav1.7 channels (specifically modification of their expression or activation) could be considered as potential targets for regulating the increased AP response of nociceptors due to hyperalgesic priming that occurs during musculoskeletal trauma.

Key membrane proteins and molecular processes that regulate AP generation in primed neurons

Because a majority of the pathways activated by distinct inflammatory mediators converge within the neuron to increase the intracellular concentrations of PKC and PKA, various membrane-bound and intracellular molecules and proteins involved in PKC and PKA signaling have been investigated for their roles in hyperalgesic priming (Ferrari et al., 2013b; Gangadharan and Kuner, 2013; Araldi et al., 2016a,b; Fang et al., 2021). In fact, numerous behavioral studies have demonstrated that activation of PKCε can consistently induce hyperalgesic priming in rats after multiple exposures to inflammatory mediators. The studies demonstrated that the generation of hyperalgesic priming in rat nociceptors depended on proteins downstream of PKC activation, such as cytoplasmic polyadenylation element-binding protein (CPEB), a protein responsible for translation of mRNA in the peripheral terminals of nociceptors. Furthermore, inhibiting downstream targets of CPEB (e.g., α calmodulin-dependent protein kinase II) or inhibiting nuclear transcription factors (e.g., via the administration of an oligodeoxynucleotide antisense) attenuated priming in rats (Aley and Levine, 1999; Reichling and Levine, 2009; Ferrari et al., 2013a,b, 2015; Araldi et al., 2016a). While we do not explicitly model molecules downstream of PKC activation that are involved in the mRNA translation of different membrane proteins, our PRCC analysis identified molecules that are precursors to PKC activation, i.e., G_αq-coupled receptors phosphorylation and PIP₂ activation, as key for increased AP firing in sensitized muscle neurons during a second inflammatory event (Figures 5C, D, black bars), further demonstrating the importance of PKC signaling in inflammation-induced neuronal response changes. However, because a behavioral pain response (i.e., paw withdrawal) was used as an indicator of the change in the neuronal response induced by hyperalgesic priming in the above-mentioned experimental studies, we still do not have electrophysiological correlates to this observed reduction in mechanical input thresholds. To develop new therapies that target the peripheral neurons, we need to understand how their response (i.e., AP firing) as well as the expression and function of key membrane proteins change during multiple inflammatory events.

Many studies that do measure electrophysiological responses of nociceptors have focused on the ability of PKC and PKA to evoke a change in the magnitude of neuronal AP firing by modifying the gating properties of key Na⁺ and K⁺ ion channels (already present on the neuronal membrane), such as Nav1.8 and Nav1.7 (Vijayaragavan et al., 2004; Wu et al., 2012; Tanaka et al., 2016), Kv1.1 (Nicol et al., 1997; King et al., 2014; D’Adamo et al., 2020), and TRPA1 (Kwan et al., 2006; Karashima et al., 2009; del Camino et al., 2010), among others, via phosphorylation. Therefore, a number of studies have investigated Nav1.8 and Nav1.7 modulation in inflammatory pain. For example, nerve growth factor (NGF)-induced thermal hyperalgesia was inhibited in Nav1.8 knockout (KO) mice compared to saline-treated mice (Kerr et al., 2001), and Nav1.8 inhibition via a pharmaceutical blocker showed reduced excitability in inflammatory and neuropathic rat models of pain (Payne et al., 2015). In addition, nociceptor-specific Nav1.7 KO mice showed an increased thermal and mechanical pain threshold and a decrease in AP response to various inflammatory mediators, including complete Freund’s adjuvant, PGE₂, and NGF (Nassar et al., 2004). However, all these studies reported changes in neuronal responses due to modification of ion channels after only one inflammatory event. Therefore, we do not know if and for how long the inflammation-induced changes in the expression and function of these channels last, because none of the studies measured the neuronal responses after a second inflammatory event.

Our analysis allowed us to characterize the neuronal AP responses after two inflammatory events in two distinct scenarios, one where we pre-exposed a group of neurons to an innocuous force (5 mN) and another where we pre-exposed them to a noxious force (100 mN), both after the administration of an inflammatory mediator. We found that, while the number of APs fired were, on average, higher for the 100 mN force input compared to 5 mN, the AP fold changes, especially after the second inflammatory event, were higher when the neurons were stimulated by an innocuous (i.e., 5 mN) force input, suggesting that the effect of hyperalgesic priming on the increase in the neurons’ excitability was higher when the initial force inputs were innocuous. Furthermore, our ability to simulate two inflammatory events led to the identification of proteins and processes key for initiating hyperalgesic priming (i.e., after the first inflammation event) and for regulating it. For example, in agreement with previous electrophysiological studies, our analysis identified Nav1.7, TRPA1, and Kv1.1 as key in driving the increase in the magnitude of AP firing after a single inflammatory event (Figures 5A, C). In addition to these three proteins, our analysis also identified Nav1.8 and Kv7.2 as well as molecular processes of phosphorylation of G_αq-coupled receptor, Nav1.8, and Nav1.7 as key for the increase in AP firing after a second inflammatory event (Figures 5B, D), which current experimental studies did not measure. Overall, our computational analyses identified five proteins that were key for the increase in AP firing to mechanical forces after two subsequent inflammatory events.

Relative contributions of model-identified key proteins/processes

While previous studies and our computational analysis have highlighted the role of individual ion channels or protein kinases in regulating pain responses during multiple inflammatory events (Nicol et al., 1997; Baker, 2005; Kwan et al., 2009; Gold and Gebhart, 2010; King et al., 2014; Nagaraja et al., 2021), we were able to compare the relative effect of modifying different proteins and kinases, one at a time, on AP firing in the same set of 50,000 simulated neurons, which is not feasible in experimental studies. For example, in a group of 2,438 distinct primed neurons, with no specific modification implemented, the average inflammation-induced increase in AP fold change in response to 5 mN force was 3.1 (SE = 0.18) and 11.1 (SE = 1.03) after the first and second inflammatory events, respectively. In the same neurons, five of the six distinct modifications we implemented (involving the proteins identified as key by our GSA), one at a time, significantly altered the average AP fold change after the second inflammatory event (Figure 7A). Of these five modifications, three increased the AP fold change by greater than 70%. Specifically, a twofold decrease in Kv1.1 expression, a twofold increase of Nav1.7 expression, and a combined modification of Kv7.2, Nav1.7, Nav1.8, and TRPA1 increased the average AP fold change by 104, 70, and 77%, with mean FC2 values of 22.6 (SE = 1.7), 18.5 (SE = 1.7), and 19.4 (SE = 1.2), respectively (Figure 7A, open, diagonally dashed, and gray bars). We saw similar results when we simulated the modifications in another group of 1,605 primed neurons using 100 mN force as input (Figure 7B). We can use this information while designing new strategies to regulate the effect of hyperalgesic priming. For example, since the modifications of both Kv1.1 and Nav1.7 individually had a significant effect on AP fold change compared to other protein modifications, we could test them as potential targets to regulate the degree of hyperalgesic priming and select the one with the least side effects.

We also simulated the neuroplastic modifications in the non-primed neuron group. We assumed that these simulations represented silent nociceptors, or non-responding sensory neurons, and aimed to determine whether there were any specific neuroplastic modifications initiated in these neurons by the first inflammatory event that could activate them to sensitize during a second inflammatory event. For example, a group of 1,800 neurons that initially (with no modifications) showed no increase in mean AP fold change after the first or second inflammatory events [Figure 7D, black bars, FC1 and FC2 of 0.50 (SE = 0.004) and 0.45 (SE = 0.004), respectively] became sensitized (or acquired nociceptive properties) after the second inflammatory event when we decreased Kv7.2 or Kv1.1 expression after the first inflammatory event. [Figure 7D, open and horizontally dashed bars, FC2 of 0.75 (SE = 0.20) and 0.74 (SE = 0.17), respectively]. Thus, in addition to identifying a panel of key proteins and molecules that could be potential targets for regulating the neuronal excitability after inflammation, our analysis also provided a one-to-one comparison of the efficacy of targeting each key protein in the same population of neurons.

Assumptions and limitations

Our computational model has several limitations arising from simplifying assumptions required to capture the complex nature of hyperalgesic priming in muscle afferent neurons. First, due to limited availability of electrophysiological data for the effects of multiple inflammatory events on AP firing in mice, we used data from rat neurons to calibrate and validate the model, which might affect model accuracy. Second, because we only represent two specific GPCR signaling pathways activated by a subset of the inflammatory mediators that can be present in an injured tissue, it is possible that proteins and molecules involved in other pathways not currently accounted for might be important for regulating hyperalgesic priming. However, as previously discussed, many inflammatory pathways are known to converge within the neuron, resulting in an increase in the concentrations of PKC and PKA (Gold and Flake, 2005; Gold and Gebhart, 2010; Gangadharan and Kuner, 2013; Pak et al., 2018), whose effects, including Epac-mediated PKC activation, we do currently model. Therefore, we could incorporate the effects of other inflammatory mediators if and when relevant. Third, in our model, we adopted many parameter values from previous computational studies developed to describe neurons from animals other than mice or from physiological tissues other than muscle (Bennett et al., 2005; Lindskog et al., 2006; Mandge and Manchanda, 2018). While we performed a validation procedure to match our computational simulations to experimental data recorded from rat neurons, we did not directly derive parameters from single ion channel current measurements in mouse neurons. This simplification could impact the accuracy of certain model parameters. Fourth, we chose a 48-h period between the two subsequent inflammatory events based on experimental evidence that shows that the reduction in the mechanical threshold due to the addition of an inflammatory mediator typically recovers in ∼24 h (Hendrich et al., 2013). However, behavioral studies have demonstrated that the effect of priming can persist for up to 1 week (Ferrari et al., 2013b,2015; Araldi et al., 2016a,b). Therefore, the mechanisms underlying the priming effect might differ based on the timing of subsequent inflammatory events in the primed neurons, which we do not currently consider but could potentially account for in future analyses. Finally, our hypotheses regarding the contributions of Kv7.2, Kv1.1, Nav1.7, Nav1.8, and TRPA1 and both Nav1.8 and Nav1.7 phosphorylation to the hyperexcitability of muscle nociceptors after two inflammatory events stem solely from simulations. These hypotheses need to be tested using experiments involving electrophysiological (patch clamp or dynamic) recordings in mice, and potentially even using human induced pluripotent stem cells, where we separately modify each protein or process and assess the effect of the modification on AP firing in response to multiple, sequential inflammatory events. Ultimately, there is always the question of translatability of the hyperalgesic mechanisms across species. Since we cannot perform in vivo investigations on human nociceptors, mouse models provide the opportunity to use genetic approaches to investigate the molecular mechanisms of nociceptive signaling. In addition, both rodents and humans are known to have the ability to sensitize following repetitive injuries and similar functional organization of the spinal cord, making them a useful substitute for human nociceptors (Fitzgerald, 2005; Toossi et al., 2021).

Conclusion

Identification of hyperalgesic priming mechanisms in muscle nociceptors is challenging due to the heterogeneity in muscle neuron subpopulations and the plethora of inflammatory mediators that can act upon them. In this study, we specifically focused on the effect of two subsequent PGE₂-induced inflammatory events on mechanical nociception in muscle tissue. The results of our computational analyses allowed us to hypothesize that (1) Kv7.2, Kv1.1, Nav1.7, Nav1.8, and TRPA1 as well as phosphorylation of G_αq-coupled receptors, Nav1.7, and Nav1.8 regulate the increase in AP firing magnitude in primed muscle nociceptors; (2) decreasing Kv1.1 expression causes the greatest increase in AP firing in primed neurons compared to individual modifications of Kv7.2, Nav1.7, Nav1.8, or TRPA1; and (3) an inflammation-induced decrease in Kv7.2 or Kv1.1 expression in silent nociceptors or non-responsive sensory neurons might cause them to become responsive and sensitize during a subsequent inflammatory event when stimulated by noxious mechanical forces.

We can use our findings to further the area of pain research in different ways. First, performing in vivo studies to experimentally test our computationally derived hypotheses could confirm the roles of specific proteins in hyperalgesic priming of nociceptors and lead to an improved understanding of chronic pain initiation in muscle tissue. In addition, we could use our findings to assess whether the same proteins that prime muscle neurons to inflammatory pain are also involved in different pathological pain scenarios, such as pain arising from other sensory processing disorders, e.g., fibromyalgia. Finally, by extending the computational model to include the kinetics of pharmaceutical drugs that act as inhibitors or enhancers for specific ion channels or receptors, we can predict their efficacy in regulating the effects of hyperalgesic priming on subsequent pain-inducing events in a dose-dependent manner.

Data availability statement

The original contributions presented in this study are included in this article/Supplementary material, further inquiries can be directed to the corresponding author.

Author contributions

SN: Conceptualization, Formal analysis, Investigation, Methodology, Software, Visualization, Writing – original draft, Writing – review and editing. ST: Conceptualization, Writing – review and editing. JR: Conceptualization, Supervision, Visualization, Writing – review and editing.

Funding

This work was supported by the Combat Casualty Care Research Program of the US Army Medical Research and Development Command (USAMRDC) under Contract No. W81XWH20C0031.

Conflict of interest

SN and ST were employed by The Henry M. Jackson Foundation for the Advancement of Military Medicine, Inc.

The remaining author declares that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Author disclaimer

The opinions and assertions contained herein are the private views of the authors and are not to be construed as official or as reflecting the views of the United States (US) Army, the US Department of Defense, or The Henry M. Jackson Foundation (HJF) for the Advancement of Military Medicine, Inc. This manuscript has been approved for public release with unlimited distribution.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fnins.2023.1254154/full#supplementary-material

References

Aley, K. O., and Levine, J. D. (1999). Role of protein kinase A in the maintenance of inflammatory pain. J. Neurosci. 19, 2181–2186. doi: 10.1523/JNEUROSCI.19-06-02181.1999

PubMed Abstract | CrossRef Full Text | Google Scholar

Amir, R., and Devor, M. (2003). Electrical excitability of the soma of sensory neurons is required for spike invasion of the soma, but not for through-conduction. Biophys. J. 84, 2181–2191. doi: 10.1016/S0006-3495(03)75024-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Araldi, D., Ferrari, L. F., and Levine, J. D. (2016a). Adenosine-A1 receptor agonist induced hyperalgesic priming type II. Pain 157, 698–709. doi: 10.1097/j.pain.0000000000000421

PubMed Abstract | CrossRef Full Text | Google Scholar

Araldi, D., Ferrari, L. F., and Levine, J. D. (2016b). Gi-protein-coupled 5-HT1b/d receptor agonist sumatriptan induces type I hyperalgesic priming. Pain 157, 1773–1782. doi: 10.1097/j.pain.0000000000000581

PubMed Abstract | CrossRef Full Text | Google Scholar

Baker, M. D. (2005). Protein kinase C mediates up-regulation of tetrodotoxin-resistant, persistent Na⁺ current in rat and mouse sensory neurones. J. Physiol. 567(Pt. 3), 851–867. doi: 10.1113/jphysiol.2005.089771

PubMed Abstract | CrossRef Full Text | Google Scholar

Bar, K. J., Natura, G., Telleria-Diaz, A., Teschner, P., Vogel, R., Vasquez, E., et al. (2004). Changes in the effect of spinal prostaglandin E2 during inflammation: Prostaglandin e (EP1-EP4) receptors in spinal nociceptive processing of input from the normal or inflamed knee joint. J. Neurosci. 24, 642–651. doi: 10.1523/JNEUROSCI.0882-03.2004

PubMed Abstract | CrossRef Full Text | Google Scholar

Bennett, M. R., Farnell, L., and Gibson, W. G. (2005). A quantitative description of the contraction of blood vessels following the release of noradrenaline from sympathetic varicosities. J. Theor. Biol. 234, 107–122. doi: 10.1016/j.jtbi.2004.11.013

PubMed Abstract | CrossRef Full Text | Google Scholar

Blair, N. T., and Bean, B. P. (2002). Roles of tetrodotoxin (TTX)-sensitive Na⁺ current, TTX-resistant Na⁺ current, and Ca²⁺ current in the action potentials of nociceptive sensory neurons. J. Neurosci. 22, 10277–10290.

Google Scholar

Cavanaugh, D. J., Lee, H., Lo, L., Shields, S. D., Zylka, M. J., Basbaum, A. I., et al. (2009). Distinct subsets of unmyelinated primary sensory fibers mediate behavioral responses to noxious thermal and mechanical stimuli. Proc. Natl. Acad. Sci. U.S.A. 106, 9075–9080. doi: 10.1073/pnas.0901507106

PubMed Abstract | CrossRef Full Text | Google Scholar

Chambers, J. D., Bornstein, J. C., Gwynne, R. M., Koussoulas, K., and Thomas, E. A. (2014). A detailed, conductance-based computer model of intrinsic sensory neurons of the gastrointestinal tract. Am. J. Physiol. Gastrointest Liver Physiol. 307, G517–G532. doi: 10.1152/ajpgi.00228.2013

PubMed Abstract | CrossRef Full Text | Google Scholar

D’Adamo, M. C., Liantonio, A., Rolland, J. F., Pessia, M., and Imbrici, P. (2020). Kv1.1 channelopathies: Pathophysiological mechanisms and therapeutic approaches. Int. J. Mol. Sci. 21:2935. doi: 10.3390/ijms21082935

PubMed Abstract | CrossRef Full Text | Google Scholar

del Camino, D., Murphy, S., Heiry, M., Barrett, L. B., Earley, T. J., Cook, C. A., et al. (2010). TRPA1 contributes to cold hypersensitivity. J Neurosci. 30, 15165–15174. doi: 10.1523/JNEUROSCI.2580-10.2010

PubMed Abstract | CrossRef Full Text | Google Scholar

Deuis, J. R., Dvorakova, L. S., and Vetter, I. (2017). Methods used to evaluate pain behaviors in rodents. Front. Mol. Neurosci. 10:284. doi: 10.3389/fnmol.2017.00284

PubMed Abstract | CrossRef Full Text | Google Scholar

Dubin, A. E., and Patapoutian, A. (2010). Nociceptors: The sensors of the pain pathway. J. Clin. Invest. 120, 3760–3772. doi: 10.1172/JCI42843

PubMed Abstract | CrossRef Full Text | Google Scholar

Fang, J., Wang, S., Zhou, J., Shao, X., Sun, H., Liang, Y., et al. (2021). Electroacupuncture regulates pain transition through inhibiting PKCepsilon and TRPV1 expression in dorsal root ganglion. Front. Neurosci. 15:685715. doi: 10.3389/fnins.2021.685715

PubMed Abstract | CrossRef Full Text | Google Scholar

Ferrari, L. F., Bogen, O., and Levine, J. D. (2013a). Role of nociceptor αCAMKII in transition from acute to chronic pain (hyperalgesic priming) in male and female rats. J. Neurosci. 33, 11002–11011. doi: 10.1523/JNEUROSCI.1785-13.2013

PubMed Abstract | CrossRef Full Text | Google Scholar

Ferrari, L. F., Bogen, O., Reichling, D. B., and Levine, J. D. (2015). Accounting for the delay in the transition from acute to chronic pain: Axonal and nuclear mechanisms. J. Neurosci. 35, 495–507. doi: 10.1523/JNEUROSCI.5147-13.2015

PubMed Abstract | CrossRef Full Text | Google Scholar

Ferrari, L. F., Levine, E., and Levine, J. D. (2013b). Role of a novel nociceptor autocrine mechanism in chronic pain. Eur. J. Neurosci. 37, 1705–1713. doi: 10.1111/ejn.12145

PubMed Abstract | CrossRef Full Text | Google Scholar

Fitzgerald, M. (2005). The development of nociceptive circuits. Nat. Rev. Neurosci. 6, 507–520. doi: 10.1038/nrn1701

PubMed Abstract | CrossRef Full Text | Google Scholar

Gangadharan, V., and Kuner, R. (2013). Pain hypersensitivity mechanisms at a glance. Dis. Model Mech. 6, 889–895. doi: 10.1242/dmm.011502

PubMed Abstract | CrossRef Full Text | Google Scholar

Gold, M. S., and Flake, N. M. (2005). Inflammation-mediated hyperexcitability of sensory neurons. Neurosignals 14, 147–157. doi: 10.1159/000087653

PubMed Abstract | CrossRef Full Text | Google Scholar

Gold, M. S., and Gebhart, G. F. (2010). Nociceptor sensitization in pain pathogenesis. Nat. Med. 16, 1248–1257. doi: 10.1038/nm.2235

PubMed Abstract | CrossRef Full Text | Google Scholar

Hameed, S. (2019). Na(v)1.7 and Na(v)1.8: Role in the pathophysiology of pain. Mol Pain 15:1744806919858801. doi: 10.1177/1744806919858801

PubMed Abstract | CrossRef Full Text | Google Scholar

Heinricher, M. M. (2016). Pain modulation and the transition from acute to chronic pain. Adv. Exp. Med. Biol. 904, 105–115. doi: 10.1007/978-94-017-7537-3_8

PubMed Abstract | CrossRef Full Text | Google Scholar

Hendrich, J., Alvarez, P., Joseph, E. K., Chen, X., Bogen, O., and Levine, J. D. (2013). Electrophysiological correlates of hyperalgesic priming in vitro and in vivo. Pain 154, 2207–2215. doi: 10.1016/j.pain.2013.07.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Hendrich, J., Alvarez, P., Joseph, E. K., Ferrari, L. F., Chen, X., and Levine, J. D. (2012). In vivo and in vitro comparison of female and male nociceptors. J. Pain 13, 1224–1231. doi: 10.1016/j.jpain.2012.09.009

PubMed Abstract | CrossRef Full Text | Google Scholar

Ho, T. C., Horn, N. A., Huynh, T., Kelava, L., and Lansman, J. B. (2012). Evidence TRPV4 contributes to mechanosensitive ion channels in mouse skeletal muscle fibers. Channels 6, 246–254. doi: 10.4161/chan.20719

PubMed Abstract | CrossRef Full Text | Google Scholar

Hodgkin, A. L., and Huxley, A. F. (1952). A quantitative description of membrane current and its application to conduction and excitation in nerve. J. Physiol. 117, 500–544. doi: 10.1113/jphysiol.1952.sp004764

PubMed Abstract | CrossRef Full Text | Google Scholar

Hucho, T., and Levine, J. D. (2007). Signaling pathways in sensitization: Toward a nociceptor cell biology. Neuron 55, 365–376. doi: 10.1016/j.neuron.2007.07.008

PubMed Abstract | CrossRef Full Text | Google Scholar

Kandasamy, R., and Price, T. J. (2015). The pharmacology of nociceptor priming. Handb. Exp. Pharmacol. 227, 15–37. doi: 10.1007/978-3-662-46450-2_2

PubMed Abstract | CrossRef Full Text | Google Scholar

Karashima, Y., Talavera, K., Everaerts, W., Janssens, A., Kwan, K. Y., Vennekens, R., et al. (2009). TRPA1 acts as a cold sensor in vitro and in vivo. Proc. Natl. Acad. Sci. U.S.A. 106, 1273–1278. doi: 10.1073/pnas.0808487106

PubMed Abstract | CrossRef Full Text | Google Scholar

Kerr, B. J., Souslova, V., McMahon, S. B., and Wood, J. N. (2001). A role for the TTX-resistant sodium channel Nav 1.8 in NGF-induced hyperalgesia, but not neuropathic pain. Neuroreport 12, 3077–3080. doi: 10.1097/00001756-200110080-00019

PubMed Abstract | CrossRef Full Text | Google Scholar

Khomula, E. V., Araldi, D., and Levine, J. D. (2019). In vitro nociceptor neuroplasticity associated with in vivo opioid-induced hyperalgesia. J. Neurosci. 39, 7061–7073. doi: 10.1523/JNEUROSCI.1191-19.2019

PubMed Abstract | CrossRef Full Text | Google Scholar

Khomula, E. V., Araldi, D., Bonet, I. J. M., and Levine, J. D. (2021). Opioid-induced hyperalgesic priming in single nociceptors. J. Neurosci. 41, 31–46. doi: 10.1523/JNEUROSCI.2160-20.2020

PubMed Abstract | CrossRef Full Text | Google Scholar

King, C. H., Lancaster, E., Salomon, D., Peles, E., and Scherer, S. S. (2014). Kv7.2 regulates the function of peripheral sensory neurons. J. Comp. Neurol. 522, 3262–3280. doi: 10.1002/cne.23595

PubMed Abstract | CrossRef Full Text | Google Scholar

Kwan, K. Y., Allchorne, A. J., Vollrath, M. A., Christensen, A. P., Zhang, D. S., Woolf, C. J., et al. (2006). TRPA1 contributes to cold, mechanical, and chemical nociception but is not essential for hair-cell transduction. Neuron 50, 277–289. doi: 10.1016/j.neuron.2006.03.042

PubMed Abstract | CrossRef Full Text | Google Scholar

Kwan, K. Y., Glazer, J. M., Corey, D. P., Rice, F. L., and Stucky, C. L. (2009). TRPA1 modulates mechanotransduction in cutaneous sensory neurons. J. Neurosci. 29, 4808–4819. doi: 10.1523/JNEUROSCI.5380-08.2009

PubMed Abstract | CrossRef Full Text | Google Scholar

Lindskog, M., Kim, M., Wikstrom, M. A., Blackwell, K. T., and Kotaleski, J. H. (2006). Transient calcium and dopamine increase PKA activity and DARPP-32 phosphorylation. PLoS Comput. Biol. 2:e119. doi: 10.1371/journal.pcbi.0020119

PubMed Abstract | CrossRef Full Text | Google Scholar

Mandge, D., and Manchanda, R. (2018). A biophysically detailed computational model of urinary bladder small DRG neuron soma. PLoS Comput. Biol. 14:e1006293. doi: 10.1371/journal.pcbi.1006293

PubMed Abstract | CrossRef Full Text | Google Scholar

Mense, S. (2010). “Functional anatomy of muscle: Muscle, nociceptors and afferent fibers,” in Muscle pain: Understanding the mechanisms, eds S. Mense and R. Gerwin (Berlin: Springer). doi: 10.1007/978-3-540-85021-2_2

CrossRef Full Text | Google Scholar

Mickle, A. D., Shepherd, A. J., and Mohapatra, D. P. (2015). Sensory TRP channels: The key transducers of nociception and pain. Prog. Mol. Biol. Transl. Sci. 131, 73–118. doi: 10.1016/bs.pmbts.2015.01.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Mitrophanov, A. Y., Rosendaal, F. R., and Reifman, J. (2015). Mechanistic modeling of the effects of acidosis on thrombin generation. Anesth. Analg. 121, 278–288. doi: 10.1213/ANE.0000000000000733

PubMed Abstract | CrossRef Full Text | Google Scholar

Nagaraja, S., Chen, L., Zhou, J., Zhao, Y., Fine, D., DiPietro, L. A., et al. (2017). Predictive analysis of mechanistic triggers and mitigation strategies for pathological scarring in skin wounds. J. Immunol. 198, 832–841. doi: 10.4049/jimmunol.1601273

PubMed Abstract | CrossRef Full Text | Google Scholar

Nagaraja, S., Queme, L. F., Hofmann, M. C., Tewari, S. G., Jankowski, M. P., and Reifman, J. (2021). In silico identification of key factors driving the response of muscle sensory neurons to noxious stimuli. Front. Neurosci. 15:719735. doi: 10.3389/fnins.2021.719735

PubMed Abstract | CrossRef Full Text | Google Scholar

Nagaraja, S., Tewari, S. G., and Reifman, J. (2023). Identification of key factors driving inflammation-induced sensitization of muscle sensory neurons. Front. Neurosci. 17:1147437. doi: 10.3389/fnins.2023.1147437

PubMed Abstract | CrossRef Full Text | Google Scholar

Nagaraja, S., Wallqvist, A., Reifman, J., and Mitrophanov, A. Y. (2014). Computational approach to characterize causative factors and molecular indicators of chronic wound inflammation. J. Immunol. 192, 1824–1834. doi: 10.4049/jimmunol.1302481

PubMed Abstract | CrossRef Full Text | Google Scholar

Nassar, M. A., Stirling, L. C., Forlani, G., Baker, M. D., Matthews, E. A., Dickenson, A. H., et al. (2004). Nociceptor-specific gene deletion reveals a major role for Nav1.7 (PN1) in acute and inflammatory pain. Proc. Natl. Acad. Sci. U.S.A. 101, 12706–12711. doi: 10.1073/pnas.0404915101

PubMed Abstract | CrossRef Full Text | Google Scholar

Nicol, G. D., Vasko, M. R., and Evans, A. R. (1997). Prostaglandins suppress an outward potassium current in embryonic rat sensory neurons. J. Neurophysiol. 77, 167–176. doi: 10.1152/jn.1997.77.1.167

PubMed Abstract | CrossRef Full Text | Google Scholar

Niederberger, E., and Geisslinger, G. (2008). The IKK-NF-κB pathway: A source for novel molecular drug targets in pain therapy? FASEB J. 22, 3432–3442. doi: 10.1096/fj.08-109355

PubMed Abstract | CrossRef Full Text | Google Scholar

Pak, D. J., Yong, R. J., Kaye, A. D., and Urman, R. D. (2018). Chronification of pain: Mechanisms, current understanding, and clinical implications. Curr. Pain Headache Rep. 22:9. doi: 10.1007/s11916-018-0666-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Parada, C. A., Reichling, D. B., and Levine, J. D. (2005). Chronic hyperalgesic priming in the rat involves a novel interaction between camp and PKCepsilon second messenger pathways. Pain 113, 185–190. doi: 10.1016/j.pain.2004.10.021

PubMed Abstract | CrossRef Full Text | Google Scholar

Payne, C. E., Brown, A. R., Theile, J. W., Loucif, A. J., Alexandrou, A. J., Fuller, M. D., et al. (2015). A novel selective and orally bioavailable Nav 1.8 channel blocker, PF-01247324, attenuates nociception and sensory neuron excitability. Br. J. Pharmacol. 172, 2654–2670. doi: 10.1111/bph.13092

PubMed Abstract | CrossRef Full Text | Google Scholar

Pozek, J. P., Beausang, D., Baratta, J. L., and Viscusi, E. R. (2016). The acute to chronic pain transition: Can chronic pain be prevented? Med. Clin. North Am. 100, 17–30. doi: 10.1016/j.mcna.2015.08.005

PubMed Abstract | CrossRef Full Text | Google Scholar

Reichling, D. B., and Levine, J. D. (2009). Critical role of nociceptor plasticity in chronic pain. Trends Neurosci. 32, 611–618. doi: 10.1016/j.tins.2009.07.007

PubMed Abstract | CrossRef Full Text | Google Scholar

Suleimanova, A., Talanov, M., Gafurov, O., Gafarov, F., Koroleva, K., Virenque, A., et al. (2020). Modeling a nociceptive neuro-immune synapse activated by ATP and 5-HT in meninges: Novel clues on transduction of chemical signals into persistent or rhythmic neuronal firing. Front. Cell Neurosci. 14:135. doi: 10.3389/fncel.2020.00135

PubMed Abstract | CrossRef Full Text | Google Scholar

Tanaka, B. S., Zhao, P., Dib-Hajj, F. B., Morisset, V., Tate, S., Waxman, S. G., et al. (2016). A gain-of-function mutation in Nav1.6 in a case of trigeminal neuralgia. Mol. Med. 22, 338–348. doi: 10.2119/molmed.2016.00131

PubMed Abstract | CrossRef Full Text | Google Scholar

Tigerholm, J., Petersson, M. E., Obreja, O., Lampert, A., Carr, R., Schmelz, M., et al. (2014). Modeling activity-dependent changes of axonal spike conduction in primary afferent C-nociceptors. J. Neurophysiol. 111, 1721–1735. doi: 10.1152/jn.00777.2012

PubMed Abstract | CrossRef Full Text | Google Scholar

Toossi, A., Bergin, B., Marefatallah, M., Parhizi, B., Tyreman, N., Everaert, D. G., et al. (2021). Comparative neuroanatomy of the lumbosacral spinal cord of the rat, cat, pig, monkey, and human. Sci. Rep. 11:1955. doi: 10.1038/s41598-021-81371-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Vijayaragavan, K., Boutjdir, M., and Chahine, M. (2004). Modulation of Nav1.7 and Nav1.8 peripheral nerve sodium channels by protein kinase a and protein kinase C. J. Neurophysiol. 91, 1556–1569. doi: 10.1152/jn.00676.2003

PubMed Abstract | CrossRef Full Text | Google Scholar

Woolf, C. J., and Ma, Q. (2007). Nociceptors–noxious stimulus detectors. Neuron 55, 353–364. doi: 10.1016/j.neuron.2007.07.016

PubMed Abstract | CrossRef Full Text | Google Scholar

Wu, D. F., Chandra, D., McMahon, T., Wang, D., Dadgar, J., Kharazia, V. N., et al. (2012). PKCepsilon phosphorylation of the sodium channel Nav1.8 increases channel function and produces mechanical hyperalgesia in mice. J. Clin. Invest. 122, 1306–1315. doi: 10.1172/JCI61934

PubMed Abstract | CrossRef Full Text | Google Scholar

Wuertz, K., Vo, N., Kletsas, D., and Boos, N. (2012). Inflammatory and catabolic signalling in intervertebral discs: The roles of NF-κB and MAP kinases. Eur. Cell Mater. 23:103–119; discussion 119–120. doi: 10.22203/ecm.v023a08

PubMed Abstract | CrossRef Full Text | Google Scholar