Identification of Key Genes Potentially Related to Intervertebral Disk Degeneration by Microarray Analysis

Abstract

Aims:

This study was designed to investigate differentially expressed genes (DEGs) in the annulus fibrosus (AF), nucleus pulposus (NP), and whole blood (WB) of intervertebral disk degeneration (IDD) patients.

Materials and Methods:

We retrieved microarray data set GSE70362, which contains the gene expression profiles of 24 AF and 24 NP samples from the Gene Expression Omnibus and identified DEGs in degenerative AF (AF-DEGs) and NP (NP-DEGs) samples compared with nondegenerative samples. We also examined gene expression profiles in WB from patients with IDD and healthy volunteers to identify DEGs in WB (WB-DEGs). We performed functional analyses on the DEGs common to AF-DEGs, NP-DEGs, and WB-DEGs. Expression of the common DEGs was partially validated by quantitative real-time-polymerase chain reaction (QRT-PCR).

Results:

In total, 846 AF-DEGs, 902 NP-DEGs, and 862 WB-DEGs were identified, and 22 DEGs were common among the three groups. Functional analyses showed that the common DEGs were enriched in 33 biological processes, 16 cellular components, 4 molecular functions, and 9 Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways; 13 of the common DEGs were included in the protein-protein interaction (PPI) network and superoxide dismutase 2 (SOD2) was identified as a hub gene in the PPI network. The QRT-PCR results for the expression of the genes protein disulfide isomerase family A member 4, FKBP prolyl isomerase 11, ectonucleotide pyrophosphatase/phosphodiesterase 4, SOD2, and actin binding LIM protein 1, were consistent with the gene chip hybridization results.

Conclusions:

This study identified key genes for future investigations of the underlying molecular mechanisms of IDD. These genes may provide future targets for the clinical treatment and diagnosis of IDD.

Introduction

Low back pain is a common public health issue that seriously affects quality of life and imposes an enormous economic burden (Andersson 1999; Luoma et al., 2000; Katz 2006; Martin et al., 2008). Extensive research exploring low back pain has focused on intervertebral disk (IVD) degeneration (Luoma et al., 2000; Carragee et al., 2006). The prevalence of IVD increases with age, and up to 90% of asymptomatic individuals aged >60 years exhibit intervertebral disk degeneration (IDD) (Brinjikji et al., 2015). Various studies have sought to identify the pathophysiology of IDD, which is thought to be multifactorial and involve genetic, mechanical, and biological factors. However, the exact mechanism underlying IDD has not been fully elucidated.

In 2015, Kazezian et al. used a gene expression array to acquire gene profiles of the annulus fibrosus (AF) and the nucleus pulposus (NP) that could be used for downstream analyses (Kazezian et al., 2015). In this study, to obtain insight into the molecular mechanism of IDD, we identified differentially expressed genes (DEGs) in the degenerative AF (AF-DEGs) and NP (NP-DEGs) compared with nondegenerative tissues based on microarray data. We also collected whole blood (WB) from eight IDD patients and eight healthy volunteers and performed microarray analysis to identify DEGs in WB (WB-DEGs). Quantitative real-time polymerase chain reaction (QRT-PCR) was conducted using WB from the same eight IDD patients and eight healthy volunteers to measure the expression levels of some DEGs common among the AF-DEGs, NP-DEGs, and WB-DEGs.

Materials and Methods

IVD microarray data

The microarray data set GSE70362 based on the Affymetrix GPL17810 platform [HG-U133_Plus_2] was downloaded from the Gene Expression Omnibus database. This data set contains gene expression data for postmortem tissues of 24 AF and 24 NP samples from 19 donors between the ages of 21 and 82 years. All samples were collected from T12/L1 to L4/L5 within <12 h after brain death was declared.

The original signal intensity data in CEL format files and classification based on the Thompson grading system (Thompson et al., 1990) were downloaded for analysis in our study. To annotate the data, we transformed the original probe IDs into gene symbols. Any probes not identified were discarded. For genes with multiple probes, only the probe with the highest fold change (FC) was retained. Within each set of 24 samples (AF or NP), 8 samples classified as Thompson grade I or Thompson grades I and II were considered nondegenerative, and 16 samples classified as Thompson grades II-V were considered degenerative.

Identification of AF-DEGs and NP-DEGs

GeneSpring GX 11.5 software (Agilent Technologies, CA) was employed to identify the DEGs. The initial data were preprocessed by using the Robust Multiarray Average procedure. Paired t-tests were performed with a threshold p-value of 0.05 and FC of 1.5 to identify DEGs between degenerative and nondegenerative AF and NP samples. For genes with multiple probes, only the probe with the highest FC was retained. The study using the IVD microarray data was performed between January 2018 and March 2018.

Human WB collection

The WB samples were collected between April 2018 and August 2018 at Sichuan Provincial Orthopedic Hospital, China. We recruited 8 patients, including 4 men and 4 women aged between 33 and 60 years with IDD confirmed by magnetic resonance imaging, and 8 volunteers, including 4 men and 4 women aged between 19 and 23 years with no clinical evidence of low back pain or sciatica. In total, 10 mL fasting WB were collected from the left median cubital vein between 7:00 and 7:30 AM from each participant between 7:00 and 7:30 AM. All WB samples were immediately incubated in a PAX gene Blood RNA tube (BD) for <72 h at −20°C and sent to Shanghai Bohao Biotechnology Co., Ltd. (Shanghai, China) for gene chip hybridization screening.

Microarray analysis of WB and identification of WB-DEGs

The microarray experiments were performed using an Agilent SurePrint G3 human gene expression microarray 8 × 60 K on an Agilent Microarray Scanner platform at Shanghai Biotechnology Co., Ltd. following the standard protocol of Agilent Technologies. The raw data obtained by the chip scan were normalized by the limma package in R Language. The Quantile algorithm was used, and the signal value was log2 normalized. Student's t-test was performed with a threshold p-value of 0.05 and FC of 1.5 to identify WB-DEGs. For genes with multiple probes, only the probe with the highest FC was retained. The complete data sets of the gene expression profiles in WB were uploaded to the GEO database and are accessible under GSE124272.

Functional analysis of the common DEGs

The enrichment analysis included Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment. The method adopted in this study was Fisher's exact test, and the clusterProfiler package from statistical environment R for Bioconductor was used to perform the enrichment analysis. The criteria selected were the number of genes that were enriched in a certain term/GO ≥2 and a p-value <0.05. The common DEGs were mapped to the Search Tool for the Retrieval of Interacting Genes (STRING website) using the default instructions. A protein-protein interaction (PPI) network of nodes with a combined score >0.15 was constructed and then visualized in Cytoscape software (version 3.7.1). Degree centrality was calculated by the plug-in CentiScaPe.

Total RNA extraction, cDNA synthesis, and QRT-PCR of common DEGs

QRT-PCR was conducted using WB to measure the expression of some DEGs common among the AF-DEGs, NP-DEGs, and WB-DEGs. Total RNA was extracted and purified using a PX Blood RNA Kit (OMEGA) following the manufacturer's instructions. Each RNA sample was reverse transcribed into cDNA using a ReverTra Ace qPCR Kit (TOYOBO, Japan), and the cDNA was subsequently used to conduct QRT-PCR. Quantification of gene expression levels was performed using a 7500 HT Sequence Detection System (ABI) with Power SYBR Green PCR Master Mix (ABI). The primers for the specific genes were designed by Primer Express (V3.0.1) and synthesized by Sangon Biotech Co. (Shanghai, China). The β-actin gene was used as an internal control, and gene expression levels were normalized to β-actin using the 2^−ΔCt method.

Statistical analysis

QRT-PCR data are presented as the mean ± SD. The statistical analysis was performed by using Student's t-tests with GraphPad Prism 6.0 (GraphPad Software, Inc., San Diego, CA). p-Values <0.05 were considered statistically significant.

Ethics statement

All procedures described in this study were authorized by the Ethics Committee of the Sichuan Provincial Orthopedic Hospital, and all patients and volunteers enrolled in the study provided written informed consent to participate in this study.

Results

Differentially expressed genes

AF-DEGs, NP-DEGs, and WB-DEGs were identified as shown in Table 1. Figure 1A-C shows hierarchical heat maps of the DEGs in AF, NP, and WB. Figure 1D shows the 22 DEGs that overlapped among the three groups and were identified as common DEGs.

FIG. 1.

Hierarchical clustering analysis of (A) AF-DEGs, (B) NP-DEGs, and (C) WB-DEGs. (D) Venn diagram illustrating the numbers of exclusive and overlapping DEGs among the AF-DEGs, NP-DEGs, and WB-DEGs. AF, annulus fibrosus; DEGs, differentially expressed genes; NP, nucleus pulposus; WB, whole blood.

Table 1.

Number of Total DEGs in the AF, NP, and WB Samples, the Proportion of Total DEGs in the Transcriptome, the Number of Upregulated DEGs, and the Number of Downregulated DEGs

	Total DEGs	Proportion of transcriptome, %	Upregulated	Downregulated
AF-DEGs	846	1.58	455	391
NP-DEGs	902	1.65	305	597
WB-DEGs	862	1.48	484	378

AF, annulus fibrosus; DEGs, differentially expressed genes; NP, nucleus pulposus; WB, whole blood.

Functional analysis of the common DEGs

Based on GO enrichment, 22 common DEGs were enriched in 33 biological processes (BP), 16 cellular components, and 4 molecular functions, only 4 of which had a p-value <0.05 (Table 2). The common DEGs were enriched in a total of nine KEGG pathways, none of which had a p-value <0.05, as shown in Table 3. Thirteen of the common DEGs were included in the PPI network, and superoxide dismutase 2 (SOD2) had the highest degree centralization calculated by the plug-in CentiScaPe (Fig. 2).

FIG. 2.

The PPI network of the common DEGs. All nodes with a combined score >0.15 are shown. Disconnected nodes in the network were excluded. The size of the node represents the degree centralization. PPI, protein-protein interaction.

Table 2.

The 4 Gene Ontology Terms for 22 Common Differentially Expressed Genes with a p-Value <0.05

GO ID	GO term	Category	p	Gene sets
GO:0007010	Cytoskeleton organization	BP	0.01	FRY, ABLIM1, SSH2, AVIL
GO:0031175	Neuron projection development	BP	0.02	FRY, ABLIM1, SSH2, AVIL
GO:0048666	Neuron development	BP	0.03	FRY, ABLIM1, SSH2, AVIL
GO:0009611	Response to wounding	BP	0.04	TREM1, ENPP4, MRVI1, TNFAIP6

The “Gene sets” column lists the gene symbols of the DEGs enriched in the corresponding GO term.

ABLIM1, actin binding LIM protein 1; BP, biological processes; ENPP4, ectonucleotide pyrophosphatase/phosphodiesterase 4; GO, Gene Ontology.

Table 3.

KEGG Pathway Enrichment Analyses for 22 Common Differentially Expressed Genes

Pathway ID	Description	p	Gene sets
hsa04060	Cytokine-cytokine receptor interaction	1.00	ACKR3
hsa04146	Peroxisome	1.00	SOD2
hsa04270	Vascular smooth muscle contraction	1.00	MRVI1
hsa04810	Regulation of actin cytoskeleton	1.00	SSH2
hsa05110	Vibrio cholerae infection	1.00	PDIA4
hsa04360	Axon guidance	1.00	ABLIM1
hsa05016	Huntington's disease	1.00	SOD2
hsa04141	Protein processing in endoplasmic reticulum	1.00	PDIA4
hsa04977	Vitamin digestion and absorption	1.00	SLC19A3

The “Gene sets” column lists gene symbols of the DEGs enriched in the corresponding KEGG pathway.

KEGG, Kyoto Encyclopedia of Genes and Genomes; PDIA4, protein disulfide isomerase family A member 4; SOD2, superoxide dismutase 2.

QRT-PCR of common DEGs

QRT-PCR was conducted using WB to verify the gene chip hybridization results of protein disulfide isomerase family A member 4 (PDIA4), FKBP prolyl isomerase 11 (FKBP11), ectonucleotide pyrophosphatase/phosphodiesterase 4 (ENPP4), SOD2, and actin binding LIM protein 1 (ABLIM1). The QRT-PCR results using WB were consistent with the gene chip hybridization results, suggesting that the microarray data correlated well with the QR-PCR data, demonstrating the reliability of the microarray results (Fig. 3A-E).

FIG. 3.

QRT-PCR results were consistent with the gene chip hybridization results. All QRT-PCR experiments were independently repeated in triplicate, and the average Ct value was used to calculate the ΔCt. The relative expression level is represented by 2^−ΔCt. P: patients with IDD; V: healthy volunteers. (A) PDIA4; (B) FKBP11; (C) ENPP4; (D) SOD2; (E) ABLIM1 (n = 8, *p < 0.05). ABLIM1, actin binding LIM protein 1; ENPP4, ectonucleotide pyrophosphatase/phosphodiesterase 4; FKBP11, FKBP prolyl isomerase 11; IDD, intervertebral disk degeneration; PDIA4, protein disulfide isomerase family A member 4; QRT-PCR, quantitative real-time polymerase chain reaction.

Expression of PDIA4, FKBP11, and ENPP4 was downregulated in WB from IDD patients compared with that in healthy volunteers, whereas expression of SOD2 in WB from IDD patients was significantly higher than that in healthy volunteers. Furthermore, the results showed that the gene expression levels of PDIA4, FKBP11, ENPP4, and SOD2 significantly differed between the patient group and volunteer group (p < 0.05). All forward and reverse primers used for PCR are listed in Table 4.

Table 4.

Sequences of Primers Used for Quantitative Real-Time Polymerase Chain Reaction

Gene	Sequence (5′ to 3′)
PDIA4	F: AATTTGAGGGTGGAGACAGAGATC
PDIA4	R: CTCAGTTTTGTGGCATGTTCTTCT
FKBP11	F: GCCAAAAGCAGGTGATTCCA
FKBP11	R: TCGCTTCTCTCCCACACACA
ABLIM1	F: TGAAACTGAGCAGGTGTGTGTTG
ABLIM1	R: CCACCTCCTCCTCTACACACAAG
ENPP4	F: TCCTTCTTGTGCCCATATCTGA
ENPP4	R: GAAGTAAAACCCTCCAAATAAAAGGAA
SOD2	F: GATCCACTGCAAGGAACAACAG
SOD2	R: GCGTGCTCCCACACATCA

F, forward; FKBP11, FKBP prolyl isomerase 11; R, reverse.

Discussion

Numerous studies have revealed that IDD is a multifactorial physiological manifestation (Adams and Roughley, 2006). The IVD consists of the AF on the outside and the NP on the inside, and these two components have different cell populations and extracellular matrix compositions. Previous studies have mostly focused on DEGs between the degenerative and nondegenerative AF or NP or whole IVD tissues, which can only be obtained through surgery. In contrast, no studies have investigated gene expression in WB from patients with IDD.

In this study, we identified 846 AF-DEGs, 902 NP-DEGs, and 862 WB-DEGs, with 22 overlapping DEGs. The expression levels of PDIA4, FKBP11, ABLIM1, ENPP4, and SOD2 in WB were validated by QRT-PCR, and these genes may be future potential targets for the clinical treatment and diagnosis of IDD.

Various studies have been conducted on the topic and are summarized in Table 5. Schubert et al. (2018) identified 5 AF markers and 6 NP markers by genome-wide expression analysis, but dysregulation of these 11 markers was not observed in WB in our study. We previously compared the gene expression profiles of degenerative AF and degenerative NP samples and identified 87 DEGs between them (Wang et al., 2018).

Table 5.

Relevant Studies on Gene Expression in Intervertebral Disk Degeneration

References	Type of samples	Species	Altered markers	Methods
Schubert et al. (2018)	Cervical disk	Human	5 AF markers 6 NP markers	Histological analysis Microarray analysis Real-time PCR
Wang et al. (2018)	AF and NP of a lumbar disk	Human	87 DEGs between degenerated AF and NP	Microarray analysis
Guo et al. (2017)	Lumbar disk	Human	35 common degeneration markers (AF/NP)	Microarray analysis
Kazezian et al. (2015)	AF of a lumbar disk	Human	17 AF degeneration markers	Microarray analysis Real-time-PCR analysis Immunohistochemistry
Zhang et al. (2010)	Peripheral blood mononuclear cells	Human	62 DEGs	Microarray analysis QRT-PCR
Sherafatian et al. (2019)	Intervertebral disk	Human	5 common DEGs among 3 expression data sets	Meta-analysis of expression data sets

Guo et al. (2017) compared the gene expression profiles of degenerative IVDs with those of nondegenerative IVDs and identified 35 DEGs that overlapped between AF-DEGs and NP-DEGs. Among these 35 common DEGs, only zinc finger protein 185 with the LIM domain was included in our common DEGs. Using microarray analysis, Kazezian et al. (2015) identified 17 molecular markers with FC >1.5 in AF degeneration, including insulin-like growth factor binding protein 3, gremlin 1, and phorbol-12-myristate-13-acetate-induced protein 1. None of these molecular markers were observed in our AF-DEGs or common DEGs.

Zhang et al. (2010) compared expression arrays of peripheral blood mononuclear cells from IDD patients with those from non-IDD individuals and identified 62 DEGs, including 39 upregulated genes and 24 downregulated genes. None of these DEGs overlap with our common DEGs. Among their 39 upregulated DEGs, annexin A3, cathelicidin antimicrobial peptide, and interleukin 1 beta were upregulated in WB from patients with IDD in our study, but cyclin B1 and myomesin (M-protein) 2 were downregulated. None of their 24 downregulated genes were observed among our WB-DEGs.

In addition, some studies identified microRNAs (miRNAs) as molecular markers for IDD. For example, a meta-analysis of three microarray studies identified five miRNAs as biomarkers for IDD (Sherafatian et al., 2019). In all of these studies, heterogeneity derived from age, gender, the level of IVD, medical circumstances, specimen collection, gene chip technology, and so on is not negligible. Further research on accurate markers for IDD is necessary.

IDD is an age-related process that is accompanied and accelerated by the accumulation of oxidative damage and inflammation (Bank et al., 1998; Sivan et al., 2006; Wang et al., 2016). Blood vessels grow into the IVD during the process of IDD (Johnson et al., 2001; Ali et al., 2008), and the resulting neovascularization exposes the avascular tissue of the IVD to high oxygen tension and subsequent oxidative injury (Nasto et al., 2013).

SOD2 plays an antioxidative role by converting superoxide radical into H₂O₂, which can be broken down into harmless H₂O and O₂ by antioxidative enzymes (Flekac et al., 2008; Chen et al., 2009) and is critical for the regulation of oxidative stress resistance (Landis and Tower, 2005). Excessive reactive oxygen species may cause protein, DNA, and membrane damage and are associated with cellular inflammatory responses by inducing the expression of cytokines and chemokines (Nel et al., 2006; Li et al., 2008; Park and Park, 2009). Thus, elimination of superoxide radical by SOD2 can be considered an anti-inflammatory process.

However, increased expression levels of many inflammatory cytokines, including tumor necrosis factor α, interferon-gamma, and interleukin-1, -2, -4, -6, -8, -10, and -17, have been found in degenerative IVD (Shamji et al., 2010), and the inflammatory response plays a critical role in the initiation and progression of IDD (Tian et al., 2013; Wang et al., 2013; Risbud and Shapiro, 2014; Molinos et al., 2015). SOD2 may be protective against IVD through antioxidant and anti-inflammatory mechanisms in the process of IDD. An age-related reduction in SOD2 expression in the NP has been observed in mice (Alvarez-Garcia et al., 2017). Data regarding the biological actions of PDIA4, FKBP11, ABLIM1, and ENPP4 in IDD are limited, and their roles in IDD require further investigation.

In summary, our study identified PDIA4, FKBP11, ENPP4, SOD2, and ABLIM1 in WB as potential biomarkers of IDD, and SOD2 was identified as a hub gene in the PPI network. These genes may hold potential for the clinical treatment and diagnosis of IDD in the future. Further research investigating the genes identified in this study is needed to determine their detailed mechanisms in IDD.

As IDD is an age-related process (Martin et al., 2002), the degree of histopathological degeneration increases with age (Boos et al., 2002). All the patients in this study were confirmed to have IDD by MRI. No study in the literature describes the lowest degree of histopathological degeneration in IDD that can be observed by MRI. Young people are more representative of people without IDD or low IDD; thus, we recruited young people as controls. With such grouping, WB-DEGs may be related to other aging processes and not only IDD. Therefore, we regarded the intersections of the WB-DEGs, AF-DEGs, and NP-DEGs as molecular markers in WB. Thus, age may be a confounding factor affecting gene expression.

In contrast, in the original literature on GSE70362, the author did not mention whether the samples were from patients with scoliosis, trauma, or other medical issues; therefore, uncertain heterogeneity may be a potential factor influencing our conclusions. With this small-sample size study, we wish to provide a preliminary demonstration of biomarkers and potential directions for follow-up research. Our results require further large-sample experimental verification and mechanism research for confirmation.

Footnotes

Acknowledgment

This research was supported by Sichuan Science and Technology Program, China (Grant No. 2018SZ0075).

Author Disclosure Statement

No competing financial interests exist.

References

Adams MA, Roughley PJ (2006) What is intervertebral disc degeneration, and what causes it? Spine (Phila Pa 1976) 31:2151-2161.

Ali

, Le Maitre

, Richardson

, et al. (2008) Connective tissue growth factor expression in human intervertebral disc: implications for angiogenesis in intervertebral disc degeneration. Biotech Histochem, 83:239-245.

Alvarez-Garcia

, Matsuzaki

, Olmer

, et al. (2017) Age-related reduction in the expression of FOXO transcription factors and correlations with intervertebral disc degeneration. J Orthop Res, 35:2682-2691.

Andersson

(1999) Epidemiological features of chronic low-back pain. Lancet, 354:581-585.

Bank

, Bayliss

, Lafeber

, et al. (1998) Ageing and zonal variation in post-translational modification of collagen in normal human articular cartilage. The age-related increase in non-enzymatic glycation affects biomechanical properties of cartilage. Biochem J 330, (Pt 1)::345-351.

Boos

, Weissbach

, Rohrbach

, et al. (2002) Classification of age-related changes in lumbar intervertebral discs: 2002 Volvo Award in basic science. Spine (Phila Pa 1976), 27:2631-2644.

Brinjikji

, Luetmer

, Comstock

, et al. (2015) Systematic literature review of imaging features of spinal degeneration in asymptomatic populations. AJNR Am J Neuroradiol, 36:811-816.

Carragee

, Alamin

, Cheng

, et al. (2006) Are first-time episodes of serious LBP associated with new MRI findings?. Spine J, 6:624-635.

Chen

, Azad

, Gibson

(2009) Superoxide is the major reactive oxygen species regulating autophagy. Cell Death Differ, 16:1040-1052.

10.

Flekac

, Skrha

, Hilgertova

, et al. (2008) Gene polymorphisms of superoxide dismutases and catalase in diabetes mellitus. BMC Med Genet, 9:30.

11.

Guo

, Zhang

, Li

, et al. (2017) Gene expression profile identifies potential biomarkers for human intervertebral disc degeneration. Mol Med Rep, 16:8665-8672.

12.

Johnson

, Evans

, Menage

, et al. (2001) Immunohistochemical detection of Schwann cells in innervated and vascularized human intervertebral discs. Spine (Phila Pa 1976), 26:2550-2557.

13.

Katz

(2006) Lumbar disc disorders and low-back pain: socioeconomic factors and consequences. J Bone Joint Surg Am 88 Suppl, 2:21-24.

14.

Kazezian

, Gawri

, Haglund

, et al. (2015) Gene expression profiling identifies interferon signalling molecules and IGFBP3 in human degenerative annulus fibrosus. Sci Rep, 5:15662.

15.

Landis

, Tower

(2005) Superoxide dismutase evolution and life span regulation. Mech Ageing Dev, 126:365-379.

16.

, Xia

, Nel

(2008) The role of oxidative stress in ambient particulate matter-induced lung diseases and its implications in the toxicity of engineered nanoparticles. Free Radic Biol Med, 44:1689-1699.

17.

Luoma

, Riihimaki

, Luukkonen

, et al. (2000) Low back pain in relation to lumbar disc degeneration. Spine (Phila Pa 1976), 25:487-492.

18.

Martin

, Deyo

, Mirza

, et al. (2008) Expenditures and health status among adults with back and neck problems. JAMA, 299:656-664.

19.

Martin

, Boxell

, Malone

(2002) Pathophysiology of lumbar disc degeneration: a review of the literature. Neurosurg Focus, 13:E1.

20.

Molinos

, Almeida

, Caldeira

, et al. (2015) Inflammation in intervertebral disc degeneration and regeneration. J R Soc Interface, 12:20150429.

21.

Nasto

, Robinson

, Ngo

, et al. (2013) Mitochondrial-derived reactive oxygen species (ROS) play a causal role in aging-related intervertebral disc degeneration. J Orthop Res, 31:1150-1157.

22.

Nel

, Xia

, Madler

, et al. (2006) Toxic potential of materials at the nanolevel. Science, 311:622-627.

23.

Park

, Park

(2009) Oxidative stress and pro-inflammatory responses induced by silica nanoparticles in vivo and in vitro. Toxicol Lett, 184:18-25.

24.

Risbud

, Shapiro

(2014) Role of cytokines in intervertebral disc degeneration: pain and disc content. Nat Rev Rheumatol, 10:44-56.

25.

Schubert

, Smink

, Arp

, et al. (2018) Quality assessment of surgical disc samples discriminates human annulus fibrosus and nucleus pulposus on tissue and molecular level. Int J Mol Sci, 19:E1761.

26.

Shamji

, Setton

, Jarvis

, et al. (2010) Proinflammatory cytokine expression profile in degenerated and herniated human intervertebral disc tissues. Arthritis Rheum, 62:1974-1982.

27.

Sherafatian

, Abdollahpour

, Ghaffarpasand

, et al. (2019) MicroRNA expression profiles, target genes and pathways in intervertebral disc degeneration; A meta-analysis of three microarray studies. World Neurosurg, 126:389-397.

28.

Sivan

, Tsitron

, Wachtel

, et al. (2006) Age-related accumulation of pentosidine in aggrecan and collagen from normal and degenerate human intervertebral discs. Biochem J, 399:29-35.

29.

Thompson

, Pearce

, Schechter

, et al. (1990) Preliminary evaluation of a scheme for grading the gross morphology of the human intervertebral disc. Spine (Phila Pa 1976), 15:411-415.

30.

Tian

, Yuan

, Fujita

, et al. (2013) Inflammatory cytokines associated with degenerative disc disease control aggrecanase-1 (ADAMTS-4) expression in nucleus pulposus cells through MAPK and NF-kappaB. Am J Pathol, 182:2310-2321.

31.

Wang

, Cai

, Shi

, et al. (2016) Aging and age related stresses: a senescence mechanism of intervertebral disc degeneration. Osteoarthritis Cartilage, 24:398-408.

32.

Wang

, Tian

, Wang

, et al. (2013) Inflammatory cytokines induce NOTCH signaling in nucleus pulposus cells: implications in intervertebral disc degeneration. J Biol Chem, 288:16761-16774.

33.

Wang

, Jiang

, Dai

, et al. (2018) Bioinformatics analysis reveals different gene expression patterns in the annulus fibrosis and nucleus pulpous during intervertebral disc degeneration. Exp Ther Med, 16:5031-5040.

34.

Zhang

, Guo

, Sun

, et al. (2010) Gene expression profiles of disc tissues and peripheral blood mononuclear cells from patients with degenerative discs. J Bone Miner Metab, 28:209-219.