Introduction

Osteoarthritis (OA) is the most common degenerative joint disease that significantly affects individuals worldwide, leading to considerable economic and property losses.1,2 OA is characterized by cartilage degeneration, joint pain, stiffness, and functional disabilities. Commonly recognized factors contributing to OA include aging, trauma, hormones, genetics, and mechanical damage to the joints.3–5 Currently, there is no effective disease-modifying medications have been identified for OA. Medications are primarily aimed at relieving symptoms in mild-to-moderate OA.6–8 While in severe OA, medications offer only short-term benefits in pain and function, necessitating joint replacement surgery as a therapeutic option.9,10

The development of OA involves a complex pathological process, including the degradation of articular cartilage, remodeling of subchondral bone, degeneration of meniscus, inflammation and fibrosis of the patellar fat pad and synovium,11–16 with the degeneration of articular cartilage being its primary manifestation.17,18 Cartilage or osteochondral damage is strongly correlated with an elevated incidence of OA.19 OA resulting from cartilage injury may advance due to the degeneration of adjacent cartilage, and catabolic factors in the cartilage surrounding the defect, such as inflammatory cytokines and biomechanical stress, may impact the progression of OA following cartilage injury.20,21 The principal features of articular cartilage surface involvement in the development of OA include chondrocyte death caused by oxidative stress and degradation of the extracellular matrix (ECM).22,23

From a biomechanical perspective, joint cartilage is composed of a large amount of dense extracellular matrix, and another characteristic of cartilage is the presence of a pericellular matrix (PCM) that surrounds the chondrocytes. The changes in the composition of these matrices leads to alterations in the mechanical environment of the cells within the cartilage matrix. If the cartilage cells are subjected to abnormal mechanical stimuli, their metabolic balance is disrupted, leading to loss of matrix and tissue degeneration, which can result in OA.24,25

Zinc finger proteins constitute the largest family of transcription factors in the human genome. They can be transcriptionally regulated by binding to specific DNA sequences and are involved in a wide range of biological processes, including cell proliferation, invasion, and metastasis.26–29 Exosomes have been demonstrated to play a crucial role in intercellular communication by transporting various biologically active substances, such as proteins, circular RNA (circRNA), microRNA (miRNA), and mRNA.30 The circ-ZNF652 mediates cell proliferation, apoptosis, and inflammation in infantile pneumonia.31 Furthermore, circ-ZNF652 has been identified as a positive regulator in the progression of certain tumors. It has a high expression level in clinical renal carcinoma tissues while its silencing leads to repressed proliferation of renal carcinoma cells.32 Additionally, exosomes contain circ-ZNF652 which contributes to hepatocellular carcinoma cell proliferation, migration, and glycolysis.33 Similarly, cell proliferation, apoptosis, and inflammation have been reported to play a role in OA pathogenesis. However, the implications of ZNF652 in the regulation of OA, and its possible role in promoting osteoarthritic chondrocyte cell death, facilitating extracellular matrix degradation, and consequently expediting OA development, remain unexplored in the current literature.

Hence, it is imperative to investigate the clinical significance and potential molecular mechanisms of ZNF652 in the context of OA. In this study, cartilage samples were collected from both OA patients and healthy individuals obtained from the GEO database. Subsequently, a screening for differentially expressed genes was performed, characterized genes associated with OA were identified using bioinformatics methods, and enrichment analysis and immune infiltration analysis were carried out. In addition, correlations were analyzed to present novel strategies and insights for the treatment of OA. Initially, 297 differentially expressed genes were identified, and following intragroup screening, three differential genes were selected, one of which was recognized as a potential biomarker for OA.

Materials and Methods Data Acquisition and Processing

We searched the GEO database using the search term “osteoarthritis” (ICD-10 code M19. 9) and narrowed the results to “Series” and “Expression profiling by array” for Homo sapiens. This yielded 100 datasets, from which we selected and downloaded six specific expression datasets including GSE10575, GSE16464, GSE117999, GSE169077, GSE178557, and GSE236924, which related to articular cartilage for our study (Table 1). Subsequently, we converted the probe matrices to gene expression matrices using the provided annotation information in the platform files. The study was approved by the Clinical Research Ethics Committee of Guangxi Medical University (Approved Number: 202008012).

Table 1 Information on the GEO Dataset Used

Data Quality Control and Data Analysis

The data sets labeled as GSE10575, GSE16464, GSE169077, GSE178557, and GSE236924 underwent standardization using the “limma” package in R software. Subsequently, data exhibiting significant differences were subjected to a log transformation (base 0.5) and then combined. PCA analysis was applied to both the pre- and post-standardization data, and validation data were batch-corrected to ensure comparability. Heatmaps of DEGs were generated using the “heatmap” package, and the data were processed using log transformation (base 0.5). A statistical significance threshold of P < 0.05 was set to discern variations in DEGs and their expressions between the experimental and control groups. Correlation heatmaps were also generated using the “corrplot” package. Furthermore, Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) analyses were concurrently conducted using the “clusterProfiler” package in R, with a threshold of p-value <0.05.

Screening Feature Genes

Three machine learning algorithms have been employed to predict disease states and identify significant feature genes. The Least Absolute Shrinkage and Selection Operator (LASSO) algorithm serves as a regression analysis tool that enhances prediction accuracy through regularization. Support Vector Machines (SVMs) have been widely utilized as a supervised machine learning technique for classification and regression, while the Recursive Feature Elimination (RFE) algorithm aids in mitigating overfitting.34 Additionally, the RF stands out as a successful integrated learning technique for large-scale classification and skewing problems, utilizing a random feature subset selection at each node.35 The LASSO regression, SVM-RFE, and Random Forest (RF) algorithms have been utilized to screen feature genes. Genes derived from these three methods are evaluated for intersections, with overlapping genes selected as the initial differential genes employing the “VennDiagram” package in the R language to construct Venn diagrams. This process includes the generation of violin and line plots showcasing the initial differential genes in each dataset.

Final Validation and Accuracy Assessment

Utilize the “pROC” package to generate a preliminary ROC curve for the differential gene and carefully evaluate the area under the curve (AUC). AUC values below 0.7 suggest lower accuracy for disease diagnosis, while values above 0.7 indicate higher accuracy. The preliminary differential genes identified were validated in the GSE117999 dataset to derive the final set for all datasets. Subsequently, violin plots and ROC curves for the final differential genes in the GSE117999 dataset were generated, and the area under the curve was computed. The analysis involved the application of relative logarithmic expression (RLE) on all datasets (p < 0.05), leading to the identification of DEGs and their associated genes. Subsequently, heatmaps and volcano plots of the DEGs were generated. Co-expression analysis was conducted using the “corrplot” package, and heatmaps depicting the associations among multiple genes were also produced.

Establishment of OA Rat Model and Collection of Articular Cartilage

All experiments were approved by the Animal Research Ethics Committee of the Guangxi Medical University. All experiments were performed following the Guangxi Medical University and national guidelines and regulations. A total of 20 8-week-old Sprague–Dawley (SD) male rats (range: 220 to 250 g) were purchased from the Beijing Vetonex Laboratory Animal Technology Co., Ltd. (Beijing, China). All rats were maintained in a specific pathogen-free (SPF) facility and randomly divided into two groups (n = 10 for each): the control group and the OA group. To create a post-traumatic OA model, the rats were subjected to an anterior cruciate ligament transaction plus medial meniscus resection (ACLT + MMx).36 Briefly, animals were anesthetized using 20% urethane (20 mL/kg) via intraperitoneal injection. After shaving and disinfection, a longitudinal incision approximately 3 cm long was made on the medial side of the knee joint, with an upward extension along the patellar ligament. The muscles were retracted, and the knee was flexed to dislocate the patella laterally for joint cavity exposure. Subsequently, both anterior and posterior cruciate ligaments and the medial collateral ligament were severed, allowing for the identification and excision of the medial meniscus. Surgical incisions were then sutured in sequence.

All rats were euthanized at week 20 postoperatively using an overdose of sodium pentobarbital (150 mg/kg via intraperitoneal injection). Five rats from each group were sacrificed for the collection of cartilage, which was immediately preserved at −80°C for subsequent analysis by a blinder (XH).

qRT-PCR

The knee joint cartilage from both groups of rats was collected. The total RNA in the cartilage was extracted using the kit method (R0027, Beyotime, China). Subsequently, the RNA was reverse-transcribed into cDNA following the HyperScripTM First-Strand cDNA Synthesis SuperMix Protocol (K1073, APExBIO, USA). Then PCR experiments were conducted using the HotStartTM Universal 2X SYBR Green qPCR Master Mix (K1170, APExBIO, USA) on an ABI Real-time PCR (7500) instrument, and the 2ΔΔCT values were calculated as per the defined method. Primer sequences are shown in Table 2.

Table 2 Primer Sequences for PCR

Western Blot

Total protein was extracted from the cartilage using RIPA buffer (strong) (K1020, APExBIO, USA) containing protease and phosphatase inhibitors. After determining the protein concentration of each sample using a BCA protein assay kit, protein samples (20 μg of total protein) underwent electrophoresis which was performed on a 10% SDS-PAGE gel prepared with the PAGE gel-fast preparation kit (EpiZyme, PG212) and transferred onto a 0.45 µm-diameter nitrocellulose membrane (MerckMillipore, IPVH00010, Germany). The membranes were blocked for 1 h at room temperature and subsequently incubated with anti-ZNF652 primary antibody (A18657, ABclonal, China, 1:1000) overnight at 4 °C. Then the secondary antibody (31460, Invitrogen, USA, 1:10000) was incubated with the blotting membrane at room temperature for 1 h. β-actin (AF7018, Affinity, USA, 1:20,000) was used as the control. Subsequently, the membranes were visualized using a chemiluminescence reagent (K1129, APExbio, USA) and analyzed quantitatively using GraphPad8 software (GraphPad, USA).

Single Gene Enrichment Analysis

We conducted GO and KEGG enrichment analysis on 2.5 single genes. Additionally, we performed GO and Gene Set Enrichment Analysis (GSEA) enrichment analysis on single gene clusters from different groups, utilizing the gene set “c2.cp.go.symbols.gmt” as a predefined gene set. To ascertain the significance of gene set correlation, we relied on the normalized enrichment score (NES), nominal p-value, and false detection rate (FDR) q-value. A NES > 1, FDR q-value < 0.25, and nominal p-value < 0.05 were considered as criteria for significant enrichment.

Single Gene Set Variant Analysis

We utilized the molecular characterization database (MSigDB, http://software.broadinstitute.org/gsea/msigdb/index.jsp) from the GSEA website (https://www.gsea-msigdb.org/gsea/index.jsp) to access the GO and KEGG gene set databases for our analysis. Subsequently, “c5.go.Hs.symbols.gmt” was chosen as the reference gene set, and gene symbol matrix files were obtained during the identification of DEGs. Following this, the gene symbol matrix files from the GO and KEGG databases were processed, and we evaluated the functional pathways based on the absolute enrichment of gene sets in each sample. Additionally, we acquired functional and pathway matrix files containing the Gene Set Variation Analysis (GSVA) scores for the functional pathways corresponding to each sample. Subsequently, we employed the “limma” package to compare the GSVA scores of OA and normal cartilage samples. The screening criteria for differential GSVA scores were set at p < 0.05 and |log FC| > 0.5.

Analysis of Immune-Related Functions

An immune-related functions analysis was conducted on the dataset of the final differential genes to determine the expression of immune functions within these genes. The CIBERSORT algorithm was utilized to assess the infiltration of immune cells in the OA, resulting in the acquisition of the immune cell infiltration matrix at a significance level of p < 0.05.

Immune Infiltrating Cell Correlation and Differential Analysis

The “ggplot2” package was utilized for correlation analysis and visualization purposes. Spearman correlation analysis was employed to investigate the relationship between the characteristic genes and differentially infiltrated immune cells. Furthermore, efforts were made to assess batch effects between different groups to examine the correlation between immune cell content and the expression of characteristic genes. The correlation coefficient (r > 0.6, p < 0.05) was applied to evaluate the correlation between the characteristic genes and immune infiltrated cells.

ceRNA Network Construction

The prediction of miRNA target genes was facilitated through the utilization of three databases: miRanda, miRDB, and TargetScan. The spongeScan database was employed to predict lncRNA-miRNA interactions, and Cytoscape was utilized for the construction of ceRNA networks.

Statistical Analysis

The bioinformatics analysis was carried out using R (R-4.4.0, Vienna University of Economics and Business, USA). qRT-PCR and Western blotting data were tested for normality using the Shapiro–Wilk test. Data were compared between the groups using an independent t-test for normally distributed variables and are presented as the mean ± standard deviation for the variables. Statistical significance was set at p < 0.05. SPSS version 20.0 (IBM Corp, Chicago, IL, USA) was used for all statistical analyses.

Results Identification of DEGs

Following the completion of batch correction and data normalization, a gene expression matrix (Figure 1A and B) was derived, revealing differential gene expression patterns associated with OA. Subsequent analysis identified 307 DEGs, comprising 183 down-regulated and 124 up-regulated genes. These findings were visually presented through the use of heatmaps and volcano maps (Figure 1C and D).

Figure 1 Initial screening of DEGs. (A and B) Batch pre-correction and post-correction for the GSE10575, GSE16464, GSE117999, GSE169077, GSE178557, and GSE236924 datasets. (C) Volcano plots of DEGs, with red color denoting up-regulated genes, green color denoting down-regulated genes, and gray color denoting non-differentiated genes. (D) Heatmap of DEGs, with highly expressed DEGs in the samples marked in red and low expressed DEGs in the samples marked in blue.

Abbreviation: DEG, differentially expressed gene.

Counting of DEGs Using Three Combined Machine Algorithms

The LASSO regression algorithm pinpointed PLIN2, DDIT3, EN2, BNIP3, PAIP2B, CLGN, GLRX, CHST7, HOXB2, PSKH1, DOCK10, and PIK3CG (Figure 2A and B). The SVM-RFE algorithm highlighted DDIT3, PSKH1, CCL19, CHST7, BNIP3, CDC42EP1, PAIP2B, EXTL2, MYOT, NR1D2, PPIC, and RRAGD (Figure 2C and D). Similarly, the Random Forest (RF) Graph Algorithm revealed thirteen characterized genes: MYL1, CEBPG, ACADL, CKM, EN2, PAIP2B, MYOT, DDIT3, RLF, TNNT1, EIF4EBP1, and SFTPB (Figure 2E and F). The genes characterized using the three methods were compared to obtain three overlapping genes: DDIT3, ZNF652, and PAIP2B. These are the initial disease-characterized genes (Figure 2G).

Figure 2 Identification of the three DEGs. (A and B) Lasso regression counting differential genes, regression graphs, and cross-validation graphs, respectively, screening 13 differential genes. (C and D) SVM machine learning methods screening feature genes, accuracy graphs, and cross-validation graphs, respectively, attaining 13 differential genes. (E and F) Gene cross-validation graphs and significance scores, respectively, screening 13 differential genes. (G) Wayne graphs of feature genes screened by the three machine learning methods to screen feature genes in a Wayne plot. Finally, three genes (ie, DDIT3, ZNF652, and PAIP2B) were screened out.

Abbreviations: DEG, differentially expressed gene; SVM, Support Vector Machine.

Intergroup Differences and Accuracy of Three Genes

In our analysis, violin plots (Figure 3A-C), total fold plots (Figure 3D), and ROC curves (Figure 3E-G) were generated to demonstrate the differential expression of the DDIT3, ZNF652, and PAIP2B genes across the GSE10575, GSE16464, GSE169077, GSE178557, and GSE236924 datasets. They reveal a decrease in DDIT3 expression and an increase in ZNF652 expression in the OA group, alongside a decrease in PAIP2B expression. The area under the ROC curves for all the three genes surpassed 0.7.

Figure 3 Differential expression and ROC curves of the three genes. (A–C) The differential expression of DDIT3, ZNF652, and PAIP2B genes in the control and OA groups. (D) The line plots of DDIT3, ZNF652, and PAIP2B genes. (E–G) The ROC curves and areas under the curves of DDIT3, ZNF652, and PAIP2B genes. *** p < 0.001.

Abbreviation: ROC, receiver operating characteristic.

ZNF652 Expression Was Elevated in the Validation Group Dataset

Upon validation of the DDIT3, ZNF652, and PAIP2B genes in the GSE117999 dataset, it was found that only ZNF652 displayed differential expression, aligning with previous findings (Figure 4A). Subsequently, the ROC curve for ZNF652 in the GSE117999 dataset (Figure 4B) indicated an area under the curve exceeding 0.7, solidifying its status as the final identified differential gene.

Figure 4 Expression of ZNF652 in the GSE117999 dataset. (A) The violin plot of ZNF652 differential expression. (B) The ROC curve and area under the curve for ZNF652. * p < 0.05.

Abbreviation: ROC, receiver operating characteristic.

Expression of Zfp652 in Cartilage of OA Rats

Upon investigation, it was determined that ZNF652 corresponds to human gene nomenclature, while Zfp652 corresponds to rat gene nomenclature. Subsequently, the mRNA and protein sequences of the two genes were obtained from the NCBI database (https://www.ncbi.nlm.nih.gov/). Then sequence comparisons were performed, revealing a high similarity in sequences between the two genes (21% and 100% of sequence similarity in mRNA and protein, respectively; Supplementary Figure 1 and 2; Supplementary Text 1 and 2). Furthermore, the qRT-PCR and Western blot results indicated that the expression of ZFN652 mRNA (1.193 ± 0.005 [95% CI 1.18 to 1.20] vs 1.000 ± 0.005 [95% CI 0.99 to 1.01], p < 0.001; Figure 5A), MMP13 mRNA (1.100 ± 0.027 [95% CI 1.05 to 1.15] vs 1.000 ± 0.019 [95% CI 0.96 to 1.04], p = 0.027; Figure 5B), and ZFN652 protein (0.981 ± 0.055 [95% CI 0.87 to 1.09] vs 0.856 ± 0.026 [95% CI 0.81 to 0.91], p = 0.012; Figure 5C and D) was elevated in the OA group than the control group, aligning with our initial prediction.

Figure 5 Expression of ZNF652 and MMP13 in cartilage of the two groups of rats. (A and B) The mRNA expression of ZNF652 and MMP13 in the cartilage of OA and control rats using qRT-PCR. (C and D) The protein expression of ZNF652 in the cartilage of the OA and control rats using Western blot. *p < 0.05; *** p < 0.001.

Abbreviation: OA, osteoarthritis.

Related Genes of ZNF652 in OA

Following the identification of the DEG, Spearman correlation analysis was conducted on the expression matrix of ZNF652 to illustrate its associated genes in a heatmap (Figure 6A). Subsequently, a visual representation depicting the correlation of the DEGs linked to ZNF652 was generated (Figure 6B).

Figure 6 Differential analysis of ZNF652 gene groupings. (A) The ZNF652 differential analysis heatmap; (B) The ZNF652 differential analysis correlation graph.

OA is Associated with Diverse Functional States and Disease Pathways in Chondrocytes

We conducted GO and KEGG enrichment analyses of the downloaded sequences using R software to discern the shared characteristics of these genes in terms of molecular functions and biological processes. Based on the results of the GO analysis, the differentially expressed genes were primarily associated with biological processes and other pathways, closely linked to functions such as elastic fiber assembly, macrophage value-added, NO anabolic regulation, and OA-associated ECM organization metabolism in cartilage (Figure 7A-E). The KEGG analysis revealed that the OA group showed enrichment in most of the disease pathways, such as the PI3K-Akt signaling pathway, MAPK signaling pathway, calcium signaling pathway, renin-angiotensin system, and Ras signaling pathway (Figure 8A-E).

Figure 7 GO analysis. (A) The GO analyzes chord chords. (B) The GO analyzes cluster plots. (C) The GO analyzes Circos circle plots. (D) The GO analyzes bar plots. (E) The GO analyzes bubble plots.

Abbreviation: GO, Gene Ontology.

Figure 8 KEGG analysis. (A) The KEGG analysis of chord chords. (B) The KEGG analysis of cluster plots. (C) The KEGG analysis of Circos circle plots. (D) The KEGG analysis of bar plots. (E) The KEGG analysis of bubble plots.

Abbreviation: KEGG, Kyoto Encyclopedia of Genes and Genomes.

ZNF652 is Closely Related to Various Pathways Associated with OA

The potential role of ZNF652 in distinguishing between diseased and normal samples was further elucidated through GSEA of the single gene associated with ZNF652, reflecting the characteristics of the regulatory pathway of ZNF652. Upon detailed scrutiny, it came to light that ZNF652 may be correlated with mitochondrial protein and gene anabolism, vascular development, glucose metabolism, transforming growth factor β (TGF-β) signaling pathway, Wnt signaling pathway, and Toll-like receptor pathway, all of which are pertinent to the ECM components of OA (Figure 9A-D). Subsequently, a comparative analysis of the activation pathways between the high and low expression groups in the single gene bound to GSVA was performed. The findings revealed that heightened levels of ZNF652 transcript in the diseased state triggered in vivo substance metabolism, while concurrently repressing metabolic processes such as growth factor binding and inhibiting the melanin synthesis pathway, the Wnt signaling pathway, and the TGF-β signaling pathway (Figure 10A and B). These alterations likely contributed to the up-regulation of ZNF652, consequently leading to cartilage modifications in OA.

Figure 9 GSEA analysis of ZNF652. (A) The enrichment function of genes in the high-expression group in ZNF652. (B) The enrichment function of genes in the low-expression group in ZNF652. (C) The enrichment pathway of genes in the high-expression group in ZNF652. (D) The enrichment pathway of genes in the low-expression group in ZNF652.

Abbreviation: GSEA, Gene Set Enrichment Analysis.

Figure 10 GSVA analysis of ZNF652. (A) The functional enrichment of genes in the high and low expression groups of genes in ZNF652. (B) The pathway enrichment of genes in the high and low expression groups of genes in ZNF652.

Abbreviation: GSVA, Gene Set Variation Analysis.

ZNF652 is Associated with Immune Cells and Their Functions in OA

The study revealed variations in the functions of immune cells, including antigen-presenting cells (APC) and mast cells, in two distinct groups associated with ZNF652 (Figure 11A). Subsequently, a correlation analysis was performed on immune cells concerning ZNF652, demonstrating that M0 macrophages, M2 macrophages, and resting mast cells exhibited correlations with ZNF652. In comparison to the control group, the number of M0 cells and M2 macrophages was reduced in OA (Figure 11B and C), showcasing a negative correlation with ZNF652. Conversely, the quantity of resting mast cells was elevated in the OA group (Figure 11D) and positively correlated with ZNF652. Additionally, Figure 11E visually illustrates the results of the correlation analysis for the three distinct cell types.

Figure 11 Analysis of ZNF652 and immune cell-related functions and immune cell correlation. (A) The functional violin plot of ZNF652 correlation with immune cells. (B) The graph of ZNF652 correlation with M0 macrophage. (C) The graph of ZNF652 correlation with M2 macrophage. (D) The graph of ZNF652 correlation with mast cell. (E) The graph of bubble plot of ZNF652 correlation with three types of immune cells. *p < 0.05; **p < 0.01.

CeRNA Network

We have established a ceRNA-based network associated with ZNF652, using the starBase and miRanda databases. The network consists of 168 nodes, including 1 marker gene, 38 miRNAs, and 130 lncRNAs, revealing intricate inter-regulatory connections. We found that hsa-miR-1229-3p, hsa-miR-1224-5p, hsa-miR-331-3p, hsa-miR-27a-5p, hsa-miR-1237-3p, hsa-miR-539-5p, hsa-let-7a-3p, hsa-miR-224-3p, hsa- miR-142-3p, hsa-miR-765, hsa-miR-618, hsa-miR-205-5p, hsa-miR-106a-5p, hsa-miR-199a-5p, hsa-miR-518d-5p, hsa-miR-570-3p, hsa-miR-590-3p, hsa-miR-576-5p, hsa-miR-135a-5p, hsa-let-7a-5p, hsa-miR-342-3p, hsa-miR-561-3p, hsa-miR-490-5p, hsa-miR-302a-3p, hsa-miR-194-5p, hsa-miR 125a-5p, hsa-miR-186-5p, hsa-miR-7-5p, hsa-miR-513a-3p, hsa-miR-130a-3p, hsa-miR-28-5p, hsa-miR-335-3p, hsa-miR-924, hsa-miR-190a-5p, hsa- miR-345-5p, hsa-miR-508-5p, hsa-miR-421, and hsa-miR-340-5p probability regulate ZNF652 (Figure 12).

Figure 12 A ceRNA networks of ZNF652 and the potential RNA regulatory pathways.

Discussion

The OA pathogenesis remains complex and not fully elucidated. Currently, the early diagnosis of OA is challenging and there is no effective disease-modifying treatment for OA. Nonetheless, several clinical studies have demonstrated an association between the development of OA and abnormal mechanical loading, aging, and specific signaling pathways.21,37 In the current study, 307 DEGs were identified as significantly associated with OA. Furthermore, after a rigorous screening and validation process employing three mechanistic learning tools, the gene ZNF652 emerged as a notable candidate which was associated with synthetic metabolism and immune response, various pathways, and immune cells and their functions in OA. The results were validated using a rat animal of OA. Our findings suggest that ZNF652 plays an important role in OA pathogenesis and can serve as an early diagnostic and therapeutic target for OA.

Elevated expression of ZNF652 has been implicated in the pathogenesis of various malignant tumors, including breast cancer, bladder cancer, hepatocellular carcinoma, and osteosarcoma.32,38–41 Furthermore, prior studies showed that the gene ZNF652 is implicated in cellular inflammatory responses.42 The expression of ZNF652 is notably elevated in LPS-induced human embryonic lung cells, consequently impeding cellular function.42 Suppression of ZNF652 promotes LPS-induced proliferation of human embryonic lung cells and mitigates the inflammatory response.42 The downregulation of ZNF652 alleviates LPS-induced injury in human embryonic lung cells by upregulating MiR-302e.42 However, the correlation between ZNF652 and OA has not been reported. On the other hand, a previous study found an association between ZNF652 and immune function alterations.39 Consistently, in the current study, we explored novel molecular detection targets through immune infiltration analysis. This connection aligns with the hypothesis that OA development is intricately linked to immune function.43

In this study, the GO and KEGG analysis revealed notable differences in biological processes, encompassing elastic fiber assembly, macrophage function, and various substance metabolic pathways. Furthermore, the DEGs in the two groups were predominantly enriched in the phosphatidylinositol 3-kinase (PI3K)/protein kinase B (PKB/AKT), mitogen-activated protein kinases (MAPK), calcium signaling, renin-angiotensin, and Ras pathways. The PI3K/Akt and MAPK signaling pathways play significant roles in cartilage alterations during OA development. The up-regulation of these pathways has a substantial impact on chondrocyte degeneration, protective biological effects, and the synthesis of cartilage matrix in OA cartilage.44–47 A previous study revealed that components associated with the PI3K/Akt signaling pathway were activated in the joints of OA rats.31 Upon inhibiting the expression of these factors, an anti-inflammatory response was observed, along with the restoration of type II collagen expression and a reduction in the expression of Matrix metalloproteinases (MMP13) and ADAMTS5. This resulted in the alleviation of clinical symptoms and a deceleration in OA progression.31 A prior study reported that the activation of the PI3K/Akt signaling pathway exacerbated chondrocyte inflammation, increased the severity of joint damage in rats, and promoted the development of OA.48 It is worth noting that the PI3K/AKT pathway can interact with the MAPK pathway to modulate OA.49 A study found that PI3K/Akt/NF-κB and MAPK pathways exert protective effects on joints in both in vivo and in vitro studies. Furthermore, the inhibition of these pathways resulted in diminished inflammatory responses and cartilage degeneration in OA mice.50 The activation of the PI3K/Akt signaling pathway promotes mitochondrial autophagy, upregulates mTOR protein expression, and attenuates OA.51 Moreover, in coordination with MAPK and other signaling pathways, the PI3K/Akt pathway can suppress M1 polarization and facilitate M2 polarization via TGF-β/SMAD, TLR/NF-κB, JAK/STAT, and other signaling pathways. This can regulate the OA inflammation, contributing to the treatment of the disease.52 Additionally, the MAPK/PI3K-AKT pathway plays a crucial role in alleviating OA by inhibiting synovial inflammation and fibrosis induced by M1 polarization of synovial macrophages.53,54 Thus, based on the findings of the current study and those of the prior studies, it is plausible to propose a regulatory interrelationship between the PI3K/AKT and MAPK pathways, macrophage polarization, and OA.

Our findings suggest that ZNF652 contributes to the immune pathogenesis of OA. As a nuclear regulator, when up-regulated ZNF652 can transmit its pathogenic effects through macrophages, potentially involving signaling pathways such as PI3K/AKT, MAPK, and macrophage polarization. Our immunofunctional results indicate an increase in mast cell number in the OA group. Numerous studies have demonstrated the involvement of mast cells in OA, wherein mast cells are abundant in synovial tissue and their presence is correlated with structural damage.55,56 Mast cells are frequently located in the joint synovial membrane, particularly in the vicinity of blood vessels and nerve endings. They are closely linked to the onset of joint pain.57 Multiple studies have documented an increase in mast cells and their degranulation in individuals with OA. A prior study found the exclusive immune cell types in the synovial tissue of OA patients, in contrast to mast cells in patients with rheumatoid arthritis.58 The activation of mast cells triggers the release of pro-inflammatory mediators, leading to synovial inflammation, bone remodeling, and cartilage damage.59

The analysis of immune cell correlation in the current study indicates differences between M0 and M2 macrophages in the OA group, both types of which were down-regulated in the OA group. Macrophages are ubiquitously present in tissues and are fundamentally involved in tissue remodeling, homeostasis, and regeneration.60 M0 macrophages are derived from monocytes and promptly migrate to the injury site and secrete inflammatory factors such as monocyte chemotactic protein-1 (MCP-1). Depending on the healing stage at the injury site, M0 polarizes to the appropriate activation pathway and, contingent upon the in vivo changes, transitions sequentially into two main distinct phenotypes, M1 macrophages and M2 macrophages. M1 macrophages are often referred to as pro-inflammatory macrophages, while M2 macrophages are termed anti-inflammatory macrophages. The activation and phenotypic transition of M0 macrophages from a pro-inflammatory (M1) to an anti-inflammatory (M2) state may significantly influence the inflammatory response and the onset of related diseases.61 The infiltration of macrophages is a prominent pathological characteristic of OA. Macrophages are recognized contributors to OA progression. Thus, the imbalance of M1 and M2 macrophages in OA is associated with the disease severity.62,63 In OA, macrophage activation results in the release of pro-inflammatory cytokines and MMPs, initiating cartilage damage and bone regrowth via catabolic and anabolic imbalances.64 The activation of pro-inflammatory M1 macrophages triggers the release of pro-inflammatory cytokines, subsequently activating chondrocytes and inducing MMP production, perpetuating a cyclical process of cartilage degradation.65 Our findings, as well as related studies, suggest a hypothesis that the activation of M1 macrophages is associated with inflammation in OA.

We also formulated the ceRNA network centering on the characterized genes. Long non-coding RNAs (lncRNAs) have the capability to compete with mRNAs for miRNA binding, ultimately regulating the expression of the ZNF652-characterized genes.66,67 The identification of miRNAs that interact directly with mRNAs, and lncRNAs that indirectly modulate mRNA expression, presents novel potential targets for the prevention or treatment of OA.

This study has some limitations. First, the human data were derived from the GEO database focusing on the cartilage of OA patients and healthy individuals. However, the sample sizes were relatively small. Moving forward, it is imperative to validate the findings of the current study in a larger dataset. Second, our results were validated in an animal (rat) model of OA but not in OA patients. Thus, there was no additional clinical evidence to confirm the diagnostic and therapeutic effect of inhibiting ZNF652 in OA. Future comparative clinical studies are needed to examine the findings of the current study. Third, although we identified an association between ZNF652 and OA, we did not comprehensively understand the biological functions of ZNF652 in OA. Therefore, more studies are required to investigate the molecular and signaling mechanisms of this gene in OA. Last, the observed differences in immune cells were limited to mast cells and M0 and M2 macrophages between the OA and the control groups. Thus, the role of other immune cells such as T cells and B cells in OA warrant further investigation.

In conclusion, this study found that ZNF652 was related to OA pathogenesis and can potentially serve as a diagnostic and therapeutic target of OA. The underlying mechanism is that ZNF652 was related to nitric oxide anabolism, macrophage proliferation, various signaling pathways, and immune cells and their functions in OA. Nevertheless, the findings need to be confirmed in clinical trials and the molecular mechanism requires further future studies.

Ethics Approval and Consent to Participate

The study was conducted with the approval of the Clinical Research Ethics Committee and Animal Research Ethics Committee of Guangxi Medical University (Approved Number: 202008012).

Acknowledgments

This study was supported by grants from National Natural Science Foundation of China (82060406, 82360429), Joint Project on Regional High-incidence Diseases Research of Natural Science Foundation of Guangxi (2022JJA141126), Advanced Innovation Teams and Xinghu Scholars Program of Guangxi Medical University, and China Postdoctoral Science Foundation (2019M650235).

Author Contributions

Yan Chen is the corresponding author of this article and has made significant contributions. All authors have participated in conceptualizing the topic, designing the experimental study, implementing, and analyzing data. They have also been involved in acquiring and interpreting graphs and figures, drafting, revising, and critically reviewing the article. Furthermore, all authors have given their final approval for the version to be published, agreed to the journal for submission, and have accepted responsibility for all aspects of the work.

Disclosure

The author(s) report no conflicts of interest in this work.

References

1. Tong L, Yu H, Huang X, et al. Current understanding of osteoarthritis pathogenesis and relevant new approaches. Bone Res. 2022;10(1):60. doi:10.1038/s41413-022-00226-9

2. Chen Y, Sun Y, Pan X, Ho K, Li G. Joint distraction attenuates osteoarthritis by reducing secondary inflammation, cartilage degeneration and subchondral bone aberrant change. Osteoarthritis Cartilage. 2015;23(10):1728–1735. doi:10.1016/j.joca.2015.05.018

3. MacKay C, Sale J, Badley EM, Jaglal SB, Davis AM. Qualitative study exploring the meaning of knee symptoms to adults ages 35-65 years. Arthritis Care Res. 2016;68(3):341–347. doi:10.1002/acr.22664

4. Prieto-Alhambra D, Judge A, Javaid MK, Cooper C, Diez-Perez A, Arden NK. Incidence and risk factors for clinically diagnosed knee, Hip and hand osteoarthritis: influences of age, gender and osteoarthritis affecting other joints. Ann Rheumatic Dis. 2014;73(9):1659–1664. doi:10.1136/annrheumdis-2013-203355

5. Wei G, Lu K, Umar M, et al. Risk of metabolic abnormalities in osteoarthritis: a new perspective to understand its pathological mechanisms. Bone Res. 2023;11(1):63. doi:10.1038/s41413-023-00301-9

6. Cho Y, Jeong S, Kim H, et al. Disease-modifying therapeutic strategies in osteoarthritis: current status and future directions. Exp Mol Med. 2021;53(11):1689–1696. doi:10.1038/s12276-021-00710-y

7. Liao Z, Umar M, Huang X, et al. Transient receptor potential vanilloid 1: a potential therapeutic target for the treatment of osteoarthritis and rheumatoid arthritis. Cell Proliferation. 2023;57(3):e13569. doi:10.1111/cpr.13569

8. Zeng D, Chen Y, Liao Z, et al. Cartilage organoids and osteoarthritis research: a narrative review. Front Bioeng Biotechnol. 2023;11:1278692. doi:10.3389/fbioe.2023.1278692

9. Carr AJ, Robertsson O, Graves S, et al. Knee replacement. Lancet. 2012;379(9823):1331–1340. doi:10.1016/S0140-6736(11)60752-6

10. Deere K, Whitehouse MR, Kunutsor SK, et al. How long do revised and multiply revised knee replacements last? An analysis of the national joint registry. Lancet Rheumatol. 2021;3(6):e438–e446. doi:10.1016/S2665-9913(21)00079-5

11. Zhu X, Chan YT, Yung PSH, Tuan RS, Jiang Y. Subchondral bone remodeling: a therapeutic target for osteoarthritis. Front Cell Dev Biol. 2020;8:607764. doi:10.3389/fcell.2020.607764

12. Battistelli M, Favero M, Burini D, et al. Morphological and ultrastructural analysis of normal, injured and osteoarthritic human knee menisci. Eur J Histochem. 2019;63(1). doi:10.4081/ejh.2019.2998

13. Emmi A, Stocco E, Boscolo-Berto R, et al. Infrapatellar fat pad-synovial membrane anatomo-functional unit: microscopic basis for Piezo1/2 mechanosensors involvement in osteoarthritis pain. Front Cell Dev Biol. 2022;10:886604. doi:10.3389/fcell.2022.886604

14. Wen CY, Chen Y, Tang HL, Yan CH, Lu WW, Chiu KY. Bone loss at subchondral plate in knee osteoarthritis patients with hypertension and type 2 diabetes mellitus. Osteoarthritis Cartilage. 2013;21(11):1716–1723. doi:10.1016/j.joca.2013.06.027

15. Chen Y, Wang T, Guan M, et al. Bone turnover and articular cartilage differences localized to subchondral cysts in knees with advanced osteoarthritis. Osteoarthritis Cartilage. 2015;23(12):2174–2183. doi:10.1016/j.joca.2015.07.012

16. Chen Y, Huang YC, Yan CH, et al. Abnormal subchondral bone remodeling and its association with articular cartilage degradation in knees of type 2 diabetes patients. Bone Res. 2017;5(1):17034. doi:10.1038/boneres.2017.34

17. Deng R, Zhao R, Zhang Z, et al. Chondrocyte membrane-coated nanoparticles promote drug retention and halt cartilage damage in rat and canine osteoarthritis. Sci Transl Med. 2024;16(735):eadh9751. doi:10.1126/scitranslmed.adh9751

18. Jeon OH, Kim C, Laberge RM, et al. Local clearance of senescent cells attenuates the development of post-traumatic osteoarthritis and creates a pro-regenerative environment. Nat Med. 2017;23(6):775–781. doi:10.1038/nm.4324

19. Tawonsawatruk T, Sriwatananukulkit O, Himakhun W, Hemstapat W. Comparison of pain behaviour and osteoarthritis progression between anterior cruciate ligament transection and osteochondral injury in rat models. Bone Joint Res. 2018;7(3):244–251. doi:10.1302/2046-3758.73.BJR-2017-0121.R2

20. Makitsubo M, Adachi N, Nakasa T, Kato T, Shimizu R, Ochi M. Differences in joint morphology between the knee and ankle affect the repair of osteochondral defects in a rabbit model. J Orthop Surg Res. 2016;11(1):110. doi:10.1186/s13018-016-0444-4

21. Chen Y, Hu Y, Yu YE, et al. Subchondral trabecular rod loss and plate thickening in the development of osteoarthritis. J Bone Mineral Res. 2018;33(2):316–327. doi:10.1002/jbmr.3313

22. Zhang M, Yang H, Wan X, et al. Prevention of injury-induced osteoarthritis in rodent temporomandibular joint by targeting chondrocyte CaSR. J Bone Miner Res. 2019;34(4):726–738. doi:10.1002/jbmr.3643

23. Legrand C, Ahmed U, Anwar A, et al. Glycation marker glucosepane increases with the progression of osteoarthritis and correlates with morphological and functional changes of cartilage in vivo. Arthritis Res Ther. 2018;20(1):131. doi:10.1186/s13075-018-1636-6

24. Danalache M, Jacobi LF, Schwitalle M, Hofmann UK. Assessment of biomechanical properties of the extracellular and pericellular matrix and their interconnection throughout the course of osteoarthritis. J Biomech. 2019;97:109409. doi:10.1016/j.jbiomech.2019.109409

25. Pettenuzzo S, Arduino A, Belluzzi E, et al. Biomechanics of chondrocytes and chondrons in healthy conditions and osteoarthritis: a review of the mechanical characterisations at the microscale. Biomedicines. 2023;11(7):1942. doi:10.3390/biomedicines11071942

26. Chen Y, Zhang B, Bao L, et al. ZMYND8 acetylation mediates HIF-dependent breast cancer progression and metastasis. J Clin Invest. 2018;128(5):1937–1955. doi:10.1172/JCI95089

27. Wu X, Zhang X, Yu L, et al. Zinc finger protein 367 promotes metastasis by inhibiting the Hippo pathway in breast cancer. Oncogene. 2020;39(12):2568–2582. doi:10.1038/s41388-020-1166-y

28. Chen L, Wu X, Xie H, et al. ZFP57 suppress proliferation of breast cancer cells through down-regulation of MEST-mediated Wnt/β-catenin signalling pathway. Cell Death Dis. 2019;10(3):169. doi:10.1038/s41419-019-1335-5

29. Liu J, Liu L, Yagüe E, et al. GGNBP2 suppresses triple-negative breast cancer aggressiveness through inhibition of IL-6/STAT3 signaling activation. Breast Cancer Res Treat. 2019;174(1):65–78. doi:10.1007/s10549-018-5052-z

30. Kahlert C, Kalluri R. Exosomes in tumor microenvironment influence cancer progression and metastasis. J Mol Med. 2013;91(4):431–437. doi:10.1007/s00109-013-1020-6

31. Liu J, Jia S, Yang Y, et al. Exercise induced meteorin-like protects chondrocytes against inflammation and pyroptosis in osteoarthritis by inhibiting PI3K/Akt/NF-κB and NLRP3/caspase-1/GSDMD signaling. Biomed Pharmacother. 2023;158:114118. doi:10.1016/j.biopha.2022.114118

32. Li Y, Zang H, Zhang X, Huang G. Exosomal Circ-ZNF652 promotes cell proliferation, migration, invasion and glycolysis in hepatocellular carcinoma via miR-29a-3p/GUCD1 axis. Cancer Manag Res. 2020;12:7739–7751. doi:10.2147/CMAR.S259424

33. Wang Y, Liu J, Ma J, et al. Exosomal circRNAs: biogenesis, effect and application in human diseases. Mol Cancer. 2019;18(1):116. doi:10.1186/s12943-019-1041-z

34. Zhang Y, Xia R, Lv M, et al. Machine-learning algorithm-based prediction of diagnostic gene biomarkers related to immune infiltration in patients with chronic obstructive pulmonary disease. Front Immunol. 2022;13:740513. doi:10.3389/fimmu.2022.740513

35. Asadi S, Roshan S, Kattan MW. Random forest swarm optimization-based for heart diseases diagnosis. J Biomed Inform. 2021;115:103690. doi:10.1016/j.jbi.2021.103690

36. Hayami T, Pickarski M, Wesolowski GA, et al. The role of subchondral bone remodeling in osteoarthritis: reduction of cartilage degeneration and prevention of osteophyte formation by alendronate in the rat anterior cruciate ligament transection model. Arthritis Rheum. 2004;50(4):1193–1206. doi:10.1002/art.20124

37. Chen Y, Zeng D, Wei G, et al. Pyroptosis in osteoarthritis: molecular mechanisms and therapeutic implications. J Inflamm Res. 2024;17:791–803. doi:10.2147/JIR.S445573

38. Lei T, Xiao B, He Y, Sun Z, Li L. High expression of ZNF652 promotes carcinogenesis and progression of breast cancer. Nan Fang Yi Ke Da Xue Xue Bao. 2020;40(12):1732–1739. doi:10.12122/j.issn.1673-4254.2020.12.06

39. Liu Y, Peng Y, Du W, et al. PD-L1-mediated immune evasion in triple-negative breast cancer is linked to the loss of ZNF652. Cell Rep. 2023;42(11):113343. doi:10.1016/j.celrep.2023.113343

40. Ke H, Zhang J, Wang F, Xiong Y. ZNF652-induced circRHOT1 promotes SMAD5 expression to modulate tumorigenic properties and nature killer cell-mediated toxicity in bladder cancer via targeting miR-3666. J Immunol Res. 2021;2021:7608178. doi:10.1155/2021/7608178

41. Cheng C, Zhang H, Dai Z, Zheng J. Circular RNA circVRK1 suppresses the proliferation, migration and invasion of osteosarcoma cells by regulating zinc finger protein ZNF652 expression via microRNA miR-337-3p. Bioengineered. 2021;12(1):5411–5427. doi:10.1080/21655979.2021.1965695

42. Wang H, Wang Q. Circ_ZNF652 regulates the proliferation and apoptosis of LPS-induced WI-38 cells via miR-302e/TLR4 axis. Transpl Immunol. 2022;74:101641. doi:10.1016/j.trim.2022.101641

43. van den Bosch MHJ, van Lent P, van der Kraan PM. Identifying effector molecules, cells, and cytokines of innate immunity in OA. Osteoarthritis Cartilage. 2020;28(5):532–543. doi:10.1016/j.joca.2020.01.016

44. Yang Q, Nanayakkara GK, Drummer C, et al. Low-intensity ultrasound-induced anti-inflammatory effects are mediated by several new mechanisms including gene induction, immunosuppressor cell promotion, and enhancement of exosome biogenesis and docking. Front Physiol. 2017;8:818. doi:10.3389/fphys.2017.00818

45. Yu F, Li M, Yuan Z, et al. Mechanism research on a bioactive resveratrol- PLA-gelatin porous nano-scaffold in promoting the repair of cartilage defect. Int J Nanomed. 2018;13:7845–7858. doi:10.2147/IJN.S181855

46. Kim S, Han S, Kim Y, et al. Tankyrase inhibition preserves osteoarthritic cartilage by coordinating cartilage matrix anabolism via effects on SOX9 PARylation. Nat Commun. 2019;10(1):4898. doi:10.1038/s41467-019-12910-2

47. He A, Ning Y, Wen Y, et al. Use of integrative epigenetic and mRNA expression analyses to identify significantly changed genes and functional pathways in osteoarthritic cartilage. Bone Joint Res. 2018;7(5):343–350. doi:10.1302/2046-3758.75.BJR-2017-0284.R1

48. Xu K, He Y, Moqbel SAA, Zhou X, Wu L, Bao J. SIRT3 ameliorates osteoarthritis via regulating chondrocyte autophagy and apoptosis through the PI3K/Akt/mTOR pathway. Int J Biol Macromol. 2021;175:351–360. doi:10.1016/j.ijbiomac.2021.02.029

49. Sun K, Luo J, Guo J, Yao X, Jing X, Guo F. The PI3K/AKT/mTOR signaling pathway in oste