The microarray data of anaphylactic patients, STEMI patients, and corresponding healthy controls were pooled from the public GEO database (https://www.ncbi.nlm.nih.gov/geo/). Datasets that satisfied the following inclusion criteria were selected: (1) gene expression profiles were based on human specimens; (2) all samples were obtained from peripheral blood; (3) the datasets included disease cases and healthy controls. GSE69063 (anaphylaxis) and GSE60993 (STEMI) datasets that meet the above criteria were used in the present study. GSE69063 contains 17 anaphylactic patients and 10 healthy controls whose blood specimens were collected 1 h after arriving at the emergency department. GSE60993 dataset consists of 7 STEMI patients and 7 healthy individuals, and peripheral blood from STEMI patients was obtained within 4 h after the attack of chest pain. Details of these datasets were summarized in Supplementary Table S1.

DEGs between case and control groups were screened by GEO2R (an official web application in NCBI that helps analyze GEO data, www.ncbi.nlm.nih.gov/geo/ge2r) (Barrett et al., 2013), with screening criteria set as |logFC (fold change)|>1 and adjusted p.value < 0.05. Probe sets without corresponding gene symbols were removed and genes with multiple probe sets were averaged. The volcano map was drawn with ggplot2 R package (https://ggplot2.tidyverse.org). TBtools software was applied to draw heatmaps to visualize the differential gene expression profiles between case and control groups (Chen et al., 2020).

Gene Ontology (GO) enrichment analysis and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis were performed with the clusterProfiler package (Yu et al., 2012). The top 20 terms of biological process (BP) and top 5 KEGG terms were visualized with bubble charts by ggplot2 R package.

The composition and distribution of immune cells in anaphylactic samples (GSE69063), STEMI samples (GSE60993), and normal samples were evaluated by CIBERSORTx (https://cibersortx.stanford.edu/) (Newman et al., 2015), which is based on a set of feature gene expression data from 22 cell subtypes (LM22 dataset). The fraction of these cell populations in diseased and healthy groups were compared using the ggpubr R package (Wilcoxon test) (https://CRAN.R-project.org/package=ggpubr).

STRING database was applied to predict PPI (https://cn.string-db.org/) (Franceschini et al., 2013) and a combined score >0.4 was considered statistically significant. Cytoscape software (Shannon et al., 2003) was used to visualize PPI networks. DEGs in the intersection of six algorithms (DMNC, Stress, MCC, Degree, Closeness, and Radiality) were determined as hub genes. The co-expression and co-localization network of hub genes was constructed and visualized by GeneMANIA (http://www.genemania.org/) (Warde-Farley et al., 2010).

Hub gene expression was confirmed in the GSE47655 (anaphylaxis) and GSE61144 (STEMI) datasets. GSE47655 contains 6 anaphylactic patients and 6 healthy controls. GSE60993 dataset consists of 7 STEMI patients and 10 healthy subjects. All samples from GSE47655 and GSE61144 were obtained from peripheral blood. Details of these datasets were shown in Supplementary Table S1. The differences in mRNA expression of hub genes between the diseased and healthy groups were compared using Student’s t-test and visualized using the ggplot2 package.

Pearson’s correlation analysis between the expression of hub genes and the abundance of immune effector cells was carried out to further explore the immunomodulatory mechanisms.

TRRUST (https://www.grnpedia.org/trrust/) (Han et al., 2018), a database that provides transcription factor-target regulatory relationships, was applied to predict transcription factors that regulate hub gene expression. The expression levels of transcription factors predicted by TRRUST were confirmed in GSE69063 (anaphylaxis) and GSE60993 (STEMI), and the differences between diseased and healthy groups were compared using t-test. ChEA3 (ChIP-X Enrichment Analysis, version 3) platform provides transcription factor enrichment analysis and integrates RNA-seq data from GTEx, TCGA, and ARCHS4, as well as CHIP-seq data from ENCODE and ReMap (Keenan et al., 2019). The hub gene list was submitted to ChEA3, and transcription factors common to those predicted by TRRUST were identified. Base on the MeanRank method, transcription factors were ranked according to their composite scores. DGIdb (Drug-gene interaction database, http://www.dgidb.org) was utilized to predict the drug-hub gene interaction (Freshour et al., 2021).

C57BL/6 mice were obtained from the Charles River Laboratories (Beijing, China). The murine model used in this study was based on an active systemic anaphylaxis model described previously (Jonsson et al., 2011) with modifications. Six-week-old female mice were sensitized subcutaneously on day 0 with bovine serum albumin (BSA, 50 μg per mouse, Sigma-Aldrich, United States) in complete Freund adjuvant (CFA, Sigma-Aldrich, United States) and boosted on day 7 and day 14 with 50 μg BSA in incomplete Freund adjuvant (IFA, Sigma-Aldrich, United States). One week after the last sensitization, mice were intravenously injected with 15 μg BSA to elicit systemic anaphylaxis. After BSA challenge, the temperature was monitored every 10 min with a rectal thermometer (TH212,China), and the severity of anaphylaxis was scored every 10 min on a scale of 0–4 based on the grading system described previously (Cloutier et al., 2018), score 0: normal; score 1: slow motions; score 2: impaired mobility, still reacting to touch; score 3: immobilized and do not react to touch; score 4: death.

When fully anesthetized with inhalant isoflurane (RWD, China), mice were immobilized in a supine position with electrodes implanted subcutaneously in the limbs. The Animal Bio Amp device (ADInstruments, Australia) and LabChart software were used to acquire and record the lead II ECGs of naïve and model mice.

One day prior to ultrasound, hair around the chest wall was carefully removed using depilatory cream (Nair, United States). Mice were anesthetized by inhalant isoflurane (RWD, China), at a concentration of 1.5%. When fully anesthetized, the mouse was positioned on the warm imaging platform ventral side up, and the medical ultrasonic gel was applied to the limb leads to generate ECG. Vevo2100 System (Visual Sonics, Canada) was used to record the transthoracic echocardiograms of naïve and model mice. The MS550D transducer was placed on the left sternal edge to obtain a parasternal short axis (PSAX) view. Left ventricular ejection fraction (EF) and fractional shortening (FS) were acquired from PSAX M-mode scans at the mid-papillary muscle level. Echocardiographic data were analyzed offline using Vevo Lab software and all measurements were averaged from 3 cardiac cycles.

Total RNA was extracted from the peripheral blood of naïve mice and model mice using TRIzol reagent (Invitrogen, United States). PrimeScript™ RT Master Mix (TaKaRa, Japan) was utilized to convert the equivalent amount (500 ng) of total RNA into cDNA. The relative mRNA expression of target genes was quantified with the SYBR master mix (Yeasen, China). All samples were normalized to the housekeeping gene β-Actin. Details of primers are listed as follows.

SPSS version 23.0 (SPSS Inc., Chicago, IL, United States) was used to perform univariate logistic analysis and construct multivariate logistic regression model with the expression profiles of hub genes as the continuous prediction variable and the physical condition (disease or not) as the categorical response variable. Statistically significant variables (p < 0.05) were included in the multivariate logistic regression model using the backward LR method for variable selection. ROC curves were plotted and the AUC value was calculated to evaluate the diagnostic efficiency of hub genes.

Statistical analysis was performed using R language (R Foundation for Statistical Computing, Vienna, Austria) and GraphPad Prism 8.0 software (GraphPad Software, San Diego, CA, United States). Data were evaluated for normal distribution using Shapiro-Wilk tests. If data were normally distributed and with similar variances, two-tailed Student’s t-test (parametric) was used to compare the differences between two groups. For normally distributed data with unequal variances, unpaired t-test with Welch’s correction was performed. For abnormally distributed data with unequal variances, Kolmogorov-Smirnov test (nonparametric) was used to determine differences between two groups. Pearson’s correlation test was used to analyze the correlation between the expression of hub genes and the abundance of immune effector cells. Quantitative data are shown as mean ± SEM. *p < 0.05, **p < 0.01, ***p < 0.001 were considered statistically significant.

The workflow of this study was shown in Supplementary Figure S1. To characterize the gene expression profile of acute anaphylaxis, we compared mRNA expression in whole blood between anaphylactic patients and healthy subjects. Venous blood was collected at 1 hour after patients’ arriving at the emergency department. A total of 265 genes were differentially expressed in the anaphylactic patients, among which 188 were upregulated and 77 were downregulated (Figure 1A). The heatmap revealed high heterogeneity of gene expression between anaphylactic patients and healthy controls (Figure 1B). To explore the pathogenesis of anaphylaxis, we performed GO and KEGG analysis, and profiling of biological processes and signaling pathways revealed the complex pathology of anaphylaxis. This disease was featured by marked activation of immune-inflammatory responses, including positive regulation of response to external stimulus, regulation of innate immune response, cytokine-mediated signaling pathways, and positive regulation of cytokine production, especially interleukin-1 (Figure 1C). In addition, the blood coagulation cascade, platelet aggregation, hemostasis, and wound healing actively participated in the development of anaphylaxis (Figure 1C). KEGG results illustrated that anaphylaxis-related DEGs were associated with transcriptional misregulation and hematopoietic cell lineage (Figure 1D).

Since anaphylaxis was characterized by marked activation of immune-inflammatory responses, the CIBERSORT deconvolution algorithm was applied to further explore the details of immune cells involved in acute anaphylaxis. As shown in Figure 2A, a wide heterogeneity of immune cell distribution was observed in anaphylactic patients compared to healthy subjects. Among the 22 immune cell subtypes, CD4⁺ memory T cells (resting), follicular helper T cells, and M0 macrophages were distributed only in anaphylactic samples but not in controls. Compared to the controls, anaphylactic specimens showed significantly increased levels of T cells gamma delta, monocytes, dendritic cells (activated) and mast cells (resting) and markedly decreased proportions of naïve B cells, CD8⁺ T cells, CD4⁺ naïve T cells, and M2 macrophages (Figure 2B).

To uncover the gene expression profiles of STEMI, we compared the mRNA expression in whole blood between STEMI patients and healthy subjects. A total of 237 STEMI-related DEGs were screened, including 199 upregulated and 38 downregulated genes, and these DEGs were visualized with the volcano plot and heatmap (Figures 3A,B). STEMI-related DEGs were primarily involved in the biological processes of innate immune response regulation, positive regulation of defense response, and positive regulation of response to external stimulus, indicating immune defense responses were essential in STEMI. Notably, myeloid leukocyte activation, positive regulation of cytokine production, leukocyte migration, and leukocyte-mediated immunity were critically involved in the development of STEMI (Figure 3C). KEGG analysis revealed that STEMI-related DEGs were associated with hematopoietic cell lineage, neutrophil extracellular trap formation, PD-L1 expression and PD-1 checkpoint pathway (Figure 3D).

As highlighted above, myeloid leukocyte activation and leukocyte-mediated immunity were essential in STEMI. To further elucidate components of immune cells key to STEMI pathology, we used CIBERSORTx to profile the composition and distribution of immune cell subtypes in STEMI and normal samples. As shown in Figure 4A, there were marked differences in the distribution of immune cells between these two groups. STEMI samples were featured by increased levels of plasma cells, M0 macrophages, dendritic cells (activated), and neutrophils, and by decreased infiltration of CD8⁺ T cells, CD4⁺ naïve T cells, NK cells (resting), and eosinophils (Figure 4B).

By extracting DEGs that were shared by anaphylaxis and STEMI, 33 common DEGs were identified (Figure 5A). The PPI network of these common DEGs contained 33 nodes and 24 edges, revealing complex interactions among these molecules (Figure 5B). These common DEGs were enriched in the secretory granules, and their molecular functions include immune receptor activity, pattern recognition receptor activity, and cytokine receptor activity (Figure 5C), suggesting immune-inflammatory regulation is essential in both anaphylaxis and STEMI. Enrichment analysis of these common DEGs demonstrated that the shared biological processes of anaphylaxis and STEMI included cytokine-mediated signaling pathway, response to reactive oxygen species, positive regulation of defense response, regulation of immune effector process, and positive regulation of NF-kappa B transcription factor activity (Figure 5C). Taking the intersection of six algorithms (DMNC, Stress, MCC, Degree, Closeness, Radiality) in cytoHubba (Figure 5D), we found 6 hub genes shared by anaphylaxis and STEMI: IL1R2, S100A12, FOS, MMP9, DUSP1, and CLEC4D (Figure 5E). The description and major functions of hub genes are listed in Supplementary Table S2. These hub genes constituted an interactive network with a co-expression rate of 96.23% and a co-localization rate of 3.77% (Figure 5F). The prominent functions of hub genes included positive regulation of defense response, response to cadmium ion, cell chemotaxis, regulation of inflammatory response, regulation of apoptotic signaling pathway, and positive regulation of DNA-binding transcription factor activity (Figure 5F).

After taking the intersection and excluding the immune cells with opposite infiltration trends in anaphylactic and STEMI groups, M0 macrophages, activated dendritic cells, CD4⁺ naïve T cells, and CD8⁺ T cells were identified as immune effector cells associated with both anaphylaxis and STEMI (Figure 6A). To further investigate the association of hub genes and these immune cell subtypes, Pearson’s correlation analysis was performed on the expression levels of hub genes and the abundance of immune cells in both anaphylactic and STEMI samples. The results illustrated that all hub genes were positively correlated with M0 macrophage (R > 0, p < 0.05), while negatively correlated (R < 0, p < 0.05) with CD4⁺ naïve T cells (Figures 6B,C ). FOS showed a significant positive correlation with activated dendritic cells, and CLEC4D had a significant negative correlation with CD8⁺ T cells (Supplementary Figure S2A, B).

To assess the robustness of hub genes identified above, we analyzed their gene expression in external human microarray datasets (GSE47655 for anaphylaxis and GSE61144 for STEMI) (Supplementary Table S1). Compared with healthy controls, 5 hub genes (IL1R2, FOS, MMP9, DUSP1, CLEC4D) were significantly upregulated in anaphylactic samples, while S100A12 showed no statistical difference (Figure 7A). In STEMI samples, the expression of all hub genes (IL1R2, S100A12, FOS, MMP9, DUSP1, CLEC4D) was markedly higher than that of controls (Figure 7B).

To further verify the findings from human microarray data, the expression levels of five hub genes (IL1R2, FOS, MMP9, DUSP1, CLEC4D) were measured in an experimentally induced murine model of acute anaphylaxis by qRT-PCR analysis. Systemic anaphylaxis was induced by intravenously injected with BSA in mice sensitized with the same antigen (Figure 8A). 10 min after the BSA challenge, mice showed obvious symptoms of anaphylaxis, including impaired mobility, lethargy, and unresponsiveness (Figure 8B). Compared with naïve mice, the anaphylactic mice exhibited ST-segment elevations on electrocardiograms (Figure 8C) and showed markedly impaired heart function (Figures 8D–F) as reflected by decreases in ejection fraction (EF; 58.68% ± 2.87% versus 27.84% ± 4.67%) and fractional shortening (FS; 30.62% ± 1.95% versus 12.53% ± 2.33%). qRT-PCR analysis of peripheral blood extracts demonstrated that the levels of IL1R2, FOS, MMP9, DUSP1, and CLEC4D in the anaphylactic mice with STEMI were significantly upregulated (Figures 8G–K), which was consistent with the human microarray results. Taken together, these results indicated that the hub genes we screened were of high reliability.

TRRUST database was applied to explore key upstream regulators for hub genes. Signal transducer and activator of transcription 3 (STAT3), cAMP-responsive element modulator (CREM), signal transducer and activator of transcription 6 (STAT6), Sp1 transcription factor (SP1), Fos proto-oncogene (FOS), serum response factor (SRF), and histone deacetylase 1 (HDAC1) were identified as transcription factors that modulated hub gene expression (Supplementary Figure S3A). The details of these transcription factors are summarized in Supplementary Table S3. With further analysis, we found that CREM and FOS were highly expressed in anaphylactic samples, while 3 transcription factors (SRF, STAT6, SP1) were significantly downregulated (Supplementary Figures S3B–H). In STEMI samples, 3 transcription factors (FOS, STAT3, SP1) were markedly upregulated, while HDAC1 was significantly downregulated (Supplementary Figures S3I–O). Moreover, ChEA3 platform was used to validate the findings from TRRUST and further rank the transcription factors according to their MeanRank composite scores. STAT3, CREM, STAT6, SP1, FOS, and SRF were identified as transcription factors of hub genes that common to those predicted by TRRUST, ranking as: FOS, STAT3, STAT6, SP1, SRF, and CREM (Supplementary Table S4). qRT-PCR analysis of peripheral blood extracts showed significantly increased expression of FOS, and markedly decreased expression of STAT6, SP1, and SRF in the murine model of anaphylaxis complicated STEMI (Figures 8L–Q), suggesting the important regulatory role of these transcription factors in anaphylaxis induced myocardial damage.

The candidate druggable hub genes were predicted via the DGIdb database. As listed in Table 1, IL1R2, FOS, MMP9, and DUSP1 showed interactions with approved drugs, suggesting the potential targeting and regulatory effects of approved drugs on these hub genes.

To predict the occurrence of anaphylaxis and STEMI, univariate logistic analysis followed by multivariate logistic regression analysis was applied to construct the diagnostic model, and receiver operating characteristic (ROC) curves were used to evaluate the diagnostic efficiency of hub genes. Univariate logistic analysis of hub genes in anaphylactic dataset GSE69063 revealed that IL1R2 (p < 0.001, AUC = 0.797, 95% CI 0.679-0.916), DUSP1 (p = 0.002, AUC = 0.719, 95% CI 0.593-0.844), MMP9 (p = 0.003, AUC = 0.706, 95% CI 0.577-0.835), and CLEC4D (p < 0.001, AUC = 0.888, 95% CI 0.803-0.974) were significantly correlated with anaphylaxis as continuous variables (Figure 9A). All statistically significant variables were then included in the multivariate logistic regression model using the backward LR method for variable selection. The results showed DUSP1 (p = 0.028) and CLEC4D (p = 0.001) were selected as significant variables and used to construct diagnostic model. The regression equation of logit (P) = −33.355 + 2.9*DUSP1+2.294*CLEC4D was established and the accuracy of this logistic model was evaluated by ROC curves. AUC values in GSE69063 and another independent anaphylactic dataset GSE47655 were 0.911 (95% CI, 0.840-0.982) and 0.938 (95% CI, 0.857-1.000), respectively (Figures 9B,C), suggesting this model was of high diagnostic efficacy in distinguishing anaphylactic patients from healthy individuals. As anaphylaxis and STEMI shared common hub genes, the diagnostic efficiency of hub genes was further evaluated in STEMI datasets. Univariate logistic analysis of hub genes in STEMI dataset GSE60993 and GSE61144 revealed that DUSP1 (p = 0.001, AUC = 0.774, 95% CI 0.649-0.898), MMP9 (p < 0.001, AUC = 0.887, 95% CI 0.794-0.980), and CLEC4D (p = 0.03, AUC = 0.649, 95% CI 0.507-0.790) were significantly correlated with STEMI (Figure 9D). Next, the efficiency of the established diagnostic model was validated in the STEMI datasets (GSE60993 and GSE61144) based on the expression profiles of DUSP1 and CLEC4D. As shown in Figure 9E, the AUC values in STEMI datasets exceeded 0.7, indicating this model also had good accuracy in discriminating STEMI. Taken together, these results suggested that two hub genes, DUSP1 and CLEC4D, held a promise for the diagnosis of anaphylaxis complicated STEMI as blood diagnostic biomarkers.