Construction of co-expression modules of HCM
A total of 14348 mRNAs were included in the WGCNA analysis. According to mRNA expression matrix, we clustered the samples and eliminated the outliers. Among the 28 HCM samples, there is no obvious outlier, which can be used for subsequent WGCNA analysis (Fig. 1A). Clinical and demographic data of 28 HCM patients were included sex, age, smoking, LAD, LVST, LVEDD, LVEF, Maxi LVWT, Maxi LVOTG and sample location (Table I). Based on mRNAs expression data and clinical data, we grouped the sample dendrogram and trait heatmap (Fig. 1B). According to the power value (β) and scale R2 value, the independence and the average connectivity of the co-expression module were defined. We checked whether the memory network approximated scale free under the selected β value. As show in figure 2B, K was negatively correlated with P (K) (R2=0.86), indicating that the selected β value could establish a scale-free gene network. As shown in Fig. 2A, the best value of the power value was 6 and the corresponding scale R2 was 0.9. Therefore, we defined the adjacency matrix with β=6 and constructed the co-expression gene module. There were 31 co-expression modules of co-expressed genes in WGCNA analysis (Fig. 2C). The number of genes in the corresponding module was shown in Table II. The gene network was visualized using a heatmap which depicted the Topological Overlap Matrix (TOM) among top 400 genes in the analysis (Fig. 2D).
Co-expression modules related to clinical traits
GSE130036 provides the following clinical data, including sex, age, smoking, LAD, LVST, LVEDD, LVEF, Maxi LVWT, Maxi LVOTG and sample location. We drew the heatmap of the correlation between gene co-expression module and clinical traits. As shown in Fig. 3, the co-expression black module was significantly associated with Maxi LVWT (maximum left ventricular wall thickness, R =0.61, p=5×10-4). We further plotted the scatter plot of gene significance for Maxi LVWT and module membership in the black module. The high correlation can be revealed that gene significance for Maxi LVWT was highly associated with module membership in black module (cor=0.74, p=2.8×10-101, Fig. 4C). Therefore, we chose the genes in black module to further analysis.
The core genes were further screened by differential expression analysis
There were total of 291 mRNAs expressed differentially between HCM samples and the normal samples (pvalue <0.05 and |logFC|>1.5). Compared with the expression of mRNAs in HCM samples, 239 mRNAs (82.13%) were upregulated in the normal samples, while the 52 mRNAs (17.87%) were downregulated (Table S1). The differences in the expression of mRNAs between the normal samples and HCM samples were hierarchical clustering analysis and visualized by hierarchical heatmaps and volcano plots (Fig. 4A, B). CES1 and CTSC were obtained from the intersection of the genes with different expression and in black module (Fig. 4D). Compared with the expression of mRNA in normal tissues, mRNA expression of CES1 and CTSC was down regulated in HCM tissues (P<0.001, Fig. 5A). To further verify the lower expression of CES1 and CTSC in HCM, we analyzed another RNA sequencing data from GEO (GSE68316). As shown in Fig. 5B, compared with the expression of mRNA in normal tissues, mRNA expression of CES1 and CTSC was also down-regulated in HCM tissues (P<0.001). Pearson correlation coefficient was used to calculate correlation between relative gene expression and Maxi LVWT. As shown in Fig. 5C, mRNA expression of CES1 and CTSC was negatively correlated with Maxi LVWT (cor=-0.48, P<0.001; cor=-0.31, P<0.001). In summary, mRNA expression of CES1 and CTSC was down-regulated in HCM and was negatively correlated with Maxi LVWT.
Clinical sample validation
In this study, 20 peripheral blood samples were obtained from HCM patients or healthy people between November 2020 and May 2021 in the Kunming Yan'an hospital (Kunming, China). We collected clinical information for all patients. The clinical features of the patients including sex, age, LAD, LVST, LVEDD, LVEF, Maxi LVWT and sample location were collected from their medical records (Table III). As shown in Fig. 6A, mRNA expression of CES1 and CTSC in HCM groups was decreased to 60.3±22.07% (P<0.001) and 82.1±25.43% (P<0.001) of the normal groups respectively. We further analyzed that mRNA expression of CES1 and CTSC was negatively correlated with Maxi LVWT in HCM patients (cor=0.54, p<0.001; cor=-0.34, p<0.001, Fig. 6B). In summary, the study confirms the abnormal low expression of CES1 and CTSC in HCM through clinical samples from HCM patients, and their expression was negatively correlated with Maxi LVWT.