Abstract Adolescent idiopathic scoliosis (AIS) is a complex three-dimensional spinal deformity that primarily affects adolescents, and its pathogenesis remains elusive. Several studies have indicated that AIS may be associated with low bone mineral density, and bone marrow mesenchymal stem cells (BMSCs) are known to play a crucial role in bone metabolism. Through high-throughput sequencing of the BMSC transcriptome, we identified 1919 differentially expressed genes and pinpointed IL6ST as a key gene influencing osteogenic differentiation in AIS. Mechanistic experiments revealed that IL6ST regulates the osteogenic differentiation of BMSCs by activating the JAK/STAT3/RANKL/OPG pathway. Moreover, knockout of IL6ST expression significantly increased the malformation rate in zebrafish models. We further explored upstream regulators of IL6ST. Using high-throughput sequencing and bioinformatics analysis, we identified CircSCAF8 as a potential upstream regulatory circRNA of IL6ST. Collectively, these findings suggest that IL6ST may be a pivotal gene in the development of AIS, deepening our understanding of its pathogenesis. Subject terms: Predictive markers, Pathogenesis, Mesenchymal stem cells __________________________________________________________________ IL6ST participates in the development of adolescent idiopathic scoliosis by regulating the osteogenic differentiation of bone marrow mesenchymal stem cells. Introduction Adolescent idiopathic scoliosis (AIS) is a complex three-dimensional deformity of the spine that mainly occurs in adolescent girls. The incidence of AIS is 1–4%, and approximately 10% of patients with AIS require medical intervention^[40]1. Currently, medical interventions for patients with AIS are limited, primarily consisting of full-time bracing and pedicle screw instrumentation correction surgery. The former can easily lead to back pain and psychological issues, while the latter involves significant surgical risks and is often associated with numerous postoperative complications^[41]2–[42]4. However, the pathogenesis of AIS remains unclear. Several hypotheses have been suggested to be related to the occurrence and development of AIS, such as abnormal bone metabolism, genetics, endocrine, and environmental factors. Among these, abnormal bone metabolism has emerged as a significant research focus in understanding the pathogenesis of AIS. Previous studies have revealed a correlation between AIS and low bone mineral density (BMD) as well as abnormalities in bone quality and strength^[43]5–[44]8, and studies have found that 27% to 38% of patients with AIS suffer from osteopenia^[45]9–[46]11. Recent studies assessing bone mechanical properties in patients with AIS have suggested that osteopenia in AIS may be the result of abnormal regulation of bone metabolism and bone modeling/remodeling^[47]12. Bone marrow mesenchymal stem cells (BMSCs) play important roles in bone metabolism. These cells can differentiate into osteoblasts, adipocytes, and chondrocytes^[48]13. The reduced osteogenic differentiation capacity of bone marrow mesenchymal stem cells in patients with AIS may be the cause of abnormal bone metabolism in these patients^[49]14. However, the underlying causes of this phenomenon and the abnormal osteogenic differentiation of BMSCs in AIS warrant further investigation. To this end, we conducted high-throughput sequencing analysis on BMSCs derived from AIS patients, aiming to identify differentially expressed mRNAs and circRNAs. Our objective was to pinpoint key genes that may be implicated in the abnormal bone development observed in AIS patients. IL6ST (gp130) is a mediator of the action of the IL-6 cytokine family. Some studies have revealed the involvement of IL6ST in the pathogenesis of various diseases. Béziat et al. reported that dominant negative mutations in IL6ST can lead to the development of hyper-IgE syndrome^[50]15. Apart from the typical symptoms of recurrent infections, pneumonia, and elevated IgE levels, a significant portion of hyper-IgE syndrome patients also exhibit bone fractures and spinal curvature, indicating that loss of IL6ST function can affect bone development^[51]16. Moreover, several studies have shown that knockout of IL6ST expression in mice can cause skeletal muscle atrophy, reduced endothelial cell generation, tumor development, and alterations in energy metabolism^[52]17–[53]20. In recent years, a series of circRNAs have been shown to play important roles in biological activities related to cancer, nervous system diseases, stem cell development, controlled pluripotency, and regulation of differentiation^[54]21–[55]23. Among them, the role of circRNAs in osteogenesis has begun to attract attention, and some studies have reported changes in circRNAs in osteogenesis through high-throughput sequencing^[56]24,[57]25. But there is a lack of relevant research on the pathogenesis of AIS. This study aimed to identify differentially expressed mRNAs that are highly associated with osteogenic differentiation of AIS BMSCs using high-throughput sequencing and to identify potential upstream regulatory circRNAs using circRNA-mRNA association analysis. Results Patients with AIS showed decreased BMD and osteogenic differentiation of BMSCs BMD measurements were collected from patients with AIS, and the results showed that the patients with AIS had lower BMD than the patients with non-AIS (Fig. [58]1A). BMD is closely related to the osteogenic differentiation of BMSCs, and BMSCs were then extracted from the blood collected from the screw path during the operation. Next, immunofluorescence staining was employed to investigate the expression of the stem cell surface marker CD73. The results revealed positive expression of CD73 in the extracted cells, confirming their identity as stem cells (Fig. [59]S1). After total mRNA and protein collection, we examined the expression levels of osteogenesis-related proteins, and the results showed that the mRNA and protein levels of Runx2 and OCN were significantly decreased in primary BMSCs from AIS (Figs. [60]1B, C and [61]S2). The differentiation potential of BMSCs was validated by Alizarin red and oil red staining, and the results showed that the primary BMSCs in the AIS group had lower osteogenic differentiation and higher adipogenic differentiation than the control BMSCs (Fig. [62]1D, E). In addition, we collected the spinal facet tissues from patients with AIS and non-AIS during the operation. By immunohistochemistry analysis, we found that the expression level of RUNX2 in the bone tissues of the patients with AIS was lower than that of the patients with non-AIS (Fig. [63]1F, G). These results suggested that patients with AIS had reduced bone mineral density, and the reduced osteogenic differentiation of BMSCs in AIS might be the cause of the low BMD. Fig. 1. BMSCs from AIS exhibit decreased osteogenic differentiation and increased adipogenic differentiation. [64]Fig. 1 [65]Open in a new tab A BMD value of patients with AIS and non-AIS, n = 15. B, C Relative protein level of OCN and RUNX2, n = 6. D, E Alizarin Red and Oil Red O staining were performed to detect the differentiation of human primary BMSCs on day 12. Scale bar: 100 μm, n = 6. F, G Representative images showing immunohistochemical staining of RUNX2. Scale bar: 100 μm. The intensity of mineral deposition was quantified by ImageJ, n = 5. Data are shown as the mean ± SD, * means p < 0.05, and ** means p < 0.01 vs. the AIS group. High-throughput sequencing identified IL6ST as a hub gene responsible for the reduced osteogenic differentiation of BMSCs in AIS Total RNA of BMSCs from patients with AIS and non-AIS was collected for sequencing analysis, and a total of 1919 differentially expressed genes were obtained, of which 1139 genes were upregulated and 880 genes were downregulated (Fig. [66]2A, Supplementary Data [67]1). Subsequently, KEGG and GO enrichment analyses of the differentially expressed genes were performed, and the results are shown (Figs. [68]2B and [69]S3). By functional annotation of genes after GO enrichment analysis, we found 34 differentially expressed genes (OBDEGs) related to osteogenic differentiation (Fig. [70]2C). Then, we conducted further analysis on these 34 genes, used the random forest method to analyze the key genes in OBDEGs, and selected the top ten genes for subsequent analysis (Fig. [71]2D). Then, we constructed a PPI based on OBDEGs using the STRING database (Fig. [72]S4), and we used BinGO to assess the biological processes enriched in these 34 OBDEGs (Fig. [73]S5). Subsequently, Cytoscape software was used to extract the hub subnetwork from the PPI network (Fig. [74]2E), and 21 hub genes were selected based on the degree value. These 21 hub genes were intercrossed with 10 key genes obtained by the random forest method to obtain 7 candidate genes (Fig. [75]2F, G). The expression of these seven genes in BMSCs from patients with AIS and non-AIS was then examined by RT-qPCR, and the results showed that JUND and SMO gene expression was increased, while BMP6 and IL6ST expression was decreased (Fig. [76]2H). Functional annotation revealed that these four genes were positively correlated with the osteogenic differentiation of BMSCs. However, the expressions of JUND and SMO were elevated, which was inconsistent with the changes of bone mineral density reduction in AIS patients, so we eliminated the JUND and SMO genes and finally selected IL6ST for further study. Western blotting detected the protein level of IL6ST in BMSCs, and the results were consistent with the RT-qPCR results (Fig. [77]2I, J). Fig. 2. Analysis of differentially expressed genes. [78]Fig. 2 [79]Open in a new tab A The volcano map shows the results of high-throughput sequencing. B KEGG enrichment analysis of the DEGs. C Ring heat map shows 34 osteoblast differentiation DEGs (OBDEGs). D Relative importance of the 34 OBDEGs calculated using the random forest analysis. E PPI network was generated based on 34 OBDEGs. Upregulated genes are marked in red and downregulated genes in blue. F Venn diagram of the intersection of hub genes obtained from the PPI and random forest. G PPI network was generated based on the 7 hub genes. Orange indicated genes upregulated, and blue indicated genes downregulated. H Relative mRNA levels of 7 hub genes, CT represents patients without AIS, n = 5. I, J Relative protein levels of IL6ST in BMSCs, n = 6. Data are shown as the mean ± SD, * means p < 0.05 vs. the control group. IL6ST positively regulated the osteogenic differentiation of BMSCs from patients with AIS To verify the effect of IL6ST knockdown on osteogenic differentiation, we designed three shRNA sequences. The control group consists of cells treated with saline. Knockdown results showed that IL6ST expression in BMSCs of patients with non-AIS was significantly downregulated by sh-RNAs #1 and #3, and sequence 1 was selected for subsequent experiments (Figs. [80]3A and [81]S6). We also used lentiviral transfection to overexpress IL6ST in AIS BMSCs (Figs. [82]3B and [83]S6). Moreover, we detected the variations in IL6ST levels during the osteogenic differentiation of BMSCs. The results indicated no significant difference in IL6ST levels between the induced differentiation group on day 10 and the uninduced group (Fig. [84]S7). We set up four groups: control, control+sh-IL6ST, AIS, and AIS+OE-IL6ST. And then, CCK-8 was used to verify the effect of knockdown and overexpression on cell viability. The results suggested that lentivirus transfection slightly damaged cell activity, but the cells were still in a state of proliferation (Fig. [85]3C). Alizarin red staining was used to verify the effect of IL6ST levels on the osteogenic differentiation of BMSCs. The results showed that knockdown of IL6ST inhibited the osteogenic differentiation of human primary BMSCs, while overexpression of IL6ST promoted the osteogenic differentiation of primary BMSCs from the AIS group (Fig. [86]3D, E). In addition, the trend of oil red staining was opposite to that of alizarin red staining (Fig. [87]3F, G). In addition, we collected the total proteins of transfected BMSCs and detected the expression levels of the osteogenesis-related proteins OCN and RUNX2 by WB. The results showed that knockdown of IL6ST reduced the expression level of osteogenesis-related proteins, while overexpression of IL6ST significantly increased the expression level of osteogenesis-related proteins (Fig. [88]3H–J). These results suggested that IL6ST could increase the osteogenic differentiation of primary BMSCs from AIS. Fig. 3. IL6ST positively regulated the osteogenic differentiation of AIS BMSCs. [89]Fig. 3 [90]Open in a new tab A The knock-down efficiency of IL6ST in BMSCs was detected by RT-qPCR, n = 6. B The overexpress efficiency of IL6ST in BMSCs, n = 6. C Cell viability was detected by CCK-8 analysis, n = 12. D, E Alizarin Red was performed to detect the osteoblast differentiation of BMSCs. Scale bar: 100 μm, n = 6. F, G Oil Red O staining was performed to detect the adipogenic differentiation of BMSCs. Scale bar: 100 μm, n = 6. H–J Relative protein level of OCN and RUNX2, n = 3. Data are shown as the mean ± SD, * means p < 0.05; ** means p < 0.01, and *** means p < 0.001 vs. the control group. # means p < 0.05 and ## means p < 0.01 vs. the AIS group. IL6ST regulates the osteogenic differentiation of BMSCs by activating the JAK/STAT3/RANKL/OPG pathway After KEGG enrichment analysis, IL6ST was predicted to be enriched upstream of the JAK/STAT3 pathway, and some studies also found that the JAK/STAT3 pathway was positively correlated with osteogenic differentiation^[91]26. We first examined the protein levels of the pathways in primary BMSCs from patients with AIS and non-AIS, and the results showed that the phosphorylation levels of JAK1 and STAT3 were significantly reduced in primary BMSCs (Fig. [92]4A, B). We subsequently knocked down IL6ST in primary BMSCs, and the results showed that the phosphorylation level of the JAK/STAT3 pathway was significantly reduced (Fig. [93]4C, D). Then, we added colivelin, a STAT3 phosphorylation activator, to treat the cells and set up control, sh-IL6ST, and sh-IL6ST+colivelin treatment groups. Immunofluorescence and RT-qPCR results for OCN and RUNX2 showed that colivelin treatment significantly increased the expression levels of osteogenesis-related proteins: OCN and RUNX2 (Figs. [94]4E, F and [95]S8). The protein levels of RANKL/OPG, an important regulatory pathway in bone metabolism, were also examined. The results showed that RANKL/OPG levels were significantly increased after IL6ST knockdown and were significantly decreased by colivelin treatment (Figs. [96]4G, H and [97]S9). Additionally, Alizarin red staining was performed to assess the impact of colivelin treatment on the impaired osteogenic differentiation caused by IL6ST knockdown. The results demonstrated a significant improvement in the reduced osteogenic differentiation upon colivelin treatment (Fig. [98]4I, J). These results suggested that IL6ST may regulate the osteogenic differentiation of BMSCs through the activation of the JAK/STAT3 pathway. Fig. 4. IL6ST regulates the differentiation of BMSCs via the JAK/STAT3/RANKL/OPG pathway. [99]Fig. 4 [100]Open in a new tab A, B Relative protein level of p-JAK1/JAK1 and p-STAT3/STAT3 in primary BMSCs, n = 6. C, D Relative protein level of p-JAK1/JAK1 and p-STAT3/STAT3 in primary BMSCs after kd-IL6ST, n = 6. E, F Representative images of OCN immunofluorescence staining in BMSCs from the control, sh-IL6ST, and sh-IL6ST+colivelin treatment groups. Scale bar: 20 μm, n = 4. G, H Relative protein level of RANKL and OPG, n = 4. I, J Alizarin Red was performed to detect the osteoblast differentiation of BMSCs. Scale bar: 100 μm, n = 6. Data are shown as the mean ± SD, * means p < 0.05, and ** means p < 0.01 vs. the control group. # means p < 0.05 and ## means p < 0.01 vs. the sh-IL6ST group. Knockout of IL6ST expression significantly increased the malformation rate of zebrafish To verify the effect of IL6ST on body development, we knocked out the expression of IL6ST in zebrafish using CRISPR‒Cas9 technology. The control group consists of zebrafish treated with saline injection. We found that the malformation rate of zebrafish was significantly increased after knockout (IL6ST^−/− 51% vs. Control 0.1%) (Fig. [101]5A, B), and the mortality rate of the knockout group was significantly higher than that of the control group (88% vs. 22%) (Fig. [102]5C). Then, we collected juvenile fish on the seventh day to extract total RNA and protein after Trizol digestion and detected the expression level of the IL6ST gene. The results showed that mRNA and protein levels of IL6ST in the knockout group were significantly decreased (Fig. [103]S10). Moreover, the mRNA and protein level of osteoblast-related proteins OCN and RUNX2 were decreased in the IL6ST knockout group (Figs. [104]5D, E and [105]S11). Alizarin red staining results showed that trunk bone development of deformed zebrafish was slower than that of the normal controls during the juvenile stage (Fig. [106]5F, G). After digestion of zebrafish soft tissues with 2% potassium hydroxide, alizarin red staining of adult zebrafish also indicates significant abnormalities in the development of zebrafish vertebrae (Fig. [107]5H). These results suggested that IL6ST knockout could cause malformation in zebrafish through abnormal osteogenesis. Fig. 5. Knocking out the IL6ST in zebrafish results in a significant increase in the deformity rate. [108]Fig. 5 [109]Open in a new tab A Representative images of zebrafish in the control group and lateral curved zebrafish after KO-IL6ST. B Incidence of the curved zebrafish phenotype in the control and IL6ST-KO groups, n = 104. C Incidence of unhatched eggs, mortality, and survival in zebrafish in the control and IL6ST-KO groups, n = 104. D, E Relative protein level of RUNX2 and OCN in zebrafish after KO-IL6ST, n = 6. F, G Alizarin Red staining of zebrafish and the number of vertebrae in zebrafish from the two groups were statistically analyzed, n = 14. H Representative images of alizarin red staining of adult zebrafish after digestion of soft tissues. Data are shown as the mean ± SD, * means p < 0.05, ** means p < 0.01, and *** means p < 0.001 vs. the control group. CircSCAF8 is a potential upstream-regulating circRNA of IL6ST After determining the abnormal expression of IL6ST, we further searched for its upstream regulators that affect the expression of IL6ST. Subsequently, we performed high-throughput sequencing of BMSCs from patients with AIS and non-AIS to detect the circRNA levels, and a total of 440 differentially expressed circRNAs (214 upregulated and 226 downregulated) were sequenced (Figs. [110]S12 and [111]6A; Supplementary Data [112]2). Then, we constructed a Sankey diagram based on the competitive endogenous RNA (ceRNA) mechanism through the RNA inter database and Cytoscape software (Fig. [113]6B). Based on these 7 hub genes, the database predicted a total of 73 differentially expressed circRNAs (Supplementary Data [114]3), six of which were potential upstream circRNAs of IL6ST: circANKRD17, circCCSER2, circFOXP1, circHIPK2, circSCAF8 and circSLAIN2. Subsequently, we performed a correlation analysis between these six circRNAs and IL6ST and found that five of them were significantly correlated with IL6ST (Fig. [115]6C). The expression of circANKRD17, circFOXP1, circSCAF8, and circSLAIN2 was detected by qPCR in primary BMSCs from patients with AIS and non-AIS (Fig. [116]6D). We then used the Coding Potential Assessment Tool (CPAT) to predict whether these circRNAs had coding ability, using the XIST sequence as a negative control and the GAPDH sequence as a positive control. The results showed that circANKRD17 and circSCAF8 had no coding ability (Fig. [117]6E), and RNA FISH was used to detect the location of these two circRNAs in BMSCs. The results showed that circANKRD17 was located in the nucleus and that circSCAF8 was located in the cytoplasm (Fig. [118]6F). The regulation of IL6ST expression by the ceRNA mechanism mainly occurs in the cytoplasm, so we finally selected circSCAF8 as our target for further study. Fig. 6. CircSCAF8 is a potential upstream-regulating circRNA of IL6ST. [119]Fig. 6 [120]Open in a new tab A The heat map of differentially expressed circRNAs in patients with and without AIS. B ceRNA networks of hub genes obtained from Fig. [121]2F. C Correlation analysis of 6 selected circRNAs and IL6ST. D Relative level of selected circRNAs in primary BMSCs, n = 5. E CPAT was used to verify the coding ability of the candidate circRNAs. The XIST transcript served as a non-coding gene control, and GAPDH served as a coding gene control. F RNA FISH showed the location of these two circRNAs in BMSCs. CircRNA probes were labeled with 488. The nuclei were stained with DAPI. Scale bar: 10 μm, n = 6. Data are shown as the mean ± SD, **** means p < 0.001 vs. the control group. Knockdown of circSCAF8 significantly inhibited the osteogenic differentiation of BMSCs by regulating the level of IL6ST We found the sequence of circSCAF8 through the Circbase database and confirmed the head-tail splicing structure (Fig. [122]7A). To further verify the regulatory effect of circSCAF8 on the osteogenic differentiation of BMSCs, we first detected the expression location of circSCAF8 and IL6ST by RNA FISH, and the results showed that circSCAF8 and IL6ST were co-expressed in the cytoplasm of BMSCs (Fig. [123]S13A). We then knocked down circSCAF8 using shRNA transfection. We designed three sh-RNAs, and two of them significantly knocked down circSCAF8 (Fig. [124]7B). Moreover, there was no significant change in SCAF8 mRNA levels (Fig. [125]S13B). Then, CCK-8 assay was used to detect the effect of knocking down circSCAF8 on BMSC activity, and the results showed that the proliferation of BMSCs in the knockdown group was lower than that in the control group (Fig. [126]7C). Subsequently, immunofluorescence was used to detect the expression of IL6ST after knocking down circSCAF8, and the results showed that knocking down circSCAF8 significantly reduced the expression level of IL6ST (Fig. [127]7D, E). In addition, the expression of the bone-related proteins OCN and Runx2 in the circSCAF8 knockdown group was significantly lower, while overexpression of IL6ST could correct their low expression (Fig. [128]7F–H). Alizarin red staining also showed that knocking down circSCAF8 reduced the osteogenic differentiation of BMSCs (Fig. [129]7I, J). All these results suggested that knocking down circSCAF8 could inhibit the osteogenic differentiation of BMSCs, and this effect was achieved by regulating the level of IL6ST (Fig. [130]8). Fig. 7. Kd-circSCAF8 significantly inhibited the osteogenic differentiation of BMSCs by regulating IL6ST. [131]Fig. 7 [132]Open in a new tab A Schematic diagram of circSCAF8 loop formation, which is formed by the partial sequence of SCAF8 gene located on chromosome 6, and the looping site is the TGTA underlined in red. B The knock-down efficiency of circSCAF8 in BMSCs was detected by qRT-PCR, n = 4. C Cell viability was detected by CCK-8 analysis, n = 24. D, E Representative images of IL6ST immunofluorescence staining in BMSCs from the control and sh-circSCAF8 treatment groups. Scale bar, 20 μm, n = 4. F–H Relative protein level of OCN and RUNX2, n = 4. I, J Alizarin Red staining was performed to detect the osteogenic differentiation of human BMSCs in control, sh-circSCAF8, and sh-circSCAF8 + OE-IL6ST treatment groups on day 15. Scale bar, 100 μm, n = 6. Data are shown as the mean ± SD, * means p < 0.05 and ** means p < 0.01 vs. the control group, # means p < 0.05 and ## means p < 0.01 vs. the sh-circSCAF8 group. Fig. 8. [133]Fig. 8 [134]Open in a new tab Schematic diagram of IL6ST in the development of AIS. Discussion AIS is a disease characterized by changes in the three-dimensional structure of the spine, and the pathogenesis is still unclear. There are several hypotheses on the pathogenesis of AIS. 1. Gene mutation theory: Some studies have identified pathogenic genes such as LBX1 and GPR326 by whole-exome sequencing or transcriptome sequencing^[135]27. 2. Neuroendocrine theory: Adolescence is the peak period of body development, and abnormal secretion of multiple hormones was found in AIS^[136]28. 3. Theory of abnormal development of bone and cartilage: AIS shows abnormal endochondral and subperiosteal osteogenesis. Our previous study found osteopenia and abnormal cartilage development in AIS^[137]28. In a cross-sectional study by Lee et al., 619 patients with AIS and 300 sex- and age-matched patients with non-AIS were enrolled, and multivariate analysis revealed that the Cobb angle was inversely and independently correlated with BMD and bone mineral content (BMC) in patients with AIS. Moreover, decreased BMD aggravated spinal deformity progression in patients with AIS^[138]29. Cheng et al. recruited a cohort of 161 AIS and control female subjects to study their bone mass and circulating Wnt16 levels. They found that both BMD and circulating Wnt16 levels were significantly reduced in AIS^[139]30. In this study, we also calculated the BMD values of AIS and patients with non-AIS, and the results showed that the BMD value of patients with AIS was lower than that of patients with non-AIS. The accumulation of bone minerals is closely related to the osteogenic/adipogenic differentiation mediated by BMSCs. Several studies have explored the role of abnormal BMSC differentiation in the pathogenesis of AIS. Li et al. found that SPRY4 was significantly downregulated in AIS MSCs, and the knockdown of SPRY4 directly decreased the osteogenic differentiation of BMSCs by regulating the phosphorylation of ERK pathway proteins^[140]31. In addition, studies have shown that CGMP-dependent protein kinase (PRKG1) is upregulated in patients with AIS, and upregulated PRKG1 is related to the CGMP-PKG signaling pathway and has been identified as a novel negative regulator of osteoblast differentiation in BMSCs^[141]32. Zhuang et al. found upregulated SMAD3 expression in AIS BMSCs, and their study demonstrated that TGF-β-activated SMAD3 inhibited the expression of Runx2 and other Runx2 target genes, thereby inhibiting the osteogenic differentiation of BMSCs^[142]33. The p38 signaling pathway has been confirmed to be dysregulated in AIS. Some studies have found that the p38 pathway is activated during the process of BMSC differentiation into cartilage induced by transforming growth factor β1, and inhibition of the p38 pathway can inhibit chondrogenesis^[143]33. In this study, we observed that the osteogenic differentiation ability of AIS BMSCs was notably lower than that of the control group. Furthermore, the expression levels of key osteogenic-related proteins, namely Runx2 and OCN, were significantly reduced in patients with AIS. Motivated by these findings, we conducted further analysis to investigate the underlying mechanisms and identify key genes that may influence the osteogenic differentiation of BMSCs in patients with AIS. By high-throughput sequencing of the transcriptome of BMSCs, 1919 differentially expressed genes were identified. To select key genes, some researchers identify the gene with the largest change in many studies, but we believe that the largest change does not mean that it is the most critical gene. GO enrichment analysis was performed to annotate the function of each differentially expressed gene, and 34 differentially expressed genes were selected that were closely related to osteogenic differentiation. Multiple bioinformatics analyses, such as random forest analysis and PPI network analysis, identified IL6ST as one of the key genes involved in the osteogenesis of AIS BMSCs. IL6ST ( gp130) is a mediator of the action of the IL-6 cytokine family. Some studies have revealed the involvement of IL6ST in the pathogenesis of various diseases. Béziat et al. reported that dominant negative mutations in IL6ST can lead to the development of hyper-IgE syndrome^[144]15. Apart from the typical symptoms of recurrent infections, pneumonia, and elevated IgE levels, a significant portion of hyper-IgE syndrome patients also exhibit bone fractures and spinal curvature, indicating that loss of IL6ST function can affect bone development^[145]16. Moreover, several studies have shown that knockout of IL6ST expression in mice can cause skeletal muscle atrophy, reduced endothelial cell generation, tumor development, and alterations in energy metabolism^[146]17–[147]20. In the process of bone metabolism, IL6ST has been found to play a role in regulating osteoblast and osteoclast function, thus impacting bone remodeling and homeostasis. And its activation can stimulate the production and activity of downstream osteoblasts^[148]34. IL6ST is essential for normal bone development; IL6ST-null mice exhibit short, malformed bones, many osteoclasts destroy newly formed bone, and few osteoblasts are present^[149]35. However, its regulation of osteogenic differentiation may be achieved by regulating the RANKL/OPG pathway^[150]36. In Standal et al., the treatment of osteoporosis with parathyroid hormone (PTH) was accompanied by increased secretion of cytokines such as IL-6, and these researchers suggested that the stimulating effect of PTH on osteogenic differentiation and bone formation required the activation of IL6ST because in mice with targeted knockout of IL6ST expression in osteocytes, stimulation with PTH in the absence of IL6ST did not result in increased osteoblast numbers and bone mineral formation^[151]37. These findings suggest that IL6ST may be an important gene regulating osteogenic differentiation and that its role in osteopenia in AIS requires further investigation. We knocked down IL6ST expression in primary BMSCs, and the results showed that knockdown of IL6ST significantly inhibited the osteogenic differentiation of BMSCs, while overexpression of IL6ST significantly promoted the osteogenic differentiation of BMSCs in AIS. Moreover, knockout of IL6ST in zebrafish significantly increased the mortality and deformity rate of zebrafish. The development of vertebral body bone was also greatly delayed, suggesting that IL6ST knockout caused malformation by delaying vertebral body bone development. KEGG and GO enrichment analysis of DEGs suggested that the possible downstream pathway of IL6ST was the JAK/STAT pathway, and existing studies have found that the activation of the JAK/STAT3 pathway can significantly promote bone formation^[152]26. Both Itoh and Li et al. found that targeted knockdown of STAT3 protein in osteoblasts or osteocytes resulted in reduced bone mass in mice due to reduced bone formation, a significant reduction in osteoblasts in the trabecular bone, and a reduced response of bone formation to mechanical stimulation^[153]38,[154]39. These results indicate that the STAT3 protein is an important protein in promoting osteogenic differentiation and osteoblast activity. Of course, in the mouse model of STAT3 hyperactivation, osteoblasts are also increased, as well as osteoclasts, suggesting that osteoblast formation is coupled to osteoclast formation, and the IL-6-dependent pathway is critical for the control of bone remodeling^[155]38,[156]40. In this study, we found that the level of phosphorylated JAK/STAT3 proteins in AIS BMSCs was significantly lower than that in the control group, and the phosphorylation of the JAK/STAT3 pathway was significantly inhibited by knockdown of IL6ST expression in primary BMSCs. In addition, the STAT3 phosphorylation agonist colivelin corrected the decreased osteogenic differentiation caused by IL6ST knockdown. This finding suggests that the IL6ST/JAK2/STAT3 pathway may be a key pathway that regulates abnormal bone development in AIS. Furthermore, the role of IL6ST in chondrogenic differentiation of BMSCs in AIS patients is also worth discussing. Nakajima et al. reported that IL6ST inhibits chondrogenic differentiation in the ATDC5 chondrocyte cell line^[157]41. Abnormal chondrogenic differentiation can lead to morphological changes in intervertebral discs and endochondral ossification. An imaging analysis of AIS patients revealed that the position and morphological variations of intervertebral discs are correlated with the severity of spinal deformity^[158]42. Yang et al. performed single-cell analysis on patients with AIS and found that the differentiation function of MSCs and chondrocyte progenitor cell populations is impaired in AIS, resulting in impaired osteoblast and chondrocyte formation^[159]43. These studies all indicate that the role of IL6ST in chondrogenic differentiation of BMSCs in AIS patients is worth exploring. To further identify the upstream regulatory factors affecting IL6ST expression, we performed whole transcriptome sequencing encompassing both mRNA and non-coding RNA profiles. Due to the critical knowledge gap in circRNA-related mechanisms underlying AIS pathogenesis, we analyzed the differentially expressed circRNAs in BMSCs from AIS and patients with non-AIS. We hope this exploration will contribute to the understanding of AIS pathogenesis. In recent years, an increasing number of studies have found a close relationship between circRNAs and bone metabolism, and the understanding of circRNAs in the regulation of BMSCs, osteoblasts, and osteoclasts has increased^[160]44. Some studies have found that circRNAs regulate multiple proteins of the bone morphogenetic protein (BMP) family. Du et al. reported that circRNA-Fgfr2 could sponge miR-133 and regulate the expression of bone morphogenetic protein 6 (BMP6)^[161]45. Xu et al. reported that circUBXN7 inhibited cell proliferation and promoted cell apoptosis by regulating the IL6ST/JAK1/STAT3 axis through sponging miR-622 in LPS-induced cell injury^[162]46. In this study, we found 440 differentially expressed circRNAs. Using miRWalk and lncBase v.2 software, we predicted the upstream miRNAs and circRNAs of 7 hub genes and constructed a ceRNA network. Six potential upstream differentially expressed circRNAs of IL6ST were then obtained. CircSCAF8 was found to be the most likely upstream circRNA of IL6ST in AIS BMSCs. Then, the expression of circSCAF8 in primary BMSCs was knocked down by shRNA transfection, and the results showed that the expression of IL6ST was significantly decreased after knockdown, and the osteogenic differentiation of BMSCs was also significantly decreased. In this study, we identified IL6ST and its upstream circSCAF8 as key genes for abnormal bone development in AIS by high-throughput sequencing. Knockdown of IL6ST or circSCAF8 expression can inhibit the osteogenic differentiation of BMSCs. In zebrafish model, we found that deletion of IL6ST could significantly increase the malformation rate and mortality of zebrafish, and the vertebral body development of deformed zebrafish was significantly later than that of normal zebrafish. All these results suggest that IL6ST may be a key gene in the development of AIS, which provides a potential direction for the future precise diagnosis and treatment of AIS. Materials and methods Patients In total, we included 45 patients with AIS between the ages of 10 and 18 years and 20 sex- and age-matched patients without AIS (Table [163]1). In the AIS group, the minimum Cobb angle of the patients was greater than 45°, and scoliosis caused by other diseases, including congenital spinal deformity, neuromuscular disease, and spinal deformity caused by metabolic and endocrine diseases, was excluded by physical examination and X-ray examination. The control group consisted of healthy individuals and patients with lumbar spine herniation or spinal fracture in our physical examination center. All patients in the control group underwent a comprehensive X-ray examination to rule out the presence of spinal deformity. Next, we assessed the lumbar vertebrae bone density data of some patients using Dual Energy X-ray Absorptiometry (DEXA). The obtained BMD value is then compared to the mean BMD of a large sample of normal individuals of the same age and gender, and a standard deviation (SD) value is calculated. This value represents the age-appropriate Z-scores, and identified samples with Z-score scores < −1 as reduced bone mass. All ethical regulations relevant to human research participants were followed. Table 1. Data of enrolled patients Items AIS control P-value Number 45 20 Age 14.2 ± 2.8 17.3 ± 3.7 0.273 Gender 18 (boy) vs 27 (girl) 12 (boy) vs 8 (girl) Bone mineral density(BMD) 0.74 ± 0.21 g/m^2 0.88 ± 0.15 g/m^2 <0.05 Z-score −1.13 ± 1.27 −0.34 ± 0.95 <0.05 Cobb angle 62.5 ± 15 [164]Open in a new tab Student’s t-test was use to tested, when p-value < 0.05, means significant difference in the two groups. In vivo experiments The Institutional Animal Care and Use Committee (IACUC) of Xiangya Hospital, Central South University, reviewed and approved all in vivo procedures and protocols. We have complied with all relevant ethical regulations for animal use. Zebrafish were bred and maintained in a flow-through system at a constant temperature of 28 °C under a 14/10 h light/dark cycle. Zebrafish embryos were obtained following standard fertilization protocols. At 1 day post-fertilization (dpf), the embryos were randomly divided into two groups: a saline-injected control group and a Cas9/RNP complexes-injected treatment group. Depending on the amount of ovulation, the number of fry in each group ranges between 50 and 100. We recorded the survival rate, hatching rate, and malformation rate of the zebrafish larvae. At 7 dpf and 30 dpf, microscopic analysis was performed to document the malformation phenotypes in both groups. Additionally, alizarin red staining was used to assess skeletal development in the zebrafish. Imaging of zebrafish and embryo samples from each treatment group was conducted using a Motic SM7 microscope (Motic, Xiamen, China). Antibodies and reagents The antibodies used in this study were IL6ST (1:200/1:1000; Abcam, Cambridge, UK), β-actin (1:1000; Cell Signaling Technology (CST), Danvers, MA, USA), RUNX2 (1:1000; CST), and OCN (1:800; CST). JAK1, p-JAK1 (1:1000; CST); STAT3, p-STAT3 (1:1000; CST), and secondary antibodies (1:10,000, China Institute of Biotechnology, Beijing, China). Alizarin red and oil red staining kits were purchased from Cyagen (HUXMA-90021; China). An OriCell BM-MSC osteogenic and adipogenic differentiation kit was purchased from Cyagen (HUXMA-90021; China). Cell Counting Kit-8 (CCK-8) was purchased from Takatsu Do (Kumamoto, Japan). Extraction of human BMSCs Screw path blood was drawn from patients undergoing surgery and centrifuged at 1000 rcf for 5 min, and the precipitated cells were treated with erythrocyte lysates. Lysed cells were resuscitated and seeded into 6-cm cell culture dishes. The cells were cultured in an incubator for 24 h, the floating cells were washed off, and the medium was changed again after 48 h. Finally, all BMSCs were cultured with 90% Ham’s F12 nutrient medium (Gibco, MA, USA) + 10% fetal bovine serum (Gibco, MA, USA) and 1% penicillin/streptomycin (Gibco). All cells were passaged twice before use. RNA extraction and sequencing Total RNA was extracted from human BMSCs using a TRIzol kit (Invitrogen, Thermo Fisher Scientific) following the manufacturer’s protocol. The quality of the RNA was assessed using an Agilent 2100 Bioanalyzer (Agilent Technologies, Palo Alto, CA, USA) and RNase-free agarose gel electrophoresis. After the extraction of total RNA, ribosomal RNA was removed to retain mRNA and non-coding (nc) RNA. Library construction and sequencing were performed on an Illumina Nova 6000 by Gene Denovo Biotechnology Co. (Guangzhou, China), generating over 12 G reads per sample. To ensure data quality, reads were further filtered using fastp^[165]47(version 0.18.0). Short reads alignment tool Bowtie2^[166]48(version 2.2.8) was used for mapping reads to ribosome RNA (rRNA) database. The rRNA-mapped reads were then removed. The remaining reads were mapped to the reference genome using HISAT2^[167]49(version 2.1.0). RNA sequencing data analysis The FPKM values for each transcript were calculated using RSEM software to analyze the RNA sequencing data. Differentially expressed lncRNAs between the AIS and control (ct) groups were identified using the deseq2 package, considering fold change ≥2 and false discovery rate (FDR, adjusted p-value) < 0.05 as statistically significant criteria. Additionally, sample gene expression correlation analysis was performed to assess the reproducibility of samples within groups. Volcano plots, heat maps, and radar plots were generated using ggplot2 and heat maps to visualize differential RNA expression profiles. Function enrichment analysis Function enrichment analyses, including Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG), were conducted to predict potential signal transduction pathways, biological functions of mRNAs in the ceRNA network, as well as biochemical metabolic pathways. GO analysis encompassed three ontologies: biological process (BP), cellular component (CC), and molecular function (MF). This allowed for mRNA annotation based on GO functional classification and analysis of the functional significance enrichment. Pathway-based analysis further facilitated understanding genes’ biological functions. The KEGG database was the primary public pathway-related resource used for pathway enrichment analysis, aiming to identify significantly enriched metabolic or signal transduction pathways in ceRNAs compared to the whole genome background. Significantly enriched GO terms were determined through a hypergeometric test with a threshold corrected p-value ≤ 0.05, while significantly differentially enriched pathways in the KEGG pathway analysis had Q values ≤ 0.05. PPI network Protein‒protein interaction (PPI) networks of DEGs were constructed from the Interactive Gene Retrieval Search Tool (STRING, [168]https://www.string-db.org/) online database. In this study, genes were identified as nodes in the network, interactions were identified as lines of the network by using STRING v10 identification, and network files were visualized using Cytoscape (v3.7.1) software to present biological interactions of core and hub genes. CircRNA-miRNA-mRNA regulatory network analysis We used the circular RNA interaction groups (CircInteractome) database ([169]https://circinteractome.nia.nih.gov/) to predict the interaction between microRNAs and DEcircRNAs. Prediction of microRNAs was performed with the online tool Venn 2.1.0 ([170]https://bioinfogp.cnb.csic.es/tools/venny/index.html) through the Venn diagram and miRNA intersection. Then, the target mRNAs of overlapping miRNAs were predicted using the miRDB database ([171]https://mirdb.org/), and the intersection of predicted mRNAs and DEmRNAs was obtained using Venny 2.1.0. Finally, the ceRNA network was constructed based on the expression of circRNAs, miRNAs, and mRNAs, and the results were visualized using Cytoscape software. Analysis of coding potential CircRNA coding potential can be analyzed by the Coding Potential Assessment Tool (CPAT) ([172]http://lilab.research.bcm.edu/cpat/). The XIST transcript was used as a negative control for the non-coding gene. GAPDH was used as a positive control for the coding gene. Real-time quantitative PCR The experimental methods used for qRT‒PCR were as described previously^[173]50. Table [174]2 lists the primers used in this study. Table 2. Primers for quantitative real-time PCR (qRT-PCR) Gene Primer sequences (5′ to 3′) IL6ST R TGTTGACGTTGCAGACTTGG F AATTGTGCCTTGGAGGAGTG JUND R TGCATTTGGATGTACAGATCG F CCCCCTCCCAAAAGAAGTAT BMP6 R GAAGGGCTGCTTGTCGTAAG F AAGAAGGCTGGCTGGAATTT SMO R TCGAGTTCTTCCATCCTGCT F TGTAACAGGTCCTCGGAAGG COL6A1 R TGGCACAAAGCGACTGGATGAAC F GTGGAGGCCGACTTCTTGTATGC ITGA11 R CGATGGGAACATTCAGGGCAGAG F ACACATTCATGGAGGCACTGGAAC IGF1 R GGCTCGCAGTAGGTAACTGG F CCCCAGAGCAATATTCCAGA SCAF8 R TGACGGCTGCACTAAAGATG F TTAGTCCCACCAGCATTTCC CircSCAF8 R CCCTGGGGATGACAAGAGTA F CCCTGGGGATGACAAGAGTA RN18s R AGAAACGGCTACCACATCCA F CCCTCCAATGGATCCTCGTT [175]Open in a new tab Immunochemistry and immunofluorescence Cells were plated on slides, fixed with 4% paraformaldehyde for 20 min, and incubated with 0.3% Triton solution for 5 min. Cells were then washed with phosphate-buffered saline, sealed with 5% bovine serum albumin for 30 min, and incubated with primary antibody overnight at 4 °C. The next day, the cells were incubated with fluorescent secondary antibodies for 1 h at room temperature. Images were acquired using a confocal microscope (Carl Zeiss AG, Jena, Germany), and the intensity of staining for target proteins was quantified using ImageJ. Alizarin red and oil red staining The cells were seeded into a 6-well plate, and differentiation was induced when the cell density reached approximately 80%. We used Human BMSCs Osteogenic/Adipogenic Induction Differentiation Kit to treat the cells according to the provided protocol for 14 days. Subsequently, cells were stained with alizarin red and oil red using staining kits according to the manufacturer’s protocol. Images were obtained using a microscope (Nikon), and the color intensity of mineral deposits was examined using ImageJ software. Cell proliferation assay BMSCs were seeded in 96-well plates at a density of 3000 cells per well according to the manufacturer’s instructions. Cell viability was measured using the CCK-8 system (Dojindo) at 0, 12, 24 and 48 h after inoculation. Absorbance was measured using a microplate reader (Tecan, Mannedorf, Switzerland). Lentivirus transfection Cells were seeded in 6-well plates at a density of 80%. After the cells attached, the culture medium was replaced, and the cells were transfected with lentiviruses at a multiplicity of infection (MOI) of 40. After 12–16 h of infection, the medium was aspirated and replaced with fresh medium without the virus. Lentiviral vectors were purchased from Obio Technology (Shanghai, China). The effective shRNA interference sequences used in this experiment were: sh-IL6ST F1: “CCGGTAAGGGATACTGG AGTGACTGCTCGAGCAGTCACTCCAGTATCCCTTATTTTTG”; R1: “AATTCA AAAATAAGGGATACTGGAGTGACTGCTCGAGCAGTCACTCCAGTATCCCTTA”. CRISPR‒Cas9 gene knockout Single guide RNA (sgRNA) for each locus was synthesized using the universal TRACR RNA template oligomer. Two gene-specific oligonucleotides were then synthesized at each site using T7 RNA polymerase (Beijing Qingke Biotechnology Co., Ltd., China). CRISPR/Cas9 ribonucleoprotein (Cas9/RNP) complexes were generated by incubating a 40 μM sgRNA mixture with 5 μM Cas9 protein (FishBio, Shanghai, China) for 5 min at room temperature prior to injection. We injected 200 to 500 pl of Cas9/RNP complexes into the yolk of zebrafish embryos at the single-cell stage. All embryos for each independent injection were obtained from mass mating of AB wild-type fish to ensure genetic diversity. After that, embryos were cultured in egg water containing methylene blue (60 μg/ml of sea salt in distilled water) at 28 °C. Statistics and reproducibility All statistical analyses were performed using GraphPad Prism 7 software (GraphPad software, San Diego, CA, USA). Data are presented as the mean ± S.D., and significant differences between groups were analyzed by unpaired Student’s t-tests. ANOVA was used to determine the significance of differences between multiple groups. P < 0.05 was considered statistically significant. Reporting summary Further information on research design is available in the [176]Nature Portfolio Reporting Summary linked to this article. Supplementary information [177]Supplementary Information^ (1.9MB, pdf) [178]42003_2025_8746_MOESM2_ESM.pdf^ (97.7KB, pdf) Description of Additional Supplementary Files [179]Supplementary Data 1^ (402.6KB, xlsx) [180]Supplementary Data 2^ (764.1KB, xlsx) [181]Supplementary Data 3^ (11.4KB, xlsx) [182]Reporting Summary^ (1.8MB, pdf) Acknowledgements