Abstract
Aims
Lumbar spinal stenosis (LSS) is a common skeletal system disease that has been partly attributed to genetic variation. However, the correlation between genetic variation and pathological changes in LSS is insufficient, and it is difficult to provide a reference for the early diagnosis and treatment of the disease.
Methods
We conducted a transcriptome-wide association study (TWAS) of spinal canal stenosis by integrating genome-wide association study summary statistics (including 661 cases and 178,065 controls) derived from Biobank Japan, and pre-computed gene expression weights of skeletal muscle and whole blood implemented in FUSION software. To verify the TWAS results, the candidate genes were furthered compared with messenger RNA (mRNA) expression profiles of LSS to screen for common genes. Finally, Metascape software was used to perform enrichment analysis of the candidate genes and common genes.
Results
TWAS identified 295 genes with permutation p-values < 0.05 for skeletal muscle and 79 genes associated for the whole blood, such as RCHY1 (PTWAS = 0.001). Those genes were enriched in 112 gene ontology (GO) terms and five Kyoto Encyclopedia of Genes and Genomes pathways, such as ‘chemical carcinogenesis - reactive oxygen species’ (LogP value = −2.139). Further comparing the TWAS significant genes with the differentially expressed genes identified by mRNA expression profiles of LSS found 18 overlapped genes, such as interleukin 15 receptor subunit alpha (IL15RA) (PTWAS = 0.040, PmRNA = 0.010). Moreover, 71 common GO terms were detected for the enrichment results of TWAS and mRNA expression profiles, such as negative regulation of cell differentiation (LogP value = −2.811).
Conclusion
This study revealed the genetic mechanism behind the pathological changes in LSS, and may provide novel insights for the early diagnosis and intervention of LSS.
Cite this article: Bone Joint Res 2023;12(6):387–396.
Article focus
-
To reveal candidate genes for lumbar spinal stenosis (LSS), and confirm the correlation between genetic variation and pathological changes in LSS.
Key messages
-
Using transcriptome-wide association study (TWAS) analysis, this study identified 374 novel susceptibility genes associated with LSS.
-
To verify the results of TWAS, the candidate genes were further compared with messenger RNA (mRNA) expression profiles of LSS and performed enrichment analysis.
-
This study found the correlation between genetic variation and the three major pathological changes of LSS.
Strengths and limitations
-
This is the first time the correlation between genetic variation and pathological changes of LSS has been confirmed via TWAS analysis.
-
TWAS analysis is a creative method that can predict gene expression in LSS, and avoid confounding from environmental differences caused by the trait that may influence expression.
-
The GWAS summary dataset and mRNA expression profiles are derived from subjects of Asian ancestry; the results of this study should be applied to other populations with caution.
Introduction
Lumbar spinal stenosis (LSS) is characterized by spinal nerve entrapment caused by bony spinal canal stenosis, and can be divided into developmental LSS and degenerative LSS according to the aetiology. The former is associated with pre-existing or congenital narrowing of the spinal canal and accounts for approximately 5% of LSS cases.1 However, most cases of LSS are attributed to the degenerative changes in the spine with ageing. The pathoanatomical phenotype of degenerative LSS is related to the spinal canal size or neural compressive element.1 The main symptoms of LSS are sciatica, numbness or weakness in buttocks or legs, and intermittent claudication or sphincter dysfunction in severe cases.2 The substantial pain and disability of LSS, which impairs ambulation, affects more than 103 million people worldwide and approximately 600,000 surgical procedures are performed in the USA each year for LSS.3
Due to the human need for upright walking, the intervertebral disc becomes the initiating factor of spinal degeneration, which further leads to the pathological changes in LSS that develop around the level of the intervertebral disc, such as buckling of the ligamentum flavum, osteophyte hypertrophy of the facet joint, and even lumbar spondylolisthesis. The severity of LSS is highly heritable, with an estimated heritability of 67% based on qualitative MRI assessment.4 A recent meta-analysis also found that genetic variation is associated with susceptibility to disc degeneration.5 The genetic heterogeneity of LSS is established at the time of parental gamete formation, independent of environmental interference in later life, and predates the appearance of LSS. It also explains the clinical phenomenon of why the degree of LSS in some patients is not always consistent with the severity of symptoms at the genetic level. Thus, the study of genetics is an ideal epidemiological tool for screening high-risk populations. It also provides a new approach for the development of biological therapies to reduce the incidence of LSS and delay the degree of spinal degeneration. Since the pathogenesis of LSS is unclear, the study of genetics may be a good tool to explore the mechanisms underlying the pathogenesis of LSS.
Previous genome-wide association studies (GWASs) have scanned a large number of genetic markers across the genome to locate the genetic variation associated with LSS.6,7 The gene AAK1 and the variants in SOX5 and CCDC26/GSDMC were correlated with the pain severity of LSS, revealing the specific genetic mechanism of individual patients’ responses to LSS.6 However, most of the genetic variants reported by GWASs are located in non-coding regions and have limited capacity for explanation at the level of gene expression.8 Compared with single nucleotide polymorphism (SNP)-based GWASs studies, transcriptome-wide association studies (TWASs) can take into account regulatory genes in non-coding regions, and better determine gene-trait associations by integrating GWASs and gene expression datasets.8 In addition, TWAS can drastically reduce the comparisons in statistical analysis and enhance the ability to detect the candidate genes of diseases with complex traits.9 In recent years, TWAS has been widely used to identify risk genes in a variety of orthopaedic diseases. For example, Qi et al10 identified 33 genes for hip osteoarthritis (OA) and 24 genes for knee OA, which provides novel clues for understanding the genetic mechanism of OA.
In the current study, by integrating the GWAS summary dataset of spinal canal stenosis derived from BioBank Japan, and precomputed gene expression weights from skeletal muscle and whole blood, we conducted a TWAS analysis to identify candidate genes associated with LSS. The significant genes were further validated by the mRNA expression profiles of LSS patients. Gene ontology (GO) enrichment analysis was performed for gene functional annotation. The findings provide novel insights into the early diagnosis of and intervention for LSS by identifying genetic variations associated with pathological changes.
Methods
GWAS summary datasets of spinal canal stenosis
A large-scale GWAS summary dataset was obtained from BioBank Japan, containing 661 diagnosed spinal canal stenosis and 178,065 controls of Asian ethnicity.11 BioBank Japan is a prospective biobank that collaboratively collected DNA and serum samples from 12 medical institutions in Japan and recruited approximately 200,000 participants, mainly of Japanese ancestry (not replicated in the UK BioBank and FinnGen). The mean age of participants at recruitment was 63 years, and 46.3% were female.11 All cases were defined by code M48 in the International Classification of Diseases, Tenth Revision.12 Genotyping was conducted using commercial platforms, such as Illumina HumanOmniExpressExome BeadChip, or a combination of the Illumina HumanOmniExpress and HumanExome BeadChip (Illumina, USA). After imputation, this study excluded variants with an imputation quality of Rsq < 0.7, resulting in 13,530,797 variants analyzed in total.11 Significant associations were considered < 5.0 × 10−8; a total of 75,230 associations across 181 traits were included in the study.11 Detailed description of sample characteristics, experimental design, and statistical analysis can be found in the published study.11
Statistical analysis
The tissue-specific TWAS analysis (including skeletal muscle and whole blood) was carried out by using Functional Summary-based Imputation software (FUSION; Gusev Lab, USA) by integrating the GWAS summary dataset and precomputed gene expression weights of different tissues. FUSION is the most commonly used TWAS method and has shown great promise in large-scale integrative omics data analysis.13,14 The imputed gene expression weights data can be viewed as a linear model based on the correlation for linkage disequilibrium (LD) among SNPs.15 Given that the skeletal muscle and whole blood were used in previous studies of LSS,16-18 the FUSION pre-computed functional weights of gene expression of the skeletal muscle and whole blood were used in our TWAS analysis. The expression weights were first calculated under different prediction models implemented in the FUSION, such as BLUP, BSLMM, LASSO, Elastic Net, and top SNPs, from the reference expression panel (GTExv8). The association testing statistics between predicted gene expression and each taxon was calculated as TWAS Z = w’Z/(w’Lw)1/2. “Z” denotes the scores of LSS, and “w”denotes the weights. “L” denotes the SNP-correlation LD matrix.15 In this study, the permutation test was run 2,000 times for each gene by Z-test within skeletal muscle and whole blood. The significant genes were defined as permutation p-value < 0.05.
The mRNA expression profiles of LSS were extracted from Gene Expression Omnibus (GEO) Datasets (GSE113212).19 Briefly, the dataset included five patients (one male, four females; age: 66 to 79 years), and the control specimen was the spine of a fresh cadaveric 30-year-old male.20 RNA sequencing was conducted on the platform of Agilent-039494 SurePrint G3 Human GE v2 8 × 60 K Microarray (Agilent, USA). The differentially expressed genes (DEGs) were screened out by using the online GEO2R tool implemented in the GEO website. The GEO2R tool provides a simple interface that allows users to perform complex R-based GEO data analysis to help identify and visualize differential gene expression.21 A p-value < 0.05 and log-fold change (| logFC |) > 1 were the criteria for differential expression. We compared the candidate genes identified by TWAS with the statistically significant DEGs detected by mRNA expression profiles of LSS to identify common candidate genes for LSS.
Metascape software (USA) was used for Gene Ontology (GO) and pathway enrichment analyses of the genes identified by TWAS, while DEGs were scanned by mRNA expression profiles.22 The enrichment results of TWAS and the mRNA expression profiles were further compared to determine the common GO terms and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways. The top GO terms in the result were selected and plotted using GO plots.
Results
TWAS analysis performed by FUSION
A total of 374 genes were identified by TWAS, including 295 for skeletal muscle, such as RCHY1 (p = 0.001), FBX025 (p = 0.001), IZUMO1 (p = 0.001), and ARNTL (p = 0.001), and 79 for whole blood, such as CRLF3 (p = 0.001), TNP03 (p = 0.001), and C2CD2l (p = 0.003) (Figure 1). Table I shows the details of the top 20 significant genes identified by TWAS. By comparing with pre-computed gene expression weights, 15 genes were expressed in both skeletal muscle and whole blood, such as FLYWCH1 (p = 0.001), DUTP6 (p = 0.004), and MBLAC1 (p = 0.005).
Fig. 1
Manhattan plot showing transcription-wide association study (TWAS)-identified genes and significantly expressed genes associated with lumbar spinal stenosis (LSS; annotated points). Each point represents a single gene, and the physical position (chromosome localization) is plotted on the x-axis, while the -log10 (p-value) of the association between gene and LSS is plotted on the y-axis.
Table I.
Top 20 genes identified by transcriptome-wide association study analysis.
| Tissue | Gene | CHR | Best GWAS ID | NSNP | TWAS.Z | P TWAS |
|---|---|---|---|---|---|---|
| Skeletal muscle | RCHY1 | 4 | rs6535523 | 427 | 3.8860 | 0.001 |
| FBXO25 | 8 | rs7822166 | 405 | −3.4367 | 0.001 | |
| IZUMO1 | 19 | rs838143 | 429 | 3.3869 | 0.001 | |
| ARNTL | 11 | rs10766064 | 550 | 3.2647 | 0.001 | |
| CNTF | 11 | rs3829218 | 412 | 3.2040 | 0.001 | |
| ATAD5 | 17 | rs216462 | 279 | −3.1973 | 0.001 | |
| FLYWCH1 | 16 | rs9930148 | 340 | 3.1917 | 0.001 | |
| PI4KA | 22 | rs4820579 | 432 | −3.1901 | 0.001 | |
| MYL12B | 18 | rs8086134 | 486 | 3.1740 | 0.002 | |
| FAM221A | 7 | rs1559012 | 439 | 3.1550 | 0.002 | |
| Whole blood | FLYWCH1 | 16 | rs9930148 | 335 | 3.3599 | 0.001 |
| CRLF3 | 17 | rs216462 | 294 | 3.3373 | 0.001 | |
| TNPO3 | 7 | rs17340646 | 359 | 3.2003 | 0.001 | |
| C2CD2L | 11 | rs643788 | 367 | −3.0000 | 0.003 | |
| MAEA | 4 | rs13117476 | 374 | 2.9289 | 0.003 | |
| DUTP6 | 1 | rs12568568 | 574 | −2.9020 | 0.004 | |
| SNAP29 | 22 | rs4820579 | 416 | 2.8489 | 0.004 | |
| MBLAC1 | 7 | rs4729565 | 265 | 2.7910 | 0.005 | |
| CNIH4 | 1 | rs6675858 | 360 | 2.7908 | 0.005 | |
| ATP6V0E2 | 7 | rs7795944 | 256 | −2.7070 | 0.007 |
-
TWAS.P and TWAS.Z values were calculated by the FUSION approach.
-
CHR, chromosome; GWAS, genome-wide association study; NSNP, number of single-nucleotide polymorphisms in the locus; TWAS, transcriptome-wide association study.
Validating the TWAS results by mRNA expression profiles of LSS
The mRNA expression profiles of LSS screened 1,702 DEGs, among which 599 were downregulated and 1,103 were upregulated. After comparing the results of the TWAS and mRNA expression profile analyses, 18 common genes were identified and are shown in Table II, such as IL15RA (PTWAS = 0.040, PDEG = 0.010; Bonferroni correction). The distribution of DEGs identified from mRNA expression profiles was visualized in the corresponding volcano plot (Figure 2).
Table II.
The top common genes identified by both transcription-wide association study analysis and mRNA expression profiles of lumbar spinal stenosis.
| Tissue | Gene | CHR | NSNP | PTWAS | PDEG | Regulation |
|---|---|---|---|---|---|---|
| Skeletal muscle | ATAD5 | 17 | 279 | 0.001 | 0.048 | Up |
| CRTAM | 11 | 516 | 0.005 | 0.018 | Down | |
| KDELR3 | 22 | 377 | 0.012 | 0.016 | Down | |
| CIDEA | 18 | 433 | 0.019 | 0.025 | Up | |
| TFCP2L1 | 2 | 384 | 0.020 | 0.005 | Down | |
| SLC16A10 | 6 | 322 | 0.024 | 0.004 | Down | |
| RIMS1 | 6 | 427 | 0.025 | 0.002 | Down | |
| NOXRED1 | 14 | 575 | 0.026 | 0.018 | Up | |
| AKAP6 | 14 | 567 | 0.028 | 0.028 | UP | |
| SLC25A10 | 17 | 328 | 0.034 | 0.046 | Up | |
| CRLF1 | 19 | 374 | 0.034 | 0.007 | Down | |
| FAM43B | 1 | 509 | 0.043 | 0.014 | Down | |
| BCR | 22 | 396 | 0.045 | 0.001 | Up | |
| ZNF331 | 19 | 482 | 0.048 | 0.043 | Up | |
| TTC23 | 15 | 545 | 0.050 | 0.001 | Up | |
| Whole blood | SPNS1 | 16 | 173 | 0.011 | 0.004 | Down |
| IL15RA | 10 | 715 | 0.040 | 0.010 | Down | |
| CAT | 11 | 605 | 0.045 | 0.003 | Up |
-
P-values were tested using Bonferroni correction.
-
DEG, differentially expressed gene; LSS, lumbar spinal stenosis; TWAS, transcriptome-wide association study.
Fig. 2
The volcano plot of messenger RNA (mRNA) expression profiles for lumbar spinal stenosis (LSS). Genes were marked in red point as differentially expressed when the following two conditions were met: p-value of < 0.05 by the moderated t-statistic, and log-fold change (|logFC|) > 1.
Gene set enrichment analysis of the TWAS results
The significant genes identified by TWAS and mRNA expression profiles of LSS were submitted to the Metascape online tool for GO and KEGG enrichment analysis, respectively. We identified 112 GO terms and five KEGG pathways enriched for TWAS results, such as negative regulation of cellular component organization (GO:0051129, LogP value = −3.249), ‘chemical carcinogenesis - reactive oxygen species’ (hsa05208, LogP value = −2.139). The significant terms were then hierarchically clustered, a subset of representative terms was selected, and these were converted into a network layout (Figure 3). The GO Chord plot shows the top GO terms enriched by TWAS results (Figure 4), and the correlations between genetic variation and pathological changes in LSS are shown in Figure 5. By integrating the DEG enrichment analysis, 71 common GO terms were identified, such as negative regulation of cell differentiation (GO:0045596). The top overlapped GO terms are listed in Table III.
Fig. 3
The network layout of top Gene Ontology (GO) terms under hierarchical clustering. In the network, each circle node represents a term, where its size is proportional to the number of input genes that fall into that term, and its colour represents its cluster identity (i.e. nodes of the same colour belong to the same cluster). Terms with a similarity score > 0.3 are linked by an edge (the thickness of the edge represents the similarity score). NOTCH4, neurogenic locus notch homolog 4.
Fig. 4
Gene Ontology (GO) Chord plot of top GO terms belonging to the biological process subontology for lumbar spinal stenosis (LSS). The genes are linked to their assigned terms via coloured ribbons. Genes are ordered according to the observed log fold change (LogFC), which is displayed in descending and ascending intensity of red and blue squares displayed next to the selected genes. TWAS, transcription-wide association study.
Fig. 5
The correlation between genetic variation and the three major pathological changes of lumbar spinal stenosis. GO, gene ontology; IL15-RA, interleukin-15 receptor α; IVDD, invertebral disc degeneration.
Table III.
Top overlapped Gene Ontology terms identified by both transcription-wide association study analysis and messenger RNA expression profiles.
| GO ID | Description | Log p-value TWAS | Log p-value mRNA expression |
|---|---|---|---|
| GO:0120035 | Regulation of plasma membrane bounded cell projection organization | −5.287 | −9.296 |
| GO:0031344 | Regulation of cell projection organization | −5.115 | −9.154 |
| GO:0010975 | Regulation of neuron projection development | −5.098 | −8.009 |
| GO:0031346 | Positive regulation of cell projection organization | −3.167 | −7.019 |
| GO:0010976 | Positive regulation of neuron projection development | −3.100 | −2.489 |
| GO:0010977 | Negative regulation of neuron projection development | −4.160 | −2.948 |
| GO:0051129 | Negative regulation of cellular component organization | −3.249 | −3.987 |
| GO:0031345 | Negative regulation of cell projection organization | −3.203 | −2.859 |
| GO:0050770 | Regulation of axonogenesis | −3.033 | −6.213 |
| GO:0045596 | Negative regulation of cell differentiation | −2.811 | −5.318 |
-
GO, Gene Ontology; mRNA, messenger RNA; TWAS, transcriptome-wide association study.
Discussion
LSS is characterized by lumbar intervertebral disc degeneration, ligamentum flavum hypertrophy, and osteophyte formation. Previous studies have identified several genome-wide significant loci associated with LSS, but the ability to interpret the relationship between the significant genes and LSS is limited.6,7 In this study, we conducted a TWAS analysis and identified genes closely associated with LSS in order to explain the relationships between these genes and the disease. To validate the results of TWAS, these genes were further compared with RNA expression profiles of LSS to screen for common genes. Our study results provide novel clues for understanding the genetic mechanism of LSS, focusing on the roles of abnormal gene transcription in the development of LSS.
The pathological evolution of all diseases includes a certain genetic variation basis, but only conducting the analysis of genetics – ignoring its correlation with pathological changes – is not conducive to early diagnosis and treatment. Starting from the pathological changes, exploring the inherent genetic variation basis for diseases can provide new diagnostic and treatment ideas. Previous GWASs did not carry out correlation analysis between LSS genetic variation and pathological changes, resulting in poor clinical application.6,7
By combining TWAS analysis and LSS mRNA expression profiling, interleukin-15 receptor α (IL-15RA) was identified on both gene and mRNA level. IL-15RA has a genetic correlation with osteophyte proliferation at the posterior margin of the vertebral body and facet joints, an important pathological change in LSS. IL-15RA is an important component of the IL-15 proinflammatory signalling that occurs in muscles after exercise.23,24 In the skeletal system, IL-15RA plays a cell-autonomous role in osteoblast function and bone mineralization.25IL-15RA SNPs are associated with susceptibility to ossification of the posterior longitudinal ligament of the spine.26 Over time, pathological ectopic bone formation in the spinal canal, including osteophytes with hyperplasia of the posterior margin of the vertebral body or facet joints, and ectopic ossification of the posterior longitudinal ligament, results in spinal stenosis of the lumbar spinal canal and progression to compressive myelopathy or radiculopathy. The hypertrophy of the ligamentum flavum is considered to be another major contributor to the development of LSS. Fibrosis is considered to be a main aetiology of the hypertrophy of the ligamentum flavum, which shows loss of elastic fibres and synthesis of collagen fibres.27,28 The gene RCHY1, an E3 ligase screened out in this study, regulates ubiquitin-mediated protein degradation, and the candidate gene FLYWCH1 can encode zinc finger protein. The aforementioned two genes are also associated with DNA damage repair, and jointly regulate the metabolism of related proteins.29,30 Yabe et al27 showed that the synthesis of proteoglycans related to protein metabolism is a key factor in the pathogenesis of ligamentum flavum hypertrophy. Thus, the genes RCHY1 and FLYWCH1 may affect the synthesis and metabolic balance of collagen fibres and elastic fibres through the regulation of DNA levels.
Intervertebral disc degeneration leads to initial relative instability and hypermobility of the intervertebral joints, which then spread to the facet joints and adjacent ligaments, subsequently resulting in hyperplasia of bony structures and ligamental hypertrophy. The gene screened in this study mainly codes for expressed Spinster Homolog 1 protein, which is a hypothetical lysosomal H+-carbohydrate transporter contributing to macroautophagy.31 He et al32 also found that hypoxia may play a protective role by enhancing macroautophagy/autophagy and retarding intervertebral disc degeneration in nucleus pulposus-derived stem cells.
The Metascape tool was used for gene functional enrichment and annotation analysis, and the data were compared with the mRNA expression profiles to screen out candidate GO terms related to LSS. Intervertebral disc degeneration (IVDD) caused by nucleus pulposus (NP) degeneration can be followed by intervertebral collapse, osteophyte hypertrophy, segmental instability, and other pathological changes, while bony spinal canal stenosis can lead to spinal cord and nerve root impingement and corresponding clinical symptoms.32 NP cells are similar to chondrocytes, which consist of collagen fibres embedded in gel-like aggrecan. The regulation of cell differentiation occurs in various pathophysiological pathways associated with LSS, especially the dedifferentiation of NP cells into fibroblast-like phenotypes during IVDD.32 Negative regulation of cell differentiation (GO:0045596) was identified by TWAS analysis and mRNA expression profile enrichment analysis in this study. In addition to cellular differentiation, IVDD is also associated with cellular component abnormalities. An increasing number of senescent cells can be observed with ageing and degeneration of intervertebral disc tissue. The time-dependent accumulation of molecular and cellular damage is the main cause of abnormal homeostasis and eventual structural and functional breakdown of intervertebral discs.33 This study found that the GO biological process negative regulation of cellular component organization (GO:0051129) is associated with LSS through negative regulation of cellular components and tissues. Senescent cells promote the functional decline in intervertebral disc tissue, which may be positively correlated with the increase in senescent markers, such as the cell cycle arrest protein p16INK4a, or the secretion of proinflammatory cytokines such as TNFα, IL-1β, IL-17, IL-6, and COX-2.34
The novelty of this study is that the correlation between genetic variation and the three major pathological changes of LSS was confirmed by TWAS for the first time. TWAS analysis is a creative method that can predict the gene expression in LSS and avoid confounding from environmental differences caused by the trait that may influence expression. In addition, the large sample size of GWAS summary dataset ensures the accuracy of our results. However, this study also has some limitations. The GWAS summary dataset and mRNA expression profiles are derived from subjects of Asian ancestry, so the findings should be applied to other populations with caution. Further TWAS analysis on other populations are needed to prove our results. In addition, the control group mRNA expression profiles exclusively comprised men; the sex differences in mRNA expression profiles may have influenced the results.
In summary, this study firstly performed TWAS to identify novel susceptibility genes associated with LSS, then integrated TWAS with mRNA expression profiles to identify common genes for LSS. We hope that the discovery of the genetic mechanism underlying LSS can be used to screen high-risk groups to achieve precise regulation of LSS through targeted biological intervention at the gene level through the development of new drugs, so as to delay degeneration of the lumbar spine, lessen clinical symptoms and disability, and reduce social and economic burdens.
References
1. Lurie J , Tomkins-Lane C . Management of lumbar spinal stenosis . BMJ . 2016 ; 352 : h6234 . Crossref , PubMed Google Scholar
2. Katz JN , Harris MB . Clinical practice. Lumbar spinal stenosis . N Engl J Med . 2008 ; 358 ( 8 ): 818 – 825 . Crossref , PubMed Google Scholar
3. Katz JN , Zimmerman ZE , Mass H , Makhni MC . Diagnosis and management of lumbar spinal stenosis: A review . JAMA . 2022 ; 327 ( 17 ): 1688 – 1699 . Crossref , PubMed Google Scholar
4. Battié MC , Ortega-Alonso A , Niemelainen R , et al. Lumbar spinal stenosis is a highly genetic condition partly mediated by disc degeneration . Arthritis Rheumatol . 2014 ; 66 ( 12 ): 3505 – 3510 . Crossref , PubMed Google Scholar
5. Wu F , Huang X , Zhang Z , Shao Z . A meta-analysis assessing the association between COL11A1 and GDF5 genetic variants and intervertebral disc degeneration susceptibility . Spine (Phila Pa 1976) . 2020 ; 45 ( 11 ): E616 – E623 . Crossref , PubMed Google Scholar
6. Suri P , Stanaway IB , Zhang Y , et al. Genome-wide association studies of low back pain and lumbar spinal disorders using electronic health record data identify a locus associated with lumbar spinal stenosis . Pain . 2021 ; 162 ( 8 ): 2263 – 2272 . Crossref , PubMed Google Scholar
7. Jiang H , Moro A , Liu Y , et al. Two GWAS-identified variants are associated with lumbar spinal stenosis and Gasdermin-C expression in Chinese population . Sci Rep . 2020 ; 10 ( 1 ): 21069 . Crossref , PubMed Google Scholar
8. Ma M , Huang D-G , Liang X , et al. Integrating transcriptome-wide association study and mRNA expression profiling identifies novel genes associated with bone mineral density . Osteoporos Int . 2019 ; 30 ( 7 ): 1521 – 1528 . Crossref , PubMed Google Scholar
9. Barfield R , Feng H , Gusev A , et al. Transcriptome-wide association studies accounting for colocalization using Egger regression . Genet Epidemiol . 2018 ; 42 ( 5 ): 418 – 433 . Crossref , PubMed Google Scholar
10. Qi X , Yu F , Wen Y , et al. Integration of transcriptome-wide association study and messenger RNA expression profile to identify genes associated with osteoarthritis . Bone Joint Res . 2020 ; 9 ( 3 ): 130 – 138 . Crossref , PubMed Google Scholar
11. Sakaue S , Kanai M , Tanigawa Y , et al. A cross-population atlas of genetic associations for 220 human phenotypes . Nat Genet . 2021 ; 53 ( 10 ): 1415 – 1424 . Crossref , PubMed Google Scholar
12. No authors listed . Code M48. International Classification of Diseases, 10th revision (ICD-10). World Health Organization . 2016 . https://icd.who.int/browse10/2016/en#/M48.0 ( date last accessed 1 June 2023 ). Google Scholar
13. Gusev A , Ko A , Shi H , et al. Integrative approaches for large-scale transcriptome-wide association studies . Nat Genet . 2016 ; 48 ( 3 ): 245 – 252 . Crossref , PubMed Google Scholar
14. Gusev A , Mancuso N , Won H , et al. Transcriptome-wide association study of schizophrenia and chromatin activity yields mechanistic disease insights . Nat Genet . 2018 ; 50 ( 4 ): 538 – 548 . Crossref , PubMed Google Scholar
15. Canela-Xandri O , Rawlik K , Tenesa A . An atlas of genetic associations in UK Biobank . Nat Genet . 2018 ; 50 ( 11 ): 1593 – 1599 . Crossref , PubMed Google Scholar
16. He N , Qi W , Zhao Y , Wang X . Relationship between severity of lumbar spinal stenosis and ligamentum flavum hypertrophy and serum inflammatory factors . Comput Math Methods Med . 2022 ; 2022 : 8799240 . Crossref , PubMed Google Scholar
17. Oh I-S , Suh D-W , Ha K-Y . Hypertrophy of the ligament flavum in degenerative lumbar stenosis associated with the increased expression of fractalkine (CX3CL1)/CX3CR1 chemokine . Connect Tissue Res . 2013 ; 54 ( 6 ): 380 – 385 . Crossref , PubMed Google Scholar
18. Weber KT , Alipui DO , Sison CP , et al. Serum levels of the proinflammatory cytokine interleukin-6 vary based on diagnoses in individuals with lumbar intervertebral disc diseases . Arthritis Res Ther . 2016 ; 18 : 3 . Crossref , PubMed Google Scholar
19. Oyadomari S , Morimoto M . The hypertrophy of the ligamentum flavum is assosciated with activation of the unfolded protein respone . Gene Expression Omnibus . 2018 . https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE113212 ( date last accessed 6 June 2023 ). Google Scholar
20. Sato N , Higashino K , Sakai T , et al. Role of epiligament in ligamentum flavum hypertrophy in patients with lumbar spinal canal stenosis:a pilot study . J Med Invest . 2018 ; 65 ( 1.2 ): 85 – 89 . Crossref , PubMed Google Scholar
21. Yuan J , Jia J , Wu T , et al. Comprehensive evaluation of differential long non-coding RNA and gene expression in patients with cartilaginous endplate degeneration of cervical vertebra . Exp Ther Med . 2020 ; 20 ( 6 ): 260 . Crossref , PubMed Google Scholar
22. Zhou Y , Zhou B , Pache L , et al. Metascape provides a biologist-oriented resource for the analysis of systems-level datasets . Nat Commun . 2019 ; 10 ( 1 ): 1523 . Crossref , PubMed Google Scholar
23. Pérez-López A , Valadés D , Vázquez Martínez C , de Cos Blanco AI , Bujan J , García-Honduvilla N . Serum IL-15 and IL-15Rα levels are decreased in lean and obese physically active humans . Scand J Med Sci Sports . 2018 ; 28 ( 3 ): 1113 – 1120 . Crossref , PubMed Google Scholar
24. Yoshida S , Fujimoto T , Takahashi T , et al. Interleukin-15 receptor subunit alpha regulates interleukin-15 localization and protein expression in skeletal muscle cells . Exp Physiol . 2022 ; 107 ( 3 ): 222 – 232 . Crossref , PubMed Google Scholar
25. Loro E , Ramaswamy G , Chandra A , et al. IL15RA is required for osteoblast function and bone mineralization . Bone . 2017 ; 103 : 20 – 30 . Crossref , PubMed Google Scholar
26. Guo Q , Lv S-Z , Wu S-W , Tian X , Li Z-Y . Association between single nucleotide polymorphism of IL15RA gene with susceptibility to ossification of the posterior longitudinal ligament of the spine . J Orthop Surg Res . 2014 ; 9 : 103 . Crossref , PubMed Google Scholar
27. Yabe Y , Hagiwara Y , Ando A , et al. Chondrogenic and fibrotic process in the ligamentum flavum of patients with lumbar spinal canal stenosis . Spine (Phila Pa 1976) . 2015 ; 40 ( 7 ): 429 – 435 . Crossref , PubMed Google Scholar
28. Sun C , Zhang H , Wang X , Liu X . Ligamentum flavum fibrosis and hypertrophy: Molecular pathways, cellular mechanisms, and future directions . FASEB J . 2020 ; 34 ( 8 ): 9854 – 9868 . Crossref , PubMed Google Scholar
29. Choi M , Choi YM , An I-S , Bae S , Jung JH , An S . E3 ligase RCHY1 negatively regulates HDAC2 . Biochem Biophys Res Commun . 2020 ; 521 ( 1 ): 37 – 41 . Crossref , PubMed Google Scholar
30. Almozyan S , Coulton J , Babaei-Jadidi R , Nateri AS . FLYWCH1, a multi-functional zinc finger protein contributes to the DNA repair pathway . Cells . 2021 ; 10 ( 4 ): 889 . Crossref , PubMed Google Scholar
31. Sasaki T , Lian S , Khan A , et al. Autolysosome biogenesis and developmental senescence are regulated by both Spns1 and v-ATPase . Autophagy . 2017 ; 13 ( 2 ): 386 – 403 . Crossref , PubMed Google Scholar
32. He R , Wang Z , Cui M , et al. HIF1A alleviates compression-induced apoptosis of nucleus pulposus derived stem cells via upregulating autophagy . Autophagy . 2021 ; 17 ( 11 ): 3338 – 3360 . Crossref , PubMed Google Scholar
33. Roh EJ , Darai A , Kyung JW , et al. Genetic therapy for intervertebral disc degeneration . Int J Mol Sci . 2021 ; 22 ( 4 ): 1579 . Crossref , PubMed Google Scholar
34. Patil P , Niedernhofer LJ , Robbins PD , Lee J , Sowa G , Vo N . Cellular senescence in intervertebral disc aging and degeneration . Curr Mol Biol Rep . 2018 ; 4 ( 4 ): 180 – 190 . Crossref , PubMed Google Scholar
Author contributions
J. Xu: Writing – original draft, Methodology, Formal analysis, Conceptualization, Resources.
H. Si: Funding acquisition, Writing – original draft, Data curation.
Y. Zeng: Methodology, Resources, Validation.
Y. Wu: Formal analysis, Methodology.
S. Zhang: Visualization, Resources, Data curation.
B. Shen: Writing – review & editing, Project administration, Supervision.
Funding statement
The authors disclose receipt of the following financial or material support for the research, authorship, and/or publication of this article: this work was supported by the National Natural Science Foundation of China (grant number 81974347 and 81802210) and the Department of Science and Technology of Sichuan Province (grant number 2021YFS0122).
ICMJE COI statement
The authors declare that they have no competing interests.
Data sharing
The GWAS summary data set is available from the BBJ (https://pheweb.jp/downloads). The datasets analyzed during the current study are available from the Gene Expression Omnibus database (https://www.ncbi.nlm.nih.gov/gds) accession number: GSE113212.
Ethical review statement
The source of the data was a publicly available database and no human participants were involved, hence ethical parameters are not applicable.
Open access funding
The open access fee for this study was funded by the National Natural Science Foundation of China (Program No.82272561 and No. 81974347) and the Foundation of the Science & Technology Department of Sichuan Province (Program No. 2022YFS0050).
Supplementary material
Tables showing the total genes identified by transcriptome-wide association study analysis and gene expression omnibus.
© 2023 Author(s) et al. This is an open-access article distributed under the terms of the Creative Commons Attribution Non-Commercial No Derivatives (CC BY-NC-ND 4.0) licence, which permits the copying and redistribution of the work only, and provided the original author and source are credited. See https://creativecommons.org/licenses/by-nc-nd/4.0/
