Molecular insights into gender-specific differences in rheumatoid arthritis: A study using high-throughput sequencing and Mendelian randomization

Abstract

Rheumatoid arthritis (RA) is a multifaceted autoimmune disorder with notable gender differences. The impact of gender-specific genetic variations on RA's pathogenesis remains unclear. This research investigates gender-specific genes in RA using Mendelian randomization (MR) and transcriptome sequencing to understand RA mechanisms, focusing on gender-specific immune responses. We used the limma package to analyze gender-differential genes in RA patients, followed by Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes analyses. Key genes were identified through least absolute shrinkage and selection operator and support vector machine recursive feature elimination, and a diagnostic nomogram was constructed and validated using the gene expression omnibus database. Immune cell and function analyses were performed using single-sample gene set enrichment analysis, and a competing endogenous ribonucleic acid (RNA) network was constructed. MR analysis was used to investigate the causal relationship between gender-differential genes and RA. Thirty differentially expressed genes were identified, with Gene Ontology analysis indicating their involvement in neutrophil-mediated killing. Six key genes (MAP7D2, AR, DNAH6, CXorf36, ORM1, and IQGAP3) were identified, and a diagnostic nomogram was developed (area under the curve: 0.956 in the training set, 0.859 in the validation set GSE68689). Furthermore, single-sample gene set enrichment analysis analysis indicated higher immune cell infiltration in female RA patients, highlighting gender's influence on immune response. The competing endogenous RNA network revealed potential RNA regulatory pathways. MR analysis found that the RETN gene has a specific role in seronegative RA patients, particularly in females. This study enhances the understanding of RA's gender-specific pathogenesis and offers a foundation for future personalized treatment and prevention strategies, aiding in developing more effective individualized treatment plans for RA patients.

Keywords: Mendelian randomization; ceRNA; gender differences; machine learning; rheumatoid arthritis; transcriptomics.

Figure 1.

Flowchart of this study. ceRNA = competing endogenous ribonucleic acid, DEGs = differentially expressed genes, GEO = gene expression omnibus, GO = Gene Ontology, KEGG = Kyoto Encyclopedia of Genes and Genomes, MR = Mendelian randomization, PBMCs = peripheral blood mononuclear cells, RA = rheumatoid arthritis.

Figure 2.

Expression levels of DEGs in female and male. (A) Volcano plot of the DEGs; (B) heatmap of the DEGs. DEGs = differentially expressed genes, FDR = false discovery rate.

Figure 3.

Mendelian randomization (MR) model. The causal estimate of MR studies must meet the following 3 assumptions: assumption 1, the selected genetic instruments must be powerfully associated with exposure (30 DEGs); assumption 2, the selected genetic instruments are unrelated to the potential confounders; assumption 3, the selected genetic instruments do not affect outcome (rheumatoid arthritis) independently of exposure (30 DEGS). DEGs = differentially expressed genes, KEGG = Kyoto Encyclopedia of Genes and Genomes.

Figure 4.

Functional analyses of the DEGs. (A) GO analysis of the DEGs; (B) KEGG analysis of the DEGs. DEGs = differentially expressed genes, GO = Gene Ontology, KEGG = Kyoto Encyclopedia of Genes and Genomes, LASSO = least absolute shrinkage and selection operator, SVM-RFE = support vector machine recursive feature elimination.

Figure 5.

Six hub DEGs were identified as diagnostic genes. (A) and (B) LASSO logistic regression, with penalty parameter tuning conducted by 10-fold cross-validation was used to select 16 genes; (C) the SVM-RFE algorithm was used to filter the 6 genes to identify; (D) hub genes obtained from the SVM-RFE and LASSO models; (E) ROC curves for the 6 hub genes from our sequencing data. DEGs: differentially expressed genes, LASSO = least absolute shrinkage and selection operator, ROC = receiver operating characteristic, SVM-RFE = support vector machine recursive feature elimination.

Figure 6.

Diagnostic model construction. (A) Construction of a nomogram model with 6 hub genes; (B) calibration curve of the model; (C) ROC curve for evaluating the diagnostic model’s performance. AUC = area under the curve, ROC = receiver operating characteristic.

Figure 7.

Validation of the diagnostic model. (A) ROC curve for validating the diagnostic model’s performance in GSE68689; (B) expression of the 6 marker genes in the validation dataset GSE68689. ROC = receiver operating characteristic.

Figure 8.

Differentially infiltrated immune cells and functions in female patients and male patients. (A) Heatmap of differential immune cells and functions; (B) correlation matrix of 16 immune cells; (C) correlation matrix of 13 immune functions; (D) the ssGSEA scores of 16 immune cells; (E) the ssGSEA scores of 13 immune functions. *P < .05, **P < .01, ***P < .001. ssGSEA = single-sample gene set enrichment analysis.

Figure 9.

Correlation among hub genes with immune cells and functions. Heatmap of correlation among 6 hub genes with immune cells and functions.

Figure 10.

ceRNA networks based on 6 hub genes. The pink circle represents the mRNAs, the orange triangle represents the miRNAs and the blue diamond represents the lncRNAs. ceRNA = competing endogenous ribonucleic acid, DEGs = differentially expressed genes, lncRNA = long non-coding RNA, mRNA= messenger RNA, RA = rheumatoid arthritis.

Figure 11.

The causality of RETN on seronegative RA risk. (A) Leave-one-out analysis. Black points depict the IVW method was used to assess the causal effect, excluding single specific variant from the analysis. The red point denotes the inverse-variance weighted estimate sing all SNPs. (B) Forest plot. The red points demonstrate the integrated estimates using all SNPs together, using IVW method. Horizontal lines represent 95% confidence intervals. (C) Scatter plot. The slope of each line denotes the estimated effect of per Mendelian randomization method. (D) Funnel plot. Vertical lines represent estimates with all SNPs. Symmetry in the funnel plot demonstrates no obvious horizontal pleiotropy. IVW = inverse-variance weighted, MR = Mendelian randomization, RA = rheumatoid arthritis, SNP = single nucleotide polymorphism.