Hypoxia-associated genes predicting future risk of myocardial infarction: a GEO database-based study

Background: Patients with unstable angina (UA) are prone to myocardial infarction (MI) after an attack, yet the altered molecular expression profile therein remains unclear. The current work aims to identify the characteristic hypoxia-related genes associated with UA/MI and to develop a predictive model of hypoxia-related genes for the progression of UA to MI.

Methods and results: Gene expression profiles were obtained from the GEO database. Then, differential expression analysis and the WGCNA method were performed to select characteristic genes related to hypoxia. Subsequently, all 10 hypoxia-related genes were screened using the Lasso regression model and a classification model was established. The area under the ROC curve of 1 shows its excellent classification performance and is confirmed on the validation set. In parallel, we construct a nomogram based on these genes, showing the risk of MI in patients with UA. Patients with UA and MI had their immunological status determined using CIBERSORT. These 10 genes were primarily linked to B cells and some inflammatory cells, according to correlation analysis.

Conclusion: Overall, GWAS identified that the CSTF2F UA/MI risk gene promotes atherosclerosis, which provides the basis for the design of innovative cardiovascular drugs by targeting CSTF2F.

Introduction

Coronary Heart Disease (CHD) conteins myocardial infarction (MI), angina pectoris, and coronary arteriosclerosis (1, 2). The incidence and mortality of CHD in developing countries, including China, are increasing year by year (3, 4). MI is one of the most serious manifestations of coronary artery disease (CAD) and the leading cause of death from non-infectious diseases worldwide (5). Atherosclerosis has been recognized as the main underlying mechanism of coronary heart disease, ultimately leading to acute coronary syndrome (ACS). ACS includes a clinical spectrum ranging from unstable angina (UA) to acute MI (6). Currently, the treatment of UA is mainly conservative drug treatment, such as anti-platelet aggregation, coronary artery dilation, lipid reduction, etc., to prevent the progress of UA, thus reducing the risk of MI (7, 8). However, the molecular mechanism from UA to MI is still poorly understood.

Ischemic/hypoxic injury is a common pathological basis for CHD, atherosclerosis, and other cardiovascular diseases. Coronary atherosclerosis can lead to narrowing of the vessel lumen, and angina pectoris caused by this rapid, transient myocardial ischemia and hypoxia is a common symptom in the chest (9), and myocardium with underdeveloped collateral circulation is susceptible to ischemic necrosis, leading to MI (10, 11). Coronary atherosclerosis is the main pathological basis of coronary heart disease, of which atherosclerotic plaque consists mainly of cholesterol and cholesteryl esters, which are highly brittle and prone to rupture or dislodge and flow with the blood, blocking small vessels or capillaries, forming thrombi, and aggravating MI (12, 13). Atherosclerotic plaque rupture has been identified as a major cause of ACS (14). Nevertheless, non-atherosclerotic factors have been identified to cause ACS (15). Therefore, further comparison of alterations in hypoxia genes in UA progression to MI must be done with association analysis of atherosclerosis-related protective factors and suppressors to determine the potential role of these genes in atherosclerosis-mediated or non-atherosclerosis-mediated ACS. Thus hypoxia-related genes may be potential targets for CAD disease progression and treatment.

The clinical occurrence of MI is often a sudden symptom. Early prevention and intervention are difficult to achieve. However, early MI has lesions such as coronary arteries, blood composition, and hormone levels. Therefore, MI can be prevented at an early stage by molecular level changes as markers (16). Recent studies in transcriptomics have been continuously applied in clinical practice, especially in routine diagnostic applications. Microarray data and RNAseq data are widely used for disease diagnosis. Numerous studies have confirmed the abnormal expression levels of some disease-related genes in the early stages of MI. Zhang et al. have then been in the process of identifying the RNA signature of coronary heart disease from the combined expression profile of lncRNA and Mrna (17). A study folded into diagnostic features byMI-related differentially expressed genes (18). We, therefore, sought to identify hypoxia-associated genes for molecular characterization from UA to MI while trying to find their association with atherosclerosis in anticipation of finding available therapeutic targets.

Materials and method

Data collection and organization

We started from GPL571, GPL6106-based platform (GSE48060, GSE61144, GSE29111, GSE97320, GSE34781) Gene Expression Omnibus (GEO) database (https://www.ncbi.nlm.nih.gov/geo/) to obtain a total of 93 samples (18 UA cases and 75 IM cases) for training set analysis. A total of 24 samples (8 UA cases and 16 MI cases) were also obtained from the GPL6884-based platform (GSE60993) for test set validation. The microarray data were quantile normalized to obtain standardization (19). The microarray probes were converted to 12,285 Ensembl gene IDs using the biomaRt software package. According to the standardized sample network connectivity Z-score < −2, the samples defined as outliers are removed (20). The batch effects of the expression profiles are processed using the operational function of the SVA package in R (21). All data processing results are shown in Supplementary Figure S1. Furthermore, total of 25 normal samples from GSE34781, GSE48060, and GSE97320 were used to compare the expression levels of the core genes in the normal, UA and MI groups. Detailed sample sizes for the GEO dataset are provided in Supplementary Table S1 and clinical information is provided in Supplementary Table S2.

Identification of genes related to hypoxia

In the training data, UA was compared with MI to identify genes with significant expression differences. We use the “limma” package R to extract differentially expressed genes (DEGs). After correction based on the false discovery rate, the cutoff value of the P value was set to 0.05 to obtain 3,569 differential genes considered as disease progression genes. Subsequently, we used “ssGSEA” to quantify the Hallmark gene set obtained from MSigDB (gsea-msigdb.org/gsea/msigdb/). Total of 4,017 differential genes related to hypoxic were identified by differential analysis. In parallel, we built a weighted gene co-expression network analysis (WGCNA) and identified hypoxia-associated modules by package “WGCNA” (22). Interactions between unique genes with hypoxic ssGSEA (23) scores were quantified by gene significance (GS), and correlations between gene expression profiles and module signature genes were indicated by module membership (MM). At a threshold of GS p < 0.05 threshold, 630 candidate genes from the “sapphire-colored module” were selected. Differential genes were analyzed by “ClusterProfiler” R (24) package to perform GO, KEGG enrichment.

Classification model construction and validation

The ten genes with the highest coefficients in the LASSO regression model were selected for the construction of the model by selecting 630 modular genes for the R package “glmet” (25). The hypoxic risk score was calculated using the corresponding coefficients of the selected features. The formula was established as follows.

$\mathrm{Risk}\mathrm{score}=\sum _{\mathrm{i}=1}^{\mathrm{n}}\mathrm{Coe}{\mathrm{f}}_{\mathrm{i}}\ast {\mathrm{x}}_{\mathrm{i}}$

where $\mathrm{Coe}{\mathrm{f}}_{\mathrm{i}}$ denotes the coefficient and ${\mathrm{x}}_{\mathrm{i}}$ denotes the expression value of each mRNA. Patients were divided into low- and high-risk groups according to the median value of risk value, and PCA downscaling analysis showed its distribution pattern of patients with different risks. ROC curves were used to assess model accuracy. A validation set was used to verify the above results.

Estimation of immune cell infiltration

The CIBERSORT (26) algorithm was used to translate normalized gene expression values into the composition of immune cells in complex samples. CIBERSORT generates an empirical P-value for the inverse convolution results using Monte Carlo sampling, which is used to compare the statistical significance of the inverse convolution results on all subsets of cells. The training set of 93 samples was included in the CIBERSORT analysis. The entire sample was included in the follow-up analysis (all P-values less than 0.05). The specific results are given in Supplementary Table S3. The ratio of immune cells in significantly enriched UA and MI samples was obtained and reported as a bar graph. We used the R package “vioplot” to compare the levels of each immune cell type. Spearman's rank correlation analysis was performed to explore the correlation between infiltrating immune cell subsets and signature genes. The R package “ggplot2” is used to visualize plots.

eQTL mapping in GTEx

SNP loci associated with UA/MI disease were obtained from the GWAS catalog database. The public eQTL database was queried in the Genotype-Tissue Expression (GTEx) project (Release V8) to determine the association between lead UA/MI GWAS variants and CSTF2T gene expression. GTEx tissue-specific eQTL were also identified for the GWAS variants using 48 different GTEx tissues sorted by NES value and eQTL p-value.

Western blot

Whole blood samples were obtained from 3 cases of UA and 3 cases of MI patients, and 3 cases of healthy volunteers in Ruijin Hospital, and kept in EDTA evacuated tubes. Subsequently, peripheral blood monocytes (PBMCs) were isolated using Ficoll-Hypaque density gradient centrifugation (TBD Science, Tianjin, China). Total protein was isolated from the isolated PBMCs in ice-cold lysis buffer using a high-intensity ultrasomic processor, according to the manufacturer's instructions. The remaining debris were removed by centrifugation, and the protein was precipitated with cold 15% Trichloroacetic for 2 h at −20 °C. After centrifugation at 20,000 × g at 4 °C for 10 min, the supernatant was discarded. The protein was dissolved in buffer, and the protein concentration was determined with a 2-D Quant kit according to the manufacturer's instructions. Equal amount of protein samples were subjected to 10% SDS-PAGE and transferred to PVDF membranes. After blocking at 37 °C for 2 h, the membranes were probed with primary antibodies at 4 °C overnight. After washing with TBST, the blots were incubated with a secondary antibody at 37 °C for 1 h. Immunoreactivity was visualized by enhanced chemiluminescence. The intensity of the bands was normalized to that of GAPDH. The statistical significance of changes in protein expression observed by Western blotting was assessed with a two-tailed t-test using GrapghPad Prism 8.0.2. The data are expressed as the mean ± standard deviation, and statistical significance was considered to be indicated by a p-value < 0.05.

Statistical analysis

All statistical analyses were performed using R software (version 4.0.3). The Wilcoxon test was used to compare the differences between the two groups. The Spearman's rank correlation analysis was used to explore the correlation between model genes and immune cells. Values of p < 0.05 were considered statistically significant.

Result

Hypoxia-related gene basis of coronary heart disease

The workflow is shown in Figure 1. Many patients with MI have rapidly progressing disease due to UA (27). We divided the patients into two groups based on two different clinical phenotypes, including 18 UA patients and 75 MI samples. Using the limma algorithm, we finally screened 3,569 DEGs between UA and MI groups based on a corrected P value < 0.05 (Figure 2A). Atherosclerosis exacerbated by ischemia/hypoxia eventually led to stable and UA and MI, and the differently enrichment of pathways between the two groups is shown in the heatmap, we then calculated the hypoxia scores of UA and MI patients according to ssGSEA (Figure 2B), and the difference analysis between high-hypoxia and low-hypoxia groups obtained 4,017 DEGs displayed in the volcano plot (Figure 2C). The 1966 intersecting genes were demonstrated by VENN plots (Figure 2D). Since our results indeed show that hypoxia is a prominent molecular feature of UA, these intersecting genes are regarded to be the molecular basis of hypoxia-related genes in UA to MI disease progression. We then performed GO and KEGG on these genes, which were highly enriched in propanoate metabolism, chronic myeloid leukemia, cell cycle, mitochondrial ribosome, H4 histone acetyltransferase complex, and cellular senescence pathways (Figure 2E), which are all related to oxygen metabolism/hypoxia (28–33).

FIGURE 1

Figure 1. The workflow of the integrative bioinformatics analyses.

FIGURE 2

Figure 2. Difference analysis showed that hypoxia was associated with the development of coronary heart disease. (A) The volcano map shows the different results between UA and MI. (B) Hypoxia ssGSEA scores were estimated in the UA and MI cohort by performing ssGSEA with a hallmark gene set. (C) The volcano map shows the different results between high-hypoxia scores and low-hypoxia scores. (D) Venn diagram shows that disease differential genes overlap with hypoxia differential genes. (E) Enrichment analysis of disease-differential genes. UA: unstable angina; MI myocardial infarction; BP: Biological Process; CC: Cellular Component; MF: Molecular Function.

Identification of hypoxia-related genes by weighted gene coexpression network analysis

To further screen for reliable genes associated with hypoxia. Based on the ssGSEA method and cancer markers from the MsigDB dataset, the hypoxia ssGSEA Zscore was calculated for the training set of patients in our dataset, and candidates associated with hypoxia were screened by WGCNA. The optimal soft threshold was determined to be 21 (Figure 3A), Dynamic Tree Cut was set to 0.2, and 7 modules were created (Figure 3B), with the turquoise-colored module having the highest correlation with UA/MI (Figure 3C), while we compared the correlation of hypoxia scores with the modules. The turquoise module possessed the highest correlation with the hypoxic phenotype (Figure 3D). From the turquoise module, 630 promising candidates were identified (Figure 3E). We performed GO, and KEGG analysis on these genes, and not surprisingly these genes were associated with Thyroid hormone synthesis, Ubiquitin mediated proteolysis, Nucleotide excision repair, Cellular senescence, and histone acetyltransferase complex (Figure 3F), which are related to oxygen metabolism/hypoxia (34–38).

FIGURE 3

Figure 3. The weight gene co-expression network analysis showed that hypoxia was associated with the development of coronary heart disease. (A) Soft threshold selection process. (B) Cluster dendrogram. Each color represents one specific co-expression module. In the colored rows below the dendrogram, the two colored rows represent the original modules and the merged modules. (C-D) The differential expression of eigengenes in UA and MI, high hypoxia and low hypoxia (FDR corrected *P < 0.05, **P < 0.01, ***P < 0.001), respectively. (E) Turquoise gene significance and membership in hypoxia network. (F) Gene ontology and KEGG pathways enrichment in turquoise modules.

Classification model construction and validation

Based on the lasso algorithm, the top 10 genes with the highest coefficients (COMMD2, CSTF2T, LNPEP, MLYCD, NEK4, NOM1, PARG, PPP3CB, PSMD5, and SKP2) were selected (Figures 4A–B). PCA analysis showed that the samples were divided into two clusters (Figure 4C). The area under the ROC curve in the training set was 0.998 (Figure 4D), demonstrating its excellent classification ability, and these results were confirmed in our test set (Figures 4E–F). We found that age was positively associated with risk score in GSE29111 (R = 0.29; P < 0.05), but not gender (Supplementary Figure S2). We also analyzed the classification ability of individual genes (Supplementary Figure S3). Furthermore, we compared Pearson's correlations between these genes and risk scores, and most of them showed high correlations, with NOM1, and LNPEP showing a strong correlation (|R|>0.6, p < 0.05) (Figure 4G). The high expression of most hypoxia-related model genes represents a low risk of MI incidence in these patients (p < 0.05) (Figure 4H).

FIGURE 4

Figure 4. Construction of disease classification model. (A-B) The lasso regression top 10 genes in training data. (C) PCA analysis in train data and (F) test data. ROC analysis in train (D) data and (E) test data. (G) Correlation analysis of risk scores and model genes. (H) Model gene expression in high and low-risk groups.

Construction of the nomogram model

The nomogram is a simple, personalized visualization tool that has been widely used in disease diagnosis and prognosis (39). We used the “rms” package in R to construct a nomogram model based on 10 hypoxia-related genes to predict the prevalence of MI in UA (Figure 5A). The column line plot c-index value is 1. The calibration curve shows that the predicted results based on the column line graph are in good agreement with the actual prognostic results (Figure 5B).

FIGURE 5

Figure 5. Construction of the nomogram. (A) Comprehensive nomogram of 10 risk genes. (B) Calibration curve.

Generation of hypoxia gene signatures

To further validate the precision of the hypoxia pattern, patients with UA/MI were classified into different genomic subtypes using a consensus clustering approach based on 10 hypoxia-associated model genes. We found the existence of two different hypoxia gene patterns (gene cluster 1 and gene cluster 2) (Figure 6A). The Sankey diagrams visualized the relationship between predicted UA and MI patients based on hypoxia-related genes and actual patients (Figure 6B). In addition, we compared the expression of these hypoxic genes between predicted UA and MI patients (Figure 6C), and the results were similar to those in the risk model, except for COMMD2, PPP3CB, and PSMB5, which were all downregulated in MI patients (p < 0.05). This again validates the accuracy of our grouping by the consensus clustering method. To further explore the relationship between these ten genes and the progression of UA/MI disease. We found that the mRNA expression of COMMD2, CSTF2T, LNPEP, MLYCD, NOM1, PARG, PSMD5, and SKP2 was not significant between the normal and UA groups, but their expressions were different between the UA and MI groups (p < 0.05) (Figure 6D). This suggests that these genes change as the disease progresses.

FIGURE 6

Figure 6. Consensus clustering of 10 important hypoxia-related genes in UA/MI. (A) Consensus matrices for k = 2. (B) Sankey diagram showing the relationship between UA, MI, and UAp, MIp modes. (C) Expression of 10 hypoxia-related genes in UAp, MIp modes. (D) Expression of 10 hypoxia-related genes among CN, UA, and MI groups.

Correlation of gene expression and immune cell infiltration

Atherosclerosis is a chronic inflammatory disease. It is characterized by complex immune interactions between resident vascular cells and specialized immune cells (40). For example, type 2 macrophages are also associated with repair and reconstruction at sites of inflammation (41). Therefore, we sought to find differences in the immune environment between UA and MI. The CIBERSORT algorithm was performed to assess the ratio of immune cell subsets in the training set. Total of 22 immune cell subpopulations are shown in Figure 7A. We found that the relative percentage of T.cells.CD4.naive was higher in the MI group compared to UA, while the relative percentages of CD4 memory resting, T cells gamma delta and relative percentages were significantly less (p < 0.05) (Figure 7B). In addition, we compared high and low-risk groups with immune cells and obtained the same results as disease groups (Figure 7C). Furthermore, we performed a correlation analysis of risk scores and immune cells. As shown in Supplementary Figure S4, among the above immune cells with differences, except for macrophage M2, which does not correlate with a risk score, the remaining 4 immune cells include T cells CD4 naive, T cells CD4 memory resting, T cells gamma delta and Eosinophils were all correlated with risk score (|R|>0.2; P < 0.05). We then did a correlation analysis of immune cells and 10 hypoxia-related genes (Figure 7D). An interesting phenomenon caught our attention, all hypoxia genes showed a negative correlation with memory B cells (|R|>0.2; p < 0.05), and we found that macrophage M2 showed a negative correlation with CSTF2F and PARG (p < 0.05), and CD8+ T cells MLYCD, CSTF2F, PSMD5 and PPP3CB showed a positive correlation (p < 0.05). Evidence for an important role of B lymphocytes in human CVD is limited. In patients with acute myocardial infarction (MI), high levels of the B-cell specific cytokines and B-cell activating factor, predict increased risk of death and recurrent MI (42). B cells are important in immune homeostasis.

FIGURE 7

Figure 7. Immune cell subsets analysis in UA/MI. (A) The proportion of immune cell subsets in UA/MI samples. (B) Differences in analysis of immune cell infiltration between UA and MI patients. (C) Differences in analysis of immune cell infiltration between high- and low-risk groups. (D) 10 hypoxia-related genes associated immune cells. NES: normalized enrichment score.

These findings indicated that the hypoxia-related model composed of 10 genes may affect the infiltration and proportion of immune cells by driving the hypoxic state of immune cells, thereby regulating the immune homeostasis of UA/MI patients.

UA/Mi-associated risk variants in the CSTF2T locus are associated with reduced CSTF2T gene expression

To further identify key genes driving UA/MI development, we intersected the genes in the model with disease risk genes and identified 7 genes that were strongly associated with disease development (Figure 8A). Notably, these 7 genes were negatively associated with the expression of atherogenic genes (CCL2, BMP4, and SELE) and positively associated with multiple atherosclerosis-protective genes (KLF2 and GCH1) (Figure 8B). Further western blot validation was performed, and the results indicated that the expression levels of seven genes (NEK4, MLYCD, LNPEP, NOM1, COMMD2, CSTF2T, PARG) related to disease and hypoxia were significantly reduced in MI patients compared with healthy volunteers, but their expression were not significant between the healthy and UA samples (Figures 8C–D). The 902 SNPs associated with UA/MI were identified by mining the GWAS catalog database to take intersections with these 7 gene SNP loci. We found that only the rs2879627 variant in the CSTF2T locus was associated with UA/MI occurrence (Figure 8E) and determined the mutation type of CSTF2T expression (Figure 8F). By mining genotype and transcriptome data in the Genotype-Tissue Expression (GTEx) database, we mapped the main GWAS risk variants associated with CSTF2 expression (rs2879627 locus) in human tissues and determined that the rs2879627 mutation was associated with CSTF2 expression in whole blood (NES < −0.1; p < 0.05) (Figure 8G).

FIGURE 8

Figure 8. Correlation between mutations at disease-associated hypoxia gene loci and their expression. (A) Venn diagram among model genes, hypoxia differential genes, and differential genes between UA and MI. (B) Association of seven disease-related hypoxia genes with atherosclerosis. (C) Western blot validation for the relative protein expression levels of seven disease-related hypoxia genes (NEK4, MLYCD, LNPEP, NOM1, COMMD2, CSTF2T, PARG) among healthy volunteers, UA and MI patients. (D) Statistical results of protein bands in C. (E) VENN plot of CAD GWAS and eQTL analysis results of seven disease-related hypoxia genes. (F) Boxplots showing expression quantitative trait locus data from the GTEx dataset for a genome-wide association study of the SNP (rs2879627) association. (G) Normalized effect sizes for single-tissue expression quantitative trait loci normalized effect sizes for rs2879627 are shown for 48 different tissues in GTEx.

Discussion

One of the major diseases that endanger human health is CAD. The most common reason for death in people with cardiovascular disease is MI, in particular (43). Therefore, rapid implementation of successful treatment requires an early, quick, and correct diagnosis and assessment of MI (44). We report for the first time that blood markers of MI from UA can be identified and a MI warning issued.

ACS is the most serious medical emergency caused by myocardial ischemia/hypoxia following coronary atherosclerotic plaque formation, with clinical manifestations of UA and acute MI (45). Hypoxia is an inevitable manifestation of atherosclerosis due to vessel wall thickening (46). Given this, we investigated the potential value of hypoxic gene signatures in this study. The hypoxic status of UA/MI was assessed by using ssGSEA and WGCNA methods, respectively, to select hypoxia-associated genes. Our discovery that differential genes based on hypoxia scores were largely the same as those of UA and MI, especially the conversion process to MI, demonstrates that hypoxia is intimately related to UA/MI. Regarding these frequent genes, GO and KEGG analyses revealed that they were enriched in propanoate metabolism, chronic myeloid leukemia, cell cycle, mitochondrial ribosome, H4 histone acetyltransferase complex, and cellular senescence pathways, which are all related to oxygen metabolism/hypoxia (28–33). These hypoxia-related pathways and functions were strengthened. In UA/MI patients, this guarantees the exclusivity and specificity of our identified gene signature. Then, to reflect the categorical gene profile of patients going from UA to MI, our WGCNA technique identified the most relevant modules for hypoxic characteristics and chose the top 10 genes with the highest weights by the LASSO regression model. The categorization algorithm was applied to forecast the risk that UA patients would advance to MI. PCA analysis showed that the classification model had different distribution patterns. The area under the ROC curve for the training set was 1 and the area under the ROC curve for the test set was 0.758, indicating that the classification model has satisfactory performance. In addition, we constructed a nomogram model based on 10 candidate hypoxia-related genes and the calibration curves showed that the decisions based on the nomogram model were accurate. We found that most of the genes in the model were significantly associated with risk scores and high and low-risk groupings in our study, and two disease patterns (cluster1 and cluster2) were identified based on 10 significant hypoxia-related genes using consensus clustering methods. These two patterns were reflected as UA prediction (UAp) and MI prediction (MIp). The Sankey diagram shows that the UAp and MIp determined based on the consensus clustering of these 10 genes are highly consistent with the actual UA, and MI patients. And the UAp and MIp groupings had the same expression levels of the 10 hypoxic genes as those of high- and low-risk patients. Thus, these results confirm the reliability of the 10 hypoxia-gene signatures for differentiating UA and MI patients.

Inflammation is a key factor in the development and progression of atherosclerosis (47). Myocardial infarction has always been considered an inflammatory disease, and the occurrence, development, and prognosis of MI are closely related to the overexpression of inflammatory cytokines (48, 49). Therefore, we tried to find the consensus and uniqueness of the immune signature of UA and MI. The percentage of immune cells in patients with UA and MI was calculated using the CIBERSORT algorithm. T cells CD4.naive was found to be higher in the MI group compared to UA, while the relative percentages of CD4 memory resting, T cell gamma delta, and macrophage M2 were significantly lower in the MI group. Recent studies have shown that the infiltrating macrophages polarization affects the expression of pro- and anti-inflammatory cytokines in the epicardial adipose tissue from UA/MI patients and that the ratio of M1/M2 macrophages is positively correlated with the severity of CAD (50). This is consistent with the phenomenon of reduced macrophage M2 in patients with MI found in our study. In addition, we found that all model genes were highly correlated with B cells. Although B cells were not previously reported to be the predominant immune cell type found in atherosclerotic lesions (51). However, they are abundant in perivascular adipose tissue, which, in addition to the spleen and bone marrow, serves as a niche for immunoglobulin (Ig) production (52). Cell production of cytokines and Ig by B-cell production of cytokines and Ig at these sites is thought to be an important regulator of inflammation in atherosclerotic lesion formation (53). B-cell depletion reduced the development of atherosclerosis in mice (53). At the same time, we learned that the hypoxia and hypoxia-inducible factor (HIF) signaling pathways are critical for B cell development and function, such as survival, proliferation, and cytokine production of B cells (54), and that improper regulation of this B cell can lead to a variety of diseases, including atherosclerosis. In patients with acute MI, high levels of the B-cell specific cytokines and B-cell activating factor, predict increased risk of death and recurrent MI (42). B cells are important in immune homeostasis. Therefore, our findings indicated that the hypoxia-related model composed of 10 hypoxia-associated genes may affect the infiltration and proportion of immune cells by driving the hypoxic state of B-cell activity, thereby regulating the immune homeostasis of UA/MI patients.

The pathological underpinning for CAD is atherosclerosis, and by crossing model genes with disease-associated hypoxia genes, we were able to identify 7 shared genes. These 7 genes were shown to be positively correlated with several atherosclerosis preventive genes and negatively correlated with the expression of atherogenic genes (CCL2, BMP4, and SELE), according to correlation analysis (KLF2 and GCH1). A large-scale genome-wide association study (GWAS) is currently being conducted has found hundreds of potential loci linked to the risk of UA and MI risk (55). However, the genetic mechanisms and causative genes of the UA/MI locus have not been fully resolved, so we sought to find correlations between UA/MI causative genes and our findings. We integrated GWAS and whole-tissue Expression Quantitative Trait Loci (eQTL) analysis to explore hypoxia-related gene expression and its disease SNP association.

Strengths and limitations

Strengths of this study focused on that we found that CAD GWAS risk variants at the CSTF2T locus were closely associated with increased CSTF2T expression in whole blood tissue, suggesting that CSTF2T may promote atherosclerosis, which may be a future direction for the treatment of atherosclerosis. However, the present study had the following limitations. First, the study was a single-centered study; therefore, the results cannot be generalized to other populations with varying demographics. Second, the less number of patients and further research is required to determine whether CSTF2T could be a better indicator for predicting the incidence of MI in UA patients. Despite these meaningful findings, we had to face some limitations. Our study is based on publicly available data. Therefore, more external cohorts are needed to validate our findings. Again, due to the difficulty of obtaining clinical information, we were unable to describe the association of risk scores with more clinical features including reinfarction, the presence of atherosclerosis.

Conclusion

In conclusion, the current study provides insight into hypoxia-related genes and establishes 10 hypoxia-related gene signatures to predict the incidence of MI in patients with UA. At the same time, a nomogram was constructed based on these genes, showing the risk of MI in patients with UA. Immune subsets analysis showed that these ten genes were mainly associated with B cells and some inflammatory cells. Furthermore, the UA/MI risk gene CSTF2F identified by GWAS promotes atherosclerosis, which provides a rationale for designing innovative cardiovascular drugs by targeting CSTF2F.

Data availability statement

The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

Ethics statement

The human study was approved by the Ethics Committee of Ruijin Hospital. Written informed consent was obtained from all participants. The patients/participants provided their written informed consent to participate in this study.

Author contributions

Data analysis of data was contributed by SL and JZ. SL wrote the first draft. Revision for important intellectual content was done by JZ, JN, and JC. Supervision was contributed by JZ, JN, and JC. All authors contributed to the article and approved the submitted version.

Funding

This work was supported by the National Natural Science Foundation of China (81570226), which was collected by Jiumei Cao.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher's note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fcvm.2023.1068782/full#supplementary-material.

References

1. Gaziano TA. Cardiovascular disease in the developing world and its cost-effective management. Circulation. (2005) 112(23):3547–53. doi: 10.1161/CIRCULATIONAHA.105.591792

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Gaziano TA. Reducing the growing burden of cardiovascular disease in the developing world. Health Aff (Millwood). (2007) 26(1):13–24. doi: 10.1377/hlthaff.26.1.13

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Hu YH, Pan XR, Liu PA, Li GW, Howard BV, Bennett PH. Coronary heart disease and diabetic retinopathy in newly diagnosed diabetes in Da qing, China: the Da qing IGT and diabetes study. Acta Diabetol. (1991) 28(2):169–73. doi: 10.1007/BF00579721

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Xu Y, Mintz GS, Tam A, McPherson JA, Iniguez A, Fajadet J, et al. Prevalence, distribution, predictors, and outcomes of patients with calcified nodules in native coronary arteries: a 3-vessel intravascular ultrasound analysis from providing regional observations to study predictors of events in the coronary tree (PROSPECT). Circulation. (2012) 126(5):537–45. doi: 10.1161/CIRCULATIONAHA.111.055004

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Mathers CD, Loncar D. Projections of global mortality and burden of disease from 2002 to 2030. PLoS Med. (2006) 3(11):e442. doi: 10.1371/journal.pmed.0030442

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Holvoet P, Klocke B, Vanhaverbeke M, Menten R, Sinnaeve P, Raitoharju E, et al. RNA-Sequencing reveals that STRN, ZNF484 and WNK1 add to the value of mitochondrial MT-COI and COX10 as markers of unstable coronary artery disease. PLoS One. (2019) 14(12):e0225621. doi: 10.1371/journal.pone.0225621

PubMed Abstract | CrossRef Full Text | Google Scholar

7. Fitzgerald DJ, Roy L, Catella F, FitzGerald GA. Platelet activation in unstable coronary disease. N Engl J Med. (1986) 315(16):983–9. doi: 10.1056/NEJM198610163151602

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Fitzgerald DJ. Platelet activation in the pathogenesis of unstable angina: importance in determining the response to plasminogen activators. Am J Cardiol. (1991) 68(7):51B–7B. doi: 10.1016/0002-9149(91)90384-w

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Gibbons RJ, Abrams J, Chatterjee K, Daley J, Deedwania PC, Douglas JS, et al. ACC/AHA 2002 guideline update for the management of patients with chronic stable angina–summary article: a report of the American college of cardiology/American heart association task force on practice guidelines (committee on the management of patients with chronic stable angina). J Am Coll Cardiol. (2003) 41(1):159–68. doi: 10.1016/s0735-1097(02)02848-6

PubMed Abstract | CrossRef Full Text | Google Scholar

10. Vasilieva EJ, Shpector AV, Raskuragev AB, Lekochmacher SS, Bespalko IA. Platelet function and plasma lipid levels in patients with stable and unstable angina pectoris. Am J Cardiol. (1991) 68(9):959–61. doi: 10.1016/0002-9149(91)90418-k

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Chakhtoura EY, Shamoon FE, Haft JI, Obiedzinski GR, Cohen AJ, Watson RM. Comparison of platelet activation in unstable and stable angina pectoris and correlation with coronary angiographic findings. Am J Cardiol. (2000) 86(8):835–9. doi: 10.1016/s0002-9149(00)01102-4

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Daugherty A, Tabas I, Rader DJ. Accelerating the pace of atherosclerosis research. Arterioscler Thromb Vasc Biol. (2015) 35(1):11–2. doi: 10.1161/ATVBAHA.114.304833

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Hokstad I, Deyab G, Wang Fagerland M, Lyberg T, Hjeltnes G, Forre O, et al. Tumor necrosis factor inhibitors are associated with reduced complement activation in spondylarthropathies: an observational study. PLoS One. (2019) 14(7):e0220079. doi: 10.1371/journal.pone.0220079

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Rioufol G, Finet G, Ginon I, Andre-Fouet X, Rossi R, Vialle E, et al. Multiple atherosclerotic plaque rupture in acute coronary syndrome: a three-vessel intravascular ultrasound study. Circulation. (2002) 106(7):804–8. doi: 10.1161/01.cir.0000025609.13806.31

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Waterbury TM, Tarantini G, Vogel B, Mehran R, Gersh BJ, Gulati R. Non-atherosclerotic causes of acute coronary syndromes. Nat Rev Cardiol. (2020) 17(4):229–41. doi: 10.1038/s41569-019-0273-3

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Song X, Zheng Y, Xue W, Li L, Shen Z, Ding X, et al. Identification of risk genes related to myocardial infarction and the construction of early SVM diagnostic model. Int J Cardiol. (2021) 328:182–90. doi: 10.1016/j.ijcard.2020.12.007

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Zhang YH, Pan X, Zeng T, Chen L, Huang T, Cai YD. Identifying the RNA signatures of coronary artery disease from combined lncRNA and mRNA expression profiles. Genomics. (2020) 112(6):4945–58. doi: 10.1016/j.ygeno.2020.09.016

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Osmak G, Baulina N, Koshkin P, Favorova O. Collapsing the list of myocardial infarction-related differentially expressed genes into a diagnostic signature. J Transl Med. (2020) 18(1):231. doi: 10.1186/s12967-020-02400-1

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Kristinsson VY. “The effect of normalization methods on the identification of differentially expressed genes in microarray data”. Högskolan i Skövde) (2007).

20. Hansen KD, Irizarry RA, Wu Z. Removing technical variability in RNA-Seq data using conditional quantile normalization. Biostatistics. (2012) 13(2):204–16. doi: 10.1093/biostatistics/kxr054

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Leek JT, Johnson WE, Parker HS, Jaffe AE, Storey JD. The sva package for removing batch effects and other unwanted variation in high-throughput experiments. Bioinformatics. (2012) 28(6):882–3. doi: 10.1093/bioinformatics/bts034

PubMed Abstract | CrossRef Full Text | Google Scholar

22. Langfelder P, Horvath S. WGCNA: an R package for weighted correlation network analysis. BMC Bioinformatics. (2008) 9:559. doi: 10.1186/1471-2105-9-559

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Barbie DA, Tamayo P, Boehm JS, Kim SY, Moody SE, Dunn IF, et al. Systematic RNA interference reveals that oncogenic KRAS-driven cancers require TBK1. Nature. (2009) 462(7269):108–12. doi: 10.1038/nature08460

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Wu T, Hu E, Xu S, Chen M, Guo P, Dai Z, et al. Clusterprofiler 4.0: a universal enrichment tool for interpreting omics data. Innovation (Camb). (2021) 2(3):100141. doi: 10.1016/j.xinn.2021.100141

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Friedman J, Hastie T, Tibshirani R. Regularization paths for generalized linear models via coordinate descent. J Stat Softw. (2010) 33(1):1–22. doi: 10.18637/jss.v033.i01

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Newman AM, Liu CL, Green MR, Gentles AJ, Feng W, Xu Y, et al. Robust enumeration of cell subsets from tissue expression profiles. Nat Methods. (2015) 12(5):453–7. doi: 10.1038/nmeth.3337

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Sheng ZQ, Li YF, Lin G, Wang Y, Lu HH. Assessment of short-term prognosis by sinus heart rate turbulence in patients with unstable angina. Exp Ther Med. (2013) 5(4):1153–6. doi: 10.3892/etm.2013.953

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Ai X, Parthun MR. The nuclear Hat1p/Hat2p complex: a molecular link between type B histone acetyltransferases and chromatin assembly. Mol Cell. (2004) 14(2):195–205. doi: 10.1016/s1097-2765(04)00184-4

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Cheloni G, Tanturli M, Tusa I, Ho DeSouza N, Shan Y, Gozzini A, et al. Targeting chronic myeloid leukemia stem cells with the hypoxia-inducible factor inhibitor acriflavine. Blood. (2017) 130(5):655–65. doi: 10.1182/blood-2016-10-745588

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Nuñez-Hernandez DM, Felix-Portillo M, Peregrino-Uriarte AB, Yepiz-Plascencia G. Cell cycle regulation and apoptosis mediated by p53 in response to hypoxia in hepatopancreas of the white shrimp litopenaeus vannamei. Chemosphere. (2018) 190:253–9. doi: 10.1016/j.chemosphere.2017.09.131

CrossRef Full Text | Google Scholar

31. Li X, Wang M, Li S, Chen Y, Wang M, Wu Z, et al. HIF-1-induced mitochondrial ribosome protein L52: a mechanism for breast cancer cellular adaptation and metastatic initiation in response to hypoxia. Theranostics. (2021) 11(15):7337–59. doi: 10.7150/thno.57804

PubMed Abstract | CrossRef Full Text | Google Scholar

32. van Vliet T, Varela-Eirin M, Wang B, Borghesan M, Brandenburg SM, Franzin R, et al. Physiological hypoxia restrains the senescence-associated secretory phenotype via AMPK-mediated mTOR suppression. Mol Cell. (2021) 81(9):2041–2052.e2046. doi: 10.1016/j.molcel.2021.03.018

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Ma S, Yeom J, Lim YH. Specific activation of hypoxia-inducible factor-2α by propionate metabolism via a β-oxidation-like pathway stimulates MUC2 production in intestinal goblet cells. Biomed Pharmacother. (2022) 155:113672. doi: 10.1016/j.biopha.2022.113672

PubMed Abstract | CrossRef Full Text | Google Scholar

34. Rezvani HR, Mahfouf W, Ali N, Chemin C, Ged C, Kim AL, et al. Hypoxia-inducible factor-1alpha regulates the expression of nucleotide excision repair proteins in keratinocytes. Nucleic Acids Res. (2010) 38(3):797–809. doi: 10.1093/nar/gkp1072

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Lin HC, Yeh CC, Chao LY, Tsai MH, Chen HH, Chuang EY, et al. The hypoxia-responsive lncRNA NDRG-OT1 promotes NDRG1 degradation via ubiquitin-mediated proteolysis in breast cancer cells. Oncotarget. (2018) 9(12):10470–82. doi: 10.18632/oncotarget.23732

PubMed Abstract | CrossRef Full Text | Google Scholar

36. Hanai S, Uchimura K, Takahashi K, Ishii T, Mitsui T, Furuya F. Hypoxia-induced thyroid hormone receptor expression regulates cell-cycle progression in renal tubule epithelial cells. Endocr J. (2021) 68(11):1309–20. doi: 10.1507/endocrj.EJ21-0245

PubMed Abstract | CrossRef Full Text | Google Scholar

37. Yang C, Liu H. Both a hypoxia-inducible EYA3 and a histone acetyltransferase p300 function as coactivators of SIX5 to mediate tumorigenesis and cancer progression. Ann Transl Med. (2022) 10(13):752. doi: 10.21037/atm-22-2663

PubMed Abstract | CrossRef Full Text | Google Scholar

38. You L, Nepovimova E, Valko M, Wu Q, Kuca K. Mycotoxins and cellular senescence: the impact of oxidative stress, hypoxia, and immunosuppression. Arch Toxicol. (2023) 97(2):393–404. doi: 10.1007/s00204-022-03423-x

PubMed Abstract | CrossRef Full Text | Google Scholar

39. Yuan HL, Zhang X, Li Y, Guan Q, Chu WW, Yu HP, et al. A nomogram for predicting risk of thromboembolism in gastric cancer patients receiving chemotherapy. Front Oncol. (2021) 11:598116. doi: 10.3389/fonc.2021.598116

PubMed Abstract | CrossRef Full Text | Google Scholar

40. Libby P. Inflammation in atherosclerosis. Arterioscler Thromb Vasc Biol. (2012) 32(9):2045–51. doi: 10.1161/ATVBAHA.108.179705

PubMed Abstract | CrossRef Full Text | Google Scholar

41. Potteaux S, Ait-Oufella H, Mallat Z. Role of splenic monocytes in atherosclerosis. Curr Opin Lipidol. (2015) 26(5):457–63. doi: 10.1097/MOL.0000000000000223

PubMed Abstract | CrossRef Full Text | Google Scholar

42. Zouggari Y, Ait-Oufella H, Bonnin P, Simon T, Sage AP, Guérin C, et al. B lymphocytes trigger monocyte mobilization and impair heart function after acute myocardial infarction. Nat Med. (2013) 19(10):1273–80. doi: 10.1038/nm.3284

PubMed Abstract | CrossRef Full Text | Google Scholar

43. Poss J, Koster J, Fuernau G, Eitel I, de Waha S, Ouarrak T, et al. Risk stratification for patients in cardiogenic shock after acute myocardial infarction. J Am Coll Cardiol. (2017) 69(15):1913–20. doi: 10.1016/j.jacc.2017.02.027

PubMed Abstract | CrossRef Full Text | Google Scholar

44. Myojo M, Ando J, Uehara M, Daimon M, Watanabe M, Komuro I. Feasibility of extracorporeal shock wave myocardial revascularization therapy for post-acute myocardial infarction patients and refractory angina pectoris patients. Int Heart J. (2017) 58(2):185–90. doi: 10.1536/ihj.16-289

PubMed Abstract | CrossRef Full Text | Google Scholar

45. Dahabreh IJ, Sheldrick RC, Paulus JK, Chung M, Varvarigou V, Jafri H, et al. Do observational studies using propensity score methods agree with randomized trials? A systematic comparison of studies on acute coronary syndromes. Eur Heart J. (2012) 33(15):1893–901. doi: 10.1093/eurheartj/ehs114

PubMed Abstract | CrossRef Full Text | Google Scholar

46. Zhang X, Zhang Y, Wang P, Zhang SY, Dong Y, Zeng G, et al. Adipocyte hypoxia-inducible factor 2alpha suppresses atherosclerosis by promoting adipose ceramide catabolism. Cell Metab. (2019) 30(5):937–951 e935. doi: 10.1016/j.cmet.2019.09.016

PubMed Abstract | CrossRef Full Text | Google Scholar

47. Klouche M, Peri G, Knabbe C, Eckstein HH, Schmid FX, Schmitz G, et al. Modified atherogenic lipoproteins induce expression of pentraxin-3 by human vascular smooth muscle cells. Atherosclerosis. (2004) 175(2):221–8. doi: 10.1016/j.atherosclerosis.2004.03.020

PubMed Abstract | CrossRef Full Text | Google Scholar

48. Orn S, Ueland T, Manhenke C, Sandanger O, Godang K, Yndestad A, et al. Increased interleukin-1beta levels are associated with left ventricular hypertrophy and remodelling following acute ST segment elevation myocardial infarction treated by primary percutaneous coronary intervention. J Intern Med. (2012) 272(3):267–76. doi: 10.1111/j.1365-2796.2012.02517.x

PubMed Abstract | CrossRef Full Text | Google Scholar

49. Van Tassell BW, Toldo S, Mezzaroma E, Abbate A. Targeting interleukin-1 in heart disease. Circulation. (2013) 128(17):1910–23. doi: 10.1161/CIRCULATIONAHA.113.003199

PubMed Abstract | CrossRef Full Text | Google Scholar

50. Hirata Y, Tabata M, Kurobe H, Motoki T, Akaike M, Nishio C, et al. Coronary atherosclerosis is associated with macrophage polarization in epicardial adipose tissue. J Am Coll Cardiol. (2011) 58(3):248–55. doi: 10.1016/j.jacc.2011.01.048

PubMed Abstract | CrossRef Full Text | Google Scholar

51. Zernecke A, Winkels H, Cochain C, Williams JW, Wolf D, Soehnlein O, et al. Meta-Analysis of leukocyte diversity in atherosclerotic mouse aortas. Circ Res. (2020) 127(3):402–26. doi: 10.1161/CIRCRESAHA.120.316903

PubMed Abstract | CrossRef Full Text | Google Scholar

52. Srikakulapu P, Hu D, Yin C, Mohanta SK, Bontha SV, Peng L, et al. Artery tertiary lymphoid organs control multilayered territorialized atherosclerosis B-cell responses in aged ApoE-/- mice. Arterioscler Thromb Vasc Biol. (2016) 36(6):1174–85. doi: 10.1161/ATVBAHA.115.306983

PubMed Abstract | CrossRef Full Text | Google Scholar

53. Ait-Oufella H, Herbin O, Bouaziz JD, Binder CJ, Uyttenhove C, Laurans L, et al. B cell depletion reduces the development of atherosclerosis in mice. J Exp Med. (2010) 207(8):1579–87. doi: 10.1084/jem.20100155

PubMed Abstract | CrossRef Full Text | Google Scholar

54. Shin DH, Lin H, Zheng H, Kim KS, Kim JY, Chun YS, et al. HIF-1alpha-mediated upregulation of TASK-2 K(+) channels augments ca(2)(+) signaling in mouse B cells under hypoxia. J Immunol. (2014) 193(10):4924–33. doi: 10.4049/jimmunol.1301829

PubMed Abstract | CrossRef Full Text | Google Scholar

55. Consortium CAD, Deloukas P, Kanoni S, Willenborg C, Farrall M, Assimes TL, et al. Large-scale association analysis identifies new risk loci for coronary artery disease. Nat Genet. (2013) 45(1):25–33. doi: 10.1038/ng.2480

PubMed Abstract | CrossRef Full Text | Google Scholar