Gene Expression Changes Precede Elevated Mechanical Sensitivity in the Mouse Intervertebral Disc Injury Model

Abstract

ABSTRACT

Background

Back pain after intervertebral disc (IVD) injury is a common clinical problem. Previous work examining early molecular changes post injury mainly used a candidate marker approach.

Methods

In this study, gene expression in the injured and intact mouse tail IVDs was determined with a nonbiased whole transcriptome approach and related to subsequent pain behavior. Mouse tail IVD injury was induced by a needle puncture. Whole murine transcriptome was determined by RNASeq. Transcriptomes of injured IVDs were compared with those of intact controls by bioinformatic methods. Mechanical allodynia was assessed by the Von Frey method.

Results

Among the 17,722 murine genes with meaningful expressions, 7242 genes were differentially expressed (P.adj ⟨ 0.01). Ontology study of upregulated genes revealed that leukocyte migration was the most enriched biological process, and network analysis showed that Tnfa had the most protein–protein interactions. The most enriched downregulated pathways were related to the pattern specification process. Mechanical allodynia persisted at the 4‐week end point.

Conclusion

The RNASeq data revealed numerous early genes that participate in inflammation and repair processes post IVD injury. Mechanical allodynia followed these gene expression changes.

Abbreviations

• 9930111H07Rik
• RIKEN cDNA 9930111H07
• Acod1
• aconitate decarboxylase 1
• Ccl2
• C-C motif chemokine ligand 2
• Cxcl
• chemokine (C-X-C motif) ligand
• D630024D03Rik
• RIKEN cDNA D630024D03 gene
• Edn2
• Endothelin 2
• Ereg
• epiregulin6
• Fabp2
• fatty acid binding protein 2
• Gm39117
• predicted gene, 39 117
• Il
• interleukin
• Kcnj6
• potassium voltage-gated channel subfamily J member 6
• Lif
• leukemia inhibitory factor
• Lrrc74b
• Leucine Rich Repeat Containing 74B
• Mmp
• matrix metallopeptidase
• Nos2
• nitric oxide synthase 2
• Rnd1
• Rho family GTPase 1
• Saa
• serum amyloid A
• Slc7a11
• solute carrier family 7 member 11
• Slpi
• secretory leukocyte peptidase inhibitor
• S100a9
• S100 calcium-binding protein A9
• Slc7a4
• Solute Carrier Family 7 Member 4
• Tnf
• Tumor necrosis factor
• Trdc
• T cell receptor delta constant

Introduction

About 20%–30% of individuals with low back pain can recall a clear injury as the origin of their pain [1]; the actual number of people who have back pain following an injury is likely higher because of recall bias: people often have difficulty accurately remembering past events, especially when it comes to minor injuries, which can lead to underreporting of potential causes of back pain. Investigating the early triggers of inflammation and determining the best timing for treatments in humans is difficult due to delayed clinic visits and limited access to biological samples. Animal models of disc injury/degeneration are valuable for addressing this knowledge gap.

The mouse intervertebral disc (IVD) injury model is informative because it consistently demonstrates molecular, histological, and biomechanical markers of disease progression [2–4]. Furthermore, the tail disc injury method offers advantages such as easy access to the tail discs and low morbidity [2]. Specific sets of inflammatory markers and changes in extracellular matrix gene expression have been well documented following injury; these markers were chosen based on experience and literature from previous studies [2, 5–7]. Recently, nonbiased methods such as whole transcriptome analysis, (e.g., RNRA sequencing), have been employed to discover new markers in response to injury. Transcriptomes of bovine and human IVD tissue have been analyzed, providing insights into the cell types present in this unique tissue and its degenerative process [8–11]. In the mouse, single-cell transcriptome analysis has been conducted during embryonic development of the intact mouse nucleus pulposus (NP) [12]. In the present study, we used bulk RNA Seq to examine the transcriptome of mouse tail IVDs 1 day after injury to gain insight into the early events post-injury. The findings may help to identify novel molecular targets and pathways to treat disc injury and slow down degeneration.

Proteoglycan loss occurs as early as 1 day post injury in the mouse tail IVD injury [13], which leads to spinal segment biomechanical dysfunction [4]. Macrophages also infiltrate the injured IVDs at 1 day post injury [14]. These morphological changes are accompanied by the elevation of a selected panel of inflammatory markers [14]. However, only a limited number of molecules were investigated in previous work. The current study specifically examined the transcriptome at an early time point to uncover previously unsuspected molecules/pathways. Discovering early changes helps determine the mechanism of transition from acute to chronic pain.

Pain-related behavior assessment is an important outcome measure because back pain related to disc degeneration is highly relevant to humans [15]. Pain behaviors have been reported in rat models of disc injury [16–19], correlating with elevated serum TNFα level [19]. Further, annular puncture with TNFα injection enhanced painful behavior [20], and inhibition of TNFα limited long-term pain (defined as 6-weeks post injury) in a rat disc degeneration model [21]. In the mouse, pain-related behavior has been studied in a lumbar disc injury model [22]. However, the behavior in response to tail IVD injury has not been examined previously, despite the widespread use of this model. Here, we used a nonbiased approach to document global gene expression changes 1 day following IVD injury, and the Von Frey method to measure mechanical allodynia over 28 days to associate injury/degeneration with pain. The National Institute of Health Research Task Force defined chronic pain in humans as pain lasting more than 6 months [23]. In rats; 6 weeks post injury was considered long-term observation [21]. In the current study, mouse pain behavior has been examined for 4 weeks and could be considered “subacute phase”

Materials and Methods

Mice

Approval for all animal experimental procedures was obtained from the Institutional Animal Care and Use Committee (IACUC) of the University of Pennsylvania, Philadelphia, PA. Briefly, 24 young adult C57BL/6 male mice (The Jackson Laboratory, Bar Harbor, ME, USA) were used. Eight mice were used for RNA extractions (4 mice underwent tail IVD injury, while 4 mice served as intact controls). Sixteen mice were used for the Von Frey assay (9 mice underwent tail IVD injuries, and 7 mice received sham injuries; 6 injured and 6 sham-injured mice were tested at day 1, 14, and 28 post injury. Three IVD-injured and 1 sham-injured mice were tested at day 7 post injury). The mice were kept in a pathogen-free environment with environmental enrichment (nestlets by Ancare, Bellmore, NY, USA). They had ad lib access to standard mouse chow and acidified water. The housing conditions included a room temperature ranging from 21.1°C to 24.4°C (equivalent to 70°F–76°F) with 30%–70% humidity, and a 12:12-h light: dark cycle.

Tail Injury Surgery

Surgery was performed under anesthesia using Ketamine (90 mg/kg) and Xylazine (10 mg/kg) subcutaneously. The skin was sterilized with Betadine, coccygeal IVDs were located, and a 26G needle was inserted into the IVD space until the needle tip reached approximately 2/3 of the disc thickness, as described previously [2]. Four consecutive levels of coccygeal (Co; 3/4, 4/5, 5/6, and 6/7) IVDs were injured with a 26 Gauge needle. We aimed to injure the AF on the needle entry side and NP, since injuring AFs on both sides caused a more severe injury [24]. For RNA extractions, the “injured” group underwent IVD injuries as described above, and “intact” control mice had no injury. For Von Frey tests, the “sham injury” mice had injury to skin and subcutaneous tissues, without penetrating the IVDs. The mice were euthanized using carbon dioxide (CO2) inhalation in accordance with American Veterinary Medical Association (AVMA) and IACUC guidelines.

Immunohistochemical Staining

The histological studies used paraffin-embedded tissues archived in our laboratory from mice on the B6 background. The histological preparation, including immunostaining, hematoxylin and eosin (H&E), and Safranin O/fast green staining, was performed by Thomas Jefferson University Sidney Kimmel Cancer Center Translational Research/Pathology Shared Resource Core Facility, with methods described previously [13, 14]. Briefly, immunohistochemical staining was performed with F4/80 (EMR1), a well-established antigen marker for mouse macrophages [25, 26]. The F4/80 antigen was examined with rabbit monoclonal antibody (Cat#: 70076, Cell Signaling Technologies, Danvers, MA). Specifically, the 4 μm paraffin sections were deparaffinized using the Shandon Varistain Gemini ES Autostainer (Ramsey, MN). Antigen retrieval was performed using 0.025% trypsin with EDTA (Cat#: 25-063-CI, Corning Life Sciences, Tewksbury, MA) for 15 min at 37°C. Immunohistochemistry staining was then performed in an intelliPATH FLX Autostainer (Biocare Medical, Pacheco, CA) with the F4/80 antibody at a 1/1000 dilution, followed by biotinylated secondary antibody (goat anti-rabbit immunoglobulin; Cat#: BA-1000, Vector Laboratories, Newark, CA) and ABC-HRP complexes (Cat#:PK6100, Vector Laboratories, Newark, CA) and signals were visualized using DAB substrate (Cat#: SK-4103, Vector Laboratories, Newark, CA). Following washing steps, slides were then counterstained with Hematoxylin, dehydrated, cleared, and coverslips applied.

Quantification of Safranin O Staining

Tissues were fixed with 4% paraformaldehyde, embedded in paraffin, sectioned to 5 μm, and stained with Safranin O as described previously [13]. The sections were digitized using the ZEISS Axionscan Z1 slide scanner (Oberkochen, Baden-Württemberg, Germany), and analyzed by ImageJ software (NIH Image). Briefly, the IVD regions were cropped and converted to HSB Stack. Red color was defined as hue color 220–255 and 0–32. The pixel number of red and the entire color spectrum (hue color 0–255) were obtained, and the percentage of pixel number of red in the entire color spectrum was calculated.

RNA Isolation

For RNA extraction, IVD tissues were collected from each injured IVD 24 h after the surgery, or from intact controls. Using 4 IVDs from 1 mouse tail yielded sufficient high-quality RNA for this study. Briefly, Co3/4, 4/5, 5/6, and Co6/7 IVDs were shaved off the cartilaginous endplate with a scalpel under a dissecting microscope (VistaVision, VWR International, Radnor, PA) and immersed in RNALater (Ambion, Foster City, CA) overnight and preserved at −80°C until extraction. Total cellular RNA was isolated by the Trizol method. Specifically, RNALater was eliminated, and the tissues were rapidly frozen with liquid Nitrogen before being transferred into Trizol (Invitrogen, Carlsbad, CA). The tissues were then homogenized using a homogenizer equipped with disposable OmniTip probes designed for hard tissue (Omni International, Kennesaw, GA). RNA was then precipitated with 70% ethanol and further purified using an RNeasy Micro Kit (Qiagen). The RNAs from all 4 discs in each mouse were pooled for RNA library preparation.

RNASeq

RNA library preparation and sequencing was performed by Azenta US (South Plainfield, NJ). Specifically, RNA samples were quantified using a Qubit 2.0 Fluorometer (ThermoFisher Scientific, Waltham, MA, USA) and RNA integrity was checked with 4200 TapeStation (Agilent Technologies, Palo Alto, CA, USA). This step was carried out by Azenta, the company that performed the RNASeq. Briefly, small amounts of samples are separated in the micro-fabricated chip channels according to their molecular weight and then examined by laser detection. The result is an electropherogram in which the amount of changed fluorescence correlates with the amount of RNA of a given size. The ratio of two ribosomal bands was calculated. The RNA integrity number (RIN) measurement is based on a machine learning algorithm that uses a capillary electrophoresis pathway and not just on the ratio of ribosomal subunits, although it is highly dependent on the ratio. RIN provides a numerical score (range 1–10) for RNA quality. A higher RIN value indicates a higher degree of RNA integrity. When the RIN values were deemed acceptable, rRNA depletion was performed to further purify mRNA before the RNASeq procedure. RNA sequencing library preparation used the NEBNext Ultra II RNA Library Preparation Kit for Illumina per the manufacturer's recommendations (NEB, Ipswich, MA, USA). Briefly, enriched RNAs were fragmented for 15 min at 94°C. First strand and second strand cDNA were subsequently synthesized. cDNA fragments were end-repaired and adenylated at their 3'ends, and universal adapters were ligated to cDNA fragments, followed by index addition and library enrichment with limited cycle PCR. Sequencing libraries were validated using the Agilent Tapestation 4200 (Agilent Technologies, Palo Alto, CA, USA), and quantified using Qubit 2.0 Fluorometer (ThermoFisher Scientific, Waltham, MA, USA) as well as by quantitative PCR (KAPA Biosystems, Wilmington, MA, USA). The sequencing libraries were then multiplexed and clustered on one flowcell. After clustering, the flowcell was loaded on the Illumina HiSeq instrument according to the manufacturer's instructions. The samples were sequenced using a 2 × 150 Pair-End (PE) configuration. Raw sequence data (.bcl files) generated from Illumina HiSeq were converted into fastq files and de-multiplexed using the Illumina bcl2fastq program version 2.20. One mismatch was allowed for index sequence identification.

Von Frey Testing to Assess Mechanical Sensitivity

A cohort of mice separate from those used for RNASeq was tested using the method described by Chaplan et al. [27] Briefly, the tail mechanical sensitivity after IVD or sham injury was evaluated in mice 1, 7, 14, and 28 days after surgery using Von Frey filaments as described previously [28]. Individual mice were placed on a wire-mesh platform (IITC Life Sciences Inc.) under a 4 × 3 × 7 cm clear acrylic cage to limit movement. Mice were habituated to the testing conditions 30 min before the first trial. During the test, a set of von Frey fibers (Stoelting Touch Test Sensory Evaluator Kit #2 to #9, ranging from 0.02- to 1 g force) was applied to the ventral surface of the tail, targeting the area below the injured IVD or sham injured area until the fibers bowed and were then held for 3 s. The threshold force required to elicit tail withdrawal (median, 50% withdrawal) was determined twice on each mouse tail with measurement intervals of 20 min.

Statistical Analysis For 50% Withdrawal Threshold

The difference in 50% withdrawal threshold was determined with a 2-factor analysis of variance (ANOVA) in repeated measures, where injured/sham was a grouping factor and time was the repeated measure. Data were expressed as mean ± standard deviation (SD). A p-value of < 0.05 was considered statistically significant. All analyses were performed using SAS statistical software (Version 9.4, SAS Institute, Cary, NC).

Data Analysis

After demultiplexing, sequence data were checked for overall quality and yield. Then, sequence reads were trimmed to remove possible adapter sequences and nucleotides with poor quality using Trimmomatic v.0.36. The trimmed reads were mapped to the reference genomes using the STAR aligner v.2.5.2b. The STAR aligner is a splice-aware aligner that detects splices junctions and incorporates.

Junctions and incorporates them to help align the entire read sequences. BAM files were generated as a result of this step. Unique gene hit counts were calculated by using feature counts from the Subread package v.1.5.2. Only unique reads within exon regions were counted. After extraction of gene hit counts, the gene-counts table was used for downstream differential expression analysis. Using DESeq2, a comparison of gene expression between the groups of samples was performed. The Wald test was used to generate p-values and Log2 fold changes, which were corrected for multiple testing using the Benjamini and Hochberg method. Genes with adjusted p-values (p.adj) < 0.01 were called as differentially expressed genes for each comparison.

Generation of Heatmap

A total of 17,722 genes had meaningful levels of expression by RNASeq. A differential gene list was generated by limiting P.adj to < 0.01, resulting in 7242 genes in this category. Then genes with no expression values detected in two or more samples were excluded, resulting in 7184 differentially expressed genes. The genes were further sorted according to log2 fold changes from high to low, and separated into upregulated genes (1905 genes; injured/intact ratio > 2; log2 > 1) and downregulated genes (1872 genes; intact/injured ratio > 2; log2 > 1). A heatmap was generated using 19 upregulated and 8 downregulated genes selected based on their highest fold changes using the R package pheatmap [29]. Rows (represent gene expression) were scaled and hierarchically clustered, and a gap was inserted among columns to differentiate the intact control and injured groups.

Gene Ontology (GO) Analysis and Protein–Protein Interaction (PPI) Network Visualization

Lists of upregulated and downregulated genes (ratio > 2; log ratio > 1) were analyzed separately for GO analysis, and biological pathways were visualized with the ggplot2 software [30]. Because of the exceptionally large numbers of genes that changed in response to injury, the gene list was further narrowed down, using ratio > 4; log2 > 2, resulting in 703 upregulated genes (injured/intact ratio > 4 fold) and 304 downregulated genes (intact/injured ration > 4 fold). The resulting gene list was imported into the STRING database (, Version: 11.5), and a PPI file was generated. The genes with a high confidence score (0.70) were imported into CytoScape software (Version: 3.9.1) [31]. Note that the confidence scores are scaled between 0 and 1, with 1 corresponding to the estimated likelihood of a given association being true. The CytoHubba plugin was used to rank nodes with the MCC (Maximal Clique Centrality) method [32]. The parameters of nodes and edges were adjusted based on node rankings and the combined interaction scores, thus visualizing PPI networks in a more intuitive and clear way [33]. Disconnected nodes in the network were not displayed, and the 25 genes with the largest number of PPI network connections were displayed in the graph.

Results

Proteoglycan Loss and Macrophage Infiltration Into the Mouse Coccygeal (Co) Intervertebral Disc (IVD)

There was rapid proteoglycan loss in the nucleus pulposus 24 h post injury (Figure 1A,B). The injured and intact IVDs were immunostained for F4/80 antigen, a classical marker for macrophages. Macrophages were easily detectable at day 1 post injury (n = 6 mice; Figure 1A′′,B′′). There were a few F4/80 positive cells in the outermost layer of the AF in the intact IVDs (Figure 1A′′,B′′). Interestingly, even in the injured discs, F4/80 positive cells were present only in the outermost layers of the annulus fibrosus.

[IMAGE OMITTED. SEE PDF]

Gene Expression Profiles of Injured and Intact IVDs were Clearly Different, and Tnfa Gene Expression had the Highest Number of Connections in the Protein–Protein Interaction PPI Network

Transcriptomes of injured and intact IVDs 24 h post injury were compared. Among the 17,722 genes with a meaningful level of expression by RNASeq, 7242 genes were differentially expressed (n = 4 mice/group; P.adj < 0.01). Among the differentially expressed genes, 1905 genes were upregulated > 2 fold in injured discs compared with controls, while 1872 genes were downregulated by > 2 fold. A heatmap was generated using 19 upregulated and 8 downregulated genes selected based on their highest fold changes, using the R package pheatmap [29]. The gene expression profiles of injured and intact IVDs are visibly different (Figure 2A). Importantly, Tnfa had the highest number of PPI network connections (Figure 2B), suggesting that regulating this gene is highly relevant to the tissue response to injury.

[IMAGE OMITTED. SEE PDF]

Leukocyte Migration and Positive Regulation of Cell Adhesion Pathways Were Enriched Among Upregulated Genes in Injured IVDs

Upregulated gene with injury (injured/intact ratio > 2; log ratio > 1; a total of 1905 genes) were inputted into the ggplot2 software. Pathways related to “leukocyte migration” and “positive regulation of cell adhesion” were enriched among upregulated genes in the injured discs (Figure 3A).

[IMAGE OMITTED. SEE PDF]

Pattern Specification and Sodium Transmembrane Transport Pathways Were Enriched Among Downregulated Genes in Injured IVDs

Downregulated genes post injury (intact/injured ratio > 2; log ratio > 1; a total of 1872 genes) were inputted into the ggplot2 software. Pathways related to “pattern specification process” and “sodium ion transmembrane transport” were enriched among downregulated genes in the injured discs (Figure 3B).

Injured Mouse Tails Develop Mechanical Allodynia

To assess mouse pain behavior post injury, the degree of mechanical allodynia under the tail IVDs was determined. Mice were randomly assigned to IVD injury (IVDs punctured) or sham injury (only skin was punctured) groups. The Von Frey test was performed at days 1, 7, 14, and 28 post injury. The 50% mechanical sensitivity threshold was significantly decreased compared with sham controls at days 14 and 28 post injury (n = 6 mice/group; p = 0.029 and 0.003, respectively). Within the IVD-injury group, the mechanical sensitivity threshold decreased significantly at days 14 and 28 compared with day 1 (n = 6 mice; p < 0.001). At day 1 post injury, there was no significant difference in mechanical sensitivity between the IVD-injured and sham groups, likely due to residual effects of anesthetics (Figure 4).

[IMAGE OMITTED. SEE PDF]

Discussion

In this study, RNASeq was performed 1 day after injury, when significant proteoglycan loss is evident [13]. The early events post injury provide mechanistic insights, and suggest that treatments to reduce inflammation and proteoglycan loss should be started immediately following injury. Among 17 722 genes with meaningful level of expression by RNASeq, 7242 genes were expressed at different levels in the injured disc model compared with intact controls. It is reassuring that genes previously used as inflammatory markers (e.g., Cxcl1, Tnfa) have been shown to upregulate in injured discs here, although most of the genes perturbed by injury have not been noted previously in similar animal disc injury models. “Leukocyte migration” was the most enriched biological process, and “pattern specification process” was the most enriched pathway among downregulated genes. TNFa had the largest number of PPI connections. Targeting the pathways and highly connected molecules may help to mitigate the degenerative cascade following disc injury. The profiles of differentially expressed genes in injured and intact discs presented here have many similarities compared with previously published work, on IVDs of mice on the DBA background, 1-week post injury [34]. There are some differences in the fold changes and biological pathways revealed by the two studies. Since both genetic background and time after injury affect gene expression profiles [35], future studies directly comparing the transcriptomes of mice on the same genetic background are indicated.

Macrophages have been identified in human degenerative IVDs, suggesting that these cells play a role in the pathophysiology of disc degeneration and back pain [36, 37]. In the mouse IVDs, macrophages enter the outermost layers of injured disc tissue early, likely directed by chemokines/cytokines released by the injured tissue. Single-cell RNASeq (scRNAseq) has been used to identify different cell types in the human and rat IVD, including NP cells, AF cells, immune cells, and progenitor cells [9, 38, 39]. In our study, in search of molecular markers of inflammation post-injury, we performed bulk RNASeq and measured the gene expression levels of the IVD tissues. Consistent with findings in humans and rats, we found “leukocyte migration” as the most enriched biological pathway by GO analysis.

There was no mechanical sensitivity by Von Frey at day 1 post injury, despite the gene expression changes and macrophage infiltration into the injured areas. This is likely due to the residual effects of anesthesia or the time needed for inflammation to spread from the IVDs to the ventral aspect of the tail. Mechanical allodynia was evident at 1–2 weeks post injury and persisted until the end of the study (28 days post injury). This is consistent with the previous observations that macrophages remained in the discs and inflammatory markers such as Cxcl1 and Adam8 remain elevated at week 4 post disc injury [14].

Genes with no expression values in two or more samples (a total of 58) were not included in the heatmap and network analysis. These genes were excluded because their expression values at baseline were very low, so the fold changes are unreliable. However, some genes (e.g., Clec4d, Il24 and Cxcr2) that were undetectable in intact discs but induced by injury may be biologically important and warrant further study. One limitation of the current study is that only young adult male mice on the C57B6 background were used. Since genetic background, sex [35], and age [40] are known to affect gene expression and behavior [41], future studies will include female mice to determine the influence of sex on the molecular markers and mouse behavior in the tail disc injury model. The advantage of the tail disc injury model is that less collateral damage occurs compared with a lumbar IVD injury procedure, which often breaches the abdominal cavity in an animal as small as a mouse. The injury model used in the present work is well established, having been the basis for numerous informative published studies of IVD pathophysiology during recent years. In the future, additional behavior analysis is warranted to confirm the relevance of the Von Frey assay to discogenic pain.

Another limitation is that the mouse tail IVD injury reflects only some aspects of the human disease since the human disc injury often results from a compression or rotation type of mechanism. The tail IVD is substantially different in load-bearing from that of the human lumbar spine. Caution must be applied in seeking to generalize the findings to other species or different types of IVD injury.

The RNASeq data revealed numerous early genes that participate in inflammation and repair processes post IVD injury. Some of these genes may play a role in the transition to chronic IVD inflammation and pain. These findings need to be compared with findings in human tissues, and the effects of altered gene expression on histology, biomechanical properties of the IVDs, and pain need to be studied further in animal models and explored in human patients.

Conclusions

With a nonbiased approach, numerous novel genes were found that could serve as markers in future mechanistic studies. The leukocyte migration pathway was the most enriched biological process, and Tnfa had the most network connections among upregulated genes. Mechanical allodynia persisted for at least 28 days, suggesting that acute injury may lead to chronic pain. Treatments targeting key pathways/molecules may be effective post-injury on disc inflammation and repair.

Clinical Significance

Novel genes highly regulated on day 1 post disc injury were identified with a nonbiased approach; they may serve as biomarkers of injury in future experiments. Enriched biological pathways and molecules with high numbers of connections may be targets for treatments post injury. Von Frey test established the relationship between disc injury/degeneration, and pain.

Acknowledgments

ZT performed mouse surgeries, collected tissues, and extracted RNA; BC performed Von Frey assays, HW, KC, BC, WO, LQ, and YZ performed data analysis. The Neurobehavior Testing Core at the University of Pennsylvania (Penn)/Institute for Translational Medicine and Therapeutics (ITMAT) and Intellectual and Developmental Disabilities Research Center (IDDRC) at the Children's Hospital of Philadelphia (CHOP)/Penn P50HD105354 assisted in IVD injury and behavior procedures. All authors made a significant contribution to the work reported, whether that is in the conception, study design, execution, acquisition of data, analysis/interpretation, or in all these areas; took part in drafting, revising, or critically reviewing the article; gave final approval of the version to be published; have agreed on the journal to which the article has been submitted; and agree to be accountable for all aspects of the work. All authors have read and approved the final submitted manuscript. We gratefully thank Dr. Martin F. Heyworth, MD, for critically editing the manuscript.

Data Availability Statement

All data are available upon request to the corresponding author (Zhang).

References

V. Nicol, C. Verdaguer, C. Daste, et al., “Chronic Low Back Pain: A Narrative Review of Recent International Guidelines for Diagnosis and Conservative Treatment,” Journal of Clinical Medicine 12 (2023): [eLocator: 1685], [DOI: https://dx.doi.org/10.3390/jcm12041685].

K. Vincent, S. Mohanty, R. Pinelli, et al., “Aging of Mouse Intervertebral Disc and Association With Back Pain,” Bone 123 (2019): 246–259.

G. E. Mosley, M. Wang, P. Nasser, et al., “Males and Females Exhibit Distinct Relationships Between Intervertebral Disc Degeneration and Pain in a Rat Model,” Scientific Reports 10 (2020): 15120–15129.