Multi-omics study on autophagic dysfunction molecular network in the pathogenesis of rheumatoid arthritis

Abstract

Background

Autophagy is associated with the development of rheumatoid arthritis (RA), but its genetic pathological mechanisms remain incompletely understood. In this study, we employed summary-data-based Mendelian randomization (SMR) and co-localization analysis to systematically investigate the relationship between autophagy-related genes and RA.

Methods

We obtained summary data on blood methylation (mQTL), gene expression (eQTL), and protein abundance (pQTL) from respective quantitative trait locus (QTL) studies. Genetic association data for RA were primarily derived from the FinnGen database, with validation performed using the UK Biobank (UKB) and GWAS Catalog databases. SMR analysis was conducted to evaluate the association between molecular characteristics of autophagy-related genes and RA. Subsequently, co-localization analysis was performed to determine whether the identified signals share the same causal genetic variants.

Results

After integrating mQTL-eQTL multi-omics data, we identified two key autophagy genes, BCL2L1 and RAF1, which may have a causal relationship with RA. Significant associations were found for BCL2L1 (cg12873919, cg13989999) and RAF1 (cg26432171) in the SMR analysis of autophagy-related mQTL, eQTL, and GWAS data (p SMR < 0.05). In the integrated mQTL-eQTL SMR analysis, cg12873919 (p SMR = 1.40E-07, OR = 0.82, 95% CI [0.76–0.88]), cg13989999 (p SMR = 1.43E-06, OR = 0.78, 95% CI [0.71–0.87]), and cg26432171 (p SMR = 9.18E-09, OR = 1.83, 95% CI [1.49–2.25]) were all significantly validated. Methylation of cg12873919 and cg13989999 in BCL2L1 was associated with increased BCL2L1 expression, consistent with their negative impact on RA risk. Conversely, the cg26432171 site in RAF1 showed a positive correlation between gene methylation and expression. In the eQTL-GWAS SMR analysis, MAPK3 expression (p SMR = 7.24E-05, OR = 0.91, 95% CI [0.87–0.95]) was negatively correlated with RA risk, a finding supported by co-localization analysis (PPH4 > 0.5), suggesting that this gene may inhibit RA pathogenesis by regulating the autophagy process. Furthermore, protein level analysis also supported the protective role of MAPK3 (p SMR = 7.53E-05, OR = 0.89, 95% CI [0.84–0.94]).

Conclusion

We identified that autophagy-related genes BCL2L1 and RAF1 may be associated with RA risk, providing strong evidence from multi-omics data. This study identifies autophagy genes related to RA, potentially offering new insights into the pathogenesis of RA.

Introduction

Autophagy is a crucial intracellular degradation pathway that maintains cellular homeostasis and regulates cellular functions by encapsulating damaged, denatured, or excess components within autophagosomes and transporting them to lysosomes for degradation [1]. Rheumatoid arthritis (RA) is a chronic inflammatory joint disease characterized by joint pain, swelling, and functional impairment, with its pathogenesis involving complex genetic, environmental, and immune factors [2]. Recent studies have revealed that autophagy plays a critical role in RA, with dysregulation of autophagy closely related to the pathophysiology of the disease [3, 4]. Autophagy not only modulates local inflammation by clearing intracellular damaged materials and inflammatory mediators but also influences disease progression by controlling apoptosis and regulating immune cell function [5]. For example, impaired autophagy in synovial cells leads to the excessive accumulation of inflammatory mediators, exacerbating local inflammation and joint destruction. Studies have shown that the autophagy-related gene ATG16L is abnormally expressed in RA patients. Mutations in the ATG16L1 gene are associated with RA susceptibility, and its expression in the synovium of RA patients is significantly reduced [6]. This decrease in expression may impair autophagosome formation and function, leading to the accumulation of inflammatory mediators in synovial cells. LC3B, an important marker of autophagosome maturation, also exhibits abnormal expression in RA patients. Research suggests that increased expression of LC3B may be related to enhanced autophagic activity in synovial tissues, which could be part of the mechanism underlying chronic inflammation in RA [7].

Recent studies have further discovered that autophagy may have a dual role in the pathogenesis of RA [8]. On the one hand, autophagy protects cells from damage by clearing damaged mitochondria and oxidative stress products, thereby inhibiting excessive activation of inflammatory responses. On the other hand, in RA patients, excessive activation of autophagy may lead to overactivation of immune cells and abnormal secretion of inflammatory mediators, exacerbating joint inflammation and damage. Therefore, abnormal expression of autophagy-related genes not only highlights the importance of autophagy in the pathogenesis of RA but also suggests that autophagy could be a potential therapeutic target for RA. Although the role of autophagy in RA is gradually being recognized, systematic studies to comprehensively explore the causal relationship between autophagy and RA are still lacking. Most existing studies rely on small-scale samples or animal models, lacking support from large-scale genetic and omics data.

Based on the above background, this study aims to systematically assess the causal effects between autophagy and RA by integrating whole-genome data with autophagy-related gene information using the SMR method (Fig. 1). We hope to provide a more comprehensive understanding of RA pathogenesis and offer strong support for the development of targeted therapeutic strategies.

Fig. 1

Study design. QTL, quantitative trait loci; RA, rheumatoid arthritis; UKB, UK Biobank; SNP, single nucleotide polymorphisms; SMR, summary-based mendelian randomization; HEIDI, heterogeneity in dependent instrument; PPH4, posterior probability of H4

Materials and methods

Data sources

The SMR analysis in this study is based on publicly available datasets. Autophagy-related genes were obtained from the Autophagy Database (http://www.autophagy.lu/), totaling 222 autophagy genes. RA GWAS data were derived from several international datasets, including the FinnGen dataset (discovery), the UKB dataset (validation), and the GWAS Catalog dataset (validation). Quantitative trait loci (QTL) reveal associations between single nucleotide polymorphisms (SNPs) and DNA methylation, gene expression, and protein abundance levels. The mQTL dataset was obtained from a meta-analysis of two European cohorts [9], including the Brisbane Systems Genetics Study (n = 614) and the Lothian Birth Cohorts (n = 1366). The eQTL dataset was obtained from the eQTLGen database [10], which includes genetic data on blood gene expression from 31,684 individuals. Genetic association data for protein levels were obtained from the pQTL study by Pietzner et al. [11], which included 10,708 Europeans. By integrating results from MR analyses at these three different levels, we identified candidate genes with potential causal relationships. Furthermore, there was no sample overlap between the exposure and outcome groups.

SMR analysis

In this study, we used the SMR (v1.3.1) tool to perform SMR analysis and HEIDI tests to evaluate the association between autophagy-related gene expression, methylation, and protein abundance with RA. By selecting the top associated cis-QTLs within ± 1000 kb upstream and downstream of the target gene with a P-value less than 5.0 × 10^− 8, we validated the advantages of the SMR method over traditional Mendelian Randomization (MR) analysis when exposure and outcome data are independent. We excluded SNPs with allele frequency differences greater than 0.2 and focused on significant causal association signals between mQTL, eQTL, and pQTL.To ensure accuracy, we set a P-value threshold of 0.05 as the statistical significance criterion. The HEIDI test was used to assess the heterogeneity of the SMR analysis results and determine whether significant differences exist between different datasets. In this study, we performed the HEIDI test and selected SNPs with a HEIDI P-value less than 0.05 as those with higher consistency between gene expression data and disease association data. For each SNP, we set the default parameters to test its consistency across different datasets and ensure the robustness and statistical significance of the results. Additionally, the developed multi-SNP SMR analysis method (SMR multi), combined with the HEIDI test, ensured that the results included in the final analysis were both significant and devoid of pleiotropy (p SMR < 0.05, p SMR multi < 0.05, and p HEIDI > 0.05), laying the foundation for subsequent co-localization and integrative analyses.

Co-localization analysis

We performed co-localization analysis using the R package coloc to identify common causal variants between autophagy gene-associated cis-QTLs (including mQTL, eQTL, pQTL) and RA. In the co-localization analysis, five different posterior probabilities were reported, corresponding to five mutually exclusive hypotheses: (H0) no trait in the region is associated with the SNP; (H1) only trait 1 is associated with the SNP; (H2) only trait 2 is associated with the SNP; (H3) both traits are associated with the SNP but with different causal variants; and (H4) both traits are associated with the SNP and share a causal variant. When selecting cis-QTLs, we filtered them based on the distance between the SNP and the target gene. We defined cis-QTLs as SNPs located within a 50 kb upstream and downstream range of the gene. These cis-QTLs were extracted from GWAS summary data and matched with genes to ensure they were within the appropriate gene expression regulatory regions. Furthermore, only those cis-QTLs that maintained consistency across multiple datasets were included to enhance the reliability of the results. Co-localization was considered successful when QTL signals with a p12 = 5 × 10^− 5 met the criterion of PPH4 > 0.5, allowing for the co-localization of signals with weaker p-values.

Statistical analysis

All statistical analyses were performed using R (v4.4.1). The R packages “ggplot2” and “ggrepel” were used for Manhattan plot generation, and “forestplot” was used for forest plot generation. The codes for generating SMR LocusPlot and SMR EffectPlot were sourced from Zhu et al. (PMID: 27019110).

Results

Autophagy genes methylation and RA

An SMR analysis was conducted to assess the relationship between autophagy gene methylation and RA (mQTL-GWAS), with full results available in Table S1. We identified 40 CpG sites significantly associated with RA, involving 23 genes (p SMR < 0.05 & p SMR multi < 0.05 & p HEIDI > 0.05, Fig. 2). The association between methylation levels at these CpG sites and RA risk was determined by SMR analysis and further validated by the HEIDI test to exclude pleiotropy. The effect estimates for different CpG sites within the same gene were not always consistent. For instance, the methylation level at CpG site cg16557858 in the ERBB2 gene was negatively correlated with RA risk (OR = 0.92, 95% CI [0.85–0.99]), while another CpG site, cg14187895, showed a positive correlation with RA risk (OR = 1.12, 95% CI [1.03–1.23]). Co-localization analysis revealed that 11 of these identified signals (corresponding to 7 genes) had strong co-localization evidence (PPH4 > 0.5, Fig. 2), including ATG16L2 (cg21806242), CTSD (cg12007048, cg15017982, cg20973931, cg22079043), EIF2AK3 (cg05707116), ERBB2 (cg05616858), RGS19 (cg10087092), TP63 (cg13518031), and TSC2 (cg00500602, cg04515572). The co-localization results for these 11 CpG sites are shown in Figure S1. Among the identified CpG sites, the site cg10087092 in RGS19 was validated in the UKB database (p SMR = 3.59E-02 & p SMR multi = 3.59E-02 & p HEIDI = 9.43E-01, OR = 0.76, 95% CI [0.58–0.98], Table S2). Additionally, the site cg13518031 in TP63 was validated in the GWAS Catalog cohort (p SMR = 1.76E-02 & p SMR multi = 1.76E-02 & p HEIDI = 2.42E-01, OR = 1.21, 95% CI [1.03–1.41], Table S3).

Fig. 2

Association between autophagy gene methylation and RA in MR analysis (Forest map). OR, odds ratio; CI, confidence interval

Autophagy genes expression levels and RA

SMR analysis was also performed to examine the relationship between autophagy gene expression levels and RA (eQTL-GWAS), with full results available in Table S4. We identified 11 autophagy-related genes associated with RA risk (p SMR < 0.05 & p SMR multi < 0.05 & p HEIDI > 0.05, Fig. 3A). The expression levels of MAPK3 (OR = 0.91, 95% CI [0.87–0.95]) and ST13 (OR = 0.84, 95% CI [0.72–0.98]) were negatively correlated with RA risk, while the expression levels of the remaining nine genes, including ATG16L2 (OR = 1.53, 95% CI [1.16–2.01]), BCL2L1 (OR = 1.49, 95% CI [1.03–2.15]), EEF2 (OR = 1.27, 95% CI [1.00-1.62]), IL24 (OR = 1.28, 95% CI [1.08–1.52]), ITGB1 (OR = 1.43, 95% CI [1.01–2.02]), KIAA0226 (OR = 1.33, 95% CI [1.06–1.67]), PPP1R15A (OR = 1.13, 95% CI [1.02–1.25]), PRKCQ (OR = 1.25, 95% CI [1.06–1.47]), and RAF1 (OR = 1.29, 95% CI [1.10–1.51]), were positively correlated with RA risk. Among these identified signals, KIAA0226 (PPH4 = 0.53) and MAPK3 (PPH4 = 0.81) had strong co-localization evidence within the co-localization region windows (Fig. 3A), with their co-localization results shown in Figure S2A. However, these findings were not validated in the UKB database or the GWAS Catalog cohort (Tables S5-6).

Fig. 3

Association between autophagy gene expression, protein expression abundance and RA in MR analysis (Forest map). A, gene expression; B, protein abundance. OR, odds ratio; CI, confidence interval

Autophagy protein abundance and RA

An SMR analysis was also conducted to assess the relationship between autophagy protein abundance and RA (pQTL-GWAS), with full results available in Table S7. We identified two autophagy-related proteins (MAPK3 and ST13) associated with RA risk (p SMR < 0.05 & p SMR multi < 0.05 & p HEIDI > 0.05, Fig. 3B). The abundance of both MAPK3 (OR = 0.89, 95% CI [0.84–0.94]) and ST13 (OR = 0.72, 95% CI [0.54–0.97]) was negatively correlated with RA risk. The chromosomal distribution of these two proteins is shown in Fig. 4A. Co-localization analysis provided strong evidence supporting the involvement of MAPK3 protein (PPH4 = 0.98), with its co-localization results shown in Figure S2B. These findings were not validated in the UKB database or the GWAS Catalog cohort (Tables S8-9).

Fig. 4

Distribution of key results obtained from integrated analysis on chromosomes (Manhattan map). The dashed line corresponds to a threshold p SMR multi = 0.05. A, Protein abundance; B, Gene methylation; C, Gene expression

Integration of methylation and expression levels of autophagy genes

Based on the key results from the SMR analysis of mQTL-GWAS and eQTL-GWAS, we suggest that the genes RAF1, ATG16L2, EEF2, and BCL2L1 may have causal associations with RA and warrant further investigation. An SMR analysis of mQTL-eQTL was performed using mQTL as the exposure and eQTL as the outcome to explore whether methylation at the identified CpG sites significantly regulates the expression of the corresponding genes, with all results available in Table S10. Further analysis revealed that the CpG sites cg12873919 (p SMR = 1.40E-07 & p SMR multi = 1.40E-07 & p HEIDI = 1.97E-01) and cg13989999 (p SMR = 1.43E-06 & p SMR multi = 1.43E-06 & p HEIDI = 2.36E-01) in BCL2L1, as well as cg26432171 (p SMR = 9.18E-09 & p SMR multi = 9.18E-09 & p HEIDI = 5.97E-02) in RAF1, were significantly validated in these analyses (Table 1). No positive cross-results were found between the mQTL-eQTL integrated analysis and the pQTL-RA SMR analysis; therefore, an eQTL-pQTL SMR analysis was not performed.

Table 1 mQTL-eQTL SMR analysis results

Integration of multi-omics evidence

After integrating multi-omics evidence, we identified two genes with primary multi-omics evidence, BCL2L1 (cg12873919, cg13989999) and RAF1 (cg26432171), that may have causal associations with RA. The chromosomal distribution of these related loci and genes is shown in Fig. 4B-C. The sites cg12873919 and cg13989999 in BCL2L1, and cg26432171 in RAF1, were significantly validated in the mQTL, eQTL, and integrated mQTL-eQTL SMR analyses. Based on OR values to determine risk associations and regulatory direction, the methylation levels of the sites cg12873919 and cg13989999 in BCL2L1 were negatively correlated with RA risk, while the gene expression levels were positively correlated with RA risk. Therefore, the methylation at sites cg12873919 and cg13989999 negatively regulates the expression level of BCL2L1. The methylation level at site cg26432171 in RAF1 was positively correlated with RA risk, as was the gene expression level. Therefore, methylation at site cg26432171 positively regulates the expression level of RAF1. In summary, the proposed mechanism might be that lower methylation levels at sites cg12873919 and cg13989999 promote the expression of BCL2L1, thereby increasing RA risk; whereas higher methylation levels at site cg26432171 promote the expression of RAF1, also increasing RA risk. The SMR LocusPlot trajectory zoom-in diagrams illustrate the genetic consistency between mQTL, eQTL, and RA GWAS within a 500 kb window for BCL2L1 (cg12873919, cg13989999) and RAF1 (cg26432171). The results indicate significant associations between SNP sites cg12873919, cg13989999 and BCL2L1 expression, as well as cg26432171 and RAF1 expression (p SMR multi < 0.05, Fig. 5). The SMR EffectPlot analysis suggests a strong correlation between methylation or expression levels at specific loci and GWAS effect sizes, indicating a potential causal relationship between methylation, gene expression, and disease risk (p SMR multi < 0.05, Fig. 6). These results further support the role of methylation and gene expression in RA-related genetic risk, particularly with high r² values in top cis-QTLs (marked in red) demonstrating robust associations.

Fig. 5

SMR LocusPlot results of BCL2L1 (cg12873919, cg13989999) and RAF1 (cg26432171). The target genes are represented in red, with p SMR multi = 0.05 as the boundary; Genes with significance within the window area are represented in blue, while the remaining genes are represented in gray; ◆, represents genes that have passed the HEIDI test (p HEIDI > 0.05); ◇, represents genes that have not passed the test. A, SMR LocusPlot results of BCL2L1; B, SMR LocusPlot results of cg12873919 and cg13989999; C, SMR LocusPlot results of RAF1; D, SMR LocusPlot results of cg26432171

Fig. 6

SMR EffectPlot results of BCL2L1 (cg12873919, cg13989999) and RAF1 (cg26432171). A-B, SMR EffectPlot results of BCL2L1; C-D, SMR EffectPlot results of RAF1

BCL2L1 and RAF1 are highly expressed in RA

Building on the previous study, we collected synovial tissue (ST) from the knees of patients with osteoarthritis (OA) and RA, isolated and cultured fibroblast like synovial cell (FLS), and performed single-cell RNA sequencing and transcriptomic analysis. Our results revealed that BCL2L1 and RAF1 were primarily expressed in FLS among various cell types (Fig. 7A). Preliminary validation showed that, compared to the OA group, BCL2L1 and RAF1 were highly expressed in both RA-ST and RA-FLS (P < 0.05, Fig. 7B). In correlation analysis, the expression level of RAF1 in the peripheral blood of RA patients showed a significant positive correlation with clinical inflammatory markers such as CRP, DAS28-CRP, and the number of swollen joints (P < 0.05, Figure S3). These findings further support the potential important role of BCL2L1 and RAF1 in the pathogenesis of RA, particularly their high expression in FLS, which may be closely related to the inflammatory response in the synovium of RA.

Fig. 7

BCL2L1 and RAF1 are highly expressed in RA. A, The single-cell sequencing results of BCL2L1 and RAF1; B, The transcriptomic sequencing results of BCL2L1 and RAF1

Discussion

In this study, we performed SMR and co-localization analyses based on multi-omics data to systematically evaluate the associations between autophagy-related mQTLs, eQTLs, and pQTLs with RA risk. We found that the autophagy-related MAPK3 gene might be associated with RA risk, while the BCL2L1 and RAF1 genes are supported by multi-omics evidence for their association with RA.

In both eQTL and pQTL analyses, we identified a significant association between MAPK3 and RA. The protein encoded by MAPK3 plays a critical role in the MAPK/ERK signaling pathway, which is involved in various cellular processes such as proliferation, differentiation, and stress response [12]. Dysregulation of this pathway has been implicated in the pathogenesis of autoimmune diseases, including RA [13]. Our findings suggest that genetic variants affecting MAPK3 expression may mediate the disease process of RA by regulating autophagy mechanisms. Specifically, under certain conditions, the ERK/MAPK pathway can promote autophagy by inhibiting the activity of Bcl-2 family proteins. However, in another context, this pathway may inhibit autophagy by activating the mTOR signaling pathway. Therefore, the mechanism by which MAPK3 influences RA may be through its ability to regulate autophagy. Our eQTL analysis showed that certain SNPs are associated with increased MAPK3 expression in immune cells, which may lead to enhanced activation of autophagy. This enhanced activation could potentially inhibit inflammatory responses, protect chondrocytes, and negatively regulate inflammatory signaling pathways, thereby exerting a protective effect in RA. Additionally, pQTL analysis provided evidence that these genetic variants also affect MAPK3 protein abundance, further supporting the role of MAPK3 in RA pathogenesis. The association between MAPK3 expression and protein levels highlights the potential of MAPK3 as a therapeutic target for RA. Targeted therapies aimed at MAPK3 or its downstream signaling components may help modulate inflammatory responses and alleviate symptoms in RA patients. Although our study suggests a protective role for MAPK3 in RA, it is also possible that, as part of the MAPK signaling pathway, MAPK3 may contribute to pro-inflammatory processes under certain conditions [14, 15]. Thus, the exact role of MAPK3 in RA may depend on specific cellular contexts and regulatory mechanisms, warranting further research to fully elucidate its role in RA pathogenesis.

The BCL2L1 gene encodes the Bcl-XL protein, which plays a key role in regulating apoptosis and maintaining cell survival [16, 17]. In RA patients, impaired autophagy may lead to an excessive response of immune cells to inflammatory factors, thereby exacerbating the synovial inflammatory response [18]. BCL2L1, by influencing the regulation of autophagy, may promote the survival of immune cells and enhance their inflammatory response. Further analysis indicates that increased expression of BCL2L1 may enhance autophagic activity, leading to excessive accumulation of inflammatory factors and thus worsening the pathological progression of RA. Additionally, BCL2L1 maintains the survival of T cells and macrophages by inhibiting the mitochondrial apoptosis pathway, potentially prolonging the lifespan of pro-inflammatory immune cells and promoting the chronicization of the local inflammatory microenvironment in the joints [19]. Our SMR analysis revealed a significant causal association between the BCL2L1 gene and RA risk. Specifically, the methylation sites cg12873919 and cg13989999 within the BCL2L1 gene exhibited significant negative regulatory effects in mQTL-eQTL analysis, with methylation levels inversely correlated with BCL2L1 expression. This negative regulation may lead to elevated BCL2L1 expression, thereby increasing RA risk. Elevated BCL2L1 expression may influence RA pathogenesis by enhancing autophagic activity. In RA patients, impaired autophagy function may lead to the excessive accumulation of inflammatory factors and abnormal local immune responses, thereby promoting disease progression [20, 21]. Therefore, BCL2L1 may influence RA pathogenesis by regulating autophagy and apoptosis. Our results suggest that BCL2L1 could be a potential therapeutic target for RA, and modulating its expression levels or autophagic function may have significant clinical implications for RA treatment.

Another key gene identified in this study is RAF1, which encodes the RAF1 protein, a critical component of the MAPK/ERK signaling pathway that regulates cell proliferation, differentiation, and stress responses [22, 23]. The study found that abnormal activation of RAF1 can interfere with autophagic flux by phosphorylating ULK1 (a key autophagy initiation protein), leading to the accumulation of damaged organelles, which in turn activates the NLRP3 inflammasome and promotes the release of inflammatory factors such as IL-1β [24]. Additionally, in the regulation of immune cell function, the RAF1-MAPK signaling pathway in Th17 cells may enhance the production of IL-17 A by regulating the transcriptional activity of RORγt, thereby exacerbating joint damage [25]. Our research demonstrated a significant positive association between the RAF1 gene and RA risk, particularly in mQTL and eQTL analyses, where the methylation site cg26432171 showed a significant positive regulatory effect, with methylation levels positively correlated with RAF1 expression. This suggests that increased methylation at the cg26432171 site may promote RAF1 expression, thereby elevating RA risk. The role of RAF1 in RA may be closely related to its function in cellular signaling. Abnormal activation of the MAPK/ERK pathway can lead to excessive activation of immune cells and abnormal secretion of inflammatory cytokines, exacerbating joint inflammation and damage [13]. The study shows that the high expression of RAF1 in the peripheral blood of RA patients is significantly positively correlated with clinical inflammatory markers such as CRP and DAS28-CRP, suggesting that RAF1 may play a promotive role in the inflammatory response of RA [26]. Elevated RAF1 expression may promote abnormal activation of this signaling pathway, further driving RA pathogenesis. Our findings underscore the critical role of RAF1 in RA pathogenesis and suggest that it could be a potential therapeutic target. Targeted therapies aimed at RAF1 or the MAPK/ERK signaling pathway may help modulate inflammatory responses and alleviate RA symptoms.

This study represents the first systematic exploration of the potential causal relationship between autophagy-related genes and RA pathogenesis through comprehensive multi-omics analysis integrating SMR and colocalization approaches. Distinct from previous single-omics or genetic association studies, our research synthesizes multidimensional data layers including gene methylation, gene expression, and protein abundance, providing more holistic insights. Particularly through SMR analysis, we not only investigated the genetic susceptibility of autophagy genes but also clarified their functional roles in RA development, establishing a solid theoretical foundation for subsequent targeted therapeutic strategies and biomarker development. Based on our research findings, the high expression of BCL2L1 and RAF1 in RA is closely related to the disease’s inflammatory response and immune cell function, making them potential biomarkers for RA. Specifically, the expression level of BCL2L1 may reflect the severity of synovial inflammation in patients and could serve as an assessment tool for disease activity. Meanwhile, RAF1 expression levels may be significantly correlated with clinical inflammation markers in patients, such as CRP and DAS28-CRP. Therefore, detecting the expression levels of BCL2L1 and RAF1 could provide important biomarkers for the early diagnosis of RA, prediction of disease progression, and clinical monitoring. These genes can be detected through blood, synovial tissue, or FLS cells, offering convenient and reliable diagnostic tools for clinical practice. Furthermore, due to the critical role of BCL2L1 and RAF1 in RA, they may not only serve as potential biomarkers but also as therapeutic targets. Treatment strategies targeting BCL2L1 could involve small molecule inhibitors or RNA interference techniques to reduce BCL2L1 expression or function, thereby improving autophagic activity in RA patients, alleviating synovial inflammation, and reducing the overactivation of immune cells. Treatment strategies targeting RAF1 could involve developing RAF1 inhibitors or blocking its downstream MAPK/ERK signaling pathway to suppress immune cell activation and proliferation, thus mitigating RA’s inflammatory response. These targeted interventions may improve the clinical symptoms of RA patients, reduce disease activity, and slow the progression of joint damage.

While existing literature has addressed aspects of this topic, our methodological integration of multi-omics perspectives and gene-disease causal analysis offers novel evidence regarding the autophagy-RA pathogenesis connection, with innovative elements in study design and dataset utilization. Notwithstanding the depth of multi-omics analysis, we acknowledge limitations in database selection. Primary databases (FinnGen, UK Biobank, GWAS Catalog) predominantly comprising European Caucasian populations may limit generalizability to other ethnic groups due to genetic and geographical diversity. Caution should be exercised when extrapolating findings to Asian or other populations without additional validation. Regarding data quality and heterogeneity, potential biases may arise from variations in sample selection criteria, processing methods, and genotyping/phenotyping accuracy across databases. While standardized processing protocols were strictly implemented, residual errors from source data characteristics remain possible. Moreover, although SMR and colocalization analyses rigorously investigated gene-RA causality, inherent limitations of public database methodologies persist. Statistical assumptions underlying SMR models and colocalization parameter settings, despite stringent controls, remain dependent on existing data frameworks and may contain unrecognized biases. While our multi-omics approach with advanced analytics has revealed novel mechanisms of autophagy-related genes in RA, we emphasize the necessity of complementary validation through diverse population studies and experimental investigations to consolidate these findings. Future research directions should prioritize cross-ethnic replications and functional studies to refine our understanding of these molecular mechanisms.

In summary, the current SMR study explored the potential causal relationships between autophagy-related gene methylation, gene expression, and protein abundance with RA and demonstrated the importance of several mitochondria-related genes and their regulation in RA pathogenesis. Future research should further validate the mechanisms of these genes and explore their clinical application value in RA patients. With the advancement of high-throughput omics technologies, integrating multi-omics data to reveal the complex pathophysiology of RA will be a key direction for future research.