Skip to main content
  • Research article
  • Open access
  • Published:

Identification of key regulators responsible for dysregulated networks in osteoarthritis by large-scale expression analysis

Abstract

Background

Osteoarthritis (OA) is a worldwide musculoskeletal disorder. However, disease-modifying therapies for OA are not available. Here, we aimed to characterize the molecular signatures of OA and to identify novel therapeutic targets and strategies to improve the treatment of OA.

Methods

We collected genome-wide transcriptome data performed on 132 OA and 74 normal human cartilage or synovium tissues from 7 independent datasets. Differential gene expression analysis and functional enrichment were performed to identify genes and pathways that were dysregulated in OA. The computational drug repurposing method was used to uncover drugs that could be repurposed to treat OA.

Results

We identified several pathways associated with the development of OA, such as extracellular matrix organization, inflammation, bone development, and ossification. By protein-protein interaction (PPI) network analysis, we prioritized several hub genes, such as JUN, CDKN1A, VEGFA, and FOXO3. Moreover, we repurposed several FDA-approved drugs, such as cardiac glycosides, that could be used in the treatment of OA.

Conclusions

We proposed that the hub genes we identified would play a role in cartilage homeostasis and could be important diagnostic and therapeutic targets. Drugs such as cardiac glycosides provided new possibilities for the treatment of OA.

Background

Osteoarthritis (OA) is the most prevalent musculoskeletal disease and a leading cause of disability in the elderly [1, 2]. OA is characterized by cartilage degeneration, the formation of osteophytes, subchondral bone remodeling, and pathological changes of the meniscus and synovitis [3]. It causes pain, joint stiffness, and disability and leads to severe economic and social burden [4]. Pain is the dominant symptom and a significant driver of clinical decision-making [5]. Treatment of OA includes alleviating pain, controlling inflammation, and slowing down tissue degradation [6]. However, there is currently no effective pharmaceutical treatment for OA that can decelerate the progression of the disease since the precise mechanisms of the pathogenesis of OA remains largely undetermined [7].

The integration of genome-scale transcriptomic profiling of different patient cohorts improves the understanding of molecular changes during OA progression and provides a scientific rationale for the development of novel treatment strategies. RNA-seq and microarray technology are widely used high-throughput genotyping methods that measure the expression of genes on a genome-wide scale with high accuracy and reliability [8,9,10]. Previously, it has been reported that during the pathogenesis of OA, the molecular signature of articular cartilage and synovial membrane underwent huge changes [11]. Expression of genes involved in the molecular matrix components, cell–matrix adhesion, and ossification had dramatic changes, which typically marked different stages of endochondral ossification or transient cartilage differentiation [11,12,13]. These changes result in remodeling of the matrix and lead to the impaired function of the tissue [14]. Notably, local low-grade inflammation of articular cartilage and synovial from patients with OA was also observed. During the progression of OA, pro-inflammatory cytokines such as IL-1, IL-6, TNF or interferon [15], were secreted from diseased tissues, which further damaged surrounding tissues and led to cartilage degradation [7]. Therapeutic strategies aimed at countering inflammation is promising but does not arrest the progressive degeneration of articular cartilage [14, 16, 17]. Thus, identifying the molecular basis of OA progression and the development of novel therapeutic strategies to improve the treatment of OA remain the key issues in OA research.

To characterize OA at the molecular level and to uncover the pathogenesis mechanisms, we collected and integrated transcriptome data of hundreds of OA and normal tissues from 7 independent studies. By comparing the gene expression pattern of OA samples with normal samples, we identified significantly differentially expressed genes (DEGs) and investigated the enriched signaling pathways of these DEGs. Using protein–protein interaction (PPI) network analysis, we discovered several hub genes that function as critical regulators in OA development and might be ideal drug targets. Furthermore, we used a computational drug repurposing method to uncover several drugs that could be repurposed to treat OA.

Methods

Literature search and data collection

We used the keywords “Osteoarthritis”, “RNA-seq or microarray”, and “Dataset” on PubMed (https://www.ncbi.nlm.nih.gov/pubmed/), GEO (http://www.ncbi.nlm.nih.gov/geo/), or ArrayExpress (https://www.ebi.ac.uk/arrayexpress/) to find relevant publications or transcriptome datasets of OA. Transcriptome data of synovial or cartilage from OA patients were included in the analysis. In the data filtering and quality control process, we selected datasets with more than three replications for each condition and measuring over 10,000 genes to ensure the coverage of genes. Altogether, we collected 132 OA samples and 74 normal samples from 7 independent datasets for further analysis.

Data normalization and removal of batch effects

For the microarray data, raw “.CEL” file was downloaded. R program (version 3.5.1; http://www.r-project.org) was used for data analysis. “oligo” [18] and “affy” [19] R packages were applied for the data processing, and then the Robust Multi-array Average (RMA) method [20] was used for background correction, normalization, and probe set summarization. For RNA-Seq data, the sequence data were aligned to the human reference genome hg38 using STAR (v2.5.3a). RSEM (v1.2.28) was applied to map aligned reads and to generate a gene count matrix by default parameters. The expression matrix was normalized using the quantile normalization method. Residual technical batch effects arising due to heterogeneous data platforms were corrected using the ComBat [21] function.

Filtering of differentially expressed genes

For microarray data, “limma” [22] was used to perform the differential gene expression analysis of each dataset. For RNA-seq data, differential expression analysis was performed using Deseq2 [23]. False discovery rate (FDR) was applied to carry out the correction of multiple testing using the Benjamini and Hochberg (BH) method. In this study, genes with |log2fold change (FC)| > 1.2 and FDR < 0.1 were selected as the threshold for differentially expressed genes (DEGs). Genes that were differentially expressed in over half of the datasets with the sample condition were selected and used for further analysis. Gene set enrichment analysis for DEGs was performed using “MAGeCKFlute” R package [24].

Integration of the PPI network and hub gene identification

The Search Tool for the Retrieval of Interacting Genes (STRING) [25] is a biological database for predicting protein interactions. The interactions between DEGs were evaluated using STRING, and gene sets with a combined score > 0.9 were defined as key DEGs. Subsequently, Cytoscape [26] (version 3.6.1; http://cytoscape.org/) was used to visualize the PPI network of the key DEGs that were identified. cytoHubba [27], a Cytoscape plugin, was used to extract the hub genes. The central elements were ranked by betweenness.

Identification of putative target genes of JUN

We searched public ChIP-seq data of JUN on the Cistrome Data Browser (http://cistrome.org/db) [28], a website that collected and integrated thousands of public ChIP-seq data of both human and mouse. We identified and downloaded 25 processed ChIP-seq data of JUN on Cistrome Data Browser. Genes with high RP scores had high likelihoods to be the target genes of JUN. We selected genes with mean RP scores > 1 as candidate target genes of JUN. For further filtering of these candidate genes, we performed co-expression analysis to identify genes that were co-expressed with JUN. Genes with absolute correlation value > 0.5 in each dataset were selected as co-expressed genes. Candidate target genes that were co-expressed in over half of the expression datasets were defined as putative target genes of JUN.

Drug repurposing analysis

Genes that were differentially expressed between OA and normal from both synovial and cartilage tissues were used for drug repurposing analysis. The probe name of the upregulated and downregulated genes were used as input on the cMap (https://portals.broadinstitute.org/cmap/) website. The “quick query” mode of the cMap algorithm was used for the drug repurposing analysis.

Statistical analysis

Statistical analyses were performed using the R software (http://www.R-project.org/). Statistical analyses gathering more than two groups were performed using ANOVA. Otherwise, for two groups, statistical analyses were performed using the unpaired t test. Multiple hypothesis testing corrections were applied where multiple hypotheses were tested and were indicated by the use of FDR.

Results

Outline of data analysis

The objective of this study was to use differentially expressed genes (DEGs) to identify dysregulated genes and pathways of OA and to uncover potential novel pharmacological strategies. To achieve this purpose, we collected hundreds of transcriptome data of patients with OA and normal tissues from 7 independent datasets (Table 1). The subsequent analyses focused on prioritizing key regulators by protein–protein interaction (PPI) network and identifying potent drugs that can be repurposed to treat OA by the computational drug repurposing method (Fig. 1).

Table 1 Information about the collected datasets
Fig. 1
figure 1

Overview of data processing step. a Transcriptomic data from 7 independent studies were selected, which included 132 OA and 74 normal samples. b Quality control and normalization for each data set. c Differential gene expression analysis of OA compared with normal was performed using “Limma” (for microarray data) or “Deseq2” (for RNA-seq data). d Strategies for integrated analysis

Molecular characteristic of OA

To get a list of OA-related DEGs, we compared gene expression profiles of OA patients with those of normal tissues. We collected transcriptome datasets from synovial or cartilage tissues of OA patients and healthy donors. To distinguish DEGs between OA and normal, we performed the differential analysis of each dataset. Among the four transcriptome datasets of synovial tissues, there were 579 DEGs from GSE1919, 2441 DEGs from GSE12021, 2431 DEGs from GSE55235, and 2423 DEGs from GSE55457. The overlap of the DEGs between each dataset was shown in Fig. 2a. Genes that were differentially expressed in two of the four datasets were selected for gene set enrichment analysis. Genes associated with the lysosome, oxidative phosphorylation, extracellular matrix organization, endopeptidase activity, skeletal system development, and collagen-containing extracellular matrix were upregulated in synovial tissues of OA patients than healthy donors (Fig. 2b). In contrast, compared with normal tissues, genes related to IL-17 signaling pathway, circadian rhythm, positive regulation of p38mapk cascade, NOD-like receptor signaling pathway, Foxo signaling pathway, cellular senescence, negative regulation of erk1 and erk2 cascade pathway, and cellular response to hypoxia were downregulated in patients with OA (Fig. 2b).

Fig. 2
figure 2

Gene ontology analysis of DEGs from synovial and cartilage tissues. a Venn plot of DEGs from synovial transcriptome datasets. b Gene ontology analysis of genes that were differentially expressed in more than two of the four synovial transcriptome datasets. “MAGeCKFlute” R package was applied to perform the enrichment analysis. c Venn plot of DEGs from cartilage transcriptome datasets. d Gene ontology analysis of genes that were differentially expressed in more than two of the three cartilage transcriptome datasets. “MAGeCKFlute” R package was applied to perform the enrichment analysis

To distinguish whether the cartilage tissues had different gene expression patterns from what we have observed in synovial tissues, we also collected three transcriptome datasets from cartilage tissues. By comparing the gene expression pattern between patients with OA and healthy donors, we obtained 6307 DEGs from E-MTAB-6266, 6037 DEGs from GSE114007, and 1578 DEGs from GSE117999. The overlap of the DEGs among these three datasets was shown in Fig. 2c. Genes that were differentially expressed in two of the three datasets were selected for further analysis. Gene ontology analysis of the DEGs revealed that genes involved in collagen-containing extracellular matrix, collagen binding, N-acetyglucosamine metabolic process, collagen fibril organization, skeletal system development, ossification, and osteoblast differentiation were upregulated in OA samples (Fig. 2d); while biological processes, including mTOR signaling pathway, Foxo signaling pathway, regulation of cell proliferation, cellular senescence, positive regulation of p38mapk cascade, and ERBB signaling pathway, were downregulated in OA samples (Fig. 2d). Although there were uniquely enriched pathways for synovial and cartilage tissues, we also found some pathways that were enriched in both synovial and cartilage tissues. Such pathways might play essential roles in the development and pathology of OA.

To better understand the differences between OA and normal on the molecular expression pattern, we focused on the genes that were differentially expressed in both synovial and cartilage tissues for further analysis. Among the 3036 DEGs of cartilage and 2620 DEGs of synovial, there were 555 genes differentially expressed in both tissues (Fig. 3a, Supplementary table 1). Next, we performed functional enrichment analysis for these 555 common DEGs. Compared with normal tissues, genes involved in glycine, serine and threonine metabolism, extracellular matrix organization, regulation of T cell activation, peptidase activator activity, bone development, ossification, and collagen-containing extracellular matrix were enriched in OA samples (Fig. 3b), while pathways related to negative regulation of erk1 and erk2 cascade, cellular senescence, response to hypoxia, circadian rhythm, positive regulation of p38mapk cascade, and HIF-1 signaling pathway were downregulated in OA samples (Fig. 3c). Genes involved in these pathways might be responsible for the cartilage homeostasis and the pathology of OA.

Fig. 3
figure 3

Gene ontology analysis of common DEGs from synovial and cartilage tissues. a Venn plot of DEGs from synovial and cartilage transcriptome datasets. b Gene ontology analysis of genes that were upregulated in both synovial and cartilage transcriptome datasets comparing OA with normal. “MAGeCKFlute” R package was applied to perform the enrichment analysis. c Gene ontology analysis of genes that were downregulated in both synovial and cartilage transcriptome datasets comparing OA with normal. “MAGeCKFlute” R package was applied to perform the enrichment analysis

Protein–protein interaction network of DEGs

To discover potential critical regulators of OA from DEGs, we performed protein–protein interaction (PPI) network analysis (Fig. 4a). We used cytoHubba, a Cytoscape plugin, to extract the hub genes in the network. The top 10 hub genes were JUN, ITGB1, VEGFA, CDKN1A, CDK2, PSMD3, CUL2, FOXO3, SGK1, and INSR (Fig.4b). These hub genes might associate with the development and pathogenesis of OA. JUN is a transcription factor that plays an essential role in regulating cellular proliferation and apoptosis. JUN is in the center of the network and interacts with a large number of DEGs in the network, suggesting the crucial role of JUN in regulating the expression of these DEGs. Integrins are major surface receptors of chondrocytes, and integrin β1 (ITGB1) has been proven to suppress chondrocyte hypertrophy and accelerate chondrocyte proliferation by in vitro studies and mouse models [34,35,36,37], suggesting that ITGB1 is involved in OA progression and OA-induced cartilage degradation. VEGF-A is the founding member of the VEGF family and is classically referred to as VEGF [38]. As reported, the expression levels of VEGF were associated with OA progression and OA-specific pathologies, such as cartilage degeneration, osteophyte formation, synovitis, and pain. Moreover, a wide range of studies suggested that inhibition of VEGF signaling reduces OA progression and benefits patients with OA [39]. CDKN1A encodes a potent cyclin-dependent kinase inhibitor and functions as a regulator of cell cycle progression at G1. Shinsuke et al. reported that CDKN1A-deficient mice were more susceptible to OA-related changes in vivo, suggesting that CDKN1A modulation might constitute a possible therapeutic strategy for OA treatment [40]. FOXO3 belongs to the forkhead family of transcription factors, which are characterized by a distinct forkhead domain [41]. FOXO3 functions as a trigger for apoptosis by regulating the expression of genes necessary for cell death [42]. Akasaki et al. reported that normal articular cartilage had a tissue-specific signature of FOXO1 and FOXO3 but not FOXO4 proteins. In OA cartilage, chondrocytes showed altered FOXO activation, which suggested FOXO might play a role in OA progression [43]. SGK1 belongs to the serine/threonine-protein kinase subfamily, which contributes to the regulation of a wide variety of physiological activities, such as membrane transports, cell growth, proliferation, and apoptosis [44]. Knockdown of SGK1 alleviates the IL-1β-induced chondrocyte anabolic and catabolic imbalance by activating autophagy in human chondrocytes [45]. Downregulation of SGK1 attenuates OA progression [46]. The biological role of these hub genes in the development and progression of OA had been studied previously, suggesting that targeting these hub genes might be an optimal novel treatment for OA. However, more experiments were needed to better elucidate the mechanism of action of these hub genes.

Fig. 4
figure 4

Protein–protein interaction network of DEGs. a STRING was used to evaluate protein interaction of DEGs between OA and normal. The interaction was visualized by Cytoscape. b Top 10 hub genes ranked by betweenness. cytoHubba was used to extract the hub genes in the network

JUN functions as a key TF that associated with the development of OA

Since JUN was in the central of the network (Fig. 4a), we focused on JUN for further analysis. We found that most genes interacting with JUN directly in the PPI network were differentially expressed between OA and normal (Fig. 5a and b), suggesting the importance of JUN in regulating the expression of these DEGs. To identify the potential target genes that were regulated by JUN, we downloaded the processed ChIP-seq data of JUN on the Cistrome Data Browser (http://cistrome.org/db) [28]. Cistrome Data Browser collected thousands of human and mouse samples with ChIP-seq data and facilitated searches for putative target genes of transcription factors by regulatory potential model [47]. Regulatory potential (RP) score for each gene was calculated, reflecting the likelihood of the transcription factor being a direct regulator of that gene, and genes with high RP scores were putative targets of the transcription factor. By analyzing 25 public ChIP-seq data of JUN, we selected 3250 genes (Supplementary table 2) with mean RP scores > 1 as candidate target genes of JUN. For further filtering of the candidate genes, we performed co-expression analysis. Among the 3250 candidate targets, 214 genes had high correlations with JUN, which were selected for further analysis (Supplementary table 3). Functional enrichment analysis of the selected 214 genes revealed that pathways including the ATF-2 transcription factor network (part of the AP-1 complex), AP-1 transcription factor network, TGF-β signaling pathway, osteoclast differentiation, and ERBB1 downstream signaling were enriched (Fig. 5c). JUN was in the AP-1 protein family, and the top 2 enriched pathways were associated with the AP-1 family network, which suggested the creditability of our data analysis method. Among these, the TGF-β signaling pathway has been reported to play a critical role in the development and progression of OA by driving chondrocytes toward hypertrophy, promoting osteoprogenitor cell differentiation into osteoblasts, mediating synthesis of cartilage-specific extracellular matrix components, and angiogenesis in subchondral bone [48]. The effects of TGF-β on the modulation of extracellular matrix components were dependent on the activation of JNK (c-Jun N-terminal Kinase), which in turn modulates the activity of c-Jun [49]. Thus, JUN might serve as a critical regulator of the development of OA by regulating the activity of TGF-β.

Fig. 5
figure 5

JUN functions as a key TF that associated with the development of OA. a A subnetwork of the first-neighbored genes with JUN. b The expression level of genes in the subnetwork in (a). Values represent mean ± s.d. *p < 0.05; **p < 0.01; ***p < 0.001; ****p < 0.0001 by Student’s t test. c Gene ontology for putative target genes of JUN. “MAGeCKFlute” R package was applied to perform the enrichment analysis

Drug repurposing for OA

Since there are no effective interventions to decelerate the progression of OA now, drug repurposing of the Food and Drug Administration (FDA)-approved therapeutic agents is a particularly attractive approach to improve the treatment of OA. Computational techniques for predictive repurposing offer a relatively efficient and authentic method of identifying testable hypotheses that may be translated into the clinic [50]. The Connectivity Map (cMap), which was established by the Broad Institute, consists of gene expression data generated by dosing of more than 1300 compounds in hundreds of cell lines [51]. cMap (https://portals.broadinstitute.org/cmap/) has been successfully used to make drug repurposing predictions for a number of disease conditions [50]. We applied cMap to explore potential drug repositioning opportunities for OA. Among the top 15 listed drugs, MG-262 ranked first (Table 2). MG-262 is a potent proteasome inhibitor that selectively and reversibly inhibits the chymotryptic activity of the proteasome [52, 53]. The proteasome inhibitors have shown anti-inflammatory activities in the animal models of arthritis, psoriasis, colitis, and other inflammatory conditions [54]. Some studies discovered that proteasome inhibitors promoted bone growth in a cell-based screen [55]. Although the biological function of MG-262 in the treatment of OA still remains an unexploited field, this analysis provided new possibilities for the treatment of OA. Of note, among the top 15 drugs, there were four drugs classified as cardiac glycosides, such as ouabain and digoxin, suggesting the potential clinical implications of cardiac glycosides in OA. Cardiac glycosides function by inhibiting the Na+/K+ ATPase (NKA) [56]. In addition, cardiac glycosides have been found to decrease inflammatory symptoms in different animal models of acute and chronic inflammation [57]. Thus, cardiac glycosides might have therapeutic benefit in the treatment of OA by countering the inflammation-induced articular cartilage degradation. Since the safety of FDA-approved drugs has been sufficiently verified, drug repurposing for approved drugs offers a less risky and more rapid potential for the investment of novel therapeutic strategies.

Table 2 Top 15 drugs predicted by cMap

Discussion

Here, we used a large-scale data integration method to characterize the molecular signatures of OA, which extended our understanding of the disease mechanisms. We also identified several essential regulators of OA, such as JUN, VEGFA, and FOXO3, which might provide a scientific rationale for the development of novel pharmacological therapies. Of note, we found JUN, a crucial dysregulated transcription factor, plays a central role in regulating the aberrant gene expression pattern in OA. Moreover, we used a computational drug repurposing method to identify potential FDA-approved drugs that can be repurposed to improve the treatment of OA.

In this study, we integrated and analyzed multiple published datasets. We were able to find genes that were consistently differentially expressed between OA and normal among different datasets and different tissues. Our analysis confirmed some findings from previous studies using genome-wide gene expression analyses. Common findings of these studies are the differential expression of genes involved in matrix-degrading enzymes (MMPs, ADAMTS), collagen organization, and inflammation [7, 58]. In our analysis, we also found increased expression of genes associated with collagen-containing extracellular matrix (DPT, GPC4, COL8A2, FBLN5, COL3A1), ossification (CDH11, COL1A1, SPARC, KAZALD1), and bone development (DYM, PTGER4, SPARC, SULF) in patients with OA, which is suggestive of active remodeling of cartilage homeostasis during OA pathogenesis. During the early stages of OA, the molecular composition and organization of the extracellular matrix are altered first [59]. The articular chondrocytes exhibit increased cell proliferation and matrix synthesis for the purpose of initiating repairing for pathological injury [59, 60]. Changes in the composition and structure of the articular cartilage further stimulate chondrocytes to produce more catabolic factors involved in cartilage degradation. Thus, the expression of genes involved in the carbohydrate metabolism and extracellular matrix components were upregulated. Moreover, chronic low-grade inflammation has also been found to contribute to the development and progression of OA [61]. During OA progression, the entire synovial joints were involved in the inflammation process [62]. Pro-inflammatory factors, such as IL-1β and TNF-α, as well as chemokines, were reported to contribute to the systemic inflammation that led to the activation of NF-κB signaling in both synovial cells and chondrocytes [63]. Based on these studies, multiple novel pharmacological strategies have emerged, including anti-inflammatory mediators (anti-IL-1 [64], anti-TNF-α [65], and anti-IL-6 [66]) and inhibition of catabolic pathways (Wnt, ADAMTS, and cathepsin K) [67]. Apart from these findings, we were able to find some essential regulators that might associate with the pathogenesis of OA, such as ITGB1, PSMD3, CUL2, and INSR. Since the association of these hub genes with the development of OA has been reported, they could serve as important diagnostic and/or therapeutic targets for OA. However, more mechanistic interpretation and clinical trial data are needed to clarify the efficacy in the treatment of OA.

Since synovial and cartilage play different roles in the bone homeostasis and pathology of OA, we compared the gene expression pattern between OA and normal using transcriptome datasets of each tissue separately. This analysis identified some tissue-specific DEGs between OA and normal. Gene sets associated with N-acetylglucosamine metabolic process, glycosaminoglycan biosynthesis and metabolism, NF-KB transcription factor activity, ERBB signaling pathway, P53 signaling pathway, and mTOR signaling pathway were uniquely differentially expressed in cartilage tissues. While biological pathways involved in lysosome, oxidative phosphorylation, endopeptidase activity, Nod-like receptor signaling pathway, and IL-17 signaling pathway were uniquely differentially expressed in synovial tissues. In the development of OA, the pathological changes of synovial and cartilage were different. The histological pattern of synovium in OA patients is characterized by synovial lining hyperplasia, increased vascularity, and sublining fibrosis [68] with the phenotypic shift of chondrocytes, including surface fibrillation, degradation of cartilage matrix, chondrocyte clusters appearance, and vascular penetration from the subchondral bone [69]. The different phenotypic changes of synovial and cartilage might explain the different OA-related DEGs between these two tissues.

Our study also identified JUN, a transcription factor, as a key regulator of these DEGs. JUN is one of the members of the Activator protein 1 (AP-1) family proteins [70]. AP-1 family proteins are basic leucine zipper (bZIP) transcription factors that consist of Jun (c-Jun, JunB, and JunD), Fos (c-Fos, FosB, Fra-1, and Fra-2), Jun dimerization partners (JDP1 and JDP2), and the closely related activating transcription factors (ATF2, LRF1/ ATF3, and B-ATF) subfamilies [71]. AP-1 family proteins are implicated in the regulation of a variety of cellular processes, including proliferation and survival, differentiation, apoptosis, cell migration, and transformation [72]. JUN has been reported to play a crucial role in regulating cell proliferation and apoptosis [73]. Ventura et al. reported that the c-Jun NH2-terminal kinase JNK signaling pathway contributes to the regulation of TGF-β-mediated biological responses [74]. TGF-β is crucial for cartilage maintenance, and the lack of TGF-β results in OA-like changes [74]. Thus, JUN might regulate the development of OA by coordinating with TGF-β signaling.

Drug discovery is a time-consuming, laborious, costly, and high-risk process. Drug repositioning is efficient, economical, and low-risk compared with the traditional drug development process. There are multiple drug repurposing methods generated, including computational approaches, biological experimental approaches, and mixed approaches. Among these, computational drug repurposing is the most powerful and the most widely used. Previously, Soul et al. constructed a data portal that provides an exploration and comparison platform for analyzed skeletal transcriptomics data [75]. This data portal also provides data repurposing mode using L1000 data, one of the cMap project resources. Using a computational drug repurposing method, we found cardiac glycosides might be repurposed in the treatment of OA. Cardiac glycosides are drugs that inhibit the Na+/K+ ATPases and are applied to treat heart failure and certain irregular heartbeats. However, recent studies reported that cardiac glycosides are a novel class of broad-spectrum senolytics for therapeutic applications in many age-related disorders [76], including osteoarthritis. Cardiac glycosides were capable of reducing the number of senescent cells, diminishing the level of local inflammation, and improved some metabolic and physical fitness parameters that decline with aging in some animal models [77, 78]. Drug repurposing can be highly attractive as a potentially cheaper and faster route to market. However, successful drug repositioning requires a deep understanding of biological mechanisms—from known overlaps of mechanisms to re-innovation of a new molecule, to find a new mechanism, dosing, route of administration, and new target. Although our drug repurposing analysis offers novel and valuable options for developing strategies to treat OA, the effectiveness of the predicted drugs in treating OA, such as cardiac glycosides, needs more systemic and detailed validation.

There are several limitations of our studies. First, osteoarthritis is typically described as a heterogeneous disease with complex pathogenesis. Different patients might have different mechanistic pathways, such that the mechanisms of OA in the elderly might be different from those after a joint injury in a younger adult or in obese individuals. The molecular signatures we discovered in this study by transcriptome data analysis might not be representative in some subtypes of OA patients. Integrated analysis of multi-omics data, including epigenetics (DNA methylation, histone post-translational modification, and/or non-coding RNA), metabolomics, and proteomics (by LC-MS) data, are needed to better elucidate the heterogeneity of OA. Such detailed and systematic analyses offer important mechanistic and potentially therapeutic insights into OA. Finally, we have highlighted several essential regulators and uncovered that cardiac glycosides might benefit the patients with OA by countering the inflammation. However, further experimental studies are needed to validate their clinical implications.

Our study integrated multiple public transcriptome data sets together, which provides more comprehensive and reliable insights into the genetic alterations associated with the disease phenotype. Moreover, we used a computational drug repurposing method to identify potent drug candidates to improve the treatment of OA.

Conclusion

In summary, we used bioinformatics analysis to identify a group of differentially expressed genes between OA and normal tissues. By protein–protein interaction network analysis, we identified several osteoarthritis-related hub genes, which might be potential diagnostic or therapeutic targets for osteoarthritis. In addition, we found JUN, a transcription factor, functions as a crucial regulator by the analysis of public ChIP-seq data and co-expression analysis. Moreover, we used a computational drug repurposing method to identify potent drugs that can be repurposed to treat OA.

Availability of data and materials

The datasets used and/or analyzed in this study are available from the corresponding author on reasonable request.

Abbreviations

OA:

Osteoarthritis

DEGs:

Differentially expressed genes

PPI:

Protein–protein interaction

RMA:

Robust Multi-array Average

ECM:

Extracellular matrix

GO:

Gene Ontology

STRING:

Search Tool for the Retrieval of Interacting Genes

RP:

Regulatory potential

References

  1. Woolf AD, Pfleger B. Burden of major musculoskeletal conditions. Bull World Health Organ. 2003;81(9):646–56.

    PubMed  PubMed Central  Google Scholar 

  2. Chou CH, Wu CC, Song IW, Chuang HP, Lu LS, Chang JH, et al. Genome-wide expression profiles of subchondral bone in osteoarthritis. Arthritis Res Ther. 2013;15(6):R190. https://doi.org/10.1186/ar4380.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  3. Katz JN. Lumbar disc disorders and low-back pain: socioeconomic factors and consequences. J Bone Joint Surg Am. 2006;88(Suppl 2):21–4.

    PubMed  Google Scholar 

  4. Shepherd C, Reese AE, Reynard LN, Loughlin J. Expression analysis of the osteoarthritis genetic susceptibility mapping to the matrix Gla protein gene MGP. Arthritis Res Ther. 2019;21(1):149. https://doi.org/10.1186/s13075-019-1934-7.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  5. Glyn-Jones S, Palmer AJ, Agricola R, Price AJ, Vincent TL, Weinans H, et al. Osteoarthritis. Lancet. 2015;386(9991):376–87. https://doi.org/10.1016/S0140-6736(14)60802-3.

    Article  PubMed  CAS  Google Scholar 

  6. Bijlsma JW, Berenbaum F, Lafeber FP. Osteoarthritis: an update with relevance for clinical practice. Lancet. 2011;377(9783):2115–26. https://doi.org/10.1016/S0140-6736(11)60243-2.

    Article  PubMed  Google Scholar 

  7. Malaise O, de Seny D. Therapeutic advances in arthritis diseases. Biochem Pharmacol. 2019;165:1–3. https://doi.org/10.1016/j.bcp.2019.04.014.

    Article  PubMed  CAS  Google Scholar 

  8. Chen G, Shi T, Shi L. Characterizing and annotating the genome using RNA-seq data. Sci China Life Sci. 2017;60(2):116–25. https://doi.org/10.1007/s11427-015-0349-4.

    Article  PubMed  CAS  Google Scholar 

  9. Nagalakshmi U, Waern K, Snyder M. RNA-Seq: a method for comprehensive transcriptome analysis. Curr Protoc Mol Biol. 2010;Chapter 4:Unit 4 11 1-13.

    PubMed  Google Scholar 

  10. Mortazavi A, Williams BA, McCue K, Schaeffer L, Wold B. Mapping and quantifying mammalian transcriptomes by RNA-Seq. Nat Methods. 2008;5(7):621–8. https://doi.org/10.1038/nmeth.1226.

    Article  PubMed  CAS  Google Scholar 

  11. Pitsillides AA, Beier F. Cartilage biology in osteoarthritis--lessons from developmental biology. Nat Rev Rheumatol. 2011;7(11):654–63. https://doi.org/10.1038/nrrheum.2011.129.

    Article  PubMed  CAS  Google Scholar 

  12. Xu YK, Ke Y, Wang B, Lin JH. The role of MCP-1-CCR2 ligand-receptor axis in chondrocyte degradation and disease progress in knee osteoarthritis. Biol Res. 2015;48(1):64. https://doi.org/10.1186/s40659-015-0057-0.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  13. Heinegard D, Saxne T. The role of the cartilage matrix in osteoarthritis. Nat Rev Rheumatol. 2011;7(1):50–6. https://doi.org/10.1038/nrrheum.2010.198.

    Article  PubMed  CAS  Google Scholar 

  14. Jaswal AP, Bandyopadhyay A. Re-examining osteoarthritis therapy from a developmental biologist's perspective. Biochem Pharmacol. 2019;165:17–23. https://doi.org/10.1016/j.bcp.2019.03.020.

    Article  PubMed  CAS  Google Scholar 

  15. Henrotin Y, Lambert C, Richette P. Importance of synovitis in osteoarthritis: evidence for the use of glycosaminoglycans against synovial inflammation. Semin Arthritis Rheum. 2014;43(5):579–87. https://doi.org/10.1016/j.semarthrit.2013.10.005.

    Article  PubMed  CAS  Google Scholar 

  16. Rannou F, Poiraudeau S. Non-pharmacological approaches for the treatment of osteoarthritis. Best Pract Res Clin Rheumatol. 2010;24(1):93–106. https://doi.org/10.1016/j.berh.2009.08.013.

    Article  PubMed  Google Scholar 

  17. Zhang W, Moskowitz RW, Nuki G, Abramson S, Altman RD, Arden N, et al. OARSI recommendations for the management of hip and knee osteoarthritis, Part II: OARSI evidence-based, expert consensus guidelines. Osteoarthritis Cartilage. 2008;16(2):137–62. https://doi.org/10.1016/j.joca.2007.12.013.

    Article  PubMed  CAS  Google Scholar 

  18. Carvalho BS, Irizarry RA. A framework for oligonucleotide microarray preprocessing. Bioinformatics. 2010;26(19):2363–7. https://doi.org/10.1093/bioinformatics/btq431.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  19. Gautier L, Cope L, Bolstad BM, Irizarry RA. affy--analysis of Affymetrix GeneChip data at the probe level. Bioinformatics. 2004;20(3):307–15.

    Article  CAS  PubMed  Google Scholar 

  20. Irizarry RA, Hobbs B, Collin F, Beazer-Barclay YD, Antonellis KJ, Scherf U, et al. Exploration, normalization, and summaries of high density oligonucleotide array probe level data. Biostatistics. 2003;4(2):249–64. https://doi.org/10.1093/biostatistics/4.2.249.

    Article  PubMed  Google Scholar 

  21. Leek JT, Johnson WE, Parker HS, Jaffe AE, Storey JD. The sva package for removing batch effects and other unwanted variation in high-throughput experiments. Bioinformatics. 2012;28(6):882–3. https://doi.org/10.1093/bioinformatics/bts034.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  22. Ritchie ME, Phipson B, Wu D, Hu Y, Law CW, Shi W, et al. limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res. 2015;43(7):e47.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  23. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15(12):550. https://doi.org/10.1186/s13059-014-0550-8.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  24. Wang B, Wang M, Zhang W, Xiao T, Chen CH, Wu A, et al. Integrative analysis of pooled CRISPR genetic screens using MAGeCKFlute. Nat Protoc. 2019;14(3):756–80. https://doi.org/10.1038/s41596-018-0113-7.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  25. Szklarczyk D, Gable AL, Lyon D, Junge A, Wyder S, Huerta-Cepas J, et al. STRING v11: protein-protein association networks with increased coverage, supporting functional discovery in genome-wide experimental datasets. Nucleic Acids Res. 2019;47(D1):D607–13. https://doi.org/10.1093/nar/gky1131.

    Article  PubMed  CAS  Google Scholar 

  26. Shannon P, Markiel A, Ozier O, Baliga NS, Wang JT, Ramage D, et al. Cytoscape: a software environment for integrated models of biomolecular interaction networks. Genome Res. 2003;13(11):2498–504. https://doi.org/10.1101/gr.1239303.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  27. Chin CH, Chen SH, Wu HH, Ho CW, Ko MT, Lin CY. cytoHubba: identifying hub objects and sub-networks from complex interactome. BMC Syst Biol. 2014;8 Suppl 4:S11.

    Article  PubMed  Google Scholar 

  28. Zheng R, Wan C, Mei S, Qin Q, Wu Q, Sun H, et al. Cistrome Data Browser: expanded datasets and new tools for gene regulatory analysis. Nucleic Acids Res. 2019;47(D1):D729–35. https://doi.org/10.1093/nar/gky1094.

    Article  PubMed  CAS  Google Scholar 

  29. Woetzel D, Huber R, Kupfer P, Pohlers D, Pfaff M, Driesch D, et al. Identification of rheumatoid arthritis and osteoarthritis patients by transcriptome-based rule set generation. Arthritis Res Ther. 2014;16(2):R84. https://doi.org/10.1186/ar4526.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  30. Huber R, Hummert C, Gausmann U, Pohlers D, Koczan D, Guthke R, et al. Identification of intra-group, inter-individual, and gene-specific variances in mRNA expression profiles in the rheumatoid arthritis synovial membrane. Arthritis Res Ther. 2008;10(4):R98. https://doi.org/10.1186/ar2485.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  31. Ungethuem U, Haeupl T, Witt H, Koczan D, Krenn V, Huber H, et al. Molecular signatures and new candidates to target the pathogenesis of rheumatoid arthritis. Physiol Genomics. 2010;42A(4):267–82. https://doi.org/10.1152/physiolgenomics.00004.2010.

    Article  PubMed  CAS  Google Scholar 

  32. Fisch KM, Gamini R, Alvarez-Garcia O, Akagi R, Saito M, Muramatsu Y, et al. Identification of transcription factors responsible for dysregulated networks in human osteoarthritis cartilage by global gene expression analysis. Osteoarthritis Cartilage. 2018;26(11):1531–8. https://doi.org/10.1016/j.joca.2018.07.012.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  33. Soul J, Dunn SL, Anand S, Serracino-Inglott F, Schwartz JM, Boot-Handford RP, et al. Stratification of knee osteoarthritis: two major patient subgroups identified by genome-wide expression analysis of articular cartilage. Ann Rheum Dis. 2018;77(3):423. https://doi.org/10.1136/annrheumdis-2017-212603.

    Article  PubMed  CAS  Google Scholar 

  34. Enomoto M, Leboy PS, Menko AS, Boettiger D. Beta 1 integrins mediate chondrocyte interaction with type I collagen, type II collagen, and fibronectin. Exp Cell Res. 1993;205(2):276–85. https://doi.org/10.1006/excr.1993.1087.

    Article  PubMed  CAS  Google Scholar 

  35. Raducanu A, Hunziker EB, Drosse I, Aszodi A. Beta1 integrin deficiency results in multiple abnormalities of the knee joint. J Biol Chem. 2009;284(35):23780–92. https://doi.org/10.1074/jbc.M109.039347.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  36. Bengtsson T, Aszodi A, Nicolae C, Hunziker EB, Lundgren-Akerlund E, Fassler R. Loss of alpha10beta1 integrin expression leads to moderate dysfunction of growth plate chondrocytes. J Cell Sci. 2005;118(Pt 5):929–36. https://doi.org/10.1242/jcs.01678.

    Article  PubMed  CAS  Google Scholar 

  37. Lian C, Wang X, Qiu X, Wu Z, Gao B, Liu L, et al. Collagen type II suppresses articular chondrocyte hypertrophy and osteoarthritis progression by promoting integrin beta1-SMAD1 interaction. Bone Res. 2019;7(1):8. https://doi.org/10.1038/s41413-019-0046-y.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  38. Ferrara N, Gerber HP, LeCouter J. The biology of VEGF and its receptors. Nat Med. 2003;9(6):669–76. https://doi.org/10.1038/nm0603-669.

    Article  PubMed  CAS  Google Scholar 

  39. Hamilton JL, Nagao M, Levine BR, Chen D, Olsen BR, Im HJ. Targeting VEGF and its receptors for the treatment of osteoarthritis and associated pain. J Bone Miner Res. 2016;31(5):911–24. https://doi.org/10.1002/jbmr.2828.

    Article  PubMed  Google Scholar 

  40. Kihara S, Hayashi S, Hashimoto S, Kanzaki N, Takayama K, Matsumoto T, et al. Cyclin-dependent kinase inhibitor-1-deficient mice are susceptible to osteoarthritis associated with enhanced inflammation. J Bone Miner Res. 2017;32(5):991–1001. https://doi.org/10.1002/jbmr.3080.

    Article  PubMed  CAS  Google Scholar 

  41. Stefanetti RJ, Voisin S, Russell A, Lamon S. Recent advances in understanding the role of FOXO3, F1000Res 7; 2018.

    Google Scholar 

  42. Ekoff M, Kaufmann T, Engstrom M, Motoyama N, Villunger A, Jonsson JI, et al. The BH3-only protein Puma plays an essential role in cytokine deprivation induced apoptosis of mast cells. Blood. 2007;110(9):3209–17. https://doi.org/10.1182/blood-2007-02-073957.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  43. Akasaki Y, Hasegawa A, Saito M, Asahara H, Iwamoto Y, Lotz MK. Dysregulated FOXO transcription factors in articular cartilage in aging and osteoarthritis. Osteoarthritis Cartilage. 2014;22(1):162–70. https://doi.org/10.1016/j.joca.2013.11.004.

    Article  PubMed  CAS  Google Scholar 

  44. Lang F, Strutz-Seebohm N, Seebohm G, Lang UE. Significance of SGK1 in the regulation of neuronal function. J Physiol. 2010;588(Pt 18):3349–54. https://doi.org/10.1113/jphysiol.2010.190926.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  45. Huang W, Cheng C, Shan WS, Ding ZF, Liu FE, Lu W, et al. Knockdown of SGK1 alleviates the IL-1beta-induced chondrocyte anabolic and catabolic imbalance by activating FoxO1-mediated autophagy in human chondrocytes. FEBS J. 2020;287(1):94–107. https://doi.org/10.1111/febs.15009.

    Article  PubMed  CAS  Google Scholar 

  46. Wang Z, Ni S, Zhang H, Fan Y, Xia L, Li N. Silencing SGK1 alleviates osteoarthritis through epigenetic regulation of CREB1 and ABCA1 expression. Life Sci. 2021;268:118733.  https://doi.org/10.1016/j.lfs.2020.118733.

  47. Tang Q, Chen Y, Meyer C, Geistlinger T, Lupien M, Wang Q, et al. A comprehensive view of nuclear receptor cancer cistromes. Cancer Res. 2011;71(22):6940–7. https://doi.org/10.1158/0008-5472.CAN-11-2091.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  48. Shen J, Li S, Chen D. TGF-beta signaling and the development of osteoarthritis, Bone Res 2; 2014.

    Google Scholar 

  49. Hocevar BA, Brown TL, Howe PH. TGF-beta induces fibronectin synthesis through a c-Jun N-terminal kinase-dependent, Smad4-independent pathway. EMBO J. 1999;18(5):1345–56. https://doi.org/10.1093/emboj/18.5.1345.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  50. Pushpakom S, Iorio F, Eyers PA, Escott KJ, Hopper S, Wells A, et al. Drug repurposing: progress, challenges and recommendations. Nat Rev Drug Discov. 2019;18(1):41–58. https://doi.org/10.1038/nrd.2018.168.

    Article  PubMed  CAS  Google Scholar 

  51. Lamb J, Crawford ED, Peck D, Modell JW, Blat IC, Wrobel MJ, et al. The Connectivity Map: using gene-expression signatures to connect small molecules, genes, and disease. Science. 2006;313(5795):1929–35. https://doi.org/10.1126/science.1132939.

    Article  PubMed  CAS  Google Scholar 

  52. Pujols L, Fernandez-Bertolin L, Fuentes-Prado M, Alobid I, Roca-Ferrer J, Agell N, et al. Proteasome inhibition reduces proliferation, collagen expression, and inflammatory cytokine production in nasal mucosa and polyp fibroblasts. J Pharmacol Exp Ther. 2012;343(1):184–97. https://doi.org/10.1124/jpet.111.190710.

    Article  PubMed  CAS  Google Scholar 

  53. Frase H, Hudak J, Lee I. Identification of the proteasome inhibitor MG262 as a potent ATP-dependent inhibitor of the Salmonella enterica serovar Typhimurium Lon protease. Biochemistry. 2006;45(27):8264–74. https://doi.org/10.1021/bi060542e.

    Article  PubMed  CAS  Google Scholar 

  54. Kisselev AF, van der Linden WA, Overkleeft HS. Proteasome inhibitors: an expanding army attacking a unique target. Chem Biol. 2012;19(1):99–115. https://doi.org/10.1016/j.chembiol.2012.01.003.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  55. Garrett IR, Chen D, Gutierrez G, Zhao M, Escobedo A, Rossini G, et al. Selective inhibitors of the osteoblast proteasome stimulate bone formation in vivo and in vitro. J Clin Invest. 2003;111(11):1771–82. https://doi.org/10.1172/JCI16198.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  56. Koren G, Soldin SJ. Cardiac glycosides. Clin Lab Med. 1987;7(3):587–606. https://doi.org/10.1016/S0272-2712(18)30733-9.

    Article  PubMed  CAS  Google Scholar 

  57. Furst R, Zundorf I, Dingermann T. New knowledge about old drugs: the anti-inflammatory properties of cardiac glycosides. Planta Med. 2017;83(12-13):977–84. https://doi.org/10.1055/s-0043-105390.

    Article  PubMed  CAS  Google Scholar 

  58. Li H, Yang HH, Sun ZG, Tang HB, Min JK. Whole-transcriptome sequencing of knee joint cartilage from osteoarthritis patients. Bone Joint Res. 2019;8(7):290–303. https://doi.org/10.1302/2046-3758.87.BJR-2018-0297.R1.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  59. Lee AS, Ellman MB, Yan D, Kroin JS, Cole BJ, van Wijnen AJ, et al. A current review of molecular mechanisms regarding osteoarthritis and pain. Gene. 2013;527(2):440–7. https://doi.org/10.1016/j.gene.2013.05.069.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  60. Im HJ, Li X, Muddasani P, Kim GH, Davis F, Rangan J, et al. Basic fibroblast growth factor accelerates matrix degradation via a neuro-endocrine pathway in human adult articular chondrocytes. J Cell Physiol. 2008;215(2):452–63. https://doi.org/10.1002/jcp.21317.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  61. Shen J, Abu-Amer Y, O'Keefe RJ, McAlinden A. Inflammation and epigenetic regulation in osteoarthritis. Connect Tissue Res. 2017;58(1):49–63. https://doi.org/10.1080/03008207.2016.1208655.

    Article  PubMed  CAS  Google Scholar 

  62. Goldring MB, Otero M. Inflammation in osteoarthritis. Curr Opin Rheumatol. 2011;23(5):471–8. https://doi.org/10.1097/BOR.0b013e328349c2b1.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  63. Chen D, Shen J, Zhao W, Wang T, Han L, Hamilton JL, et al. Osteoarthritis: toward a comprehensive understanding of pathological mechanism. Bone Res. 2017;5(1):16044. https://doi.org/10.1038/boneres.2016.44.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  64. Singer JW, Fleischman A, Al-Fayoumi S, Mascarenhas JO, Yu Q, Agarwal A. Inhibition of interleukin-1 receptor-associated kinase 1 (IRAK1) as a therapeutic strategy. Oncotarget. 2018;9(70):33416–39. https://doi.org/10.18632/oncotarget.26058.

    Article  PubMed  PubMed Central  Google Scholar 

  65. Grunke M, Schulze-Koops H. Successful treatment of inflammatory knee osteoarthritis with tumour necrosis factor blockade. Ann Rheum Dis. 2006;65(4):555–6. https://doi.org/10.1136/ard.2006.053272.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  66. Brenn D, Richter F, Schaible HG. Sensitization of unmyelinated sensory fibers of the joint nerve to mechanical stimuli by interleukin-6 in the rat: an inflammatory mechanism of joint pain. Arthritis Rheum. 2007;56(1):351–9. https://doi.org/10.1002/art.22282.

    Article  PubMed  CAS  Google Scholar 

  67. Pap T, Bertrand J. Syndecans in cartilage breakdown and synovial inflammation. Nat Rev Rheumatol. 2013;9(1):43–55. https://doi.org/10.1038/nrrheum.2012.178.

    Article  PubMed  CAS  Google Scholar 

  68. Mathiessen A, Conaghan PG. Synovitis in osteoarthritis: current understanding with therapeutic implications. Arthritis Res Ther. 2017;19(1):18. https://doi.org/10.1186/s13075-017-1229-9.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  69. Goldring MB. Articular cartilage degradation in osteoarthritis. HSS J. 2012;8(1):7–9. https://doi.org/10.1007/s11420-011-9250-z.

    Article  PubMed  PubMed Central  Google Scholar 

  70. Karin M, Liu Z, Zandi E. AP-1 function and regulation. Curr Opin Cell Biol. 1997;9(2):240–6. https://doi.org/10.1016/S0955-0674(97)80068-3.

    Article  PubMed  CAS  Google Scholar 

  71. Vesely PW, Staber PB, Hoefler G, Kenner L. Translational regulation mechanisms of AP-1 proteins. Mutat Res. 2009;682(1):7–12. https://doi.org/10.1016/j.mrrev.2009.01.001.

    Article  PubMed  CAS  Google Scholar 

  72. Shaulian E, Karin M. AP-1 in cell proliferation and survival. Oncogene. 2001;20(19):2390–400. https://doi.org/10.1038/sj.onc.1204383.

    Article  PubMed  CAS  Google Scholar 

  73. Kappelmann M, Bosserhoff A, Kuphal S. AP-1/c-Jun transcription factors: regulation and function in malignant melanoma. Eur J Cell Biol. 2014;93(1-2):76–81. https://doi.org/10.1016/j.ejcb.2013.10.003.

    Article  PubMed  CAS  Google Scholar 

  74. Ventura JJ, Kennedy NJ, Flavell RA, Davis RJ. JNK regulates autocrine expression of TGF-beta1. Mol Cell. 2004;15(2):269–78. https://doi.org/10.1016/j.molcel.2004.06.007.

    Article  PubMed  CAS  Google Scholar 

  75. Soul J, Hardingham TE, Boot-Handford RP, Schwartz JM. SkeletalVis: an exploration and meta-analysis data portal of cross-species skeletal transcriptomics data. Bioinformatics. 2019;35(13):2283–90. https://doi.org/10.1093/bioinformatics/bty947.

    Article  PubMed  CAS  Google Scholar 

  76. Martin N, Soriani O, Bernard D. Cardiac glycosides as senolytic compounds. Trends Mol Med. 2020;26(3):243–5. https://doi.org/10.1016/j.molmed.2020.01.001.

    Article  PubMed  CAS  Google Scholar 

  77. He S, Sharpless NE. Senescence in health and disease. Cell. 2017;169(6):1000–11. https://doi.org/10.1016/j.cell.2017.05.015.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  78. Paez-Ribes M, Gonzalez-Gualda E, Doherty GJ, Munoz-Espin D. Targeting senescent cells in translational medicine. EMBO Mol Med. 2019;11(12):e10234. https://doi.org/10.15252/emmm.201810234.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

Download references

Acknowledgments

The authors wish to thank the Gene Expression Omnibus (GEO) and ArrayExpress database for making their processed data publicly available.

Funding

This research was partially supported by the National Natural Science Foundation of China (82002282).

Author information

Authors and Affiliations

Authors

Contributions

SS, FW, YL, and FL were involved in the study design. SS, FW, ZZ, RT, YL, and MJ participated in the data collection. SS and FW performed all the data analyses. SS, FW, FL, and YL drafted the manuscript. All authors contributed and approved to the discussion and revision of the manuscript.

Corresponding authors

Correspondence to Fan Liu or Yake Liu.

Ethics declarations

Ethics approval and consent to participate

Not applicable

Consent for publication

All the authors approved the manuscript for publication.

Competing interests

The authors declare that they have no conflicts of interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1: Sup.table 1

: Differentially expressed genes (DEGs) between OA and normal.

Additional file 2: Sup.table 2

: RP scores of public ChIP-seq data of JUN.

Additional file 3: Sup.table 3

: Putative targets of JUN.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Shi, S., Wan, F., Zhou, Z. et al. Identification of key regulators responsible for dysregulated networks in osteoarthritis by large-scale expression analysis. J Orthop Surg Res 16, 259 (2021). https://doi.org/10.1186/s13018-021-02402-9

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s13018-021-02402-9

Keywords