Distinct transcriptional responses of mouse sensory neurons in models of human chronic pain conditions

Background: Sensory neurons play an essential role in almost all pain conditions, and have recently been classified into distinct subsets on the basis of their transcriptomes. Here we have analysed alterations in dorsal root ganglia (DRG) gene expression using microarrays in mouse models related to human chronic pain. Methods: Six different pain models were studied in male C57BL/6J mice: (1) bone cancer pain using cancer cell injection in the intramedullary space of the femur; (2) neuropathic pain using partial sciatic nerve ligation; (3) osteoarthritis pain using mechanical joint loading; (4) chemotherapy-induced pain with oxaliplatin; (5) chronic muscle pain using hyperalgesic priming; and (6) inflammatory pain using intraplantar complete Freund’s adjuvant. Microarray analyses were performed using RNA isolated from dorsal root ganglia and compared to sham/vehicle treated controls. Results: Differentially expressed genes (DEGs) were identified. Known and previously unreported genes were found to be dysregulated in each pain model. The transcriptomic profiles for each model were compared and expression profiles of DEGs within subsets of DRG neuronal populations were analysed to determine whether specific neuronal subsets could be linked to each of the pain models. Conclusions: Each pain model exhibits a unique set of altered transcripts implying distinct cellular responses to different painful stimuli. No simple direct link between genetically distinct sets of neurons and particular pain models could be discerned.

Introduction Chronic pain is a major clinical problem affecting roughly 20% of the population with more than 6% suffering debilitating levels of pain 1,2 . Despite the huge clinical burden, little progress has been made in developing more effective analgesic agents. Dorsal root ganglion (DRG) sensory neurons are particularly exciting targets for drug development because of their essential role in driving pain sensations in the central nervous system. With the wider availability of high-throughput RNA-seq, there has been a major effort to define DRG transcriptomes, both at the level of whole ganglia 3 and single neurons 4,5 . This has allowed a new classification of sets of sensory neurons based on genetic identity, rather than on the rate of action potential propagation (A fibres and C fibres). As well as neurons, non-neuronal cells such as glia 6 and immune system cells are found within sensory ganglia, with leukocytes alone making up 5-10% of cells in DRG 7,8 .
We wondered if the newly defined sensory neuron subsets were differentially activated in different pain models. We therefore performed microarray analyses of RNAs isolated from the DRG of male mice subjected to six different interventions that mimic common chronic pain in human patients. DRG were dissected after the development of peak pain behaviours in each case. The six different chronic pain models were: cancer cell injection in the intramedullary space of the femur as a model for bone cancer pain; partial sciatic nerve ligation (PSL) for neuropathic pain; mechanical joint loading for osteoarthritis pain; oxaliplatin-induced painful neuropathy for chemotherapyinduced pain; hyperalgesic priming model for chronic muscle pain and intraplantar complete Freund's adjuvant (CFA) for inflammatory pain. The rodent models for PSL 9 , CFA 10 , chemotherapy induced pain with oxaliplatin 11 , bone cancer pain 12 and hyperalgesic priming 13-15 are well characterized. The mechanical joint loading model robustly recapitulates behavioural signs and symptoms of osteoarthritis 16 .
By comparing altered levels of gene expression in six distinct chronic pain models we hoped to gain insights into the neuronal subsets that are particularly important in different pain conditions. We also could determine genes and related pathways that are common or unique to each condition. Here we show that there is no simple relationship between transcriptional changes in the newly defined subsets of sensory neurons and particular painful insults. However, we have identified several known and previously unreported genes altered in expression in each pain condition. Relevant pathways and upstream regulators were also analysed using Ingenuity Pathway Analysis including cytokines, transcription factors (TFs), G-protein coupled receptors (GCPRs) and ion channels. The results from these experiments highlight the complexity of sensory neuron cell types involved in nociceptive responses. Although altered DRG gene expression may not necessarily play a causal role in pain, some dysregulated genes are potential analgesic drug targets.

Methods
All experiments were performed in accordance with the UK Animals (Scientific Procedures) Act 1986 with prior approval under a Home Office project licence (PPL 70/7382). Mice were kept on a 12-h light/dark cycle and provided with food and water ad libitum. All animals were acclimatized for 1 week to the facility before the start of the experiment. Mice were housed in individually ventilated cages (Techniplast GM500 Mouse IVC Green line) containing Lignocel bedding with a maximum of 5 adult mice per cage. All tests and surgeries were conducted using adult male C57BL/6J mice supplied by Charles River and Envigo (specific numbers detailed in the figure legends and ages detailed below for each model). Sample size for each behavioural model was calculated using G*Power (Ver. 3.1.9.2) for a power of 0.8 17 . Surgical procedures were performed by trained researchers and under aseptic conditions. All efforts were made to ameliorate any suffering of the animals with surgery performed under anaesthesia. Mice were euthanized by CO 2 asphyxiation followed by cervical dislocation. microCT (uCT) for bone cancer Femurs were post-fixed in 4% PFA (Sigma-Aldrich) for 24h and kept in 70% ethanol (VWR) till scanning with microCT. Images were acquired with Skyscan software at 6.41um/pixel with 0.6 degrees rotation steps and 2 frame averaging. Image data was reconstructed with NRecon software. A CT-analyzer was used to select a 1mm volume of interest (VOI) region, starting at 0.6mm from the growth plate. Bone mineral density was determined by plotting attenuation coefficients against a standard curve determined with two phantoms with known density. Representative images were binarized with ImageJ (1.52c) and 3D Viewer extension (Java 3D 1.6 with 3D Viewer plugin 4.0.2) used for 3D reconstructions of femur scans.
Partial sciatic nerve ligation (PSL) model Surgical procedures were performed under isoflurane anaesthesia (2-3%). A partial nerve injury in 14 week old mice was induced by tying a tight ligature with 6-0 silk suture around approximately 1/3 to 1/2 the diameter of the sciatic nerve, similar to the approach described in rats 9 . Ipsilateral and contralateral L3-L5 DRGs were dissected prior to RNA isolation on postsurgery day 16.

Mechanical joint loading (MJL) model
In vivo loading: Osteoarthritis was induced by subjecting 12 week old mice to a two week loading regimen using an electronic materials testing machine (Bose 3100). Throughout the loading episodes, mice were anaesthetised with isoflurane 3.5% (±0.5%) with the right tibia positioned vertically between two custommade cups which fixate the knee and ankle joint in deep flexion 20 . The loading regimen was the same as previously described 16 .
Loading was repeated three times per week for two consecutive weeks. The non-loaded control group consisted of age-and cage-matched mice that were not subjected to a loading regimen but received isoflurane anaesthesia for the same duration as loaded mice. Ipsilateral L3-L5 DRGs were dissected from loaded and non-loaded mice after three weeks.

Behavioural tests
For behavioural experiments, animals were acclimatized to the equipment for at least 2 days prior to testing. Observers who performed behavioural experiments were blind to the test/sham groups. The cancer induced bone pain model used a limb score assessment, as previously described 19 . A limb score of 2 (i.e. significant limping) was used as the checkpoint for ending the experiment with DRGs and ipsilateral femurs collected.
Mechanical sensitivity was measured using the up-down method 24 and static weight bearing 18 assays. Thermal nociceptive thresholds were determined by measuring paw withdrawal latency using the Hargreaves' apparatus 25 . The response to noxious cold was measured using the cold plantar assay 26 . Differential gene expression (DGE) and splice variant analysis Differential gene expression (DGE) and splice variant analysis of transcriptome array data was performed using Transcriptome Analysis Console (TAC Ver 4.0, Thermo Fisher). Normalization was carried out according to the manufacturer's default protocol converting probe cell intensity (CEL files) into signal data (CHP files). An ANOVA (eBayes) cut-off p-value of <0.05 was used to identify DEGs, with the fold change filter then used to sort genes. Genes without a curated gene symbol associated with the Affymetrix probe set were excluded from downstream DGE analyses. Differentially expressed genes (DEGs) were correlated with known expression in DRG neuronal subtypes 4 . Splice variant analysis was performed using the Gene + Exon function in TAC which uses a gene normalized intensity value (ratio of probe set intensity to gene expression level). The cut off was set at Splicing Index of 2, splicing ANOVA with p<0.05, gene expressed in both conditions, probe selection region (PSR)/ junction expressed in at least one condition and gene containing at least one PSR. Non-coding genes were excluded from the splicing analyses. An exon event score indicated consensus alternative splicing event score where a score of 1 represents a high likelihood of splicing event. Heatmaps, volcano plots and statistical analyses were generated using GraphPad Prism 7.

RNA extraction and microarrays
Ingenuity Pathway Analysis (IPA) Ingenuity Pathway Analysis (IPA, Qiagen) for upstream regulators was performed according to the manufacturer's instructions. IPA primarily utilizes published results of knockdown or knockout studies of genes characterized as growth factor, kinase/phosphatase, transcription/translation regulator, cytokine, transmembrane receptor, enzyme, ion channel, transporter, G-protein coupled receptor, peptidase, and microRNA. This analysis examines target genes from the array dataset, compares their direction of change to what is expected from the literature to predict likely upstream regulators using overlap p-value (Fishers Exact Test, p<0.01) and an activation z-score corrected for bias. The top upstream regulators with an activation z-score were examined for each pain pathway. The activation z-score was used to make a prediction regarding activation or inhibition of the upstream regulator.

Results
We modelled six different chronic pain conditions in mice and used behavioural assays to confirm the development of mechanical or thermal hypersensitivity, followed by transcriptomic analyses of whole dorsal root ganglia tissue. Significant DGE between naïve and injured conditions and alternative splicing was assessed for each model. DGE values for every gene on the mouse transcriptome array and alternative splicing changes are available in Supplementary Table 1.

Cancer-induced bone pain
Carcinoma cells were introduced into the mouse femur and pain behaviour measured daily to determine the limb use score (4: normal use of limb, 3: slight limping, 2: clear limping; see Methods). Mice were culled when a limb score of 2 was reached, which ranged from Day 8-16 post surgery ( Figure 1A). Weight bearing tests confirmed a significant reduction in use of the affected limb ( Figure 1B). Bone mineral density (BMD) was evaluated as an additional indicator of the cancer phenotype and was shown to be significantly decreased in cancer mice compared to sham controls ( Figure 1C, D). This is also shown in Figure 1E using a micro-CT reconstruction.
Dorsal root ganglia (ipsilateral L2-L4) were isolated from test mice and paired sham controls and RNA expression levels analysed using microarrays. The volcano plot ( Figure 1F) shows that all significantly dysregulated genes (p<0.05) were differentially expressed between ±10-fold with only 1 non-coding gene (Gm25931) falling outside of this range (fold change of -10.75).
Interestingly, very few of the differentially expressed DRG genes in the cancer induced bone pain model were shared with other models ( Figure 1G). A heat map was generated for genes that were significantly up or down regulated by 2 fold or more ( Figure 1H). Forty-five genes were found to be upregulated >2 fold and 15 genes downregulated. Dysregulated genes were coding and non-coding, with many genes annotated as 'predicted' (Supplementary Table 1). Excluding these 'predicted' genes, the top 5 most upregulated genes were Nts, miR6898, Snord49b, Sprr1a and Gal ( Figure 1H). The cancer induced bone pain model induced significant changes in 23 microRNAs in the >1.5 fold range, indicating the potential role of miRNAs in cancer pain modulation (Supplementary Figure 2A).

Partial sciatic nerve ligation (PSL) model
In the PSL neuropathic pain model mice developed mechanical allodynia from Day 7 post surgery, as characterized by a significant reduction in the 50% paw withdrawal threshold (PWT) in the injured (ipsilateral) compared to the non-injured (contralateral) paw ( Figure 2A). Pain persisted for 15 days post-surgery, with ipsilateral and contralateral L3-L5 DRGs isolated on day 16 and RNA analysed by microarrays.
The volcano plot ( Figure 2B) shows that all significantly dysregulated genes (p<0.05) were differentially expressed between ±10-fold with the exception of Nts, Gpr151 and Npy (respective fold changes of 13.21, 13.88 and 66.18). Similar to the cancer model, very few PSL stimulated genes were shared with other pain models ( Figure 2C). The PSL model stimulated changes in the largest cohort of DRG genes in the 2 fold change range (97 genes) compared to the other 5 models ( Figure 2D). The top 5 most upregulated genes, Npy, Gpr151, Nts, Gal and Atf3, have all been previously reported as being upregulated in neuropathic pain 29-32 . Eighty-six genes were identified as being upregulated more than 2 fold with several of these being noncoding RNAs (Supplementary Table 1). Fourteen microRNAs were upregulated >1.5-fold and 5 were downregulated <-1.5-fold (p<0.05, ANOVA) (Supplementary Figure 2B).

Mechanical joint loading (MJL) model
Osteoarthritis was induced through a two-week mechanical joint loading protocol after which animals developed significant mechanical hypersensitivity and altered weight bearing ( Figure 3A). Ipsilateral L3-L5 DRGs from loaded and nonloaded mice were extracted 3 weeks after the loading regimen, RNA extracted and analysed for DEGs. The volcano plot ( Figure 3B) shows that all significantly dysregulated genes (p<0.05) were differentially expressed between ±7-fold. Only fourteen genes were dysregulated >2 fold ( Figure 3D) with only 1-2 of these genes also differentially expressed in the other pain models ( Figure 3C). The top three upregulated genes (Ppbp, DLEU2_5 and Scarna9) have no previous association with osteoarthritic pain.

Oxaliplatin-induced painful neuropathy
Oxaliplatin is a chemotherapeutic drug known to induce mechanical and cold allodynia 33 . Mice exposed to oxaliplatin were found to have a reduced cold pain threshold in the cold plantar test ( Figure 4A). Twenty-eight days following the first oxaliplatin treatment, ipsilateral L3-L5 DRGs were extracted from test and vehicle treated mice and analysed for differential gene expression. The volcano plot ( Figure 4B) shows that all significantly dysregulated genes (p<0.05) were differentially expressed between ±10-fold. Oxaliplatin treatment resulted in 58 dysregulated genes in the >2 fold range, which in contrast to other models, the majority (36 genes) being down-rather than up-regulated ( Figure 4D). Out of the 58 DEGs, 7 were shared with the CFA model, 3 with the PSL model and 2 with the mechanical joint loading model ( Figure 4C).

Chronic muscle pain (CMP)
A hyperalgesic priming model was used to study chronic muscle pain, whereby the gastrocnemius muscle was injected twice, four days apart, with 3% carrageenan. Prolonged mechanical hypersensitivity resulted from the hyperalgesic priming ( Figure 5A). On Day 29, ipsilateral L3-L5 DRGs were extracted from test and vehicle treated mice and analysed for DGE. The volcano plot ( Figure 5B) shows that all significantly dysregulated genes (p<0.05) were differentially expressed between ±5-fold. Interestingly, only ten genes were dysregulated >2 fold ( Figure 5D) with Snord49b shared with the cancer-induced bone pain model and Gm19860 shared with the mechanical joint loading model ( Figure 5C).

CFA inflammatory pain model
Intraplantar injection of CFA was used to induce a robust inflammatory pain phenotype, with significant thermal hypersensitivity observed thirty-six hours post CFA injection ( Figure 6A). Ipsilateral L3-L5 DRGs were extracted 2 days post CFA and RNA analysed for differential gene expression. The volcano plot ( Figure 6B) shows that all significantly dysregulated genes (p<0.05) were differentially expressed between ±9-fold. Seventy genes were differentially expressed >2 fold with the majority of genes (52) being upregulated ( Figure 6D). Once again, very few dysregulated genes were shared between CFA and the other models ( Figure 6C).     Analysis of upstream regulators of neuronal DEGs in chronic pain models Next, we analysed pathways and upstream regulators implicated in the four pain conditions with the highest number of DEGs using the Ingenuity Pathway Analysis program. The top upstream regulators in the bone cancer pain model are shown in Figure 8A (p-value of overlap <0.05, ANOVA) and includes cytokine IFNG (activation z-score 1.9); growth factor LEP (activation z-score 0.4); and transcription regulator NFE2L2 (activation z-score 1.9). The analysis was unable to determine whether these regulators were activated or inhibited. The top upstream regulators in the PSL model are depicted in Figure 8B (p-value of overlap<0.05, ANOVA), and includes growth factors TGFB1 (activation z-score 2.2, predicted to be activated) and   BMP7 (activation z-score -1.9); cytokines TNF, IL-6, IL-1-beta (activation z-scores > 2, all predicted to be activated in PSL); kinase MTOR (activation z-score 2.3, predicted to be activated); and the transcription regulator GATA4 (activation z-score -2.8). In the oxaliplatin model ( Figure 8C) the top upstream regulators (p<0.05, ANOVA) involved included cytokines IFNG, IFNL1, PRL (activation z scores< -2, all predicted to be inhibited) and IL1RN (activation z score >2, predicted to be activated); G-protein coupled receptor ACKR2 (activation z score 3, predicted to be activated); and transcription regulators TRIM24, SIRT1, STAT3 (activation z scores > 2, all predicted to be activated) and STAT1 (activation z score <-3, predicted to be inhibited). Finally, for CFA the top upstream regulators (p<0.05, ANOVA) included cytokines IFNG, IL-1B and IFN-B1 (activation z scores>1.7, all predicted to be activated); transcription regulator STAT1 (activation z score 1.8, predicted to be activated) and TRIM24 (activation z score -1.47, predicted to be inhibited) ( Figure 8D). We observed interesting patterns of activation compared to inhibition of upstream regulators between pain models. For example, STAT1 is activated in oxaliplatin, but inhibited in the CFA model; IFNG is activated in CFA and PSL, but inhibited in the oxaliplatin model.

Discussion
DRG sensory neurons are a heterogeneous population that have recently been classified according to their gene expression profiles 4 . By modelling six chronic pain conditions in mice we aimed to (1) identify dysregulated DRG genes in each model; (2) map dysregulated genes to one or more of the 11 neuronal subsets; and (3) determine whether different pain models shared similar dysregulated genes and/or neuronal subsets. Overall our analysis highlighted how dissimilar the different models were in respect to their transcriptomic profiles with the majority of dysregulated genes not shared between models (see Venn diagrams in Figure 1-Figure 6). This implies that chronic pain can arise from a diverse set of mechanisms involving different genes and pathways. Although transcriptional changes may not be causative in the pain induction, the present study suggests distinct cell and molecular mechanisms are involved with associated transcriptional correlates in the pain models studied. Although some key genes were shared between models, such as Neurotensin and Galanin being upregulated in cancer induced bone pain and PSL neuropathic pain models, the majority of dysregulated genes were exclusive to each pain model.
We also largely failed to link narrow subsets of DRG neurons to particular chronic pain models using microarray analyses (see Figure 7). Our aim was to identify dysregulated genes from the microarrays and map these genes to the neuronal populations previously outlined using single cell RNAseq data 4 . However, this approach was limited by the fact that the published RNAseq data is derived from naïve mice (i.e. not in pain) and potentially the gene expression patterns across the neuronal subsets may dynamically change in the pain models. This is certainly true for genes such as Galanin which is markedly upregulated and expressed in ~40-50% of DRG neurons following axotomy 32 and BDNF, which is upregulated in medium and large diameter DRG neurons in neuropathic pain 35 . Furthermore, genes that have no DRG expression in the Usoskin et al. dataset (i.e. zero sequencing reads) were not included in our transcriptomic analyses, but could be expressed following injury. For example, Neuropeptide Y is dramatically upregulated in DRG neurons following peripheral nerve injury, but is not expressed in basal conditions 29 . Thirdly, many genes are broadly expressed in several DRG populations, meaning that restricted critical subsets are difficult to highlight.
Recent advances in deep sequencing have dramatically increased the number of known transcribed genes 36 . The Mouse Transcriptome Array that we used gives a comprehensive coverage of the mouse transcriptome, with >214,000 transcripts represented. Genes that are represented include both coding and non-coding genes (such an lncRNAs, small RNAs and expressed pseudogenes). Although such genes are now known to be expressed, gene annotation and understanding of gene function is lagging behind. This is highlighted in our transcriptomic analyses where a high number of predicted genes (typically with the 'Gm' prefix) showed differential expression. A major challenge is to determine the function of these types of genes and to ascertain the contribution, if any, to the pain phenotypes observed. Similarly, although the dense probe coverage on the array gives the opportunity to assess differential splicing patterns (Supplementary Figure 3), experimental verification is needed to prove that the changes in splicing patterns are relevant.
Given the unique sets of underlying genes, we attempted to identify potential common upstream regulators using IPA analyses. We identified several common cytokines as upstream regulators of DEGs in our pain models consistent with the role of cytokines in generation and maintenance of pain 37,38 , pointing to the interaction of DRG neurons and immune cells as a key point of focus. T cells are known to infiltrate DRGs in models of neuropathic pain 39-41 , however whether the cytokine related downstream DEGs identified in our pain models are being induced by immune cells remains to be investigated.
In summary, we present a systematic study of the transcriptomic changes in six different chronic mouse pain models as a resource for the pain community. Although gene expression changes may not be causally linked to the pain phenotypes observed, the transcriptomic profiles reported here highlight the diversity of cell and molecular mechanisms apparent in different pain conditions. Potential therapeutic gene targets for the distinct chronic pain conditions studied may also be identified after further mechanistic studies.

Data availability
Microarray data has been deposited at Gene Expression Omnibus Array Express for public use with reference number E-MTAB-6864. Table 1), replicates for Figure 1E (3D reconstruction images) and raw values in a GraphPad Prism file for Figure 1A-D, Figure 2A, Figure 3A, Figure 4A, Figure 5A and Figure 6A for all replicates are deposited in OSF: https://doi. org/10.17605/OSF.IO/DG7Z3 42 Data are available under the terms of the Creative Commons Zero "No rights reserved" data waiver (CC0 1.0 Public domain dedication).

Competing interests
No competing interests were disclosed.

Supplementary material
Supplementary Splice variant analysis was performed using the Gene + Exon function in TAC which uses a gene normalized intensity value (ratio of probe set intensity to gene expression level). Click here to access the data.

Supplementary File 1: Supplementary figures.
Click here to access the data.  (exon SI >2 or <-2, Exon p-val <0.05). Exon event score calculated by TAC software (see methods) where 1 represents a high likelihood for the spliced event. Exon events were defined in TAC as: Intron retention-when a sequence is spliced out as an intron or remains in the mature mRNA transcript; Cassette Exon (skipped exon) -when one exon is spliced out of the primary transcript together with its flanking introns; Alternative 3`acceptor site when two or more splice sites are recognized at the 5` end of an exon, an alternative 3` splice junction is used changing the 5` boundary of the downstream exon; Alternative 5` donor site when two or more splice sites are recognized at the 3` end of an exon, an alternative 5` splice junction is used changing the 3` boundary of the upstream exon. The paper discovers differentially expressed genes in mouse primary sensory neurons (DRG) in a variety of pain models that are believed to be models of human pain states. Looking in the primary neurons the authors find that, although pain levels appear similar and significant in all the models studied, the changes in gene expression for each of the pain states is different. The paper is of great interest, and I have enjoyed trying to grapple with its contents. This is a way of saying that I like it, but I wanted to clarify a couple of things, that have raised a flag.

Open Peer Review
The discussion is written as though no gene expression data have ever been obtained before in pain models, and I agree with David Wynick that a more helpful approach would require that some context is offered in the discussion, particularly as this group of authors are generating a benchmark. For example, with reference to the findings on local inflammation, in Figure 2D I do not see a single sodium channel gene affected, and yet Na+ channel genes are reported by others to be dysregulated, for example NaV1.3 expression is reported to double in DRG in inflammatory circumstances . Do the authors think this is surprising and what is the reader to make of this? In fact, is Na+ channel gene expression in DRG neurons affected at all by any of the induced pain states? What has happened to other potentially important pain related ion channels such as KCNQ and HCN channels?
Is one correct deduction then that while there may indeed be stereotyped gene expression changes taking place (given that somehow the maladaptive changes in the nervous system exist for all these pain scenarios) they must be in second and higher order neurons, involving glutamate gated ion channel subunits and so on? Secondly, as it seems inescapable to conclude that the membrane behaviour of primary sensory neurons must be somehow 'abnormal' in pain states driven from the periphery, is the correct conclusion that this is precipitated without any necessity of modification of ion channel expression, completely by a combination of already latent functionality and post-translational modification.
I would like a little more clarification on Figure 7. I think the authors are looking a for a pattern of gene dysregulation that might conspicuously 'fit' the already discovered baseline gene expression profiles of sub-populations of neurons. I think I understand that fitting the profiles against each subtype (previously described) reveals nothing useful, but I would be happy with a bit more explanation for this (bottom page 10 and Figure 7), as I want to rid myself of (I think !) the erroneous impression that single cell analysis is going on. 1 Place the current dataset in the context of the previous published studies that have used transcriptional profiling on the DRG after different pain models (and in some cases comparing two different models), and whether the previous findings and genes identified are consistent with the current data. Identifying genes in the DRG on the basis of their fold-change after a pain model does not mean that those genes are necessarily involved in mediating pain. Indeed, many of the targets of the current analgesics used in clinical practice do not show >two-fold change in expression in mouse models. It is quite possible that many of the genes that do show large changes in each model play roles that relate to apoptosis, survival and/or regeneration, rather than pain transmission and nociception. The main finding of this paper provides further evidence that it may not be possible to identify targets in the final common nociceptive pathway that will target all causes of pain and would be amendable to drug discovery. A more targeted and mechanistic approach is likely to be necessary to focus more on specific pathways that are activated and thus on the individual diseases that cause the pain.

If applicable, is the statistical analysis and its interpretation appropriate? Yes
Are all the source data underlying the results available to ensure full reproducibility? Yes

Are the conclusions drawn adequately supported by the results? Yes
No competing interests were disclosed.

Competing Interests:
I have read this submission. I believe that I have an appropriate level of expertise to confirm that it is of an acceptable scientific standard.