Identification of Key Genes as Early Warning Signals of Acute Myocardial Infarction Based on Weighted Gene Correlation Network Analysis and Dynamic Network Biomarker Algorithm

Purpose: The specific mechanisms and biomarkersunderlying the progression of stable coronary artery disease (CAD) to acute myocardial infarction (AMI) remain unclear. The current study aims to explore novel gene biomarkers associated with CAD progression by analyzing the transcriptomic sequencing data of peripheral blood monocytes in different stages of CAD.

Material and Methods: A total of 24 age- and sex- matched patients at different CAD stages who received coronary angiography were enrolled, which included 8 patients with normal coronary angiography, 8 patients with angiographic intermediate lesion, and 8 patients with AMI. The RNA from peripheral blood monocytes was extracted and transcriptome sequenced to analyze the gene expression and the differentially expressed genes (DEG). A Gene Oncology (GO) enrichment analysis was performed to analyze the biological function of genes. Weighted gene correlation network analysis (WGCNA) was performed to classify genes into several gene modules with similar expression profiles, and correlation analysis was carried out to explore the association of each gene module with a clinical trait. The dynamic network biomarker (DNB) algorithm was used to calculate the key genes that promote disease progression. Finally, the overlapping genes between different analytic methods were explored.

Results: WGCNA analysis identified a total of nine gene modules, of which two modules have the highest positive association with CAD stages. GO enrichment analysis indicated that the biological function of genes in these two gene modules was closely related to inflammatory response, which included T-cell activation, cell response to inflammatory stimuli, lymphocyte activation, cytokine production, and the apoptotic signaling pathway. DNB analysis identified a total of 103 genes that may play key roles in the progression of atherosclerosis plaque. The overlapping genes between DEG/WGCAN and DNB analysis identified the following 13 genes that may play key roles in the progression of atherosclerosis disease: SGPP2, DAZAP2, INSIG1, CD82, OLR1, ARL6IP1, LIMS1, CCL5, CDK7, HBP1, PLAU, SELENOS, and DNAJB6.

Conclusions: The current study identified a total of 13 genes that may play key roles in the progression of atherosclerotic plaque and provides new insights for early warning biomarkers and underlying mechanisms underlying the progression of CAD.

Introduction

Coronary artery disease (CAD) remains the leading cause of disease burden worldwide (1). Based on clinical presentation, myocardium injury biomarkers, electrocardiography characteristics, and the extent of myocardium injury, CAD is generally classified as stable CAD and acute coronary syndrome (ACS), which included unstable angina and acute myocardial infarction (AMI). Stable CAD is primarily caused by lumen stenosis caused by atherosclerotic plaque and subsequent oxygen demand-supply mismatch, while ACS is the rupture of vulnerable plaque and subsequent occlusive thrombosis formation and myocardial necrosis (2).

However, the exact mechanism underlying the formation, progression, and rupture of plaque is unclear. The current well-established mechanism of this process involves lipid-driven inflammation (3). Monocyte-derived macrophages are one type of the key inflammatory cells within the plaque and participated in each stage of atherosclerosis plaque formation (4). The rapid development of high-throughput omics technology, such as genomics, transcriptomics, proteomics, and metabolomics, provides new insights into the mechanisms and biomarkers of CAD. For instance, peripheral RNA expression differed significantly between CAD patients and normal control (5) and within CAD patients (6).

There has been rich evidence in studies on biomarkers for diagnosis and risk stratification of coronary heart disease. The validated biomarkers are related to different pathophysiological processes of coronary heart disease, including myocardial injury, altered myocardial stress, inflammation, and vascular endothelial dysfunction. For example, cardiac troponin T (cTnT) and cardiac troponin I (cTnI) demonstrate myocardial tissue specificity (7) and are released into the blood when myocardial tissue suffers damage, leading to an elevated concentration level in peripheral blood (8). Likely, a heart-type fatty acid-binding protein is released into peripheral blood during AMI, and studies confirm its prognostic value for patients with suspected ACS and negative cardiac troponin test results (9, 10). The natriuretic peptide family is a set of typical biomarkers associated with myocardial stress. In 2017, B-type natriuretic peptide (BNP) and N-terminal pro-B-type natriuretic peptide (NT-proBNP) were recommended for the diagnosis, evaluation, and management of heart failure patients by AHA (11). Also, a recent study shows that BNP and NT-proBNP can guide risk stratification in patients with coronary heart disease (12). The soluble suppressor of tumorgenicity 2 (sST2) is a typical biomarker associated with inflammation and is proved to be an independent risk factor for long-term all-cause death in a stable coronary disease cohort (13). Endothelin, converted from its relatively stable primer big endothelin-1, has a strong constrictive effect on the coronary arteries, leading to endothelial dysfunction (14, 15). Research evidence shows that big endothelin–1 has the ability to predict long–term prognosis for both stable coronary heart disease and acute myocardial infarction (16, 17).

However, previous studies mainly focused on the biomarkers associated with the occurrence of CAD, while the biomarkers associated with the progression of CAD remain lacking.

Weighted correlation network analysis (WGCNA), proposed by scholars Zhang and Horvath, is an efficient and accurate bioinformatics method for analyzing microarray data (18). WGCNA analysis divides genes into several modules based on the similarity of gene expression profile and identifies the gene module and corresponding hub genes that are highly correlated with the clinical trait of interest.

WGCNA methods have been successfully applied to identify hub genes in many diseases, including cardiovascular disease (19), cancer (20), psychological disease (21), etc.

However, the primary objective of WGCNA, like most traditional bioinformatics methods, is to distinguish a disease state from a normal state or to diagnose the disease state rather than the “predisease” state (22). In other words, WGCNA may fail to accurately predict the early onset of disease before its development. To overcome this limitation, Chinese scholars Chen et al. proposed the dynamic network biomarker (DNB) theory based on the dynamic features of molecules within the biological system (23). Compared with traditional bioinformatics methods such as WGCNA, the DNB method mainly aims to diagnose the “predisease state,” which can detect the early warning signals and achieve early diagnosis before the onset of disease. The DNB method has been successfully used to detect the critical transition point and key molecules in many pathological process (22, 24–26). We hypothesize that the combined use of both WGCNA and the DNB method will contribute to the identification of robust biomarkers to predict disease deterioration.

The current study aims to explore novel gene biomarkers as early warning signals of AMI by analyzing the gene expression profile of peripheral blood monocytes in different stages of CAD based on transcriptomic sequencing technology, combined with both WGCNA and DNB methods.

Materials and Methods

Study Design and Participants

The study design is shown in Figure 1. We collected blood samples from patients who agreed to provide blood samples from June 2011 to March 2015 at Fuwai Hospital. Eligible patients had symptoms indicating CAD and were receiving elective coronary angiography. We excluded patients with rheumatic heart disease, organic heart diseases and cardiomyopathy, severe liver and renal dysfunction, severe infectious diseases, malignant tumors, immune system and connective tissue diseases, and metabolic diseases, including hyperthyroidism and Cushing syndrome. Fasting blood samples and demographic information were collected at Fuwai Hospital after obtaining informed consent. Specifically, blood sample of patients in the AMI group were collected when they were rehospitalized due to ischemia symptoms. Sampling was performed prior to angiography in the fasting status. We randomly selected 8 patients with stable coronary heart disease and angiography–confirmed intermediate lesion (IML) (visual stenosis 50%–70%), as well as 8 age– and sex–matched patients with normal coronary angiography (NCA) and patients with AMI. AMI was diagnosed according to the third universal diagnosis of acute myocardial infarction, with cardiac biomarker (primarily cTn) elevation above 99% upper reference limit with ischemia symptoms, new–onset ischemic ST–T segment change or left bundle branch block, and Q wave formation. Patients in the IML group had objective evidence of ischemia such as chest pain and other myocardial ischemia symptoms or positive exercise stress test results. In recent years, the coronary physiology evaluation method fractional flow reserve (FFR) is the “gold standard” for identifying lesions with physiological significance. Recent studies have also confirmed that a high proportion of intermediate lesions, despite diameter stenosis <80%, may cause myocardial ischemia as detected by FFR. In addition, patients with intermediate lesions varied significantly in their long–term prognosis despite a similar degree of lesion stenosis. Therefore, the current study enrolled patients with intermediate lesions to represent the disease stage of “stable coronary heart disease.” The study protocol complied with the Declaration of Helsinki and was approved by the ethics committee of Fuwai Hospital (No. 2012–431).

FIGURE 1

Figure 1 Study flow chart. A total of 24 age– and sex–matched coronary artery disease (CAD) patients at different disease stages were enrolled, which included 8 patients with normal coronary angiography (NCA), 8 patients with intermediated coronary lesion (IML), and 8 patients with acute myocardial infarction (AMI). Peripheral monocytes were isolated and RNA was extracted. Transcriptomic sequencing was performed, and various analytical methods were used to analyze gene expression profile across different groups. The overlapping genes between different analytical methods were identified as key genes in the progression of CAD. CAD, coronary artery disease, NCA, normal coronary angiography, IML, intermediate coronary lesion, AMI, acute myocardial infarction, DEG, differential expressed genes, WGCNA, weighted gene correlation network analysis, DNB, dynamic network biomarker.

Peripheral Blood Monocyte Isolation, RNA Isolation, and Sequencing

Fasting blood samples were collected in the morning before the angiography procedure. Leukocytes were isolated by centrifugation, and monocytes were isolated by using EasySep™ Human CD14 Positive Selection Kit (#18058, Stemcell, USA) following the manufacturer’s instructions. RNA was extracted from leukocytes by using TRIzol^® Reagent (#15596018, Invitrogen, USA) and frozen at −80°C. The purity of RNA was assessed by NanoPhotometer Spectrophotometers (IMPLEN, USA), and the quantity and quality of RNA were assessed by Agilent 2100 Bioanalyzer and Agilent RNA 6000 Nano Assay Kits.

A total of 3 µg of qualified RNA was used to construct the cDNA library. mRNA was purified from total RNA using poly–T oligo–attached magnetic beads and fragmented with NEBNext First Strand Synthesis Reaction Buffer (5×). The sequencing was performed by Annoroad Gene Technology Corporation (Beijing) by using the Illumina Novaseq S2 platform (Illumina, USA) with the PE–150 module.

Data Processing

Basic data processing included three major steps: data quality control, transcriptome mapping, and gene expression analysis. Raw data were filtered by the Cutadapt software to remove low–quality sequencing reads and generate high–quality data (clean reads). The human reference database and annotation files were downloaded from the ENSEMBL database (version Homo_sapiens.GRCh38.91.chr). The clean reads were then aligned to the reference genome by using HISAT2 v2.1.0 software.

Fragments per kilobase per million mapped fragments (FPKM) were calculated to assess gene expression level and used for subsequent further analysis. The algorithm for calculating FPKM is as follows:

$$FPKM=\frac{{10}^3*F}{NL/{10}^6}$$

where F is the number of reads mapped to the gene, N is the total number of mapped reads or fragments, and L is the gene length.

Differential Expression Analysis and GO Enrichment Analysis

Differentially expressed gene (DEG) analysis was performed by using DESeq2 R packages, and genes with fold change ≥1.5 or fold change ≤0.67 and an adjusted p–value <0.05 were identified as differentially expressed genes. The volcano plot was used to visualize differentially expressed genes between groups and was generated by the ggplot packages. GO enrichment analysis was used to investigate the biological function of differentially expressed genes, which was performed and visualized on the Metascape website (https://metascape.org/gp/index.html#/main/step1) (27). The parameters were set as follows: min overlap of 3, p–value cutoff of 0.01, and minimum enrichment factor of 1.5. The top 20 enriched pathways were selected for visualization.

WGCNA Analysis

The WGCNA analysis was performed by using the R package “WGCNA” (28), which mainly included three steps: construction of a gene coexpression network, identify gene modules, and gene module–clinical trait correlation analysis. The expression matrix was constructed with the original FPKM. We included genes with FPKM greater than zero in more than 8 samples. The soft–threshold power was chosen by the pickSoftThreshold function. The gene network was constructed by the blockwiseModules function, and the parameters were set as follows: minModuleSize = 30, reassignThreshold = 0, and mergeCutHeight = 0.25. The corresponding eigengenes were obtained to summarize the expression profile of each gene module. The association between each gene module and clinical trait was performed by using the Pearson’s correlation method, which correlate the module eigengenes with each clinical trait. Genes with absolute gene modulemembership > 0.8 and genetraitsignificance > 0.2 were identified as hub genes (29).

DNB Analysis

A detailed theoretical foundation and computational algorithm for DNB analysis were described previously (23). Based on the DNB theory, the progression of the disease is not smooth but abrupt (Figure 2), and there is a critical transition point after which the system shifts abruptly from one state to another. According to this transition point, the disease progression process can be divided into three stages: “normal state,” which is a relatively stable state where the disease undergoes gradual and slow change, “predisease state,” which is the limit of the normal state just before the transition to the disease state, and the “disease state,” which is another relatively stable state and is usually irreversible to the normal state. Based on the DNB theory, there exists a group of molecules (genes, proteins or metabolites, etc.) in the predisease state, which can be used for predicting disease. This group of molecules characterizes the dynamic features of the underlying system and are termed as DNB and satisfies the following three criteria:

FIGURE 2

Figure 2 Illustration of DNB theory. Based on the DNB theory, there is a critical transition point (tipping point) during the progression of the disease. The disease progression process can be divided into three stages according to the tipping point: the normal state, the predisease state, and the disease state. The abbreviations are the same as above.

1. The average standard deviations (SDs) of DNB molecules increase drastically.

2. The average Pearson’s correlation coefficients (PCCs) of DNB molecules drastically increase in the predisease state.

3. The average PCCs of molecules between DNB molecules and non–DNB molecules drastically decrease.

The above three criteria can which can be represented as a composite index (CI):

$$CI=\frac{S{D}_d\text{x}PC{C}_d}{PC{C}_O}$$

whereas SD_d is the average standard deviation (SD) for molecules inside the DNB module, PCC_d is the average Pearson’s correlation coefficient (PCC) in absolute value for molecules inside the module, and PCC_o is the average PCC in absolute value for molecules between DNB and non–DNB. The CI is expected to increase abruptly and significantly before the critical transition to the disease state and can serve as an early warning signal.

Statistical Analysis

SPSS 26.0 was used for statistical analysis. Continuous variables were expressed as mean ± standard deviation, and comparisons between groups were performed by analysis of variance. Categorical variables were expressed as frequency (percentage) using the chi–square test or Fisher’s exact test for comparison between groups. The difference was considered statistically significant when the bilateral P value was less than 0.05.

Results

Clinical Characteristics of the Recruited Patients

The baseline characteristics of the 8 IML patients and 8 gender– and–age–matched NCA and AMI patients are demonstrated in Table 1. Among the 24 patients, 17 (70.83%) were men, with a median age of 64 years. Additionally, no significant difference was found in terms of age, gender, and body mass index (BMI) between the different groups. As expected, CAD patients had a significantly higher rate of family history of CAD and hyperlipidemia, and more tended to be current smokers, but the difference was not statistically significant.

TABLE 1

Table 1 Baseline characteristics across groups.

Differential Gene Expression Analysis of Different Phases of CAD Patients

Gene expression profiles differed significantly across groups. As shown in Figure 3, there were a total of 192 DEGs between the NCA and IML groups, 2,269 DEGs between the AMI and NCA groups, and 385 DEGs between the AMI and IML groups. To analyze the potential biological function of DEGs, the online analytic tool Metascape was used, and differential genes were uploaded on the website. The top enriched pathways of the DEGs between the IML group and the NCA group were associated with cell adhesion, cell migration, and cell response to inflammation, the top enriched pathways of the DEGs between the AML group and the NCA group were associated with cellular response to stress and cellular cytokine production, and the top enriched pathways of the DEGs between the AML group and the IML group were associated with mRNA metabolic process, cellular signal transduction, and cellular response to growth factor stimulus (Supplementary Figures S1–S3).

FIGURE 3

Figure 3 Volcano plot of differentially expressed genes between the IML and NCA groups (A), AMI and NCA groups (B), as well as AMI and IML group (C). The abbreviations are the same as above.

WGCNA Analysis

The soft threshold of 8 was used to construct the scale–free network (Figure 4). Nine modules were identified based on average hierarchical clustering and dynamic tree clipping (Figure 5). The number of genes in each gene module is shown in Table 2.

FIGURE 4

Figure 4 Determination of the soft threshold for the WGCNA analysis. Analysis of the scale–free index for various soft–threshold powers (A) and mean connectivity (B) for various soft–threshold powers. The abbreviations are the same as above.

FIGURE 5

Figure 5 Clustering dendrogram of all differentially expressed genes based on dissimilarity, with an assigned color for each gene module.

TABLE 2

Table 2 The number of genes and median MM and GS in each gene module.

To identify the genes associated with the progression of atherosclerosis, we evaluated the association between each gene module and clinical trait by calculating the module significance (MS) for each module–trait correlation (Figure 6). Disease stage is a trichotomous variable that indicates NCA, IML, and AMI. The pink module had the highest positive association with disease stage (r² = 0.64, p = 8e−04), followed by the turquoise module (r² = 0.55, p = 0.005).

FIGURE 6

Figure 6 Identification of the gene module associated with clinical traits. Heatmap of the correlation between the module eigengenes and clinical traits of CAD patients. Disease stage is a trichotomous variable (NCA, IML, and AMI). The number in each cell is the correlation coefficient (the corresponding p–value). The abbreviations are the same as above.

Genes within each gene module, as well as the corresponding module membership and gene significance, are shown in Supplementary Table S1. The top 3 GO–BP–enriched pathways in pink module are endoplasmic reticulum organization, regulation of alpha–beta T cell, and peptidyl–serine phosphorylation (Figure not shown). The top enriched pathways in the turquoise module are associated with cell response to inflammatory stimuli, lymphocyte activation, cytokine production, and the apoptotic signaling pathway (Supplementary Figure S4), which has been reported to be closely associated with the atherosclerotic process.

The scatter plot for gene significance and module membership for the above gene modules are shown in Supplementary Figure S5.

DNB Genes

The average DNB score for each group is shown in Figure 7. The average DNB score of the IML group is greater than that of the NCA and AMI groups, which indicates that the critical transition point of atherosclerosis is at the IML stage. This is also consistent with the pathophysiological process of atherosclerosis progression in clinical practice.

FIGURE 7

Figure 7 The average DNB score and corresponding parameters for the three groups. x–axis represents the three groups, and y–axis represents the average DNB score for each group.

There are a total of 103 DNB genes (Supplementary Table S2). The top 3 GO–BP enrichment pathways included regulation of activation of Janus kinase activity, regulation of cell adhesion, and negative regulation of proteolysis (Supplementary Figure S6).

Genes Potentially Associated with CAD progression

To identify robust biomarkers that may serve as early warning signals of AMI, we analyzed the common genes (overlapping genes) between differentially expressed gene sets and DNB gene sets and identified the SGPP2 genes that may play key roles in the process of atherosclerosis progression (Figure 8).

FIGURE 8

Figure 8 Venn diagram showing the common genes between DEG and DNB (A), as well as DNB and WGCNA analyses (B). The abbreviations are the same as above.

We also analyzed the overlapping genes between DNB gene sets and WGCNA gene sets and identified the following 12 genes: DAZAP2, INSIG1, CD82, OLR1, ARL6IP1, LIMS1, CCL5, CDK7, HBP1, PLAU, SELENOS, and DNAJB6. Gene expression levels are expressed in the violin plot in Figure 9.

FIGURE 9

Figure 9 The violin plot showing the gene expression level of the 13 overlapping genes between DEG/WGCNA and DNB analyses. The abbreviation.

Discussion

CAD remains the leading cause of disease burden worldwide. Elevated LDL–c levels are the well–established principal risk factor for the onset and progression of CAD. However, a substantial number of patients still suffer from adverse cardiovascular events despite the reduction of LDL–c level to 40 mg/dl or less (30, 31), suggesting additional “residual cardiovascular risks” are yet to be explored.

In the microenvironment of coronary atherosclerotic plaques, the main type of immune cell is the macrophage, which is differentiated from peripheral blood monocyte. In atherogenic conditions, circulation monocytes sense changes in the microenvironment and adopt specific gene expression in response. These monocytes are then recruited to the subintimal layers of the artery wall and differentiated into macrophages. Therefore, studying the gene expression of peripheral blood monocytes in different stages of CAD through RNA–sequencing will help to reveal the underlying mechanism.

Disease progression is a complex nonlinear process, which is not necessarily “smooth” but “abrupt” (23). There is usually a drastic change during disease progression, which causes the critical transition from the normal state to the disease state (Figure 2).

Traditional bioinformatics methods primarily compare the molecular characteristics between disease and normal state and poorly detect the predisease state, or accurately predict the onset/deterioration of disease before its occurrence, due to the similarity of molecular characteristics between the pre–disease state and the normal state. To overcome this limitation, a novel computational algorithm, dynamic network biomarkers, was proposed to detect the critical transition point in the complex biological process, such as the progression of atherosclerosis in our study. The DNB theory was based on solid nonlinear dynamic system theory and was reported to serve as an early–warning signal prior to disease deterioration.

Several of these genes have been reported to be closely associated with the onset or progression of atherosclerosis, further demonstrating the robustness of our conclusion.

C–C motif chemokine ligand 5 (CCL5) encodes a 68–amino acid chemokine, which functions as a chemoattractant for blood immune cells and the natural ligand for the chemokine receptor chemokine (C–C motif) receptor 5 (CCR5). CCL5 is involved in a wide range of inflammatory processes, including advanced atherosclerosis and myocardial reperfusion injury. CCL5 has been detected in atherosclerosis plaque (32). CCR5 deficiency reduces the development of diet–induced atherosclerosis in mice (33). Inhibition of CCL5 reduces myocardial reperfusion in atherosclerosis mice (34).

Oxidized low–density lipoprotein receptor 1 (OLR1) encodes a low–density lipoprotein receptor, which binds, internalizes, and degrades oxidized low–density lipoprotein. The previous meta–analysis demonstrated a significant association between OLR1 gene polymorphisms and CAD risk (35). OLR1 promotes endothelial dysfunction by inducing pro–atherogenic signaling via the endothelial uptake of oxidized LDL (oxLDL), which contributes to the initiation, progression, and destabilization of atheromatous plaques (36). In addition to its expression in endothelial cells, OLR1 is also expressed in immune cells, including macrophages, lymphocytes, and neutrophils, further implicating this receptor in multiple aspects of atherosclerotic plaque formation. In conclusion, OLR1 holds promise as a novel diagnostic and therapeutic target for atherosclerosis and CHD.

Plasminogen activator, urokinase (PLAU) encodes a secreted serine protease that converts plasminogen to plasmin. PLAU has been implicated in a broad spectrum of biological and pathological processes, including chemotaxis, cell adhesion, migration and growth, fibrinolysis, proteolysis, angiogenesis, inflammation, and neointima formation (37). PLAU has also been reported to be associated with atherosclerosis plaque formation and AMI: macrophage–specific overexpression of the PLAU gene accelerated atherosclerosis, coronary artery occlusions, and premature death in ApoE^−/− mice (38). PLAU has already been reported to be closely associated with AMI: an SNP rs4065 of the PLAU gene is associated with AMI risk in the Chinese Han population (39). Taken together, the above evidence suggests the PLAU gene as a novel therapeutic target for the treatment of atherosclerosis.

Several genes have been reported to be involved in metabolic or inflammatory processes and, therefore, may potentially be involved in the pathophysiological process of atherosclerosis. However, this hypothesis still needs to be demonstrated in future studies.

Sphingosine–1–phosphate phosphatase 2 (SGPP2) is differentially expressed between both the AMI and NCA groups, as well as the AMI and IML groups. SGPP2 encodes sphingosine–1–phosphate phosphatase 2, which can degrade sphingosine 1–phosphate (S1P) to produce sphingosine. SGPP2 is expressed in human umbilical vein endothelial cells and neutrophils (40). The results of previous studies suggest that the SGPP2 gene may be involved in the onset and progress of atherosclerosis. The SGPP2 gene affects the endothelial barrier function via altering the expression of interleukin 1–β (IL–1β) (40) in endothelial cells. IL–1β is a proinflammatory cytokine that can induce endothelial cell inflammation and destroy the endothelial barrier function (41). In addition, SGPP2 is involved in the inflammatory response process: SGPP2 knockout mice showed reduced expression of proinflammatory factors, including IL–6, TNF–α, and IL–1β, molecular mechanisms, and inhibition of inflammation–induced signal transducer and activator of transcription 3 (STAT–3) signal pathway activation related. The activation of the STAT–3 signaling pathway plays an important role in the process of macrophage inflammation and polarization (42). Endothelial barrier dysfunction, inflammatory cell infiltration, and release of inflammatory factors are the key factors leading to plaque progression and rupture.

Insulin–induced gene 1 (INSIG1) gene encodes an endoplasmic reticulum membrane protein that regulates cholesterol metabolism, lipogenesis, and glucose homeostasis in various tissues (43). INSIG1 acts as a negative regulator of cholesterol biosynthesis by mediating the retention of the SCAP–SREBP complex in the endoplasmic reticulum, thereby blocking the processing of sterol regulatory element–binding proteins (SREBPs) (44). INSIG1 gene single nucleotide polymorphisms were associated with coronary heart disease risk in the Chinese Han population (45). Knockdown of INSIG1 resulted in a significant reduction of cholesterol efflux to HDL (46). INSIG1 variation may contribute to statin–induced changes in plasma TG in a sex–specific manner (47).

CD82 encodes a membrane glycoprotein, a member of the transmembrane 4 superfamily. The primary research area of CD82 is tumors, and CD82 has been recognized as a tumor metastasis suppressor gene (48). CD82 inhibits pathological angiogenesis. Endothelial cells CD82 knockout enhanced the migration and invasion capabilities of endothelial cells (49). CD82 also plays a key role in the regulation of endothelial–monocyte interactions, which include monocyte recruitment and migration (50). The above studies suggest that CD82 may affect atherosclerotic plaque formation, which requires validation in future studies.

Selenoprotein S (SELENOS) encodes a transmembrane protein, which is involved in the degradation process of misfolded proteins in the ER and may also have a role in inflammation control. SELENOS has been reported to be associated with the risk of diabetes: genetic polymorphisms of SELENOS genes are associated with diabetes risk in the Chinese population (51, 52). The serum source of SELENOS is primarily from hepatocytes, and the serum level of SELENOS was associated with the risk of DM and its macrovascular complications (53). Given that SELENOS is closely associated with inflammation, oxidative stress, as well as glucose metabolism (53), the above evidence indicates SELENOS played a key role in the pathophysiology process of AS, which requires future validation.

The biological function and their roles in atherosclerosis progression for the other genes, including DAZAP2, ARL6ip, CDK7, LIMS, HBP1, and DNAJB6, remain to be investigated. The current study has several limitations: Firstly, sample size of the current study was limited, and 8 samples were tested using RNA–seq for each group. Setting up biological replicates is necessary to eliminate errors in sequencing studies, and sequencing technique or statistical tools cannot fully eliminate biological variability. One recent study recommended no fewer than 6 biological replicates should be included in a single group in RNA–seq study (54). Also, to compensate for this deficiency, the DNB algorithm was applied. One advantage of the DNB algorithm is that it can identify the critical transition stage even based on small samples of high–throughput data (55). The first publication proposed and demonstrated the validity of the DNB algorithm by detecting the “early warning signals” prior to acute lung injury in carbonyl chloride inhalation–induced acute lung injury, based on 2–5 samples of lung tissue transcriptomic data at each sampling period (23). Secondly, since the current study was a retrospective study, future studies are required to validate whether the key genes in our study can predict adverse events in a prospective cohort, which is an ongoing project. Thirdly, since the current study was based on high–throughput sequencing data, quantitative validation of the gene expression level is required based on the qPCR method. Finally, basic studies are needed to explore the biological functions of the key genes identified in our study.

Conclusion and Future Prospects

In conclusion, based on the peripheral blood mononuclear cell transcriptome sequencing data from patients at different disease progression stages of coronary artery disease, combined with traditional DEG analysis, WGCNA analysis, and novel DNB methods, the current study identified a total of 13 genes that may play key roles involved in the progression of atherosclerotic plaque and provides new insights for early warning biomarkers and potential underlying mechanisms underlying the progression of CAD.

Data Availability Statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: https://www.ncbi.nlm.nih.gov/geo/, GSE166780.

Ethics Statement

The studies involving human participants were reviewed and approved by the ethics committee of Fuwai Hospital. The patients/participants provided their written informed consent to participate in this study.

Author Contributions

CS and ZQ carried out the sequence data analysis, participated in the sequence alignment, and drafted the manuscript. LC and JG created the method of the DNB analysis. SY, XB, CW, RZ, LJ, and QL participated in the design of the study and performed the statistical analysis. RF and KD conceived of the study and participated in its design and coordination and helped to draft the manuscript. All authors read and approved the final manuscript.

Funding

This work was supported by the CAMS Innovation Fund for Medical Sciences (CIFMS) (Grant No. 2021–I2M–1–008) awarded to KD and was supported by The Natural Science Foundation of China (No. 81870277). This work was also supported by Special Fund for Science and Technology Innovation Strategy of Guangdong Province (Nos. 2021B0909050004, 2021B0909060002); Major Key Project of PCL (No.PCL2021A12); JST Moonshot R&D Project (No.JPMJMS2021).

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.

Acknowledgments

We acknowledged all the staffs for data collection and patients’ follow–up.

Supplementary Material

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

References

1. Roth GA, Mensah GA, Johnson CO, Addolorato G, Ammirati E, Baddour LM, et al. Global Burden of Cardiovascular Diseases and Risk Factors 1990–2019: Update From the Gbd 2019 Study. J Am Coll Cardiol (2020) 76(25):2982–3021. doi: 10.1016/j.jacc.2020.11.010

CrossRef Full Text | Google Scholar

2. Bentzon JF, Otsuka F, Virmani R, Falk E. Mechanisms of Plaque Formation and Rupture. Circ Res (2014) 114(12):1852–66. doi: 10.1161/circresaha.114.302721

CrossRef Full Text | Google Scholar

3. Soehnlein O, Libby P. Targeting Inflammation in Atherosclerosis – From Experimental Insights to the Clinic. Nat Rev Drug Discov (2021) 20(8):589–610. doi: 10.1038/s41573–021–00198–1

CrossRef Full Text | Google Scholar

4. Moore KJ, Tabas I. Macrophages in the Pathogenesis of Atherosclerosis. Cell (2011) 145(3):341–55. doi: 10.1016/j.cell.2011.04.005

CrossRef Full Text | Google Scholar

5. Liao J, Wang J, Liu Y, Li J, Duan L. Transcriptome Sequencing of LncRNA, Mirna, Mrna and Interaction Network Constructing in Coronary Heart Disease. BMC Med Genomics (2019) 12(1):124. doi: 10.1186/s12920–019–0570–z

CrossRef Full Text | Google Scholar

6. Peng XY, Wang Y, Hu H, Zhang XJ, Li Q. Identification of the Molecular Subgroups in Coronary Artery Disease by Gene Expression Profiles. J Cell Physiol (2019) 234 (9): 16540–8. doi: 10.1002/jcp.28324

CrossRef Full Text | Google Scholar

7. Mohammed AA, Januzzi JL Jr. Clinical Applications of Highly Sensitive Troponin Assays. Cardiol Rev (2010) 18(1):12–9. doi: 10.1097/CRD.0b013e3181c42f96

CrossRef Full Text | Google Scholar

8. Apple FS, Collinson PO. Analytical Characteristics of High–Sensitivity Cardiac Troponin Assays. Clin Chem (2012) 58(1):54–61. doi: 10.1373/clinchem.2011.165795

CrossRef Full Text | Google Scholar

9. Viswanathan K, Kilcullen N, Morrell C, Thistlethwaite SJ, Sivananthan MU, Hassan TB, et al. Heart–Type Fatty Acid–Binding Protein Predicts Long–Term Mortality and Re–Infarction in Consecutive Patients With Suspected Acute Coronary Syndrome Who Are Troponin–Negative. J Am Coll Cardiol (2010) 55(23):2590–8. doi: 10.1016/j.jacc.2009.12.062

CrossRef Full Text | Google Scholar

10. Moon MG, Yoon CH, Lee K, Kang SH, Youn TJ, Chae IH. Evaluation of Heart–Type Fatty Acid–Binding Protein in Early Diagnosis of Acute Myocardial Infarction. J Korean Med Sci (2021) 36(8):e61. doi: 10.3346/jkms.2021.36.e61

CrossRef Full Text | Google Scholar

11. Chow SL, Maisel AS, Anand I, Bozkurt B, de Boer RA, Felker GM, et al. Role of Biomarkers for the Prevention, Assessment, and Management of Heart Failure: A Scientific Statement From the American Heart Association. Circulation (2017) 135(22):e1054–e91. doi: 10.1161/cir.0000000000000490

CrossRef Full Text | Google Scholar

12. Kragelund C, Grønning B, Køber L, Hildebrandt P, Steffensen R. N–Terminal Pro–B–Type Natriuretic Peptide and Long–Term Mortality in Stable Coronary Heart Disease. New Engl J Med (2005) 352(7):666–75. doi: 10.1056/NEJMoa042330

CrossRef Full Text | Google Scholar

13. Dieplinger B, Egger M, Haltmayer M, Kleber ME, Scharnagl H, Silbernagel G, et al. Increased Soluble St2 Predicts Long–Term Mortality in Patients With Stable Coronary Artery Disease: Results From the Ludwigshafen Risk and Cardiovascular Health Study. Clin Chem (2014) 60(3):530–40. doi: 10.1373/clinchem.2013.209858

CrossRef Full Text | Google Scholar

14. Papassotiriou J, Morgenthaler NG, Struck J, Alonso C, Bergmann A. Immunoluminometric Assay for Measurement of the C–Terminal Endothelin–1 Precursor Fragment in Human Plasma. Clin Chem (2006) 52(6):1144–51. doi: 10.1373/clinchem.2005.065581

CrossRef Full Text | Google Scholar

15. Kolettis TM, Barton M, Langleben D, Matsumura Y. Endothelin in Coronary Artery Disease and Myocardial Infarction. Cardiol Rev (2013) 21(5):249–56. doi: 10.1097/CRD.0b013e318283f65a

CrossRef Full Text | Google Scholar

16. Khan SQ, Dhillon O, Struck J, Quinn P, Morgenthaler NG, Squire IB, et al. C–Terminal Pro–Endothelin–1 Offers Additional Prognostic Information in Patients After Acute Myocardial Infarction: Leicester Acute Myocardial Infarction Peptide (Lamp) Study. Am Heart J (2007) 154(4):736–42. doi: 10.1016/j.ahj.2007.06.016

CrossRef Full Text | Google Scholar

17. Zhou BY, Guo YL, Wu NQ, Zhu CG, Gao Y, Qing P, et al. Plasma Big Endothelin–1 Levels at Admission and Future Cardiovascular Outcomes: A Cohort Study in Patients With Stable Coronary Artery Disease. Int J Cardiol (2017) 230:76–9. doi: 10.1016/j.ijcard.2016.12.082

CrossRef Full Text | Google Scholar

18. Zhang B, Horvath S. A General Framework for Weighted Gene Co–Expression Network Analysis. Stat Appl Genet Mol Biol (2005) 4:Article17. doi: 10.2202/1544–6115.1128

CrossRef Full Text | Google Scholar

19. Niu X, Zhang J, Zhang L, Hou Y, Pu S, Chu A, et al. Weighted Gene Co–Expression Network Analysis Identifies Critical Genes in the Development of Heart Failure After Acute Myocardial Infarction. Front Genet (2019) 10:1214. doi: 10.3389/fgene.2019.01214

CrossRef Full Text | Google Scholar

20. Giulietti M, Occhipinti G, Righetti A, Bracci M, Conti A, Ruzzo A, et al. Emerging Biomarkers in Bladder Cancer Identified by Network Analysis of Transcriptomic Data. Front Oncol (2018) 8:450. doi: 10.3389/fonc.2018.00450

CrossRef Full Text | Google Scholar

21. Zhang ZQ, Wu WW, Chen JD, Zhang GY, Lin JY, Wu YK, et al. Weighted Gene Coexpression Network Analysis Reveals Essential Genes and Pathways in Bipolar Disorder. Front Psychiatry (2021) 12:553305. doi: 10.3389/fpsyt.2021.553305

CrossRef Full Text | Google Scholar

22. Liu R, Wang X, Aihara K, Chen L. Early Diagnosis of Complex Diseases by Molecular Biomarkers, Network Biomarkers, and Dynamical Network Biomarkers. Med Res Rev (2014) 34(3):455–78. doi: 10.1002/med.21293

CrossRef Full Text | Google Scholar

23. Chen L, Liu R, Liu ZP, Li M, Aihara K. Detecting Early–Warning Signals for Sudden Deterioration of Complex Diseases by Dynamical Network Biomarkers. Sci Rep (2012) 2:342. doi: 10.1038/srep00342

CrossRef Full Text | Google Scholar

24. Yang Y, Lin X, Lu X, Luo G, Zeng T, Tang J, et al. Interferon–Microrna Signalling Drives Liver Precancerous Lesion Formation and Hepatocarcinogenesis. Gut (2016) 65(7):1186–201. doi: 10.1136/gutjnl–2015–310318

CrossRef Full Text | Google Scholar

25. Yang B, Li M, Tang W, Liu W, Zhang S, Chen L, et al. Dynamic Network Biomarker Indicates Pulmonary Metastasis at the Tipping Point of Hepatocellular Carcinoma. Nat Commun (2018) 9(1):678. doi: 10.1038/s41467–018–03024–2

CrossRef Full Text | Google Scholar

26. Zhang C, Zhang H, Ge J, Mi T, Cui X, Tu F, et al. Landscape Dynamic Network Biomarker Analysis Reveals the Tipping Point of Transcriptome Reprogramming to Prevent Skin Photodamage. J Mol Cell Biol (2021) 13 (11): 822–33. doi: 10.1093/jmcb/mjab060

CrossRef Full Text | Google Scholar

27. Zhou Y, Zhou B, Pache L, Chang M, Khodabakhshi AH, Tanaseichuk O, et al. Metascape Provides a Biologist–Oriented Resource for the Analysis of Systems–Level Datasets. Nat Commun (2019) 10(1):1523. doi: 10.1038/s41467–019–09234–6

CrossRef Full Text | Google Scholar

28. Langfelder P, Horvath S. Wgcna: An R Package for Weighted Correlation Network Analysis. BMC Bioinf (2008) 9:559. doi: 10.1186/1471–2105–9–559

CrossRef Full Text | Google Scholar

29. Liu Y, Gu HY, Zhu J, Niu YM, Zhang C, Guo GL. Identification of Hub Genes and Key Pathways Associated With Bipolar Disorder Based on Weighted Gene Co–Expression Network Analysis. Front Physiol (2019) 10:1081. doi: 10.3389/fphys.2019.01081

CrossRef Full Text | Google Scholar

30. Sabatine MS, Giugliano RP, Keech AC, Honarpour N, Wiviott SD, Murphy SA, et al. Evolocumab and Clinical Outcomes in Patients With Cardiovascular Disease. New Engl J Med (2017) 376(18):1713–22. doi: 10.1056/NEJMoa1615664

CrossRef Full Text | Google Scholar

31. Schwartz GG, Steg PG, Szarek M, Bhatt DL, Bittner VA, Diaz R, et al. Alirocumab and Cardiovascular Outcomes After Acute Coronary Syndrome. New Engl J Med (2018) 379(22):2097–107. doi: 10.1056/NEJMoa1801174

CrossRef Full Text | Google Scholar

32. Veillard NR, Kwak B, Pelli G, Mulhaupt F, James RW, Proudfoot AE, et al. Antagonism of Rantes Receptors Reduces Atherosclerotic Plaque Formation in Mice. Circ Res (2004) 94(2):253–61. doi: 10.1161/01.Res.0000109793.17591.4e

CrossRef Full Text | Google Scholar

33. Braunersreuther V, Zernecke A, Arnaud C, Liehn EA, Steffens S, Shagdarsuren E, et al. Ccr5 But Not Ccr1 Deficiency Reduces Development of Diet–Induced Atherosclerosis in Mice. Arterioscler Thromb Vasc Biol (2007) 27(2):373–9. doi: 10.1161/01.ATV.0000253886.44609.ae

CrossRef Full Text | Google Scholar

34. Braunersreuther V, Pellieux C, Pelli G, Burger F, Steffens S, Montessuit C, et al. Chemokine Ccl5/Rantes Inhibition Reduces Myocardial Reperfusion Injury in Atherosclerotic Mice. J Mol Cell Cardiol (2010) 48(4):789–98. doi: 10.1016/j.yjmcc.2009.07.029

CrossRef Full Text | Google Scholar

35. Salehipour P, Rezagholizadeh F, Mahdiannasser M, Kazerani R, Modarressi MH. Association of Olr1 Gene Polymorphisms With the Risk of Coronary Artery Disease: A Systematic Review and Meta–Analysis. Heart Lung (2021) 50(2):334–43. doi: 10.1016/j.hrtlng.2021.01.015

CrossRef Full Text | Google Scholar

36. Akhmedov A, Sawamura T, Chen CH, Kraler S, Vdovenko D, Lüscher TF. Lectin–Like Oxidized Low–Density Lipoprotein Receptor–1 (Lox–1): A Crucial Driver of Atherosclerotic Cardiovascular Disease. Eur Heart J (2021) 42(18):1797–807. doi: 10.1093/eurheartj/ehaa770

CrossRef Full Text | Google Scholar

37. Chavakis T, Kanse SM, May AE, Preissner KT. Haemostatic Factors Occupy New Territory: The Role of the Urokinase Receptor System and Kininogen in Inflammation. Biochem Soc Trans (2002) 30(2):168–73. doi: 10.1042/bst0300168

CrossRef Full Text | Google Scholar

38. Cozen AE, Moriwaki H, Kremen M, DeYoung MB, Dichek HL, Slezicki KI, et al. Macrophage–Targeted Overexpression of Urokinase Causes Accelerated Atherosclerosis, Coronary Artery Occlusions, and Premature Death. Circulation (2004) 109(17):2129–35. doi: 10.1161/01.Cir.0000127369.24127.03

CrossRef Full Text | Google Scholar

39. Xu J, Li W, Bao X, Ding H, Chen J, Zhang W, et al. Association of Putative Functional Variants in the Plau Gene and the Plaur Gene With Myocardial Infarction. Clin Sci (Lond) (2010) 119(8):353–9. doi: 10.1042/cs20100151

CrossRef Full Text | Google Scholar

40. Mechtcheriakova D, Wlachos A, Sobanov J, Kopp T, Reuschel R, Bornancin F, et al. Sphingosine 1–Phosphate Phosphatase 2 Is Induced During Inflammatory Responses. Cell Signal (2007) 19(4):748–60. doi: 10.1016/j.cellsig.2006.09.004

CrossRef Full Text | Google Scholar

41. Béguin EP, van den Eshof BL, Hoogendijk AJ, Nota B, Mertens K, Meijer AB, et al. Integrated Proteomic Analysis of Tumor Necrosis Factor Α and Interleukin 1β–Induced Endothelial Inflammation. J Proteomics (2019) 192:89–101. doi: 10.1016/j.jprot.2018.08.011

CrossRef Full Text | Google Scholar

42. Chen Q, Lv J, Yang W, Xu B, Wang Z, Yu Z, et al. Targeted Inhibition of Stat3 as a Potential Treatment Strategy for Atherosclerosis. Theranostics (2019) 9(22):6424–42. doi: 10.7150/thno.35528

CrossRef Full Text | Google Scholar

43. Dong XY, Tang SQ. Insulin–Induced Gene: A New Regulator in Lipid Metabolism. Peptides (2010) 31(11):2145–50. doi: 10.1016/j.peptides.2010.07.020

CrossRef Full Text | Google Scholar

44. Yang T, Espenshade PJ, Wright ME, Yabe D, Gong Y, Aebersold R, et al. Crucial Step in Cholesterol Homeostasis: Sterols Promote Binding of Scap to Insig– a Membrane Protein That Facilitates Retention of Srebps in Er. Cell (2002) 110(4):489–500. doi: 10.1016/s0092–8674(02)00872–3

CrossRef Full Text | Google Scholar

45. Liu X, Li Y, Wang L, Zhao Q, Lu X, Huang J, et al. The Insig1 Gene, Not the Insig2 Gene, Associated With Coronary Heart Disease: Tagsnps and Haplotype–Based Association Study. The Beijing Atherosclerosis Study. Thromb Haemost (2008) 100(5):886–92. doi: 10.1160/TH08-01-0050

CrossRef Full Text | Google Scholar

46. Orekhov AN, Pushkarsky T, Oishi Y, Nikiforov NG, Zhelankin AV, Dubrovsky L, et al. Hdl Activates Expression of Genes Stimulating Cholesterol Efflux in Human Monocyte–Derived Macrophages. Exp Mol Pathol (2018) 105(2):202–7. doi: 10.1016/j.yexmp.2018.08.003

CrossRef Full Text | Google Scholar

47. Theusch E, Kim K, Stevens K, Smith JD, Chen YI, Rotter JI, et al. Statin–Induced Expression Change of Insig1 in Lymphoblastoid Cell Lines Correlates With Plasma Triglyceride Statin Response in a Sex–Specific Manner. Pharmacogen J (2017) 17(3):222–9. doi: 10.1038/tpj.2016.12

CrossRef Full Text | Google Scholar

48. Feng J, Huang C, Wren JD, Wang DW, Yan J, Zhang J, et al. Tetraspanin Cd82: A Suppressor of Solid Tumors and a Modulator of Membrane Heterogeneity. Cancer Metastasis Rev (2015) 34(4):619–33. doi: 10.1007/s10555–015–9585–x

CrossRef Full Text | Google Scholar

49. Akin I, Nienaber CA. "Obesity Paradox" in Coronary Artery Disease. World J Cardiol (2015) 7(10):603–8. doi: 10.4330/wjc.v7.i10.603

CrossRef Full Text | Google Scholar

50. Yeung L, Hickey MJ, Wright MD. The Many and Varied Roles of Tetraspanins in Immune Cell Recruitment and Migration. Front Immunol (2018) 9:1644. doi: 10.3389/fimmu.2018.01644

CrossRef Full Text | Google Scholar

51. Li F, Mao A, Fu X, She Y, Wei X. Correlation Between Seps1 Gene Polymorphism and Type 2 Diabetes Mellitus: A Preliminary Study. J Clin Lab Anal (2019) 33(8):e22967. doi: 10.1002/jcla.22967

CrossRef Full Text | Google Scholar

52. Zhao L, Zheng YY, Chen Y, Ma YT, Yang YN, Li XM, et al. Association of Genetic Polymorphisms of Sels With Type 2 Diabetes in a Chinese Population. Biosci Rep (2018) 38(6): BSR20181696. doi: 10.1042/bsr20181696

CrossRef Full Text | Google Scholar

53. Yu SS, Men LL, Wu JL, Huang LW, Xing Q, Yao JJ, et al. The Source of Circulating Selenoprotein S and Its Association With Type 2 Diabetes Mellitus and Atherosclerosis: A Preliminary Study. Cardiovasc Diabetol (2016) 15:70. doi: 10.1186/s12933–016–0388–3

CrossRef Full Text | Google Scholar

54. Schurch NJ, Schofield P, Gierliński M, Cole C, Sherstnev A, Singh V, et al. How Many Biological Replicates Are Needed in an Rna–Seq Experiment and Which Differential Expression Tool Should You Use? RNA (New York NY) (2016) 22(6):839–51. doi: 10.1261/rna.053959.115

CrossRef Full Text | Google Scholar

55. Aihara K, Liu R, Koizumi K, Liu X, Chen L. Dynamical Network Biomarkers: Theory and Applications. Gene (2022) 808:145997. doi: 10.1016/j.gene.2021.145997

CrossRef Full Text | Google Scholar