Corrigendum: Identification of immune-associated genes in diagnosing osteoarthritis with metabolic syndrome by integrated bioinformatics analysis and machine learning
- 1The Graduate School, Heilongjiang University of Chinese Medicine, Harbin, China
- 2The Third Affiliated Hospital of Heilongjiang University of Chinese Medicine, Harbin, China
- 3Teaching and Research Section of Orthopedics and Traumatology, Heilongjiang University of Chinese Medicine, Harbin, China
- 4The Second Department of Orthopedics and Traumatology, The Second Affiliated Hospital of Heilongjiang University of Chinese Medicine, Harbin, China
- 5The Bone Injury Teaching Laboratory, Heilongjiang University of Chinese Medicine, Harbin, China
Background: In the pathogenesis of osteoarthritis (OA) and metabolic syndrome (MetS), the immune system plays a particularly important role. The purpose of this study was to find key diagnostic candidate genes in OA patients who also had metabolic syndrome.
Methods: We searched the Gene Expression Omnibus (GEO) database for three OA and one MetS dataset. Limma, weighted gene co-expression network analysis (WGCNA), and machine learning algorithms were used to identify and analyze the immune genes associated with OA and MetS. They were evaluated using nomograms and receiver operating characteristic (ROC) curves, and finally, immune cells dysregulated in OA were investigated using immune infiltration analysis.
Results: After Limma analysis, the integrated OA dataset yielded 2263 DEGs, and the MetS dataset yielded the most relevant module containing 691 genes after WGCNA, with a total of 82 intersections between the two. The immune-related genes were mostly enriched in the enrichment analysis, and the immune infiltration analysis revealed an imbalance in multiple immune cells. Further machine learning screening yielded eight core genes that were evaluated by nomogram and diagnostic value and found to have a high diagnostic value (area under the curve from 0.82 to 0.96).
Conclusion: Eight immune-related core genes were identified (FZD7, IRAK3, KDELR3, PHC2, RHOB, RNF170, SOX13, and ZKSCAN4), and a nomogram for the diagnosis of OA and MetS was established. This research could lead to the identification of potential peripheral blood diagnostic candidate genes for MetS patients who also suffer from OA.
1 Introduction
Osteoarthritis (OA) is one of the common degenerative diseases of the musculoskeletal system in orthopedics, often involving one or more joints, and has a high prevalence worldwide (1), and according to WHO, there is a 10% chance of OA occurring in people over 60 years of age worldwide (2). So far, the causative factors and pathogenesis of OA are not clear, but most scholars believe that its pathogenesis is a combination of several factors, such as obesity, aging, trauma, excessive joint strain, metabolic disorders, inflammation, and genetics (3).
Metabolic syndrome (MetS) is a metabolic disorder characterized by abdominal obesity, dyslipidemia, hyperglycemia, and hypertension (4). According to the National Institute for Health and Nutrition Examination Survey, around one in five adults in the United States suffer from MetS (5). Most patients with metabolic syndrome are associated with varying degrees of obesity, which increases stress on the entire musculoskeletal system and poses a higher risk for the development of osteoarthritis and various other musculoskeletal disorders. Osteoarthritis and metabolic syndrome share common mechanisms of inflammation, oxidative stress, and metabolic dysfunction in their etiology (6). Traditionally, osteoarthritis is a non-inflammatory disease affected by trauma or metabolic dysregulation, and age-related joint degeneration is thought to be a causal factor in the development of the disease (7, 8). However, a growing body of evidence suggests that low-grade inflammation may be a key factor driving the pathogenesis of OA (9). Inflammation and metabolic disorders play a very important role in the progression of osteoarthritis (10, 11). Risk factors, including diabetes, hypertension, and hyperlipidemia, are largely involved in osteoarthritis through the release of inflammatory and adipokines that accelerate the progression of osteoarthritis by driving articular cartilage degeneration and bone marrow lesions (12). Several theories describe how MetS risk factors affect the progression of OA, such as high blood pressure can cause subchondral ischemia, abnormal lipids can cause lipid deposition in chondrocytes, and high blood glucose can cause oxidative stress and low inflammation, eventually leading to cartilage destruction (13). In addition, central obese patients also have abnormal leptin and lipocalin levels, which can further aggravate the development of OA (6, 14). To date, metabolic syndrome-associated osteoarthritis (MetS-OA) has been well characterized as a distinct phenotype of osteoarthritis (15, 16), and the goal of this model is to investigate the relationship between inflammatory response and metabolic disorders in order to improve the concept of treating MetS and to aid in the reduction of inflammatory response in OA patients. There is some evidence to support this, but the mechanism of action in metabolic syndrome-related osteoarthritis is still being investigated. In recent years, comprehensive bioinformatics analysis has been used to identify novel genes associated with various diseases that can serve as diagnostic and prognostic biomarkers. However, the common diagnostic and interlinked genes of OA and MetS are not known. Therefore, this study used a bioinformatics approach to screen for biomarkers associated with immune infiltration in both, which could help identify immune-related potential diagnostic markers for OA in patients with MetS.
2 Materials and methods
2.1 Microarray data
The datasets used were all from the NCBI Gene Expression Omnibus (GEO) database (17), and the OA datasets were GSE169077, GSE55457, and GSE55235 (18); while the MetS dataset was GSE98895 (19). Figure 1 depicts the study flowchart.
Figure 1 Study flowchart, GSE, gene expression omnibus series; WGCNA, weighted gene co-expression network analysis; Limma, linear models for microarray data; DEGs, differentially expressed genes.
2.2 Data processing and differential gene screening
The following steps were all carried out using the R software (https://www.bioconductor.org/), starting with “affy” to calibrate, log-transform, and normalize the three OA datasets, then merging them and using “SVA” to remove batch effects, with p-values set to 0.05 and |log2 fold change (FC)|>1.5, and “limma” for differential gene screening. Table 1 presents detailed dataset information, including the microarray platform, sample groups, and numbers.
2.3 Weighted gene co-expression network analysis
WGCNA (20) was used to investigate the gene modules most closely linked to MetS. First, the top 50% of genes with the highest median absolute deviation (MAD) were filtered. The expression matrix was then filtered to remove ineligible data. Third, use a “soft” threshold power (β) for co-expression of similarity to calculate adjacency. Then, using dynamic tree cuts and hierarchical clustering, a topological overlap matrix (TOM) was created to group genes into modules by random colors, and a gene dendrogram was constructed using a TOM-based measure of phase dissimilarity and a minimum gene cluster size (n=100). Fifth, for the next step of the study, the dissimilarity of the module genes was calculated, and the average linkage hierarchy clustering was performed. It is finally visualized.
2.4 Functional enrichment analysis
Sangerbox (http://www.sangerbox.com/tool) was used for Gene Ontology (GO) (21) and Kyoto Encyclopedia of Genes and Genomes (KEGG) (22) enrichment analysis.
2.5 Machine learning
Machine learning algorithms are used to screen the core genes for OA diagnosis. “glmnet” (23) was used for LASSO (24) regression, and “randomForest” (25) was used for RF (26) analysis. The two intersecting genes were used as the core genes for OA diagnosis.
2.6 Nomogram construction and evaluation of recipient operating characteristics
Nomogram was created for the screened genes using the “rms” (27) R package, and its value in the clinic was determined by the area under the curve (AUC) and 95% CI using ROC. When the AUC exceeds 0.7, it is considered to have diagnostic value.
2.7 Immune infiltration analysis
The CIBERSORT algorithm is used to determine the proportion of immune cells in cells or tissues (28). The bar graphs show the proportion of each type of immune cell in various samples, and the “corrplot” (29) R package is used to generate a heat map of the correlation between 22 immune cells. The vioplot was used to visualize the differences between the OA and normal immune cell groups.
2.8 Statistical analysis
R software version 4.2.2 and Sangerbox were used for all statistical analyses. To compare normally distributed continuous variables between two groups, the Independent Student’s t-test was used, and the Wilcoxon rank sum test was used to analyze non-normally distributed variables. Every statistical test was two-sided. The statistical significance level was set at p-Value <0.05.
3 Results
3.1 Differentially expressed genes
The integrated OA dataset yielded a total of 2263 DEGs after LIMMA analysis, with 1341 up-regulated genes and 922 down-regulated genes. Figures 2A, B show the heat map and volcano map generated from the above data. Furthermore, the MetS dataset yielded 1449 DEGs, including 605 up-regulated genes and 844 down-regulated genes (Figures 3A, B).
Figure 2 Heatmap and valcano plot for the DEGs identified from the integrated OA dataset. (A) Each row shows the DEGs, and each column refers to one of the samples of OA cases or controls. The red and blue represent DEGs with upregulated and downregulated gene expression, respectively. (B) Red and green plot triangles represent DEGs with upregulated and downregulated gene expression, respectively. OA, Osteoarthritis: DEGs, differentially expressed genes.
Figure 3 Identification of DEGs via Limma and module genes via WGCNA in MetS. (A) The heatmap displays the top 50 upregulated and downregulated DECs identified from MetS dataset. Each row represents the intersection of genes, and each column represents one of MetS cases or controls. Red and blue represent upregulated and downregulated ge expression. (B) The volcano plot shows all DEGs, of which red and green triangles refer to significant DEGs. (C, D) b= 14 is selected as the soft threshold with the combined analysis of scale independence and average connectivity. (E) Clustering dendrogram of the MetS and control samples. (F) Gene co-expression modules represented by different colors under the gene tree (G) Heatmap of eigengene adjacency. (H) Heatmap of the association between modules and Mers. The yellow module is shown to be correlated significantly with MetS. Numbers at the top and bottom brackets represent the correlation coefficient and p-value, respectively. (I) Correlation plot between module membership and gene significance of genes included in the yellow module WGCNA, weighted gene co-expression network analysis, Limma, linear models for microarray data; DEGs, differentially expressed genes.
3.2 Selection of key modules
In MetS, use WGCNA to filter critical modules. When the soft threshold =14, the scale-free network performs best (Figures 3C, D). Figure 3E depicts the clustering tree graph for both groups and the ten randomly colored gene modules obtained (Figures 3F, G). The yellow module (691 genes) had the highest MetS correlation (Figure 3H) and can be used as a key module in the following analysis. Figure 3I shows the results of a correlation analysis for the genes in the yellow module, which showed a positive correlation (r=0.62).
3.3 Functional enrichment analysis of metabolic syndrome
The intersection of 1449 DEGs and 691 modular genes yielded a total of 108 genes (Figure 4A). The “PI3K-Akt signaling pathway” and “glycerophospholipid metabolism” were primarily enriched in KEGG analysis (Figure 4B). According to GO analysis, the main enrichment in biological process (BP) terms was in “B cell activation involved in immune response” and “cellular response to DNA damage stimuli” (Figure 4C). It is primarily found in “cell membrane,” “nuclear lumen,” and “vesicles” in terms of cellular component (CC) ontology (Figure 4D). It was primarily enriched in “molecular function regulators,” “enzyme regulator activity,” and “signaling receptor binding” in molecular function (MF) analysis (Figure 4E).
Figure 4 Enrichment analysis of the intersection of genes in MetS. (A) Venn diagram shows that 108 genes are identified from the intersection of DEGs via Limma and yellow module genes via WGCNA. (B) KEGG pathway analysis of the intersection of genes. Different colors represent various significant pathways and related enriched genes. (C–E) GO analysis of the Intersection of genes, including biological process, cellular component, and molecular function, respectively. The y-axis represents different GO terms, the x-axis represents gene ratio enriched in relative GO terms, the circle size refers to gene numbers, and the color represents p-value. MetS, Metabolic syndrome: KEGG, Kyoto Encyclopedia of Genes and Genomes GO, Gene Ontology: WGCNA, weighted gene co-expression network analysis: Limma, linear models for microarray datac DEGs, differentially expressed genes.
3.4 Enrichment analysis of osteoarthritis with metabolic syndrome
To further explore whether key MetS-related genes might be associated with the pathogenesis of OA, 82 genes were identified by the Venn diagram from the intersection of OA DEGs and MetS key module genes (Figure 5A). According to KEGG analysis, 82 genes were primarily enriched in “inflammatory mediator regulation of TRP channels,” “cellular senescence,” and “MAPK signaling pathway,” all of which are closely related to the immune system (Figure 5B). GO analysis revealed that they were primarily enriched in “immune response,” “cell death,” and “immune system process” (BP); “integral component of endoplasmic reticulum membrane,” “stereocilium mem-brane” and “immunological synapse,” (CC); and “catalytic activity,” “ATP binding,” and “sequence-specific DNA binding” (MF) (Figures 5C–E).
Figure 5 Enrichment analysis of common genes from OA with MetS. (A) Venn diagram shows that 82 common genes are identified from the intersection of genes in OA using Limma and MetS using WGCNA. (B) KEGG analysis of 82 common genes. (C–E) GO analysis (hiological process, cellular component, and molecular function) of 82 common genes. OA. Osteoarthritis: MetS, Metabolic syndrome, WGCNA, weighted gene co-expression network analysis.
3.5 Core genes screening using machine learning
The Lasso regression and RF algorithms were used to identify core genes and create relevant nomograms for ROC analysis. Lasso regression screened 10 candidate genes (Figures 6A, B), and the RF algorithm identified 30 most important genes (Figures 6C, D), and the two were taken to intersect (Figure 6E), resulting in the identification of eight genes (FZD7, IRAK3, KDELR3, PHC2, RHOB, RNF170, SOX13, and ZKSCAN4).
Figure 6 Machine learning in screening candidate diagnostic biomarkers for OA with MetS (A, B) Biomarkers screening in the Lasso model. The number of genes (n-10) corresponding to the lowest point of the curve is the most suitable for OA with MetS diamos. (C, D) The random forest algorithm shows the error in OA; control group and genes are ranked based on the importance score. (E) Venn diagram shows that eight candidate diagnostic genes are identified via the above two algorithms OA. Osteoarthritis: MetS Metabolic syndrome.
3.6 Determining diagnostic value
We created a nomogram (Figure 7A) and plotted ROC curves based on the eight candidate genes to assess the diagnostic value of each gene. The calculated AUCs and 95% confidence intervals were as follows: FZD7 (AUC 0.86, CI 0.96–0.76), IRAK3 (AUC 0.92, CI 0.99–0.84), KDELR3 (AUC 0.94, CI 1.00–0.87), PHC2 (AUC 0.89, CI 0.99–0.80), RHOB (AUC 0.96, CI 1.00–0.90), RNF170 (AUC 0.82, CI 0.94–0.70), SOX13 (AUC 0.83, CI 0.94–0.72), ZKSCAN4 (AUC 0.90, CI 0.99–0.81) (Figures 7B–I). The findings indicated that the acquired genes had a high value for the diagnosis of OA in combination with MetS.
Figure 7 Nogrim comtraction and the diagnostic value evaluation (A) The visible nomogram for diagnosing OA with MetS (B–I) The ROC curve of each candidate gene (RHOR, KDELR3, TRAKA, PHC2, ZKSCANA, SOX13, RNF 170 and FZD7) and nomogram show the significant OA with MetS diagnostic value. OA. Osteoarthritis, MetS, Metabolic syndrome; AUC, area under the curve.
3.7 Immune infiltration analysis
We discovered that genes associated with MetS can also play a role in OA, primarily in immune regulation. An in-depth examination of the nomogram and ROC revealed that they could be used as potential biomarkers for the diagnosis of OA, which was confirmed by immune infiltration.
The percentage of 22 immune cells in each sample is shown in the bar graph for both datasets (Figure 8A). Voltammograms showed higher levels of B cells memory, Macrophages M0, Dendritic cells activated, Mast cells resting and Eosinophils in OA patients, and lower levels of T cells CD4 memory resting, Macrophages M2, Mast cells activated and Neutrophils (Figure 8B). Correlation analysis of the 22 immune cell types showed that NK cells resting was positively correlated with Neutrophils (r=0.62), T cells CD4 naive was positively correlated with T cells CD4 memory activated (r=0.52), Mast cells activated was negatively correlated with Mast cells resting (r=-0.67), and Macrophages M2 was negatively correlated with T cells gamma delta (r=-0.52) (Figure 8C). In conclusion, OA patients have varying degrees of multiple immune cell infiltration that may be potential regulatory points for the treatment of OA.
Figure 8 Immune cell infiltration analysis between OA and control. (A) The proportion of 22 kinds of immune cells in different samples visualized from the burplot. (B) Comparison regarding the proportion of 22 kinds of immune cells between OA and control groups visualized by the vioplot. (C) Correlation of 22 immune cell type compositions, *,p < 0.05, **,p < 0.01, ***,p < 0.001. Both horizontal and vertical axes demonstrate immune cell subtypes OA, Osteoarthritis.
4 Discussion
OA is one of the most significant causes of disability in middle-aged and older adults and places a significant burden on public health (30). Several biomarkers for OA have been identified in recent studies, including osteoclast protein (OPN), cartilage oligomeric matrix protein (COMP), cartilage acidic protein 1, and CRTAC1 (31, 32). Studies that combine these two diseases are relatively rare. There were currently no markers for diagnosing OA using nomogram and machine learning methods, so we combined both bioinformatics analysis and machine learning and evaluated their diagnostic value using nomogram and ROC. Notably, we identified eight key immune-related candidate genes (FZD7, IRAK3, KDELR3, PHC2, RHOB, RNF170, SOX13, and ZKSCAN4) and developed a nomogram for the diagnosis of OA in MetS patients.
We were able to approximate the odds of developing OA in MetS patients by simply testing the peripheral blood of MetS patients with known expression of core genes, given that the MetS samples used in this study were all from peripheral blood. Because of its simplicity and effectiveness, peripheral blood testing is used in the diagnosis of many diseases (33, 34). Following that, we will develop a more refined model that can accurately reflect gene expression and assign it a numerical value to make it more accurate for diagnosis (35). When significant changes in target indicators occur in MetS patients, early monitoring and intervention can be performed, which is more valuable for diagnosing MetS combined with OA.
Frizzled 7 (FZD7) belongs to the G protein-coupled receptor family and is a 7-channel transmembrane cell surface receptor and a key factor in the Wnt signaling pathway (36). Frizzled proteins have a universal structure: a cysteine structural domain (CRD) that fills the extracellular space, followed by a structural domain containing seven putative transmembrane fragments (37). CRD can interact with Wnts and form a complex that reduces β-catenin phosphorylation (38). The FZD family proteins mediate the Wnt/β-catenin signaling pathway and contribute to the overall pathological progression of OA, especially in the remodeling of cartilage (39, 40). Exos derived from hADSC have been shown to contain MiR-376c-3p, which has the ability to target WNT3 or WNT9a and inhibit the Wnt/β-catenin signaling pathway, reducing synovial fibrosis and chondrocyte degradation in OA patients. In conclusion, FZD7 can regulate cell proliferation, differentiation, migration, and tumorigenesis development by activating downstream signaling pathways through binding to Wnt ligands (41, 42).
IRAK3, an interleukin-1 receptor-associated kinase (IRAK) family member, is primarily expressed in immune cells such as macrophages, monocytes, dendritic cells, and epithelial cells (43). As an important regulatory protein of the TLRs/IL-1R pathway, it is involved in the inhibition of Toll-like receptor signaling. In addition, IRAK1 and IRAK2 can dissociate from the myosin complex, and IRAK3 inhibits it, rendering IRAK1 and IRAK2 unable to interact with TRAF6 and thus inhibiting NF-B pathway activation, and reduces inflammation by reducing the production of pro-inflammatory factors (44). Of course, another possibility is that the MyD88-IRAK puposome interacts with IRAK3, which stimulates NF-B activation via complex formation, where both MEKK3 and TRAF6 can form complexes with IRAK3, which is also dependent on the NF-B pathway, to promote anti-inflammatory-related expression (45). An experiment discovered that miR-33b-3p regulated IRAK3 and that it was effective in alleviating IL-1-induced apoptosis and inflammation. Overall, this appears to be a promising OA target (46).
KDELR3, a member of the KDEL family, is in charge of encoding a protein related to the endoplasmic reticulum (ER) (47). The ER influences not only the synthesis and transport of lipids and steroids but also the activity of hormones and glucose metabolism. As a result, KDELR3 expression varies significantly between individuals, with higher levels of expression in non-atherosclerotic tissues. In addition, several studies have found that KDELR3 is significantly aberrantly expressed in a variety of malignancies, including prostate adenocarcinoma and hepatocellular carcinoma (48, 49). As a result, we hypothesize that KDELR3 plays a critical role in OA patients.
The HD1 structural domain of PHC2 is close to its FCS zinc finger structural domain, which is capable of interacting with the catalytic ring structural domain of RING1B. RING1B, a core component of Polycomb Repressive Complex 1 (PRC1), can monoubiquitinate histone H2A (H2AK119ub1) and plays a catalytic role at lysine 119. It is a key regulator of the estrogen receptor alpha (ER) transcriptional process (50), through which estrogen regulates biological processes like reproductive maturation, energy homeostasis, and skeletal growth by binding to ERα (51). Most importantly, estrogen action on ERα can increase miR-140 expression levels and decrease the expression levels of matrix metalloproteinase-13 (MMP-13) in human articular chondrocytes (52). It has also been shown that PHC2 can bind to the Vcam1 locus and act to reduce systemic immunodeficiency (53). It is reasonable to believe that PHC2 can have an impact on the progression of OA.
RHOB is a Rho GTPase with a C-terminus and an N-terminus, the latter of which contains a G domain, also known as the RhoA-like domain (54). Because chromosome 2 retains RHOB’s genetic information and lacks the alternative RHOB precursor mRNA, only one RHOB protein sequence can be translated (55). RHOB can be induced by transforming growth factor beta (TGF) signaling and is immune to p53 regulation, allowing for a timely response to non-genotoxic (polysaccharide, hypoxia, inflammatory factors) and genotoxic stimuli (radiation) (56). RHOB deficiency inhibited pathological angiogenesis in ischemic retina patients, implying that RHOB alleviates symptoms by promoting the formation of lymphatic vessels following injury or in the presence of inflammation (57). Interestingly, RHOB also activates IL-1β (interleukin 1β), LPS (lipopolysaccharide), and TNFα (tumor necrosis factor α) in an inflammatory setting (54), contrary to previous observations. However, its mechanism of action in OA is unclear and remains to be investigated. In general, RHOB can regulate a variety of cellular processes, including vesicle trafficking, apoptosis, DNA repair, angiogenesis, proliferation, migration, and invasion (58).
RNF170 is a novel E3 ubiquitin ligase that mediates the ubiquitination of the endoplasmic reticulum calcium channel sarcoglycan 1,4,5 trisphosphate receptor leading to the degradation of the sarcoglycan 1,4,5 trisphosphate receptor via the proteasome pathway and subsequently affecting the calcium flow and content in the endoplasmic reticulum lumen and cytoplasm. By specifically targeting TLR receptor molecules, RNF170 mediates polyubiquitination of lysine at position 766 on the TIR domain of TLR3 (59), resulting in the degradation of TLR3 receptor molecules via the proteasome pathway and reducing the effect on TLR3 downstream signaling. activation of the TLR3 downstream signaling pathway, thereby inhibiting the production of inflammatory factors and type I interferons (60). A search of the GEO database revealed that RNF170 protein expression was significantly increased in human monocytes after infection with human immunodeficiency virus (HIV) (61), and when the body’s fibroblasts underwent an immune response to cytomegalovirus infection, the expression level of RNF170 protein in fibroblasts showed a trend of first increasing and then decreasing (62), and in general, RNF170 is extremely important in natural immunity.
SOX13 is a member of the SOX protein family, and its coding gene is located at 1q31.3-32.1, which is found in a wide range of cells and tissues. SOX13 contains 604 amino acids, including three specific functional regions, namely the leucine zipper region, the glutamine-rich sequence region, and the HMG functional region (63). Its expression in three embryonic cell lineages suggests that it may be involved in various developmental processes, with expression detected in the neural tube, developing brain, kidney, liver, mesenchymal, and chondrogenic progenitor cells, and has a significant impact on sex determination, neurogenesis and endochondral bone formation (64). Innate lymphocytes, including type 1, type 2 (ILC2) and type 3 (ILC3) subpopulations, are increasingly recognized as key regulators of tissue homeostasis and inflammation through the release of cytokines (65), with NCR ILC3s producing mainly IL-17A, which is regulated by DNA binding inhibitor 2 (ID2), transcription factor 1, and SRY (sex-determining region Y)-box transcription factor 13 (Sox13) (66). Fida et al. simultaneously tested 297 cases of primary biliary sclerosis (PBC), 22 cases of autoimmune cholangitis, 29 cases of autoimmune hepatitis, and 90 patients with T1DM for SOX13-Ab and found a high rate of positivity (67). In addition, Sox13, as an autoimmune antigen, can also modulate T-cell specificity, and in summary, it is reasonable to suspect that SOX13 can modulate the inflammatory response in OA.
ZKSCAN4, a zinc finger protein family member, localizes to chromosome 6p21-p22.1 and regulates genomic stability, stem cell generation, and telomere elongation (68, 69). It is found in a variety of tissues, including the cervix, testes, and trachea, as well as the kidneys, adrenal glands, mouth, skin, lungs, brain, spleen, uterus, liver, intestines, and muscles (70). Several authors examined European and Asian populations using genome-wide association studies (GWAS) of genes and discovered significant differences in some genes between healthy people and rheumatoid arthritis patients, especially ZKSCAN4, ABCF1, BTN3A3, BRD2, FLOT1, HLA-DMA, HLA-G, HLA-F, HSP90AB1, IER3, and TUBB (71). ZKSCAN4 improves the HDM2 promoter region’s association with YY1 by modifying the HDM2 chromatin structure to promote its expression (72). HDM2 is a ubiquitin protein ligase that not only promotes the degradation but also the transcription of the tumor suppressor p53. Under hypoxia, p53 promotes cell survival by promoting the production of cellular metabolic energy in OA, inhibiting reactive oxygen species (ROS) production and ROS detoxification (73). Although there are no functional studies to validate ZKSCAN4, it can be speculated that it is tightly linked to OA.
Previous research has revealed that immune and inflammatory responses are present in all phases of OA. The position of inflammatory factors and cellular infiltration in the progression of OA, as well as the various manifestations of OA in different immune settings, provide a theoretical foundation for further research into the relationship between immunity and OA. Many diseases are now recognized as being influenced by the immune microenvironment, and in OA patients, both cartilage damage and repair processes involve immune cells (74). It has been established that B cells and macrophages, as the primary immune cells, are involved in cartilage damage and repair, in addition to NK cells, T cells, and DCs (75). In addition, mast cells (MC), T cells, and macrophages can be found in large numbers in the synovial tissue of OA patients. A unique pattern of immune infiltration has recently been identified, characterized by increased polarization of CD4+ T cells to activated Th1 cells and increased secretion of immunomodulatory cytokines (76). This has similarities to our findings, with OA patients having higher levels of CD8+ T cells, B-cell memory, activated dendritic cells, M0 macrophages, and eosinophils, and lower levels of Mast cells activated, T cells CD4 memory activated, M2 macrophages and Neutrophils. Above all, understanding inflammatory signaling mechanisms is critical for OA diagnosis and treatment.
5 Conclusion
Based on bioinformatics analysis and machine learning, we systematically identified eight related candidate genes (FZD7, IRAK3, KDELR3, PHC2, RHOB, RNF170, SOX13, and ZKSCAN4) and provided a template for the diagnosis of OA combined with MetS. We also noticed that the immune system of MetS patients with OA is out of balance, that the percentage of immune cells can be affected by the immune microenvironment, and that the screened genes could be used for clinical diagnosis and treatment.
6 Limitation
There are some restrictions on our research. First, despite pooling three OA datasets, the samples are still tiny, and because of the small sample size, the diagnostic value of the column line graphs is quite high. Additionally, we wanted to choose a different sample to verify the diagnostic results. The eligible ones weren’t accessible, and the Mets dataset was already small. As a result, we were unable to confirm the results. After that, the findings should be verified in a bigger study with a larger sample size. Second, even though the interactions between candidate hub genes and dysregulated immune cells warrant further study, the eight candidate hub genes are primarily concentrated in regulatory immune pathways.
Data availability statement
The original contributions presented in the study are included in the article/Supplementary Material. Further inquiries can be directed to the corresponding author.
Author contributions
All authors listed have made a substantial, direct, and intellectual contribution to the work, and approved it for publication.
Funding
The State Administration of Traditional Chinese Medicine Construction Project of Inheritance Studio of National Famous Traditional Chinese Medicine Experts (No. 2021-2024). Youth Qihuang Scholar Support Project of the State Administration of Troditional Chinese Medicine (No. 2021-2023). Key Laboratory of Basic and Clinical Research of Osteonecrosis in Heilongjiang Province.
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fimmu.2023.1134412/full#supplementary-material
References
1. 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 Rheum Dis (2014) 73(9):1659–64. doi: 10.1136/annrheumdis-2013-203355
2. Disease GBD, Injury I, Prevalence C. Global, regional, and national incidence, prevalence, and years lived with disability for 354 diseases and injuries for 195 countries and territories, 1990-2017: A systematic analysis for the global burden of disease study 2017. Lancet (London England) (2018) 392(10159):1789–858. doi: 10.1016/S0140-6736(18)32279-7
3. Blagojevic M, Jinks C, Jeffery A, Jordan KP. Risk factors for onset of osteoarthritis of the knee in older adults: A systematic review and meta-analysis. Osteoarthritis Cartilage (2010) 18(1):24–33. doi: 10.1016/j.joca.2009.08.010
4. Fahed G, Aoun L, Bou Zerdan M, Allam S, Bou Zerdan M, Bouferraa Y, et al. Metabolic syndrome: Updates on pathophysiology and management in 2021. Int J Mol Sci (2022) 23(2):786. doi: 10.3390/ijms23020786
5. Beltran-Sanchez H, Harhay MO, Harhay MM, McElligott S. Prevalence and trends of metabolic syndrome in the adult U.S. population, 1999-2010. J Am Coll Cardiol (2013) 62(8):697–703. doi: 10.1016/j.jacc.2013.05.064
6. Zhuo Q, Yang W, Chen J, Wang Y. Metabolic syndrome meets osteoarthritis. Nat Rev Rheumatol (2012) 8(12):729–37. doi: 10.1038/nrrheum.2012.135
7. 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:16044. doi: 10.1038/boneres.2016.44
8. Courties A, Gualillo O, Berenbaum F, Sellam J. Metabolic stress-induced joint inflammation and osteoarthritis. Osteoarthritis Cartilage (2015) 23(11):1955–65. doi: 10.1016/j.joca.2015.05.016
9. Terkawi MA, Ebata T, Yokota S, Takahashi D, Endo T, Matsumae G, et al. Low-grade inflammation in the pathogenesis of osteoarthritis: Cellular and molecular mechanisms and strategies for future therapeutic intervention. Biomedicines (2022) 10(5):1109. doi: 10.3390/biomedicines10051109
10. Zhang H, Cai D, Bai X. Macrophages regulate the progression of osteoarthritis. Osteoarthritis Cartilage (2020) 28(5):555–61. doi: 10.1016/j.joca.2020.01.007
11. Zheng L, Zhang Z, Sheng P, Mobasheri A. The role of metabolism in chondrocyte dysfunction and the progression of osteoarthritis. Ageing Res Rev (2021) 66:101249. doi: 10.1016/j.arr.2020.101249
12. Ramirez-Perez S, Reyes-Perez IV, Martinez-Fernandez DE, Hernandez-Palma LA, Bhattaram P. Targeting inflammasome-dependent mechanisms as an emerging pharmacological approach for osteoarthritis therapy. iScience (2022) 25(12):105548. doi: 10.1016/j.isci.2022.105548
13. Nemet M, Blazin T, Milutinovic S, Cebovic T, Stanojevic D, Zvekic Svorcan J. Association between metabolic syndrome, its components, and knee osteoarthritis in premenopausal and menopausal women: A pilot study. Cureus (2022) 14(7):e26726. doi: 10.7759/cureus.26726
14. Le Clanche S, Bonnefont-Rousselot D, Sari-Ali E, Rannou F, Borderie D. Inter-relations between osteoarthritis and metabolic syndrome: A common link? Biochimie (2016) 121:238–52. doi: 10.1016/j.biochi.2015.12.008
15. Courties A, Sellam J, Berenbaum F. Metabolic syndrome-associated osteoarthritis. Curr Opin Rheumatol (2017) 29(2):214–22. doi: 10.1097/BOR.0000000000000373
16. Jiang X, Zhong R, Dai W, Huang H, Yu Q, Zhang JA, et al. Exploring diagnostic biomarkers and comorbid pathogenesis for osteoarthritis and metabolic syndrome Via bioinformatics approach. Int J Gen Med (2021) 14:6201–13. doi: 10.2147/IJGM.S325561
17. Barrett T, Wilhite SE, Ledoux P, Evangelista C, Kim IF, Tomashevsky M, et al. Ncbi geo: Archive for functional genomics data sets–update. Nucleic Acids Res (2013) 41(Database issue):D991–5. doi: 10.1093/nar/gks1193
18. 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. doi: 10.1186/ar4526
19. D'Amore S, Hardfeldt J, Cariello M, Graziano G, Copetti M, Di Tullio G, et al. Identification of mir-9-5p as direct regulator of Abca1 and hdl-driven reverse cholesterol transport in circulating Cd14+ cells of patients with metabolic syndrome. Cardiovasc Res (2018) 114(8):1154–64. doi: 10.1093/cvr/cvy077
20. Langfelder P, Horvath S. Wgcna: An r package for weighted correlation network analysis. BMC Bioinf (2008) 9:559. doi: 10.1186/1471-2105-9-559
21. The Gene Ontology C. The gene ontology resource: 20 years and still going strong. Nucleic Acids Res (2019) 47(D1):D330–D8. doi: 10.1093/nar/gky1055
22. Kanehisa M, Goto S. Kegg: Kyoto encyclopedia of genes and genomes. Nucleic Acids Res (2000) 28(1):27–30. doi: 10.1093/nar/28.1.27
23. Yang C, Delcher C, Shenkman E, Ranka S. Machine learning approaches for predicting high cost high need patient expenditures in health care. BioMed Eng Online (2018) 17(Suppl 1):131. doi: 10.1186/s12938-018-0568-3
24. Ellis K, Kerr J, Godbole S, Lanckriet G, Wing D, Marshall S. A random forest classifier for the prediction of energy expenditure and type of physical activity from wrist and hip accelerometers. Physiol Meas (2014) 35(11):2191–203. doi: 10.1088/0967-3334/35/11/2191
25. Zhang M, Zhu K, Pu H, Wang Z, Zhao H, Zhang J, et al. An immune-related signature predicts survival in patients with lung adenocarcinoma. Front Oncol (2019) 9:1314. doi: 10.3389/fonc.2019.01314
26. Alderden J, Pepper GA, Wilson A, Whitney JD, Richardson S, Butcher R, et al. Predicting pressure injury in critical care patients: A machine-learning model. Am J Crit Care (2018) 27(6):461–8. doi: 10.4037/ajcc2018525
27. Pan X, Jin X, Wang J, Hu Q, Dai B. Placenta inflammation is closely associated with gestational diabetes mellitus. Am J Transl Res (2021) 13(5):4068–79.
28. Newman AM, Liu CL, Green MR, Gentles AJ, Feng W, Xu Y, et al. Robust enumeration of cell subsets from tissue expression profiles. Nat Methods (2015) 12(5):453–7. doi: 10.1038/nmeth.3337
29. Hu K. Become competent within one day in generating boxplots and violin plots for a novice without prior r experience. Methods Protoc (2020) 3(4):64. doi: 10.3390/mps3040064
30. Pereira D, Peleteiro B, Araujo J, Branco J, Santos RA, Ramos E. The effect of osteoarthritis definition on prevalence and incidence estimates: A systematic review. Osteoarthritis Cartilage (2011) 19(11):1270–85. doi: 10.1016/j.joca.2011.08.009
31. Silverwood V, Blagojevic-Bucknall M, Jinks C, Jordan JL, Protheroe J, Jordan KP. Current evidence on risk factors for knee osteoarthritis in older adults: A systematic review and meta-analysis. Osteoarthritis Cartilage (2015) 23(4):507–15. doi: 10.1016/j.joca.2014.11.019
32. Watt FE. Osteoarthritis biomarkers: Year in review. Osteoarthritis Cartilage (2018) 26(3):312–8. doi: 10.1016/j.joca.2017.10.016
33. Ohe Y, Fushida S, Yamaguchi T, Kinoshita J, Saito H, Okamoto K, et al. Peripheral blood platelet-lymphocyte ratio is good predictor of chemosensitivity and prognosis in gastric cancer patients. Cancer Manag Res (2020) 12:1303–11. doi: 10.2147/CMAR.S241069
34. Rao C, Jin J, Lu J, Wang C, Wu Z, Zhu Z, et al. A multielement prognostic nomogram based on a peripheral blood test, conventional mri and clinical factors for glioblastoma. Front Neurol (2022) 13:822735. doi: 10.3389/fneur.2022.822735
35. OH B, Gransar H, Callister T, Shaw LJ, Schulman-Marcus J, Stuijfzand WJ, et al. Development and validation of a simple-to-Use nomogram for predicting 5-, 10-, and 15-year survival in asymptomatic adults undergoing coronary artery calcium scoring. JACC Cardiovasc Imaging (2018) 11(3):450–8. doi: 10.1016/j.jcmg.2017.03.018
36. MacDonald BT, Tamai K, He X. Wnt/Beta-catenin signaling: Components, mechanisms, and diseases. Dev Cell (2009) 17(1):9–26. doi: 10.1016/j.devcel.2009.06.016
37. Wang Y, Chang H, Rattner A, Nathans J. Frizzled receptors in development and disease. Curr Top Dev Biol (2016) 117:113–39. doi: 10.1016/bs.ctdb.2015.11.028
38. Schulte G. International union of basic and clinical pharmacology. lxxx. the class frizzled receptors. Pharmacol Rev (2010) 62(4):632–67. doi: 10.1124/pr.110.002931
39. Lories RJ, Corr M, Lane NE. To wnt or not to wnt: The bone and joint health dilemma. Nat Rev Rheumatol (2013) 9(6):328–39. doi: 10.1038/nrrheum.2013.25
40. Mariani E, Pulsatelli L, Facchini A. Signaling pathways in cartilage repair. Int J Mol Sci (2014) 15(5):8667–98. doi: 10.3390/ijms15058667
41. Fernandez A, Huggins IJ, Perna L, Brafman D, Lu D, Yao S, et al. The wnt receptor Fzd7 is required for maintenance of the pluripotent state in human embryonic stem cells. Proc Natl Acad Sci U.S.A. (2014) 111(4):1409–14. doi: 10.1073/pnas.1323697111
42. Simonetti M, Agarwal N, Stosser S, Bali KK, Karaulanov E, Kamble R, et al. Wnt-fzd signaling sensitizes peripheral sensory neurons Via distinct noncanonical pathways. Neuron (2014) 83(1):104–21. doi: 10.1016/j.neuron.2014.05.037
43. Miyata M, Lee JY, Susuki-Miyata S, Wang WY, Xu H, Kai H, et al. Glucocorticoids suppress inflammation Via the upregulation of negative regulator irak-m. Nat Commun (2015) 6:6062. doi: 10.1038/ncomms7062
44. Hubbard LL, Moore BB. Irak-m regulation and function in host defense and immune homeostasis. Infect Dis Rep (2010) 2(1):e9. doi: 10.4081/idr.2010.e9
45. Freihat LA, Wheeler JI, Wong A, Turek I, Manallack DT, Irving HR. Irak3 modulates downstream innate immune signalling through its guanylate cyclase activity. Sci Rep (2019) 9(1):15468. doi: 10.1038/s41598-019-51913-3
46. Tao T, Zhang Y, Wei H, Heng K. Downregulation of Irak3 by mir-33b-3p relieves chondrocyte inflammation and apoptosis in an in vitro osteoarthritis model. Biosci Biotechnol Biochem (2021) 85(3):545–52. doi: 10.1093/bbb/zbaa105
47. Trychta KA, Back S, Henderson MJ, Harvey BK. Kdel receptors are differentially regulated to maintain the er proteome under calcium deficiency. Cell Rep (2018) 25(7):1829–40 e6. doi: 10.1016/j.celrep.2018.10.055
48. Bai Y, Qi W, Liu L, Zhang J, Pang L, Gan T, et al. Identification of seven-gene hypoxia signature for predicting overall survival of hepatocellular carcinoma. Front Genet (2021) 12:637418. doi: 10.3389/fgene.2021.637418
49. Xu Z, Xu L, Liu L, Li H, Jin J, Peng M, et al. A glycolysis-related five-gene signature predicts biochemical recurrence-free survival in patients with prostate adenocarcinoma. Front Oncol (2021) 11:625452. doi: 10.3389/fonc.2021.625452
50. Zhang Y, Chan HL, Garcia-Martinez L, Karl DL, Weich N, Slingerland JM, et al. Estrogen induces dynamic eralpha and Ring1b recruitment to control gene and enhancer activities in luminal breast cancer. Sci Adv (2020) 6(23):eaaz7249. doi: 10.1126/sciadv.aaz7249
51. Mauvais-Jarvis F, Clegg DJ, Hevener AL. The role of estrogens in control of energy balance and glucose homeostasis. Endocr Rev (2013) 34(3):309–38. doi: 10.1210/er.2012-1055
52. Liang Y, Duan L, Xiong J, Zhu W, Liu Q, Wang D, et al. E2 regulates mmp-13 Via targeting mir-140 in il-1beta-Induced extracellular matrix degradation in human chondrocytes. Arthritis Res Ther (2016) 18(1):105. doi: 10.1186/s13075-016-0997-y
53. Bae J, Choi SP, Isono K, Lee JY, Park SW, Choi CY, et al. Phc2 controls hematopoietic stem and progenitor cell mobilization from bone marrow by repressing Vcam1 expression. Nat Commun (2019) 10(1):3496. doi: 10.1038/s41467-019-11386-4
54. Eckenstaler R, Hauke M, Benndorf RA. A current overview of rhoa, rhob, and rhoc functions in vascular biology and pathology. Biochem Pharmacol (2022) 206:115321. doi: 10.1016/j.bcp.2022.115321
55. Schaefer A, Reinhard NR, Hordijk PL. Toward understanding rhogtpase specificity: Structure, function and local activation. Small GTPases (2014) 5(2):6. doi: 10.4161/21541248.2014.968004
56. Wang XH, Wang Y, Diao F, Lu J. Rhob is involved in lipopolysaccharide-induced inflammation in mouse in vivo and in vitro. J Physiol Biochem (2013) 69(2):189–97. doi: 10.1007/s13105-012-0201-z
57. Gerald D, Adini I, Shechter S, Perruzzi C, Varnau J, Hopkins B, et al. Rhob controls coordination of adult angiogenesis and lymphangiogenesis following injury by regulating Vezf1-mediated transcription. Nat Commun (2013) 4:2824. doi: 10.1038/ncomms3824
58. Gutierrez E, Cahatol I, Bailey CAR, Lafargue A, Zhang N, Song Y, et al. Regulation of rhob gene expression during tumorigenesis and aging process and its potential applications in these processes. Cancers (Basel) (2019) 11(6):818. doi: 10.3390/cancers11060818
59. Song XQ. Negative regulation of Tlr3-triggered innate immune response by E3 ligase Rnf170 [D]. Peking Union Medical College (2018).
60. Oshiumi H, Matsumoto M, Funami K, Akazawa T, Seya T. Ticam-1, an adaptor molecule that participates in toll-like receptor 3-mediated interferon-beta induction. Nat Immunol (2003) 4(2):161–7. doi: 10.1038/ni886
61. Tilton JC, Johnson AJ, Luskin MR, Manion MM, Yang J, Adelsberger JW, et al. Diminished production of monocyte proinflammatory cytokines during human immunodeficiency virus viremia is mediated by type I interferons. J Virol (2006) 80(23):11486–97. doi: 10.1128/JVI.00324-06
62. Hertel L, Mocarski ES. Global analysis of host cell gene expression late during cytomegalovirus infection reveals extensive dysregulation of cell cycle gene expression and induction of pseudomitosis independent of Us28 function. J Virol (2004) 78(21):11988–2011. doi: 10.1128/JVI.78.21.11988-12011.2004
63. Mertin S, McDowall SG, Harley VR. The DNA-binding specificity of Sox9 and other sox proteins. Nucleic Acids Res (1999) 27(5):1359–64. doi: 10.1093/nar/27.5.1359
64. Wang Y, Ristevski S, Harley VR. Sox13 exhibits a distinct spatial and temporal expression pattern during chondrogenesis, neurogenesis, and limb development. J Histochem Cytochem (2006) 54(12):1327–33. doi: 10.1369/jhc.6A6923.2006
65. Vivier E, Artis D, Colonna M, Diefenbach A, Di Santo JP, Eberl G, et al. Innate lymphoid cells: 10 years on. Cell (2018) 174(5):1054–66. doi: 10.1016/j.cell.2018.07.017
66. Malhotra N, Narayan K, Cho OH, Sylvia KE, Yin C, Melichar H, et al. A network of high-mobility group box transcription factors programs innate interleukin-17 production. Immunity (2013) 38(4):681–93. doi: 10.1016/j.immuni.2013.01.010
67. Fida S, Myers MA, Whittingham S, Rowley MJ, Ozaki S, Mackay IR. Autoantibodies to the transcriptional factor Sox13 in primary biliary cirrhosis compared with other diseases. J Autoimmun (2002) 19(4):251–7. doi: 10.1006/jaut.2002.0622
68. Lee K, Gollahon LS. Zscan4 interacts directly with human Rap1 in cancer cells regardless of telomerase status. Cancer Biol Ther (2014) 15(8):1094–105. doi: 10.4161/cbt.29220
69. Zalzman M, Falco G, Sharova LV, Nishiyama A, Thomas M, Lee SL, et al. Zscan4 regulates telomere elongation and genomic stability in Es cells. Nature (2010) 464(7290):858–63. doi: 10.1038/nature08882
70. Ecker K, Lorenz A, Wolf F, Ploner C, Bock G, Duncan T, et al. A ras recruitment screen identifies Zkscan4 as a glucocorticoid receptor-interacting protein. J Mol Endocrinol (2009) 42(2):105–17. doi: 10.1677/JME-08-0087
71. Zhu H, Xia W, Mo XB, Lin X, Qiu YH, Yi NJ, et al. Gene-based genome-wide association analysis in European and Asian populations identified novel genes for rheumatoid arthritis. PloS One (2016) 11(11):e0167212. doi: 10.1371/journal.pone.0167212
72. Wei SY, Shih YT, Wu HY, Wang WL, Lee PL, Lee CI, et al. Endothelial yin yang 1 phosphorylation at S118 induces atherosclerosis under flow. Circ Res (2021) 129(12):1158–74. doi: 10.1161/CIRCRESAHA.121.319296
73. Humpton TJ, Vousden KH. Regulation of cellular metabolism and hypoxia by P53. Cold Spring Harb Perspect Med (2016) 6(7):a026146. doi: 10.1101/cshperspect.a026146
74. Li M, Yin H, Yan Z, Li H, Wu J, Wang Y, et al. The immune microenvironment in cartilage injury and repair. Acta Biomater (2022) 140:23–42. doi: 10.1016/j.actbio.2021.12.006
75. de Lange-Brokaar BJ, Ioan-Facsinay A, van Osch GJ, Zuurmond AM, Schoones J, Toes RE, et al. Synovial inflammation, immune cells and their cytokines in osteoarthritis: A review. Osteoarthritis Cartilage (2012) 20(12):1484–99. doi: 10.1016/j.joca.2012.08.027
76. Rosshirt N, Hagmann S, Tripel E, Gotterbarm T, Kirsch J, Zeifang F, et al. A predominant Th1 polarization is present in synovial fluid of end-stage osteoarthritic knee joints: Analysis of peripheral blood, synovial fluid and synovial membrane. Clin Exp Immunol (2019) 195(3):395–406. doi: 10.1111/cei.13230
Keywords: differentially expressed genes, osteoarthritis, metabolic syndrome, machine learning, immune infiltration
Citation: Li J, Wang G, Xv X, Li Z, Shen Y, Zhang C and Zhang X (2023) Identification of immune-associated genes in diagnosing osteoarthritis with metabolic syndrome by integrated bioinformatics analysis and machine learning. Front. Immunol. 14:1134412. doi: 10.3389/fimmu.2023.1134412
Received: 30 December 2022; Accepted: 27 March 2023;
Published: 17 April 2023.
Edited by:
Michael V. Volin, Midwestern University, United StatesCopyright © 2023 Li, Wang, Xv, Li, Shen, Zhang and Zhang. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Xilin Xv, ZWFzdG9waEBzaW5hLmNvbQ==; Xiaofeng Zhang, enhmZW5naGxqQDE2My5jb20=
†These authors have contributed equally to this work and share first authorship