A construction and comprehensive analysis of the immune-related core ceRNA network and infiltrating immune cells in peripheral arterial occlusive disease

Background: Peripheral arterial occlusive disease (PAOD) is a peripheral artery disorder that increases with age and often leads to an elevated risk of cardiovascular events. The purposes of this study were to explore the underlying competing endogenous RNA (ceRNA)-related mechanism of PAOD and identify the corresponding immune cell infiltration patterns.

Methods: An available gene expression profile (GSE57691 datasets) was downloaded from the GEO database. Differentially expressed (DE) mRNAs and lncRNAs were screened between 9 PAOD and 10 control samples. Then, the lncRNA-miRNA-mRNA ceRNA network was constructed on the basis of the interactions generated from the miRcode, TargetScan, miRDB, and miRTarBase databases. The functional enrichment and protein–protein interaction analyses of mRNAs in the ceRNA network were performed. Immune-related core mRNAs were screened out through the Venn method. The compositional patterns of the 22 types of immune cell fraction in PAOD were estimated through the CIBERSORT algorithm. The final ceRNA network and immune infiltration were validated using clinical tissue samples. Finally, the correlation between immune cells and mRNAs in the final ceRNA network was analyzed.

Results: Totally, 67 DE_lncRNAs and 1197 DE_mRNAs were identified, of which 130 DE_mRNAs (91 downregulated and 39 upregulated) were lncRNA-related. The gene ontology enrichment analysis showed that those down- and upregulated genes were involved in dephosphorylation and regulation of translation, respectively. The final immune-related core ceRNA network included one lncRNA (LINC00221), two miRNAs (miR-17-5p and miR-20b-5p), and one mRNA (CREB1). Meanwhile, we found that monocytes and M1 macrophages were the main immune cell subpopulations in PAOD. After verification, these predictions were consistent with experimental results. Moreover, CREB1 was positively correlated with naive B cells (R = 0.55, p = 0.035) and monocytes (R = 0.52, p = 0.049) and negatively correlated with M1 macrophages (R = −0.72, p = 0.004), resting mast cells (R = −0.66, p = 0.009), memory B cells (R = −0.55, p = 0.035), and plasma cells (R = −0.52, p = 0.047).

Conclusion: In general, we proposed that the immune-related core ceRNA network (LINC00221, miR-17-5p, miR-20b-5p, and CREB1) and infiltrating immune cells (monocytes and M1 macrophages) could help further explore the molecular mechanisms of PAOD.

Introduction

Peripheral arterial occlusive disease (PAOD) is an atherosclerotic condition involving non-cardiac and non-cerebral arteries. Nowadays, it has developed into a widespread disease with more than 200 million people affected worldwide and has become the third most common cause of death from cardiovascular disease (Fowkes et al., 2013). The importance of it is growing by virtue of its increasing incidence. Patients with PAOD often suffer from chronic limb ischemia that results in intermittent claudication, resting pain, disability, and even death. As PAOD is one common manifestation of systemic atherosclerosis, it is essential to study peripheral atherosclerosis for exploring the potential pathogenesis and progression of PAOD and effective therapeutic targets.

In recent years, more and more studies have shown that atherosclerosis, as a chronic inflammatory disease, is significantly associated with the infiltration of immune cells such as neutrophils, macrophages, T cells, and B cells into the inner layer of the vessel wall (Hansson and Hermansson, 2011; Baptista et al., 2018). In atherosclerosis, hypercholesterolemia leads to the accumulation of plasma low-density lipoprotein (LDL) in the artery wall, which stimulates and recruits monocytes and elicits local inflammation. Then, monocytes are infiltrated to differentiate locally into macrophages, and the lipid metabolic disorders and efferocytosis of macrophages are reduced, leading to irreversible inflammation (Tabas and Bornfeldt, 2016). Macrophages polarized to M1 and M2 exert pro-inflammatory and anti-inflammatory effects, respectively (Moore and Tabas, 2011). T cells account for about 40% of the total number of immune cells in human atherosclerotic lesions. Among them, regulatory T (Treg) cells produce transforming growth factor β, which inhibits the proliferation of T-helper type 1 (Th1) and T-helper type 17 (Th17) cells (Fantini et al., 2006). Th1 cells and natural killer (NK) cells secrete pro-inflammatory factors, which destroy collagen fibers and promote the transformation of atherosclerotic plaques to vulnerable phenotypes (Konkel et al., 2017). Th17 cells are a subtype of T cells, which can promote the formation of thick collagen fibers and contribute to the stability of plaques (Gisterå et al., 2013). Dendritic cells (DCs) can make the innate and adaptive immune responses to act as important modulators in atherosclerosis (Cybulsky et al., 2016). Several B-lymphocyte subsets contribute to the inflammatory process of atherosclerosis through cellular and humoral responses (Tsiantoulas et al., 2014). However, in PAOD, the landscape of immune cell infiltration has not been fully elucidated. Moreover, the relationship between immune-related genes and immune cells in PAOD is largely unknown.

Noncoding RNAs (ncRNAs) regulate gene expression at transcriptional and posttranscriptional levels without coding proteins in the transcriptome, including microRNAs (miRNAs) and long noncoding RNAs (lncRNAs). miRNAs are a class of highly conserved single-stranded noncoding small RNAs, which contain approximately 19–25 nucleotides and have post-transcriptional regulatory activity. lncRNAs are defined as a type of ncRNAs that are longer than 200 nucleotides in length, with multilevel regulatory functions in gene expression, such as transcription, translation, and epigenetics. Accumulating evidence has shown that functional ncRNAs play an important role in the pathogenesis and progression of many diseases, such as cancer, digestive system diseases, and cardiovascular diseases. In recent years, a competing endogenous RNA (ceRNA) network hypothesis has been proposed (Salmena et al., 2011). In the ceRNA network, lncRNAs can serve as endogenous molecular sponges for miRNAs to regulate the expression of messenger RNAs (mRNAs) indirectly. In this way, the function of ncRNAs can be linked to the function of mRNAs that encode proteins. Given their complexity, the dysregulation of the lncRNA-miRNA-mRNA network is closely related to the pathogenesis and progression of many human diseases, such as cardiovascular diseases. For instance, Ye et al. (2019) found that lncRNA MIAT upregulates CD47 expression by sponging miR-149-5p and inhibits efferocytosis in advanced atherosclerosis. Yang et al. (2021) discovered that the lncRNA XIST serves as a ceRNA and promotes atherosclerosis by increasing miR-599-mediated expression of TLR4. Nevertheless, few data-based studies have been conducted to analyze the relationship between the immune-related ceRNA regulatory network and infiltrating immune cells in PAOD.

In this study, we compared differentially expressed (DE) mRNAs and lncRNAs between 9 PAOD and 10 control samples downloaded from the Gene Expression Omnibus (GEO) database. Then, the target miRNAs of DE_lncRNAs and DE_mRNAs of target miRNAs were predicted. Subsequently, protein–protein interaction (PPI) analysis among the predicted DE_mRNAs was conducted, and hub DE_mRNAs were identified by the Cytoscape’s cytoHubba plugin. The overlapping genes between the hub DE_mRNAs and immune-related genes were identified as core mRNAs to construct the potential immune-related core ceRNA regulatory network of PAOD. Meanwhile, the CIBERSORT method was used to analyze the different patterns of immune cell infiltration in PAOD. The immune-related core ceRNA network and immune infiltration were validated using clinical tissue samples. Finally, we investigated and visualized the correlation between the core mRNAs and infiltrating immune cells in an effort to better understand the molecular immune mechanism during the progression of PAOD.

Materials and methods

Data acquisition

In this study, the microarray dataset GSE57691 (Biros et al., 2015) that assessed the relative gene expression in human abdominal aortic aneurysm (AAA) and PAOD, and GSE137580 that studied the global miRNAs in atherosclerotic models that oxidative LDL treated human aortic endothelial cells (HAEC) were obtained from the Gene Expression Omnibus database (GEO, http://www.ncbi.nlm.nih.gov/geo) (Barrett et al., 2013). The specimens of GSE57691 were obtained from 20 patients with small AAA, 29 patients with large AAA, 9 PAOD patients, and 10 organ donors. Then, the data on PAOD and normal artery were picked out for further analysis. The test platforms of GSE57691 and GSE137580 were GPL10558 Illumina HumanHT-12 V4.0 expression beadchip and GPL24741 Agilent-070156 Human_miRNA_V21.0_Microarray 046064 (gene name version), respectively. The GSE137580 data were used to validate the relative expression level of miRNAs in the potential ceRNA network. Additionally, immune-related genes were downloaded from the Immunology Database and Analysis Portal (ImmPort) database (http://www.immport.org/) (Bhattacharya et al., 2014).

Differential expression analysis

First, the probe sets were converted into corresponding gene symbols according to the platform profile with annotation information. If multiple probe sets correspond to the same gene, their mean value was calculated by R software (version 4.1.0). Then, based on the gene annotation information included in the ENSEMBL database (https://asia.ensembl.org/Homo_sapiens/Info/Index), the expression profile dataset was divided into lncRNA and mRNA groups. The linear models for the microarray data (LIMMA) package of R software were utilized to normalize raw data and perform DE_RNA analysis between PAOD and normal artery groups (Ritchie et al., 2015). p-values were adjusted by the Benjamini–Hochberg (BH) false discovery rate (FDR) method (Glickman et al., 2014). The cut-off value of DE_RNAs was set as adj. p-value < 0.01 and |fold change (FC)| > 1.5 (Meng et al., 2019). The heatmap of the DE_lncRNAs and the volcano plot of all RNAs were constructed for data visualization by the “pheatmap” (https://CRAN.R-project.org/package=pheatmap) and “ggplot2” packages in R software, respectively (Ginestet, 2011).

Prediction of lncRNA-miRNA-mRNA interactions

The highly conserved miRNA families of the miRcode database (http://www.mircode.org/) were applied to predict interactions between DE_lncRNAs and potential miRNAs (Jeggari et al., 2012). Subsequently, the TargetScan (http://www.targetscan.org/) (Agarwal et al., 2015), miRDB (http://www.mirdb.org/) (Chen and Wang, 2020), and miRTarBase (http://mirtarbase.mbc.nctu.edu.tw/) (Huang et al., 2020) databases were utilized to forecast miRNA–mRNA pairs. Only those genes that concurrently existed in all three databases were considered as candidate targets of mRNAs for further analysis.

Venn method

The Venn method was employed to analyze overlapping genes. Intersections between DE_lncRNA predicted mRNAs and DE_mRNAs, as well as immune-related genes and lncRNA-related DE_mRNAs (the intersections between DE_lncRNA predicted mRNAs and DE_mRNAs), were calculated using an online tool. (http://bioinformatics.psb.ugent.be/webtools/Venn/). Moreover, the tool was also used to identify the hub genes.

Functional enrichment and protein–protein interaction (PPI) analysis

First, the lncRNA-related DE_mRNAs were divided into expression upregulated and downregulated groups. Then, the “clusterProfiler” package (Yu et al., 2012) in R software was used to perform gene ontology (GO) (Gene Ontology Consortium, 2006) and Kyoto Encyclopedia of Genes and Genomes (KEGG) (Kanehisa et al., 2010) enrichment analyses. GO enrichment analysis included three categories: cellular component (CC), molecular function (MF), and biological process (BP). The “ggplot2” package in R software was utilized to draw the bubble charts for visualization of the results of the GO and KEGG enrichment analyses. p-value < 0.05 was considered significantly enriched when screening. Subsequently, the online Search Tool for the Retrieval of Interacting Genes (STRING, https://string-db.org/) database (Szklarczyk et al., 2019) was used to determine the relationship between the lncRNA-related DE_mRNAs. The minimum required interaction score was 0.7, and the disconnected nodes in the network were hidden. Then, Cytoscape software (version 3.8.2) was used to develop the PPI network. Furthermore, the Maximal Clique Centrality (MCC), Degree, and Maximum Neighborhood Component (MNC) algorithms in the cytoHubba plugin were used to screen out the top 10 mRNAs in the PPI network. The overlapping genes obtained by the three different aforementioned algorithms were identified as hub genes. Finally, we identified the overlapping genes between the hub genes before and the immune-related genes as immune-related core genes.

Construction of the immune-related core ceRNA network

First, the interactions between lncRNAs, miRNAs, and immune-related core mRNAs were confirmed, as described in the aforementioned item. Second, the immune-related lncRNA-miRNA-mRNA ceRNA network was developed, and the “ggalluvial” package (http://corybrunson.github.io/ggalluvial/) in R software was used to draw a Sankey diagram for data visualization. Subsequently, a correlation analysis between lncRNA and immune-related core genes in the ceRNA network was performed by the “LIMMA” package in R software. Ultimately, the relative expression levels of miRNAs in immune-related core ceRNA were validated in the GSE137580 dataset and visualized by GraphPad Prism (version 7.0). p-value < 0.05 was considered statistically significant.

Estimation of immune cell infiltration

CIBERSORT (https://cibersort.stanford.edu/) is a versatile analytical tool that uses gene expression data to quantify the cell fractions from complex tissues and has been confirmed by flow cytometry (Newman et al., 2015). To analyze the proportion of 22 infiltrating immune cells in atherosclerotic plaques of PAOD patients and normal controls, the mRNA expression data were uploaded to the CIBERSORT platform. Only samples that had a CIBERSORT algorithm output of p-value < 0.05 were filtered out, and the immune cell infiltration matrix was obtained for further analysis. Histograms and heatmaps were drawn to show the rate of immune cell infiltration in different samples. Subsequently, the Wilcoxon rank-sum test was performed to assess the differential composition of infiltrating immune cells between PAOD patients and controls. Results were visualized by the “pheatmap” and “vioplot” (https://github.com/TomKellyGenetics/vioplot) packages in R software. Furthermore, Pearson’s correlation analysis was adopted to explore the correlation among 22 immune cell subtypes. A correlation heatmap was drawn by the “corrplot” package (https://github.com/taiyun/corrplot) in R software to visualize the correlation analysis results. Finally, the “ggstatsplot” package (https://github.com/IndrajeetPatil/ggstatsplot) was used to perform correlation analysis on the immune-related core DE_mRNAs and infiltrating immune cells, and the “ggplot2” package was used to visualize the results in R software.

Real-time quantitative PCR

From June 2021 to February 2022, we recruited five patients with PAOD who underwent femoral endarterectomy and five organ donors who donated normal iliac arteries in the First Affiliated Hospital of Anhui Medical University. All tissue specimens were divided into two parts and frozen in liquid nitrogen immediately when they were isolated. One part of these samples was pretreated, and total RNA was extracted with TRIzol reagent (Invitrogen Life Technologies, United States) for performing real-time quantitative polymerase chain reaction (RT-qPCR). Spectrophotometry was used to measure the purity of RNA. RNAs were reverse transcribed into complementary DNAs using the Bestar qPCR RT Kit (DBI, Germany), following the instructions of the manufacturer. Subsequently, complementary DNA was amplified by RT-qPCR using an Applied Biosystems SYBR Green mix kit (ABI, United States). GAPDH was used as an internal reference for lncRNAs and mRNAs, while U6 was used as a reference for miRNAs. Primer sequences were obtained from PrimerBank and miRprimer2 databases (Supplementary Table S1). The reactions were measured on the ABI 7900HT Real-Time PCR system (ABI, United States), and the 2^−ΔΔCT method was used for analysis.

Hematoxylin and eosin (H&E) and immunofluorescence staining

The other part of those samples was demineralized after fixation in 4% paraformaldehyde and then embedded in paraffin. All samples were cut into 4-μm slices for further staining. H&E staining was performed to assess the atherosclerotic lesions. The protein expression levels of CD68 and iNOS were analyzed by immunofluorescence staining. Antibodies (i.e., CD68 and iNOS) were purchased from Proteintech (Chicago, IL, United States). All procedures were conducted according to the recommendations of the manufacturer. Images were observed with a fluorescence microscope (Leica DMI6000B, Germany) and analyzed by ImageJ software.

Statistical analysis

The data are presented as the mean ± standard deviation (SD). SPSS (SPSS Inc., Chicago, IL, United States, version: 19.0) was used to conduct statistical analysis. Student’s t-test was used for comparisons between two groups. A p-value of less than 0.05 was considered statistically significant. All experiments were performed at least three times.

Results

Identification of DE mRNAs and lncRNAs

In order to clarify the process of this research and make it easier for readers to read, a schematic representation is provided in Figure 1. Raw data were downloaded from the GSE57691 dataset in the GEO database. In total, RNA-seq data from 9 PAOD and 10 normal tissues were analyzed using criteria of |fold change (FC)| > 1.5 and adj. p-value < 0.01. A total of 2,142 lncRNAs and 18,219 mRNAs were re-annotated according to the platform profile with annotation information. Moreover, we identified 1264 DE RNAs, including 67 lncRNAs (27 downregulated and 40 upregulated) and 1197 mRNAs (752 down-regulated and 445 up-regulated) meeting the thresholds. Then, the corresponding volcano plot and heatmap of DE_lncRNAs are shown in Figure 2.

FIGURE 1

FIGURE 1. Schematic representation of our analytic process.

FIGURE 2

FIGURE 2. Differentially expressed (DE) mRNAs and lncRNAs in tissues between PAOD patients and normal controls. (A) Volcano plot of DE_mRNAs. (B) Heatmap of DE_lncRNAs.

Functional enrichment analysis of lncRNA-related DE_mRNAs

In order to construct the ceRNA network, DE_lncRNAs were further analyzed. The miRcode database was employed to predict potential DE_lncRNA-targeted miRNAs. Then, potential miRNA–mRNA pairs were analyzed using the miRTarBase, TargetScan, and miRDB databases. A total of 34 miRNAs were identified as DE_lncRNA-predicted miRNAs, and 130 mRNAs (91 downregulated and 39 upregulated) were identified as DE_lncRNA-predicted mRNAs. Subsequently, the Venn method was used to analyze the overlapping genes between DE_mRNAs and DE_lncRNA-predicted mRNAs. As shown in Figure 3A, the intersection contains all of the DE_lncRNA-predicted mRNAs. These mRNAs were also called lncRNA-related DE_mRNAs. To determine the potential mechanisms of lncRNA-related DE_mRNAs, these mRNAs were divided into down- and upregulated groups for further GO and KEGG enrichment analyses (Figures 4A–D). A biological process analysis showed that mRNAs in the downregulated group were significantly enriched in dephosphorylation, regulation of autophagy, and inositol phosphate catabolic processes, while mRNAs in the upregulated group were mostly enriched in regulation of translation, positive regulation of the cellular catabolic process, and positive regulation of fat cell differentiation. A cellular component analysis showed that mRNAs in the downregulated group were significantly enriched in microtubules, while mRNAs in the upregulated group were mostly enriched in the endoplasmic reticulum–Golgi intermediate compartment. A molecular function analysis showed that mRNAs in the downregulated group were significantly enriched in ubiquitin-like protein ligase binding and phosphatase binding, while mRNAs in the upregulated group were mostly enriched in Rho GTPase binding and vinculin binding. The KEGG pathway enrichment analysis showed mRNAs in the downregulated group were significantly enriched in the PI3K-Akt signaling pathway, cellular senescence, and regulation of the actin cytoskeleton. However, mRNAs in the upregulated group were mostly enriched in human cytomegalovirus infection and the cGMP-PKG signaling pathway.

FIGURE 3

FIGURE 3. Identification of lncRNA-related DE-mRNAs (A) and hub genes (B). The details of intersection in the Venn diagram are listed as follows. Construction of the immune-related core ceRNA network in PAOD (C). MCC, maximal clique centrality; MNC, maximum neighborhood component.

FIGURE 4

FIGURE 4. Functional enrichment analysis of lncRNA-related DE_mRNAs. (A,B) GO enrichment analysis of downregulated and upregulated lncRNA-related DE-mRNAs. (C,D) KEGG enrichment analysis of downregulated and upregulated lncRNA-related DE-mRNAs.

PPI network construction and hub gene identification

The PPI network of lncRNA-related DE_mRNAs containing 130 nodes and 61 edges was constructed based on the STRING online database and visualized by Cytoscape software (Figure 5). Subsequently, Cytoscape’s plugin cytoHubba was used to identify the top 10 genes based on three commonly used classification methods (MCC, Degree, and MNC) (Supplementary Table S1). By overlapping these genes, seven hub genes (YWHAH, PTEN, PIK3R1, FOXO3, YWHAE, CREB1, and FOXO1) were consequently identified, as shown in Figure 3B.

FIGURE 5

FIGURE 5. PPI network analysis. The PPI network consisting of 130 nodes and 61 edges was visualized in Cytoscape. Red and green notes symbolize the upregulated and downregulated genes. The size of each node is positively correlated with its degree value.

Construction of the immune-related core ceRNA network

To construct the immune-related core ceRNA network, the Venn method was used to analyze the intersection between lncRNA-related DE_mRNAs and immune-related genes obtained from the ImmPort database. Consequently, a total of seven genes (TNFAIP3, PLXND1, CREB1, PIK3R1, HMGB1, PLAU, and KRAS) were identified as lncRNA-related immune DE_mRNAs (Figure 3A). Subsequently, two genes (PIK3R1 and CREB1) were identified as immune-related core genes by overlapping the lncRNA-related immune DE_mRNAs and the seven hub genes acquired in Cytoscape software. After that, immune-related core DE_mRNAs and their paired miRNAs and lncRNAs were chosen to develop the ceRNA regulatory network. In total, the immune-related core ceRNA network contained eight lncRNAs, three miRNAs, and two mRNAs (Figure 3C). Then, the correlation between the expression of immune-related core DE_mRNAs and their paired lncRNAs was analyzed. The results illustrated that the expression of CREB1 was positively correlated with LINC00221 (R = 0.861, p < 0.001) and MEG3 (R = 0.492, p = 0.033) (Figures 6A,B). Similarly, PIK3R1 was positively correlated with that of HSPB2-C11orf52 (R = 0.508, p = 0.026), (Figure 6C). Additionally, the relative expression levels of miRNAs in the potential immune-related core ceRNA network were validated in the GSE137580 dataset and visualized by GraphPad Prism 7. As shown in Figures 6D–F, compared with the negative control group, the relative expression levels of miR-107, miR-20b-5p, and miR-17-5p were all low in the atherosclerosis group, due to the expression trend of miR-20b-5p and miR-17-5p being consistent with that of prediction, while that of miR-107 was opposite to that of prediction. Moreover, the correlation ship between the expression of LINC00221 and CREB1 was very close. Hence, the final potential immune-related core ceRNA network in this study contained one lncRNA (LINC00221), two miRNAs (miR-20b-5p and miR-17-5p), and one mRNA (CREB1) in Table 2.

FIGURE 6

FIGURE 6. Correlation between the expression of immune-related core DE_mRNAs and their paired lncRNAs (A–C). Validation of the miRNAs of the immune-related core ceRNA network in GSE137580 (D–F). The results with statistically significant differences (p < 0.05) are shown.

Composition of infiltrating immune cells

The composition of 22 infiltrating immune cells in atherosclerosis tissues of PAOD patients and normal controls was estimated using the CIBERSORT algorithm (Figures 7A,B). Since the output p-values of GSM1386842, GSM1386843, GSM1386849, and GSE1386836 were greater than 0.05, they were excluded for further analysis. The distribution of 22 immune cell types in each sample varied significantly, among which M2 macrophages accounted for the largest proportion. The relationships among 22 immune cells are presented in Figure 7C. Monocytes were negatively correlated with M1 macrophages (R = −0.75). Activated mast cells were positively correlated with eosinophils (R = 0.75) and activated dendritic cells (R = 0.74). Activated memory CD4 T cells were positively correlated with naive CD4 T cells (R = 0.75) and activated dendritic cells (R = 0.71). Naive CD4 T cells were positively correlated with activated dendritic cells (R = 0.98). Memory B cells were positively correlated with plasma cells (R = 0.88). Other immune cell subpopulations were weakly to moderately correlated. The violin plot of the immune cell infiltration difference showed that, compared with the normal control sample, two types of immune cells, monocytes and M1 macrophages, were differentially expressed. Monocytes were upregulated, while M1 macrophages were downregulated in atherosclerosis tissues of PAOD patients (Figure 7D).

FIGURE 7

FIGURE 7. Composition of infiltrating immune cells assessed using the CIBERSORT algorithm in PAOD tissues. (A) Distribution of immune cell infiltration in each sample. (B) Heatmap of immune cell types. (C) Correlation among infiltrating immune cells. (D) Violin plot of infiltrating immune cells.

Biological experiments

To validate the immune-related core ceRNA network in PAOD, RT-qPCR was used to detect the expression levels of the core genes. As it is shown in Figures 8A–D, compared with the normal control group, the expression of LINC00221 and CREB1 in the PAOD group was increased, while the expression of miR-20b-5p and miR-17-5p was decreased (all p-values < 0.05). Intimal structure disorder, lipid infiltration, and smooth muscle cell atrophy after endarterectomy of the femoral artery are found in Figure 8E in comparison with normal arterial intima. Subsequently, immunofluorescence staining of monocyte- and macrophage-associated molecules CD68 (monocytes) and iNOS (M1) showed that monocyte and M1 macrophage infiltration were significantly increased and decreased, respectively, in the PAOD group (all p-values < 0.05). The results of positive area analysis of CD68 and iNOS are shown in Figures 8F,G. These results are consistent with the bioinformatics results we predicted.

FIGURE 8

FIGURE 8. RT-qPCR validation of the immune-related core ceRNA network and evaluation of monocytes and M1 macrophages in tissue samples. Quantification of the relative expression levels of LINC00221 (A), miR-17-5p (B), miR-20b-5p (C), and CREB1 (D) using RT-qPCR. (E) Representative pictures of H&E and immunofluorescence staining (magnification, ×20). Quantification of molecules in CD68 (monocyte) Fig: (F) and iNOS (M1 macrophage) Fig: (G) positive areas.

Correlation analysis

Finally, the correlation between CREB1 and infiltrating immune cells was estimated. In this analysis, the Wilcoxon test was adopted and significantly correlated pairs with a p-value < 0.05. As shown in Figure 9, CREB1 was positively correlated with naive B cells (R = 0.55, p = 0.035) and monocytes (R = 0.52, p = 0.049) and negatively correlated with M1 macrophages (R = −0.72, p = 0.004), resting mast cells (R = −0.66, p = 0.009), memory B cells (R = −0.55, p = 0.035), and plasma cells (R = −0.52, p = 0.047).

FIGURE 9

FIGURE 9. Correlation between CREB1 and infiltrating immune cells. The size of each dot is positively correlated with the strength of the correlation between genes and immune cells.

Discussion

PAOD is a peripheral artery disease increasing with age and often leads to distal limb ischemia that results in reduced quality of life and death. Despite improvements in surgical, interventional, and pharmacological therapy of PAOD, the progression of atherosclerosis is still not prevented or curbed efficiently. Hence, in order to improve the treatment of PAOD, it is of vital significance to explore the underlying molecular mechanisms. For the past few years, the hypothesis of the ceRNA network has greatly raised the interest of researchers, and it makes the link between ncRNAs and mRNAs. On the basis of the ceRNA hypothesis, lncRNAs can act as miRNA sponges to regulate the expression of mRNAs. In this study, for exploring the potential pathogenesis of PAOD, we constructed the immune-related ceRNA regulatory network of PAOD based on the GSE57691 microarray dataset. After comprehensive analysis and careful verification, the final immune-related core ceRNA network, including one lncRNA (LINC00221), two miRNAs (miR-17-5p and miR-20b-5p), and one mRNA (CREB1), was identified. The interaction between LINC00221, miR-17-5p, and miR-20b-5p may regulate peripheral atherosclerosis via CREB1.

As post-transcriptional regulators, miRNAs play an important role in influencing the expression of downstream target genes, which are involved in a variety of physiological and pathological processes (Valencia-Sanchez et al., 2006; Bartel, 2009). During recent years, the dysregulation of miRNA expression has been found to be associated with the pathogenesis and progression of many diseases, including atherosclerosis. In fact, several miRNAs have been found to be associated with the pathogenesis of atherosclerosis. Tan et al. (2019) demonstrated that inhibiting the expression of miR-17-5p can suppress inflammation and reduce lipid accumulation in atherosclerosis. An et al. (2019) reported that the expression of miR-17-5p is significantly decreased, and lncRNA SNHG16 can promote proliferation and inflammatory response of macrophages through the miR-17-5p/NF-κB signaling pathway in patients with atherosclerosis. Shen et al. (2019) identified that the expression of circRNA0044073 was upregulated, and the expression of miR-107 was downregulated in atherosclerotic blood cells. Moreover, circRNA-0044073 can suppress the levels of miR-107 via a sponge mechanism and increase the proliferation and invasion of cells in atherosclerosis. As for miR-20b-5p, the reports related to its role are mainly in various cancers, but there are few reports on atherosclerosis. More research studies are needed on the role of miR-20b-5p in atherosclerosis in the future.

The CREB1 gene encodes cyclic-AMP response-binding protein 1, which is a transcription factor that is a member of the leucine zipper family of DNA-binding proteins. CREB1 has been shown to be involved in both positive and negative regulation of atherosclerosis. On the one hand, CREB1 upregulation is observed in the vessels isolated from normal mice compared to atherosclerotic mice (Schauer et al., 2010). On the other hand, CREB1 can activate the pro-inflammatory cytokine IL-17 that is directly responsible for macrophage accumulation and the ensuing inflammation in the atherosclerotic plaque in mice (Kotla et al., 2013). Being equally ambiguous is the role of CREB1 in the endothelial dynamic balance. A large body of data seems to suggest that the deletion of CREB1 in endothelial cells may result in an enhanced inflammatory response and barrier dysfunction (Chava et al., 2012; Xiong et al., 2020). On the contrary, CREB1 can promote leukocyte adhesion by directly binding to human umbilical vein endothelial cell ICAM1 and activating its transcription (Hadad et al., 2011). Similarly, CREB1 interacts with BAF47 (BRG1-associated factor 47) and recruits BAF47 to the proximal neogenin 1 promoter, leading to neogenin 1 trans-activation that contributes to endothelial dysfunction (Li et al., 2022). These discrepant roles of CREB1 may allude to its coupling to various signaling pathways targeting either the stimulation or suppression of progression in atherosclerosis. Spatiotemporally controlled CREB1 transgenic animal models should be employed in future studies for delineating the implied role of CREB1 in atherosclerosis.

Additionally, in this study, we described the infiltrating immune cells in peripheral arterial plaques, analyzed the differences in the abundance of immune cells between the groups, and estimated the correlation between the immune-related core genes (CREB1) and infiltrating immune cells. The results demonstrated that monocytes and macrophages were the main immune cell subpopulations in atherosclerosis of peripheral artery tissues. In addition, atherosclerotic plaques had increased infiltration of monocytes while having decreased proportions of M1 macrophages compared to normal artery tissues. The result that M1 macrophages are found at higher levels in normal tissues appears to be in contrast with the current knowledge of proatherogenic cells. The immune cell number is a relative percentage among each group. Atherosclerotic tissues have more immune cell infiltration, so there is a relatively lower content of M1 macrophages in these tissues than in the control tissues. Owing to an increased infiltration of immune cells, the diminution in percentage may not manifest a decreased number of M1 macrophages. If the infiltration score is readjusted, the absolute percentage of M1 macrophages in the PAOD group may not be lower than that in the normal group.

Monocytes are derived from bone marrow-derived progenitor cells, and the early stage of monocytes development may be regulated by the content of cellular cholesterol in a manner that can affect atherogenesis. Recent insight suggests that the key initiating step of incipient atherogenesis in human and animal models indicates subendothelial accumulation of apolipoprotein B-containing lipoproteins (apoB-LPs) (Williams and Tabas, 1995). The pivotal early inflammatory response to accumulated apoB-LPs is the activation of surficial endothelial cells, which leads to the recruitment of circulating monocytes (Mestas and Ley, 2008). Activated endothelial cells secrete chemokines and interact with cognate chemokine receptors on monocytes in a manner that promotes monocytes’ migration into the intima, where they differentiate into macrophages and phagocytize lipoproteins, leading to foam cell formation. Importantly, atherogenesis can be prevented or retarded in mouse models of atherosclerosis through preventing monocyte infiltration by blocking chemokines or their receptors (Mestas and Ley, 2008). The functions of macrophages within plaques are shaped largely by external stimuli such as intracellular energy metabolism (Van den Bossche et al., 2017), gut microbiota metabolites (Wang et al., 2011), and genetic and epigenetic factors including ncRNAs (Erbilgin et al., 2013; Amit et al., 2016). Traditionally, macrophages are classified into pro-inflammatory M1 macrophages (activated by lipopolysaccharide and interferon-γ) and anti-inflammatory M2 macrophages (induced by interleukin-4 and interleukin-13) (Johnson and Newby, 2009; Rath et al., 2014). In general, M1 macrophages perform processes that promote atherosclerosis progression, whereas M2 macrophages carry out functions that can restrain plaque progression and facilitate plaque regression (Peled and Fisher, 2014). M1 macrophages, through secreting cytokines, proteases, and other factors, increase the cellular expansion of lesions and cause changes in plaque morphology that result in plaque rupture and acute lumenal thrombosis. Two key changes in plaque morphology promoted by M1 macrophages are plaque necrosis and fibrous cap thinning. A previous study has confirmed that M1 macrophages can secrete matrix metalloproteinases (MMPs), such as MMP9 and MMP2, which may contribute to plaque rupture, and another study shows that MMPs co-localize with M1 macrophages in advanced plaques (Huang et al., 2012). In this study, we found that the expression of the CREB1 gene was positively correlated with monocytes (R = 0.52, p = 0.049) and negatively correlated with M1 macrophages (R = −0.72, p = 0.004). Potentially, the crosstalk between the ceRNA network and the infiltration of immune cells plays a crucial part in regulating atherosclerosis progression. Further research studies are required to clarify the complex interactions between these genes and immune cells.

However, this article has some limitations. First, only one dataset (GSE57691) that is from the west and lacks ethnic diversity was used as a data source. Moreover, the analysis of immune cell infiltration was based on the CIBERSORT algorithm, the immune cell types of which were not comprehensive. Finally, the underlying regulatory mechanisms of the ceRNA network and immune cells were not elucidated clear enough, and more functional biological experiments with larger sample sizes are needed to further verify this in the future.

Conclusion

Taken together, in this study, the immune-related ceRNA network, by the composition of one lncRNA (LINC00221), two miRNAs (miR-17-5p and miR-20b-5p), and one mRNA (CREB1), was first constructed in PAOD. Afterward, the immune cell infiltration analysis was performed to estimate the abundance and differences of different immune cells. The final results show that monocytes and M1 macrophages were considered to be important immune cells associated with PAOD formation. Moreover, the expression of the CREB1 gene was positively correlated with monocytes (R = 0.52, p = 0.049) and negatively correlated with M1 macrophages (R = −0.72, p = 0.004). These findings provide new insights into the pathogenesis and progression of PAOD and novel potential therapeutic targets. Perhaps in the future, new drugs can be developed for these novel potential therapeutic targets to delay the progression of PAOD and improve the long-term patency rate of vascular lumen in PAOD patients undergoing surgery or interventional therapy.

Data availability statement

The original contributions presented in the study are included in the article/Supplementary Material; further inquiries can be directed to the corresponding author.

Ethics statement

The studies involving human participants were reviewed and approved by the Medical Ethics Committee of Anhui Medical University (No. 20210015).

Author contributions

HO and YL contributed to the conception and design of the study. ZC and JX performed data analysis for the study. BZ contributed to software analysis. JL conducted the biological experiments. ZC and JX prepared the original manuscript. HO and YL reviewed and edited the manuscript. All authors contributed to the manuscript revision, read, and approved the final version.

Funding

This research was supported by the Natural Science Foundation of Anhui Province of China (Grant No. 2108085QH308) and the National Natural Science Foundation of China (Grant No. 81874063).

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/fgene.2022.951537/full#supplementary-material

References

Agarwal, V., Bell, G. W., Nam, J. W., and Bartel, D. P. (2015). Predicting effective microRNA target sites in mammalian mRNAs. Elife 4. e05005. doi:10.7554/eLife.05005 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Amit, I., Winter, D. R., and Jung, S. (2016). The role of the local environment and epigenetics in shaping macrophage identity and their effect on tissue homeostasis. Nat. Immunol. 17, 18–25. doi:10.1038/ni.3325 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

An, J. H., Chen, Z. Y., Ma, Q. L., Wang, H. J., Zhang, J. Q., and Shi, F. W. (2019). LncRNA SNHG16 promoted proliferation and inflammatory response of macrophages through miR-17-5p/NF-κB signaling pathway in patients with atherosclerosis. Eur. Rev. Med. Pharmacol. Sci. 23, 8665–8677. doi:10.26355/eurrev_201910_19184 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Baptista, D., Mach, F., and Brandt, K. J. (2018). Follicular regulatory T cell in atherosclerosis. J. Leukoc. Biol. 104, 925–930. doi:10.1002/JLB.MR1117-469R | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Barrett, T., Wilhite, S. E., Ledoux, P., Evangelista, C., Kim, I. F., Tomashevsky, M., et al. (2013). NCBI geo: Archive for functional genomics data sets-update. Nucleic Acids Res. 41, D991–D995. doi:10.1093/nar/gks1193 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Bartel, D. P. (2009). MicroRNAs: Target recognition and regulatory functions. Cell 136, 215–233. doi:10.1016/j.cell.2009.01.002 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Bhattacharya, S., Andorf, S., Gomes, L., Dunn, P., Schaefer, H., Pontius, J., et al. (2014). ImmPort: Disseminating data to the public for the future of immunology. Immunol. Res. 58, 234–239. doi:10.1007/s12026-014-8516-1 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Biros, E., Gäbel, G., Moran, C. S., Schreurs, C., Lindeman, J. H., Walker, P. J., et al. (2015). Differential gene expression in human abdominal aortic aneurysm and aortic occlusive disease. Oncotarget 6, 12984–12996. doi:10.18632/oncotarget.3848 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Chava, K. R., Tauseef, M., Sharma, T., and Mehta, D. (2012). Cyclic AMP response element-binding protein prevents endothelial permeability increase through transcriptional controlling p190RhoGAP expression. Blood 119, 308–319. doi:10.1182/blood-2011-02-339473 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, Y., and Wang, X. (2020). miRDB: an online database for prediction of functional microRNA targets. Nucleic Acids Res. 48, D127–D131. doi:10.1093/nar/gkz757 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Cybulsky, M. I., Cheong, C., and Robbins, C. S. (2016). Macrophages and dendritic cells: Partners in atherogenesis. Circ. Res. 118, 637–652. doi:10.1161/CIRCRESAHA.115.306542 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Erbilgin, A., Civelek, M., Romanoski, C. E., Pan, C., Hagopian, R., Berliner, J. A., et al. (2013). Identification of CAD candidate genes in GWAS loci and their expression in vascular cells. J. Lipid Res. 54, 1894–1905. doi:10.1194/jlr.M037085 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Fantini, M. C., Becker, C., Tubbe, I., Nikolaev, A., Lehr, H. A., Galle, P., et al. (2006). Transforming growth factor beta induced FoxP3+ regulatory T cells suppress Th1 mediated experimental colitis. Gut 55, 671–680. doi:10.1136/gut.2005.072801 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Fowkes, F. G., Rudan, D., Rudan, I., Aboyans, V., Denenberg, J. O., McDermott, M. M., et al. (2013). Comparison of global estimates of prevalence and risk factors for peripheral artery disease in 2000 and 2010: A systematic review and analysis. Lancet 382, 1329–1340. doi:10.1016/S0140-6736(13)61249-0 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Gene Ontology Consortium (2006). The gene ontology (GO) project in 2006. Nucleic Acids Res. 34, D322–D326. doi:10.1093/nar/gkj021 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Ginestet, C. (2011). ggplot2: Elegant graphics for data analysis. J. R. Stat. Soc. Ser. A 174, 245–246. doi:10.1111/j.1467-985x.2010.00676_9.x |

CrossRef Full Text | Google Scholar

Gisterå, A., Robertson, A. K., Andersson, J., Ketelhuth, D. F., Ovchinnikova, O., Nilsson, S. K., et al. (2013). Transforming growth factor-β signaling in T cells promotes stabilization of atherosclerotic plaques through an interleukin-17-dependent pathway. Sci. Transl. Med. 5, 196ra100. doi:10.1126/scitranslmed.3006133 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Glickman, M. E., Rao, S. R., and Schultz, M. R. (2014). False discovery rate control is a recommended alternative to Bonferroni-type adjustments in health studies. J. Clin. Epidemiol. 67, 850–857. doi:10.1016/j.jclinepi.2014.03.012 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Hadad, N., Tuval, L., Elgazar-Carmom, V., Levy, R., and Levy, R. (2011). Endothelial ICAM-1 protein induction is regulated by cytosolic phospholipase A2α via both NF-κB and CREB transcription factors. J. Immunol. 186, 1816–1827. doi:10.4049/jimmunol.1000193 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Hansson, G. K., and Hermansson, A. (2011). The immune system in atherosclerosis. Nat. Immunol. 12, 204–212. doi:10.1038/ni.2001 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Huang, W. C., Sala-Newby, G. B., Susana, A., Johnson, J. L., and Newby, A. C. (2012). Classical macrophage activation up-regulates several matrix metalloproteinases through mitogen activated protein kinases and nuclear factor-κB. PLoS ONE 7, e42507. doi:10.1371/journal.pone.0042507 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Huang, H. Y., Lin, Y. C., Li, J., Huang, K. Y., Shrestha, S., Hong, H. C., et al. (2020). miRTarBase 2020: updates to the experimentally validated microRNA-target interaction database. Nucleic Acids Res. 48, D148–D154. doi:10.1093/nar/gkz896 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Jeggari, A., Marks, D. S., and Larsson, E. (2012). miRcode: a map of putative microRNA target sites in the long non-coding transcriptome. Bioinformatics 28, 2062–2063. doi:10.1093/bioinformatics/bts344 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Johnson, J. L., and Newby, A. C. (2009). Macrophage heterogeneity in atherosclerotic plaques. Curr. Opin. Lipidol. 20, 370–378. doi:10.1097/MOL.0b013e3283309848 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Kanehisa, M., Goto, S., Furumichi, M., Tanabe, M., and Hirakawa, M. (2010). KEGG for representation and analysis of molecular networks involving diseases and drugs. Nucleic Acids Res. 38, D355–D360. doi:10.1093/nar/gkp896 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Konkel, J. E., Zhang, D., Zanvit, P., Chia, C., Zangarle-Murray, T., Jin, W., et al. (2017). Transforming growth factor-β signaling in regulatory T cells controls T helper-17 cells and tissue-specific immune responses. Immunity 46, 660–674. doi:10.1016/j.immuni.2017.03.015 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Kotla, S., Singh, N. K., Heckle, M. R., Tigyi, G. J., and Rao, G. N. (2013). The transcription factor CREB enhances interleukin-17A production and inflammation in a mouse model of atherosclerosis. Sci. Signal. 6, ra83. doi:10.1126/scisignal.2004214 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, N., Liu, H., Xue, Y., Chen, J., Kong, X., and Zhang, Y. (2022). Upregulation of neogenin-1 by a CREB1-BAF47 complex in vascular endothelial cells is implicated in atherogenesis. Front. Cell Dev. Biol. 10, 803029. doi:10.3389/fcell.2022.803029 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Meng, L. B., Shan, M. J., Qiu, Y., Qi, R., Yu, Z. M., Guo, P., et al. (2019). TPM2 as a potential predictive biomarker for atherosclerosis. Aging (Albany NY) 11, 6960–6982. doi:10.18632/aging.102231 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Mestas, J., and Ley, K. (2008). Monocyte-endothelial cell interactions in the development of atherosclerosis. Trends cardiovasc. Med. 18, 228–232. doi:10.1016/j.tcm.2008.11.004 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Moore, K. J., and Tabas, I. (2011). Macrophages in the pathogenesis of atherosclerosis. Cell 145, 341–355. doi:10.1016/j.cell.2011.04.005 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Newman, A. M., Liu, C. L., Green, M. R., Gentles, A. J., Feng, W., Xu, Y., et al. (2015). Robust enumeration of cell subsets from tissue expression profiles. Nat. Methods 12, 453–457. doi:10.1038/nmeth.3337 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Peled, M., and Fisher, E. A. (2014). Dynamic aspects of macrophage polarization during atherosclerosis progression and regression. Front. Immunol. 5, 579. doi:10.3389/fimmu.2014.00579 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Rath, M., Müller, I., Kropf, P., Closs, E. I., and Munder, M. (2014). Metabolism via arginase or nitric oxide synthase: Two competing arginine pathways in macrophages. Front. Immunol. 5, 532. doi:10.3389/fimmu.2014.00532 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Ritchie, M. E., Phipson, B., Wu, D., Hu, Y., Law, C. W., Shi, W., et al. (2015). Limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res. 43, e47. doi:10.1093/nar/gkv007 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Salmena, L., Poliseno, L., Tay, Y., Kats, L., and Pandolfi, P. P. (2011). A ceRNA hypothesis: The rosetta stone of a hidden RNA language. Cell 146, 353–358. doi:10.1016/j.cell.2011.07.014 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Schauer, I. E., Knaub, L. A., Lloyd, M., Watson, P. A., Gliwa, C., Lewis, K. E., et al. (2010). CREB downregulation in vascular disease: A common response to cardiovascular risk. Arterioscler. Thromb. Vasc. Biol. 30, 733–741. doi:10.1161/ATVBAHA.109.199133 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Shen, L., Hu, Y., Lou, J., Yin, S., Wang, W., Wang, Y., et al. (2019). CircRNA-0044073 is upregulated in atherosclerosis and increases the proliferation and invasion of cells by targeting miR-107. Mol. Med. Rep. 19, 3923–3932. doi:10.3892/mmr.2019.10011 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Szklarczyk, D., Gable, A. L., Lyon, D., Junge, A., Wyder, S., Huerta-Cepas, J., et al. (2019). STRING v11: Protein-protein association networks with increased coverage, supporting functional discovery in genome-wide experimental datasets. Nucleic Acids Res. 47, D607–D613. doi:10.1093/nar/gky1131 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Tabas, I., and Bornfeldt, K. E. (2016). Macrophage phenotype and function in different stages of atherosclerosis. Circ. Res. 118, 653–667. doi:10.1161/CIRCRESAHA.115.306256 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Tan, L., Liu, L., Jiang, Z., and Hao, X. (2019). Inhibition of microRNA-17-5p reduces the inflammation and lipid accumulation, and up-regulates ATP-binding cassette transporterA1 in atherosclerosis. J. Pharmacol. Sci. 139, 280–288. doi:10.1016/j.jphs.2018.11.012 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Tsiantoulas, D., Diehl, C. J., Witztum, J. L., and Binder, C. J. (2014). B cells and humoral immunity in atherosclerosis. Circ. Res. 114, 1743–1756. doi:10.1161/CIRCRESAHA.113.301145 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Valencia-Sanchez, M. A., Liu, J., Hannon, G. J., and Parker, R. (2006). Control of translation and mRNA degradation by miRNAs and siRNAs. Genes Dev. 20, 515–524. doi:10.1101/gad.1399806 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Van den Bossche, J., O'Neill, L. A., and Menon, D. (2017). Macrophage immunometabolism: Where are we (going). Trends Immunol. 38, 395–406. doi:10.1016/j.it.2017.03.001 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, Z., Klipfell, E., Bennett, B. J., Koeth, R., Levison, B. S., Dugar, B., et al. (2011). Gut flora metabolism of phosphatidylcholine promotes cardiovascular disease. Nature 472, 57–63. doi:10.1038/nature09922 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Williams, K. J., and Tabas, I. (1995). The response-to-retention hypothesis of early atherogenesis. Arterioscler. Thromb. Vasc. Biol. 15, 551–561. doi:10.1161/01.atv.15.5.551 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Xiong, S., Hong, Z., Huang, L. S., Tsukasaki, Y., Nepal, S., Di, A., et al. (2020). IL-1β suppression of VE-cadherin transcription underlies sepsis-induced inflammatory lung injury. J. Clin. Invest. 130, 3684–3698. doi:10.1172/JCI136908 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Yang, K., Xue, Y., and Gao, X. (2021). LncRNA XIST promotes atherosclerosis by regulating miR-599/TLR4 Axis. Inflammation 44, 965–973. doi:10.1007/s10753-020-01391-x | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Ye, Z. M., Yang, S., Xia, Y. P., Hu, R. T., Chen, S., Li, B. W., et al. (2019). LncRNA MIAT sponges miR-149-5p to inhibit efferocytosis in advanced atherosclerosis through CD47 upregulation. Cell Death Dis. 10, 138. doi:10.1038/s41419-019-1409-4 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Yu, G., Wang, L. G., Han, Y., and He, Q. Y. (2012). clusterProfiler: an R package for comparing biological themes among gene clusters. OMICS 16, 284–287. doi:10.1089/omi.2011.0118 | |

PubMed Abstract | CrossRef Full Text | Google Scholar

Glossary

AAA abdominal aortic aneurysm

apoB-LPs apolipoprotein B-containing lipoproteins

BP biological process

CC cellular component

CD cluster of differentiation

ceRNA competing endogenous RNA

CREB1 cAMP-responsive element-binding protein 1

DCs dendritic cells

DE differentially expressed

FC fold change

FDR false discovery rate

FOXO1 forkhead box O1

FOXO3 forkhead box O3

GEO Gene Expression Omnibus

GO gene ontology

HAEC human aortic endothelial cell

HMGB1 high mobility group box 1

IL interleukin

IMMPORT Immunology Database and Analysis Portal

iNOS inducible nitric oxide synthase

KEGG Kyoto Encyclopedia of Genes and Genomes

KRAS KRAS proto-oncogene GTPase

LDL low-density lipoprotein

LIMMA linear models for microarray data

lncRNA long noncoding RNA

MCC maximal clique centrality

MEG3 maternally expressed 3

MF molecular function

miRNA microRNA

MMPs matrix metalloproteinases

MNC maximum neighborhood component

mRNA messenger RNA

ncRNA noncoding RNA

NK natural killer

PAOD peripheral arterial occlusive disease

PIK3R1 phosphoinositide-3-kinase regulatory subunit 1

PLAU plasminogen activator urokinase

PLXND1 plexin D1

PPI protein–protein interaction

PTEN phosphatase and tensin homolog

RNA ribonucleic acid

RT-qPCR real-time quantitative polymerase chain reaction

STRING Search Tool for the Retrieval of Interacting Genes

Th1 T-helper type 1

Th17 T-helper type 17

TNFAIP3 TNF alpha-induced protein 3

Treg regulatory T

YWHAE tyrosine 3-monooxygenase/tryptophan 5-monooxygenase activation protein epsilon

YWHAH tyrosine 3-monooxygenase/tryptophan 5-monooxygenase activation protein eta