Comprehensive Analysis of RNA Modifications Related Genes in the Diagnosis and Subtype Classification of Dilated Cardiomyopathy

Introduction

Dilated cardiomyopathy (DCM) is a common primary cardiomyopathy marked by left ventricular (LV) dilation and systolic impairment.¹ The clinical manifestation of DCM ranges from asymptomatic conditions to heart failure, ventricular arrhythmias, and sudden cardiac death. The causes of DCM are diverse, encompassing several hereditary and non-genetic factors. Non-genetic etiologies encompass infections, autoimmune disorders, myocarditis, and exposure to alcohol, chemicals, and poisons. Genetic mutations affecting genes related to cytoskeletal, sarcomeric components, and mitochondrial and nuclear envelope proteins comprise around 30–50% of patients.^2,3 These genetic abnormalities are consistently identified during the screening and clinical diagnosis of familial DCM.⁴ Patients with DCM due to non-genetic causes may potentially be affected by their genetic makeup.⁵ Consequently, a comprehensive comprehension of the correlations between unique DCM phenotypic attributes and genetic and molecular traits, along with advancements in early screening, would enhance the diagnosis, prevention, and treatment of DCM.

Growing evidence has shown the correlation between DCM clinical characteristics and epigenetic regulation, alongside genetic variations.^6,7 RNA modifications are regarded as vital posttranscriptional regulators of gene expression programs, significantly influencing numerous developmental processes.⁸ In recent years, RNA modifications have garnered attention due to the advent of RNA sequencing techniques for mapping modification sites, novel methods for accurate detection and quantification of modifications, and the comprehensive characterization of RNA modification effectors: “writers” “erasers” and “readers”⁹.

More than 170 distinct chemical modifications have been reported on RNAs, including N6-methyladenosine (m⁶A), N1-methyladenosine (m¹A), 5-methylcytosine (m⁵C), N7-methylguanosine (m⁷G), and alternative polyadenylation (APA) alterations. These RNA modifications are reversible and dynamic, typically facilitated by “writers” “erasers” and “readers” and play a role in different physiological and pathological processes. For instance, the m⁶A demethylase FTO was downregulated in a diabetic cardiomyopathy mouse model, and its overexpression enhanced cardiac function by mitigating myocyte hypertrophy and myocardial fibrosis.¹⁰ The m⁶A “reader” protein YTHDC1 facilitates proper contractile function, and its absence leads to the onset of DCM (characterized by significant LV enlargement and severe systolic dysfunction) in mice due to the aberrant splicing of Titin (TTN).¹¹ Moreover, mutations in TTN are recognized to cause DCM, present in roughly 18% of sporadic instances and 25% of familial instances.¹² Levels of m6A in mRNA from human heart failure samples were found to be elevated compared to control samples.¹³ Additionally, interactions may exist among these RNA modifications. For instance, the prominent m⁶A “readers” YTHDF1-3 and YTHDC1 can also function as “readers” that directly attach to m¹A in RNA.¹⁴ Nevertheless, contemporary research has predominantly concentrated on the individual genes associated with m⁶A modification in DCM, with limited evidence regarding the functions of various RNA modifications in DCM as a whole.

This study conducted an integrated investigation of multiple RNA modification-related genes in DCM. Genes associated with RNA modification exhibiting unique expression patterns in DCM were analyzed, leading to the identification of biomarkers for DCM to develop a diagnostic model. Furthermore, we endeavored to categorize DCM patients into different populations exhibiting varying genetic and molecular traits based on RNA modification-related genes. Ultimately, we employed clinical data to validate the association between the hub RNA modification-related genes and DCM. This work elucidates the correlations between DCM phenotypic characteristics and epigenetic regulation.

Methods

Data Acquisition

The microarray datasets GSE141910 and GSE120895 were obtained from the Gene Expression Omnibus (GEO) database. The GSE141910 dataset served as the training set, consisting of 332 LV free-wall tissue samples from 166 non-failing healthy donors (control) and 166 DCM patients, which were collected for analysis. The data in GSE141910 was produced using the GPL16791 Illumina HiSeq 2500 platform. The GSE120895 dataset included 55 tissue samples from 47 patients with DCM and 8 persons with normal LV systolic function (controls). The data in GSE120895 was produced using the GPL570 Affymetrix Human Genome U133 Plus 2.0 Array platform. The GSE120895 dataset used as an external validation set.

RNA Modification-Related Genes

A total of 63 RNA modification-related genes were identified from prior research,^15,16 and the expression levels of these genes were extracted from the GSE141910 dataset.

Differential Expression Analysis

Differentially expressed RNA modification-related genes (DERMGs) between the DCM and control groups (CTRL) were identified using the “Limma” package with a cut-off value of FDR < 0.05. Similarly, differential expression genes between two molecule clusters were also screened using the “Limma” package with cut-off values of FDR < 0.05 and |log2FC| > 0.263.

Function Enrichment Analyses

Enrichment analysis for biological process keywords in Gene Ontology (GO) annotations and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways was performed using DAVID (version 6.8) with a threshold of FDR < 0.05. Furthermore, pathways that were differentially enriched between two molecular clusters were examined using gene set variation analysis (GSVA). MSigDB (https://www.gsea-msigdb.org/gsea/msigdb/index.jsp), a pre-defined set of genes, used as an enrichment reference. Each pathway was measured utilizing GSVA (version 1.36.3), followed by the identification of differential pathways between two molecular clusters using an FDR < 0.05.

Establishment and Evaluation of Diagnostic Model

To identify genes linked to DCM, the DERMGs were examined using univariate logistic regression analysis provided in the “rms” package (version 6.3–0), in which DERMGs with P < 0.05 were selected. Next, the LASSO algorithm provided in the “lars” package (version 1.2) was further used to screen optimal DERMGs. Finally, a diagnostic classifier was established based on these optimal DERMGs using the Support Vector Machine (SVM) method provided in the “e1071” package (version 1.6–8). The predictive performance of the diagnostic classifier was evaluated by the receiver operator characteristic (ROC) curve, which was plotted using the “pROC” package (version 1.12.1).

Consensus Clustering Analysis

The best DERMGs were utilized to categorize 166 DCM patients into distinct molecular clusters with ConsensusClusterPlus (version 1.54.0). The RNA modification scores for each sample were quantified via GSVA (version 1.36.3) and compared among molecular clusters using the Kruskal−Wallis test.

Immune Infiltration Status

The CIBERSORT algorithm (https://cibersort.stanford.edu/index.php) was utilized to assess the proportions of 22 invading immune cells in tissue samples. The immune and stromal scores of tissue samples were inferred using the “ESTIMATE” package (http://127.0.0.1:29606/library/estimate/html/estimateScore.html).

Real-Time Quantitative PCR (RT-qPCR)

The study of our clinical samples adheres to the Declaration of Helsinki. Six myocardial samples were collected from Shaanxi Provincial People’s Hospital, comprising three normal samples (Control group, CTRL) from donors without DCM, and three samples from patients diagnosed with DCM (DCM group). The Ethical Committee of Shanxi Provincial People’s Hospital sanctioned this study, and the corresponding patient supplied written informed permission. Total RNA was extracted from myocardial samples with the Total RNA Kit I (Omega Bio-tec, Inc, USA). Total RNA was transcribed to cDNA on ice utilizing the PrimeScript™ RT reagent kit (TaKaRa, Dalian, China). Quantitative PCR (qPCR) was conducted with the Green PCR Master Mix System (Thermo Fisher Scientific). GAPDH served as an internal reference. qPCR was performed using an HT7900 Real-Time PCR System (Applied Biosystems). The relative fold change was assessed with the 2^−ΔΔCT method. The key pivotal genes in myocardial samples from Mouse DCM model were confirmed also using q-PCR. The primer sequences employed are detailed in Supplementary Table S1.

Construction of Animal Models with DCM

Eight male C57/BL mice, six weeks old and weighing 20.24 ± 2.33 g, were acquired from the Laboratory Animal Center of Xi’an Jiaotong University in China. The mice were arbitrarily assigned to two groups: a control (CTRL) group and a DCM group. Four mice in each group. Each mouse in the DCM group received a cumulative dose of 4 mg/kg doxorubicin (DOX, Sigma, St. Louis, USA) through four consecutive intraperitoneal injections over 28 days, with an injection occurring once every seven days. The CTRL was administered injections of equivalent volumes of 0.9% saline at identical time intervals. Transthoracic echocardiography was conducted in the fifth week to assess cardiac function during deep anesthesia with pentobarbital sodium (50 mg/kg, Sigma, USA). At the fifth week, cardiac tissue samples were obtained for q-PCR analysis, and all animals were euthanized by exposure to a high concentration of isoflurane. An echocardiography was performed before sacrifice, and body weight was documented during modeling.

Echocardiography

A dose of 30 mg/kg pentobarbital was administered intraperitoneally to anesthetize the mice. The KOLO SiliconWave 30 ultra-high-resolution ultrasound imaging instrument for small animals (Suzhou, China) was utilized for echocardiography. Utilizing long- and short-axis images of the left ventricle, the following parameters were assessed: left ventricular ejection fraction (LVEF), fractional shortening (FS), and left ventricular internal diameters during diastole (LVIDd) and systole (LVIDs). Data from three cardiac cycles were evaluated to evaluate cardiovascular function.

Correlation Between Hub Genes and LV Parameters in Echocardiography

The LV parameters (EF%, FS%, and LVIDs) were examined using the Pearson algorithm, and the findings were shown with the R package “ggplot2”.

Statistical Analysis

All analyses were conducted utilizing R version 3.6.1. Correlation analysis was conducted utilizing the “cor” function (http://77.66.12.57/R-help/cor.test.html), and correlation heatmaps were generated using pheatmap (version 1.0.8). The proportions of immune cell infiltration between the DCM and CTRL groups, as well as between the two molecular clusters, were analyzed using the Kruskal−Wallis test. The Kruskal−Wallis test was employed to evaluate immune and stromal scores between two molecular clusters. The expression of co-stimulatory molecules and human leukocyte antigen (HLA) family genes between two molecular clusters was evaluated using Student’s t-test. A P-value of less than 0.05 was deemed statistically significant.

Results

Expression Pattern of RNA Modification-Related Genes Changed in DCM

The primary analysis processes are displayed in Figure 1.

Figure 1 The workflow of this study.

A differential expression analysis was performed to ascertain the expression pattern of RNA modification-related genes between DCM and CTRL samples. Of the identified genes, 26 exhibited significantly altered expression patterns in the DCM samples relative to the CTRL, with m⁵C genes comprising fifty percent of this subset (Figure 2A). The expression of most genes was markedly downregulated in DCM, including the m⁶A genes YTHDC2, YTHDF3, and CBLL1. In contrast, seven genes were markedly upregulated in DCM, including two APA genes (PABPN1 and CSTF3) and five m⁵C genes (DNMT3B, TET1, NSUN6, UNG and ZBTB4). A correlation heatmap indicated that the expression of the m⁵C gene ZBTB38 was positively correlated with the APA genes PCF11 and CPSF2 but negatively correlated with the APA gene PABPN1. Such correlations were also observed for the m⁶A genes YTHDC2 and YTHDF3. Moreover, the expression of the m⁶A genes YTHDC2, YTHDF3, and CBLL1 showed a strong positive correlation with the m⁵C genes ZBTB38 and ZBTB33 (Figure 2B). These DERMGs were significantly enriched in mRNA methylation, mRNA polyadenylation, the C-5 methylation of cytosine, and the mRNA surveillance pathway (Figure 2C and D).

Figure 2 RNA modification-related genes in dilated cardiomyopathy. (A) Boxplot showing the expression levels of differentially expressed RNA modification-related genes in dilated cardiomyopathy (DCM) vs Control. (B) Correlation heatmap showing the correlations among the 26 differentially expressed RNA modification-related genes; The significantly enriched biological processes (C) and KEGG pathways (D) for differentially expressed RNA modification-related genes. Compared with CTRL group, 0.01<*P<0.05, 0.005<**P<0.01, ***P<0.005.

RNA Modification-Related Genes-Based Diagnostic Model

A univariate logistic regression analysis was conducted on 26 DERMGs to identify genes associated with DCM, all of which shown substantial correlations with the condition (Figure 3A). LASSO analysis was further conducted to identify the most valuable genes, and an optimal gene signature with 13 DERMGs was ultimately determined (Figure 3B). Among these, the expression of TRMT61B, MBD2, YTHDC2, NOP2, TRMT10C, WDR4, and CPSF2 were found to decrease in DCM, while the expression of CSTF3, ZBTB4, UNG, NSUN6, TET1, and DNMT3B were increased in DCM (Figure 3C). An SVM diagnostic classifier was then developed utilizing these 13 genes. The ROC curve demonstrated that this gene profile exhibited strong predictive performance for DCM, with an area under the curve (AUC) of 0.980 (Figure 3D). The efficacy of our diagnostic classifier was corroborated using an external dataset, GSE120895, which demonstrated commendable prediction ability for DCM, yielding an AUC of 0.902 (Figure 3E–F).

Figure 3 RNA modification-related genes-based diagnostic model. (A) Forest plot showing the 26 differentially expressed RNA modification-related genes significantly associated with dilated cardiomyopathy in univariate Cox regression analysis. (B) LASSO coefficient distribution of RNA modification-related genes and 10-fold cross-validated likelihood deviance of the LASSO coefficient for parameter selection. (C) Heatmap showing the expression pattern of 13 optimal diagnostic genes in GSE141910 dataset. (D) ROC curve showing the diagnostic performance of the diagnostic model. (E) Heatmap showing the expression pattern of 13 optimal diagnostic genes in GSE120895 dataset. (F) ROC curve showing the diagnostic performance of the diagnostic model.

Immune Infiltration Status of Left Ventricular Tissue in DCM Differed from the Control

The state of immunological infiltration was assessed utilizing CIBERSORT. Among 22 immune cell types, naive B cells, plasma B cells, resting memory CD4⁺ T cells, M2 macrophages, and activated mast cells exhibited a comparatively significant infiltration fraction in the left ventricular free-wall tissue. Furthermore, 11 immune cells exhibited notable disparities in infiltration quantity between the DCM and CTRL samples. DCM samples exhibited a significantly elevated presence of naive B cells, M1 macrophages, activated mast cells, and neutrophils in comparison to the CTRL samples (Figure 4A). We subsequently examined the associations of the 13 optimum diagnostic genes with invading immune cells that varied between the DCM and CTRL groups. The expression levels of TRMT61B, MBD2, YTHDC2, and CPSF2 exhibited a robust positive connection with the infiltration abundance of plasma B cells and eosinophils. It remained adversely linked with the amount of infiltrating CD8 T cells and regulatory T cells (Tregs) (Figure 4B).

Figure 4 Associations of RNA modification-related genes with immune infiltration. (A) Boxplot showing the difference in infiltration abundance of 22 immune cells between dilated cardiomyopathy (DCM) and control samples; (B) Correlation heatmap showing the correlations between RNA modification-related genes expression and infiltration abundance of differential immune cells. Compared with CTRL group, 0.01<*P<0.05, 0.005<**P<0.01, ***P<0.005.

RNA Modification-Related Genes-Based Molecule Clusters in DCM

Given the significance of the 13 optimum diagnostic genes in DCM, a consensus clustering analysis was performed based on the expression of these genes, categorizing 166 DCM patients into two molecular clusters: 116 samples in Cluster 1 and 50 samples in Cluster 2 (Figure 5A). Cluster 1 exhibited a markedly elevated RNA modification score compared to Cluster 2 (Figure 5B). The immune infiltration state varied between the two molecular clusters. Cluster 1 exhibited a markedly elevated stromal score (Figure 5C) and a significantly low immune score (Figure 5D) compared to Cluster 2. A total of 14 immune cells showed significant differences between Cluster 1 and 2 (Figure 5E).

Figure 5 Consensus clustering analysis. (A) Consensus clustering grouped dilated cardiomyopathy samples into two clusters based on expression of RNA modification-related genes; Boxplots showing the difference in RNA modification score (B), stromal score (C) and immune score (D) between two clusters; (E) Boxplot showing the difference in infiltration abundance of 22 immune cells between two clusters. Compared with Cluster 1 group, 0.01<*P<0.05, ***P<0.005.

Cluster 1 exhibited a greater infiltrating abundance of naive B cells, plasma B cells, and resting memory CD4⁺ T cells, whereas Cluster 2 demonstrated a larger infiltrating abundance of M2 macrophages and monocytes. Furthermore, these two molecular clusters exhibited distinct expression patterns for HLA family genes (Figure 6A) and several co-stimulatory molecules (Figure 6B). For example, multiple HLA genes (HLA-DRB1, HLA-DPB1, and HLA-C) were highly expressed in Cluster 2 compared to Cluster 1.

Figure 6 Expression of HLA family genes and co-stimulatory molecules. Boxplots showing the difference in expression of HLA family genes (A) and co-stimulatory molecules (B) between two clusters. Compared with Cluster 1 group, 0.01<*P<0.05, 0.005<**P<0.01, ***P<0.005.

Gene Expression Patterns and Pathways Between Two Clusters

To validate the two molecular clusters associated with RNA modification-related genes, their respective gene expression patterns and pathways were subsequently examined. As anticipated, the two clusters of molecules associated with RNA modification exhibited significantly distinct gene expression patterns, encompassing 742 differently expressed genes (Figure 7A and B). These genes were primarily associated with the inflammatory response, ion transmembrane transport mechanisms (such as chloride and sodium ions), neuroactive ligand-receptor interactions, and calcium signaling pathways (Figure 7C and D). Additionally, GSVA identified many pathways that were differentially enriched between the two molecular clusters. Cluster 2 exhibited considerable activation of various metabolism-related pathways, including those for arachidonic acid, pyrimidine, glutathione, glycine, serine, and threonine metabolism. Conversely, primary immunodeficiency, the renin-angiotensin system, and neuroactive ligand-receptor interaction pathways were activated in Cluster 1 (Figure 7E).

Figure 7 Gene expression and pathways between two clusters. Heatmap (A) and volcano plot (B) showing the different expression pattern of genes between two clusters; Bar charts showing the significantly enriched biological process (C) and KEGG pathways (D) for differentially expressed genes; (E) The pathways differentially enriched between two clusters in gene set variation analysis.

Q-PCR Validated the Relationship of 13 RNA Modification-Related Genes’ Levels Between Control and DCM

RT-qPCR results indicated that, in comparison to the CTRL group, the expression levels of six optimal RNA modification-related genes (CSTF3, ZBTB4, UNG, NSUN6, TET1 and DNMT3B) among the 13 optimal genes were considerably elevated in the DCM group (P<0.05). The seven optimum genes (TRMT61B, MBD2, YTHDC2, NOP2, TRMT10C, WDR4 and CPSF2) in the DCM group were significantly downregulated (P<0.05). TET1, DNMT3B, MBD2 and TRMT61B exhibited the most pronounced variations between the two groups (P<0.01) (Figure 8). Consequently, the subsequent murine DCM model was employed to ascertain the correlation between the four genes (TET1, DNMT3B, MBD2, and TRMT61B) and cardiac function in DCM. The validation outcomes affirmed the dependability of this diagnostic paradigm.

Figure 8 Verification of 13 optimal RNA modification-related genes by qRT-PCR in Control and DCM group. Compared with CTRL group, 0.01<*P<0.05, 0.005<**P<0.01.

Experimental Validations Of four Hub Genes Expression in DCM Mice

In the modeling phase, the DCM group exhibited significantly reduced body weight compared to the CTRL group (P<0.05) (Figure 9A). The qPCR analysis of cardiac tissue samples revealed increased mRNA expression levels of TET1 and DNMT3B, alongside decreased mRNA expression levels of MBD2 and TRMT61B (P<0.05) (Figure 9C). Echocardiography demonstrated a reduction in EF% and FS% (P < 0.05) in the DCM group relative to the CTRL group, alongside a significantly elevated LVIDs (P < 0.05) (Figure 9B and D). We investigated the correlation between cardiac function and four hub genes (TET1, DNMT3B, MBD2 and TRMT61B). The number of TET1 PCR cycles correlated positively with EF% (R = 0.8891, P = 0.0031) and FS% (R = 0.8707, P = 0.0007) but negatively with LVIDs (R = −0.8756, P =0.0006). The number of DNMT3B PCR cycles correlated positively with EF% (R = 0.8883, P=0.0032) and FS% (R = 0.8384, P= 0.0093) and negatively with LVIDs (R = −0.7506, P =0.0319). The number of MBD2 PCR cycles correlated negatively with EF% (R = −0.9715, P <0.0001), FS% (R = −0.9405, P = 0.0005), and positively with LVIDs (R = 0.9328, P = 0.0007). The number of TRMT61B PCR cycles correlated negatively with EF% (R = −0.9527, P =0.0003), FS% (R = −0.9664, P < 0.0001), and positively with LVIDs (R = 0.9511, P = 0.0003) (Figure 9E).

Figure 9 The relationship between TET1, DNMT3B, MBD2 and TRMT61B mRNA levels, and DCM in mice model. (A) The body weight in the Control and DCM group. (B and D) The echocardiography features in the Control and TAAD groups. (C) TET1, DNMT3B, MBD2 and TRMT61B mRNA levels in the CTRL and DCM group. (E) The correlations between four optimal hub genes TET1, DNMT3B, MBD2 and TRMT61B mRNA levels, and cardiovascular functional parameters in the CTRL and DCM group. Compared with the CTRL group, 0.01<*P<0.05, 0.005<**P<0.01, ***P<0.005.

Discussion

Despite growing data indicating the role of RNA modifications in various human diseases, pertinent studies on DCM are still absent. This study revealed that several RNA modification-related genes were dysregulated in DCM, including m⁵C regulators. The bulk of the dysregulated genes, including m⁶A regulators, were significantly downregulated in DCM. Research has shown that m⁶A alteration is essential for cardiomyocyte contractile performance and cardiac remodeling, both of which are vital in heart failure resulting from DCM.^17,18 This validated the deregulation of m⁶A regulators in DCM. However, fewer research have reported on the role of m⁵C regulators and other RNA modification regulators in DCM or the cardiac remodeling process. This analysis revealed substantial connections among dysregulated RNA modification-related genes, indicating potential interactions among various RNA modifications.

From these dysregulated RNA modification-related genes, 13 genes were identified as diagnostic biomarkers for DCM, including TRMT61B, MBD2, YTHDC2, NOP2, TRMT10C, WDR4, CPSF2, CSTF3, ZBTB4, UNG, NSUN6, TET1, and DNMT3B. The SVM diagnostic model constructed utilizing these 13 genes demonstrated strong predictive ability for DCM. Mitochondria are essential for structural and functional cardiac remodeling, and their malfunction is associated with heart failure, ventricular hypertrophy, and other cardiac disorders.¹⁹ TRMT61B is mostly localized in the mitochondria and encodes a mitochondria-specific tRNA methyltransferase that catalyzes the formation of m¹A at position 58.²⁰ Furthermore, TRMT10C is a regulator responsible for m1A in mitochondrial tRNA and mitochondrial ND5 mRNA,²¹ and its recessive mutations result in defects in mitochondrial RNA processing.²² DNMT3B is a major methyltransferase expressed in human and mouse hearts, and its cardiac-specific deletion contributes to myocardial thinning and the progression of severe systolic insufficiency.²³ In addition, the increased expression of DNMT3B and decreased expression of MBD2 have been reported to be involved in epigenetic alterations in ischaemic cardiomyopathy.²⁴ The potential function of other genes has not been reported in DCM but in various other diseases, such as the predictive value of CPSF2 in papillary thyroid carcinoma²⁵ and CPSF3 in non-small cell lung cancer,²⁶ as well as YTHDC2 in acute myocardial infarction.²⁷

Analysis of immune infiltration revealed a comparatively significant quantity of naive B cells, M1 macrophages, activated mast cells, and neutrophils in DCM samples. Tang et al discovered that the anomalous distribution of transitional B cells and B1 cells contributed to the etiology of idiopathic DCM. They observed that the proportion of B1 cells was closely linked with DCM severity.²⁸ Cardiac macrophages are known to modulate adaptive myocardial remodeling, and the reduction of CCR2 macrophages in a DCM mouse model has been associated with increased mortality and impaired ventricular remodeling.²⁹ Cardiac mast cells are involved in myocardial fibrosis and are considered the core of poor adverse myocardial remodeling.^30,31 The neutrophil/lymphocyte ratio correlates with the severity of DCM and serves as a predictive diagnostic for the condition.^32,33 Additionally, we identified numerous genes (eg, MBD2) that had substantial correlations with infiltrating immune cells in DCM, corroborating prior research. The ablation of MBD2 inhibited TH17 differentiation.³⁴ In experimental autoimmune myocarditis, MBD2 functioned as a target gene, facilitating pathogenesis and CD4+ T cell immunometabolic dysfunction.³⁵ Cluster 1 exhibited a superior RNA modification score and a worse immune score compared to Cluster 2. Co-stimulatory molecules are important for T-cell activation,³⁶ and the expression of multiple co-stimulatory molecules has been associated with the development of acute myocarditis and direct myocardial damage in DCM.³⁷

We found that various co-stimulatory molecules, including as CD40 and CD27, were dysregulated between the two clusters. Human leukocyte antigens (HLA) are essential genetic markers linked to numerous illnesses, including DCM. The HLA-DR3 antigen has been documented to confer a protective effect in idiopathic DCM,³⁸ while HLA-DRB1 and HLA-DQB1 gene polymorphism is related to susceptibility to DCM and myocarditis.^39,40 We evaluated the expression of HLA genes between the two molecular clusters and discovered that multiple HLA genes (HLA-DRB1 and HLA-DPB1) were substantially elevated in Cluster 2 compared to Cluster 1. Furthermore, disparities in routes were noted between the two clusters. Cluster 2 exhibited considerable activation of various metabolism-related pathways, whereas Cluster 1 demonstrated activation of primary immunodeficiency and renin-angiotensin system pathways. During heart failure, the renin-angiotensin system (RAS) is activated and is recognized as a significant therapeutic target for the condition.⁴¹ Consequently, RAS gene polymorphism has been regarded as a potential modulator of the phenotypes associated with hypertrophic cardiomyopathy and DCM.⁴² Taken together, these findings revealed the significant differences between the two RNA modification- molecular clusters.

In our prior research, we integrated bulk and single-cell RNA sequencing data to discover hub genes of fibroblasts in DCM.⁴³ In the rat DCM model, we also found that the immune microenvironment around cardiac fibroblasts in pathological tissue was disordered, and RNA modification could lead to immune microenvironment disorder. Therefore, in order to further study the pathogenesis of DCM, based on the work, we analyzed the application value of RNA modifications related genes in the diagnosis and subtype classification of DCM in this study, and confirmed that RNA-modification-related genes (TET1, DNMT3B, MBD2 and TRMT61B) were most closely associated with decreased cardiac function in DCM. Nevertheless, our research has some limitations. First, because to the considerable heterogeneity exhibited by DCM, further research necessitating a bigger sample size and clinical data will be required. Secondly, there is a lack of linkage among the subtypes, diagnostic genes, and clinical aspects, necessitating more inquiry. The association between the subtypes of RNA modification-related genes and the progression of DCM, along with their functions and processes in the pathological process of DCM, will be further examined.

Conclusions

This work is the inaugural investigation into the systematic role of RNA modification-related genes in DCM. Alongside the established gene signature for DCM, two molecular subtypes of DCM, based on RNA modification-related genes, were discerned; these subtypes exhibited significant variations across multiple dimensions. This study had numerous drawbacks. Initially, over 170 varieties of chemical modifications were recognized on RNAs, and this analysis encompassed just genes associated with many prevalent RNA modifications. The potential functions of genes associated with other RNA modifications in DCM, such as pseudouridine and 2’-O-methylation, were overlooked, despite being prevalent RNA modifications. Secondly, this was a preliminary exploratory study, and several studies should further corroborate the primary findings. The expression and clinical significance of the 13 critical genes, together with the functions of infiltrating immune cells in DCM, are crucial. Furthermore, we successfully produced two distinct RNA modification-related molecular clusters, which exhibited significant differences across various dimensions, demonstrating the viability of RNA modification-based molecular clusters. Nevertheless, it was necessary to examine whether clinical pathological traits varied between the two groups, given the absence of clinical data in the dataset. Dysregulation in the expression patterns of many RNA modification-related genes was noted, leading to the establishment of a diagnostic model based on 13 critical genes for DCM. Two RNA modification molecular clusters with distinct clinical characteristics were subsequently found for DCM. Collectively, this study provides novel insights into the functions of RNA modification in DCM, which may contribute to advances in DCM diagnosis and treatment.