Bioinformatics Identification of Candidate Biomarkers in Endomyocardial Biopsy and Peripheral Blood for Cardiac Allograft Rejection.

BACKGROUND Cardiac allograft rejection is still a crucial barrier to achieving satisfactory outcomes after surgery. In this study, we propose to find candidate biomarkers from endomyocardial biopsy (EMB) and peripheral blood (PB) samples for efficient diagnosis and treatment of cardiac allograft rejection. MATERIAL AND METHODS Microarray datasets were obtained from the Gene Expression Omnibus (GEO). Differentially expressed genes (DEGs) of cardiac allograft rejection patients and control subjects from EMB and PB samples were screened using the online tool GEO2R. Gene ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis of all samples' DEGs were performed with the DAVID online tool. Protein-protein interaction (PPI) networks were constructed and visualized using Cytoscape and the top 10 hub genes were selected. Finally, the most highly enriched GO and KEGG pathways of the top 10 hub genes were determined. RESULTS A total of 57 502 genes from EMB samples and 131 624 genes from PB samples were identified. Gene characteristics and enrichment analysis indicated that both EMB and PB samples contained DEGs involved in antigen presentation, immune cells activation, inflammatory process, and cellular injuries. In EMB samples, there were some DEGs related to heart tissue injury and cardiac malfunction. Moreover, DEGs that regulates hypoxia-induced factors and erythrocyte function in response of ischemia and hypoxia stress were present in PB samples but were absent in EMB samples. CONCLUSIONS The screened differentially expressed genes (DEGs) from EMB and PB samples of patients with cardiac graft rejection are potential candidate biomarkers of diagnosis and treatment.

Background

Heart transplantation is one of the most effective options to treat patients with end-stage heart disease. However, acute immune rejection is still an intractable obstacle for prolonging survival and promoting quality of life after surgery [1]. To solve this urgent issue, the first and most crucial step is to find biomarkers that can be used in rapid and accurate diagnosis for early-stage immune rejection disorder. Additionally, finding efficient therapeutic targets for alleviating rejection injury before the occurrence irreversible graft failure and multiple organ failure is also important.

In clinical practice, endomyocardial biopsy (EMB) is the criterion standard to identify acute cardiac rejection [2]. However, drawbacks of EMB have been realized since this technology has been initially applied decades ago. On the one hand, EMB is an invasive procedure that inevitably causes complications, including cardiac tamponade, arrhythmias, severe heart internal-structure injury, bleeding, and thrombosis [3]. On the other hand, EMB is mainly based on the histopathologic interpretation, which is relatively insensitive to detect early rejection [4]. To find a better way to diagnose and intervene in heart transplant rejection, in recent years, researchers have been searching for find biomarkers that can improve the accuracy of EMB and make it more suitable for early detection of rejection [5,6]. For example, it was found that a set of microRNAs are expressed differently in EMB samples from patients with heart transplant rejection in earlier stage [7]. Meanwhile, other researchers are committed to the development of minimally invasive and non-invasive cardiac allograft rejection monitoring, including traditional imaging diagnosis, intramyocardial electrocardiogram, molecular imaging, and testing biomarkers in peripheral blood (PB). Among the minimally-invasive and non-invasive methods currently used, testing of PB samples has been regarded as the most repeatable, convenient, and economical method for diagnosis of rejection [8]. It can determine whether the type of rejection is acute cellular rejection (ACR), antibody-mediated rejection (AMR), or mixed. And novel biomarkers such as cell-free DNA in peripheral blood enable this test more sensible for diagnosis [4]. However, limitations still exist, as the peripheral blood test has poor ability to distinguish rejection and infection from immunological stability [8,9], so further investigations to identify biomarkers of cardiac allograft rejection from EMB and PB are warranted.

In this study, we analyzed microarray data of rejection and non-rejection heart transplantation recipients from EMB and PB samples by bioinformatics methods, which were downloaded from the Gene Expression Omnibus (GEO) database. Then, we screened out the differentially expressed genes (DEGs) from the 2 types of samples, displayed their significant Gene Ontology (GO) functional annotation analysis and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment, and constructed the protein–protein interaction (PPI) networks. Through further analysis, the hub genes of EMB and PB samples and their specific bio-function in cardiac graft rejection was discovered. By performing in-depth analysis and comparing characteristics of the DEGs and pivotal roles they likely play in cardiac graft rejection, we selected the candidate biomarkers of the 2 sample types that possess huge potential value to improve future cardiac allograft rejection diagnosis and treatment.

Material and Methods

Microarray Source

By searching public resources in the Gene Expression Omnibus (GEO) database (http://www.ncbi.nlm.nih.gov/geo/), the expression profiling of human endomyocardial biopsies (EMB) samples of heart transplantation recipients with or without rejection were included in this study, which were GSE2596, GSE4470, and GSE4315. The profiling was conducted on the GPL1053 (NIAID – Hsaband) and GPL96 ([HG-U133A] Affymetrix Human Genome U133A Array) platforms. As for rejection and non-rejection data from peripheral blood (PB) samples, GSE4315, GSE5967, GSE33970, and GSE87301 were included, and the platforms of these profiles were GPL96 and GPL570 ([HG-U133_Plus_2] Affymetrix Human Genome U133 Plus 2.0 Array). All of the data are from open-access sources available online, and did not involve any human or animal experiments.

Data Processing and DEG Identification

The raw microarray data including EMB and PB samples were obtained from the GEO database. After filtering out chronic rejection cases, data of both samples were analyzed separately by the online tool GEO2R (http://www.ncbi.nlm.nih.gov/geo/geo2r) to screen out different genes expressed between rejection and non-rejection samples. We downloaded data files analyzed by GEO2R and we used the OmicStudio tools (https://www.omicstudio.cn/tool) to show the volcano plots of genes expressed in EMB and PB samples from rejection and non-rejection groups (Figure 1). Among these genes, those that meet the standard of P value <0.05 and |logFC| ≥1.0 were regarded as significant differentially expressed genes (DEGs). According to the absolute value of logFC ranked from high to low, the top 10 DEGs from both samples were selected (Tables 1, 2).

Figure 1

Differential expression of genes from endomyocardial biopsy and peripheral blood samples

(A) Volcano plot of 57 502 genes expressed differentially in endomyocardial biopsy of patients with and without cardiac rejection. (B) Volcano plot of 131 624 genes expressed differentially in peripheral blood of patients with and without cardiac rejection. Red plots represent significantly upregulated genes (logFC ≥1, P value <0.05), blue plots are significantly downregulated genes (logFC ≤−1, P value <0.05), while the gray ones demonstrate genes with non-differential expression. Volcano plot graphs were created using the OmicStudio tools available at https://www.omicstudio.cn/tool.

Table 1

Top 10 differentially expressed genes (DEGs) from endomyocardial biopsy samples.

DEGs	log FC	p-Value
IGHG3	3.2153	1.56E-08
SCYA5	3.091389	1.46E-11
CCL18	3.05464	7.80E-07
NPPB	2.60278392	0.034928
CSPG2	2.535519	7.17E-07
CD2	2.438863	1.10E-10
UBD	2.228056	2.01E-06
COPG	2.21468	1.30E-07
IGL@	2.1198	1.25E-05
ZAP70	2.105274	1.16E-08

The criteria for DEGs screening are |logFC| ≥1 and P value <0.05.

Table 2

Top 10 differentially expressed genes (DEGs) from peripheral blood samples.

DEGs	log FC	p-Value
IFI27	2.43246027	0.016955
HLA-DQB1	−2.2071461	0.00089723
CA1	2.1195333	0.00011024
HLA-DQA1	−2.1127547	0.00500552
FECH	2.0857495	0.00020633
BPGM	2.0250519	0.00062333
XK	1.7945742	0.00260448
HEMGN	1.7774263	0.037994
IFI27	1.7642483	0.04744193
RNF182	1.7447703	0.00845521

The criteria for DEGs screening are |logFC| ≥1 and P value <0.05.

Enrichment Analyses of DEGs

DEGs expression profiles were submitted to the DAVID online tool (Database for Annotation, Visualization and Integrated Discovery, version 6.8, http://david.ncfcrf.gov/) for analysis. To assess the significance of the biological function of DEGs, GO and KEGG pathway annotation were performed for EMB and PB samples. P<0.05 and gene counts ≥10 were set as the threshold for EMB samples. In the case of PB samples, P<0.05 and gene counts ≥5 were considered statistically significant (Tables 3, 4).

Table 3

Enrichment of Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis of significance for differentially expressed genes (DEGs) from endomyocardial biopsy samples.

Category	Term	Description	Count	p-Value
GOTERM_MF	GO: 0005515	Protein binding	57	0.01841913
GOTERM_MF	GO: 0042803	Protein homodimerization activity	10	0.013700626
GOTERM_MF	GO: 0042802	Identical protein binding	12	0.00173391
GOTERM_MF	GO: 0005102	Receptor binding	11	1.49E-05
GOTERM_CC	GO: 0016021	Integral component of membrane	39	0.00321182
GOTERM_CC	GO: 0005794	Golgi apparatus	10	0.026519452
GOTERM_CC	GO: 0005886	Plasma membrane	34	0.001736587
GOTERM_CC	GO: 0070062	Extracellular exosome	28	3.47E-04
GOTERM_CC	GO: 0005887	Integral component of plasma membrane	21	1.46E-05
GOTERM_CC	GO: 0016020	Membrane	21	0.004716471
GOTERM_CC	GO: 0005576	Extracellular region	19	7.84E-04
GOTERM_CC	GO: 0005615	Extracellular space	18	2.82E-04
GOTERM_CC	GO: 0006955	Immune response	17	6.07E-10
GOTERM_CC	GO: 0000139	Golgi membrane	11	6.63E-04
GOTERM_BP	GO: 0045087	Innate immune response	13	2.72E-06
GOTERM_BP	GO: 0009986	Cell surface	13	1.41E-05
GOTERM_BP	GO: 0007165	Signal transduction	12	0.042676364
GOTERM_BP	GO: 0050852	T cell receptor signaling pathway	10	8.77E-08
GOTERM_BP	GO: 0006954	Inflammatory response	10	1.79E-04
KEGG_PATHWAY	hsa04612	Antigen processing and presentation	12	4.60E-11
KEGG_PATHWAY	hsa05168	Herpes simplex infection	12	5.45E-07
KEGG_PATHWAY	hsa04514	Cell adhesion molecules (CAMs)	10	4.23E-06
KEGG_PATHWAY	hsa05164	Influenza A	10	2.21E-05

The criteria are count ≥10 and P value <0.05. BP – biological process; CC – cellular component; MF – molecular function.

Table 4

Enrichment of Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis of significance for differentially expressed genes (DEGs) from peripheral blood samples.

Category	Term	Description	Count	p-Value
GOTERM_MF_	GO: 0005515	Protein binding	114	1.90E-06
GOTERM_MF_	GO: 0044822	poly(A) RNA binding	19	0.022148293
GOTERM_MF_	GO: 0042803	Protein homodimerization activity	16	0.004263137
GOTERM_CC_	GO: 0005737	Cytoplasm	66	0.006773839
GOTERM_CC_	GO: 0005829	Cytosol	55	2.23E-05
GOTERM_CC_	GO: 0016020	Membrane	33	0.008967662
GOTERM_CC_	GO: 0005622	Intracellular	23	0.007774602
GOTERM_CC_	GO: 0005741	Mitochondrial outer membrane	6	0.013942611
GOTERM_CC_	GO: 0005840	Ribosome	6	0.021196797
GOTERM_CC_	GO: 0045121	Membrane raft	6	0.046903966
GOTERM_BP_	GO: 0035556	Intracellular signal transduction	11	0.006858197
GOTERM_BP_	GO: 0043066	Negative regulation of apoptotic process	10	0.03717623
GOTERM_BP_	GO: 0000184	Nuclear-transcribed mRNA catabolic process, nonsense-mediated decay	6	0.006570574
GOTERM_BP_	GO: 0006413	Translational initiation	6	0.011674787
GOTERM_BP_	GO: 0000209	Protein polyubiquitination	6	0.03611137
GOTERM_BP_	GO: 0006661	Phosphatidylinositol biosynthetic process	5	0.002580151
GOTERM_BP_	GO: 0006614	SRP-dependent cotranslational protein targeting to membrane	5	0.014114136
GOTERM_BP_	GO: 0019083	Viral transcription	5	0.025108741

The criteria are count ≥5 and P value <0.05. BP – biological process; CC – cellular component; MF – molecular function.

PPI Network Construction and Hub Genes Identification

Protein–protein interaction (PPI) networks of DEGs for EMB and PB samples were constructed and visualized with Cytoscape software (version 3.8.2). A confidence score cutoff >0.4 was set as the standard of significance for PPI pairs (Figure 2). The top 10 hub genes of each PPI network were selected and entered into CytoHubba, a plugin in Cytoscape software, according to the degree of nodes rank (≥10) from high to low (Figure 3). Then, we analyzed GO and KEGG pathway enrichment of the top 10 hub genes (Tables 5, 6).

Figure 2

PPI Networks of significant DEGs in endomyocardial biopsy and peripheral blood sample constructed with Cytoscape 3.8.2 software

(A) Significant DEGs PPI networks of myocardium. (B) Significant DEGs PPI networks of peripheral blood. Nodes represent genes. Edges represent the interaction among proteins. the structure of proteins can be checked in nodes. Colors of edges represent the interaction between the proteins proved by different kinds of evidence. PPI – protein–protein interaction.

Figure 3

Networks of top 10 hub genes with connective degree of ≥10 achieved with CytoHubba, a plugin in Cytoscape

(A) Ten hub genes of DEGs from endomyocardial biopsy samples. (B) Ten hub genes of DEGs from peripheral blood. The depth of color of each gene reflects its degree of connection.

Table 5

Top 10 hub genes from endomyocardial biopsy samples enriched in Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG).

Category	Description	Genes from Top 10 Hub
GOTERM_MF	Protein binding	ZAP70, CD74, TYROBP, CCL5, B2M, PSMB9
GOTERM_CC	Integral component of membrane	HLA-DQA1, ZAP70, CD74, TYROBP, HLA-B, HLA-C, HLA-DRA
GOTERM_CC	Plasma membrane	HLA-DQA1, ZAP70, CD74, TYROBP, HLA-B, HLA-C, HLA-DRA, B2M
GOTERM_CC	Membrane	HLA-DQA1, CD74, HLA-B, HLA-C, B2M
GOTERM_CC	Cell surface	CD74, TYROBP, HLA-B, HLA-C, HLA-DRA
GOTERM_CC	Extracellular exosome	HLA-B, HLA-C, HLA-DRA, B2M, PSMB9
GOTERM_BP	Immune response	HLA-DQA1, CD74, ZAP70, HLA-B, HLA-C, CCL5, B2M
GOTERM_BP	T cell receptor signaling pathway	HLA-DQA1, ZAP70, HLA-DRA, PSMB9
KEGG_PATHWAY	Antigen processing and presentation	HLA-DQA1, CD74, HLA-B, HLA-C, HLA-DRA, B2M
KEGG_PATHWAY	Herpes simplex infection	HLA-DQA1, CD74, HLA-B, HLA-C, CCL5, HLA-DRA
KEGG_PATHWAY	Cell adhesion molecules (CAMs)	HLA-DQA1, HLA-B, HLA-C, HLA-DRA

color indicates upregulated genes and indicates downregulated genes. BP – biological process; CC – cellular component; MF – molecular function.

Table 6

Top 10 hub genes from peripheral blood samples enriched in Gene Ontology (GO).

Category	Description	Genes from Top 10 Hub
GOTERM_MF	Protein binding
GOTERM_CC	Cytoplasm
GOTERM_CC	Cytosol

color indicates upregulated genes and indicates downregulated genes. CC – cellular component; MF – molecular function.

Results

Identification of the DEGs

Genes expressed in EMB and PB samples of patients with and without rejection after heart transplantation were obtained by submitting microarray data to the GEO2R analyzer, resulting in 57 502 genes from EMB samples and 131 624 genes from PB samples. The results were demonstrated via volcano plot (Figure 1) where red squares represented upregulated DEGs (53 genes for EMB and 164 for PB) (logFC ≥1, P<0.05), blue triangles represent downregulated DEGs (88 genes for EMB and 118 for PB) (logFC ≤−1, P<0.05), and the gray circles represent non-differential expression genes (57 361 genes for EMB and 131 342 for PB) (−1

Biofunctional Enrichment Analysis of DEGs

The GO and KEGG pathway enrichment analyses of DEGs were performed using the DAVID online tool. In EMB samples, there were 19 GO and 4 KEGG pathway categories mostly enriched. For GO, cellular component (CC) was the most significant category, including protein binding, protein homodimerization activity, identical protein binding, receptor binding, integral component of membrane, Golgi apparatus, plasma membrane, extracellular exosome, integral component of plasma membrane, membrane, extracellular region, extracellular space, immune response, and Golgi membrane. Biological process (BP) is the second category, consisting of innate immune response, cell surface, signal transduction, T cell receptor signaling pathway, and inflammatory response. For molecular function (MF), protein binding, protein homodimerization activity, identical protein binding, and receptor binding were enriched. Four categories of KEGG pathway in EMB were Antigen processing and presentation, Herpes simplex infection, Cell adhesion molecules (CAMs), and Influenza A (Table 3). However, in PB samples, DEGs were only significantly enriched in GO, which contained 4 categories of CC (cytosol, membrane, mitochondrial outer membrane, and cytoplasm), 4 categories of BP (protein targeting to Golgi, Nuclear-transcribed mRNA catabolic process, Odontogenesis of dentin-containing tooth, and phosphatidylinositol biosynthetic process), and 3 categories of MF (protein binding, ubiquitin conjugating enzyme binding, and protein homodimerization activity) (Table 4).

Construction of the PPI networks and Hub genes analysis

DEGs for EMB and PB samples were imported independently into Cytoscape software (version 3.8.2) for protein–protein interaction analysis. After setting a confidence score cutoff >0.4 as the threshold value and removing proteins without obvious connection with majority members, PPI networks for each group were constructed and visualized using Cytoscape software (Figure 2). Next, the top 10 hub genes from PPI networks of EMB and PB samples were identified with degrees ≥10, which contained HLA-DQA1, CD74, TYROBP, HLA-B, HLA-C, CCL5, HLA-DRA, B2M, and PSMB9 in EMB samples and SLC4A1, NEDD4L, PIKFYVE, AHSP, SIAH2, ALAS2, RNF14, SNCA, UBA6, and GYPB in PB samples (Figure 3). The top 10 hub genes from EMB were obviously enriched in GO, including protein binding, integral component of membrane, plasma membrane, membrane, immune response, cell surface, extracellular exosome, and T cell receptor signaling pathway. In the KEGG pathway, antigen processing and presentation, herpes simplex infection, and cell adhesion molecules (CAMs) were included. In PB samples, the most enriched GO by top 10 hub genes was protein binding (MF) then followed by cytoplasm, cytosol, intracellular (CC), protein polyubiquitination B, intracellular signal transduction, phosphatidylinositol biosynthetic process, and negative regulation of apoptotic process (BP) (Tables 5, 6).

Discussion

Diagnosing cardiac allograft rejection is mostly based on EMB and PB testing. EMB possesses relatively higher accuracy, but it requires special equipment and skilled physicians, and it is associated with some complications in patients [10]. On the contrary, PB is the most convenient minimally invasive test for cardiac allograft rejection diagnosis and can be performed repeatedly, but its results can be easily influenced by infection or other pathological conditions [11]. Therefore, comprehensively considering the research and analysis of gene expression changes in EMB and PB can potentially improve the diagnosis and treatment of allogeneic cardiac rejection.

In this study, the first step completed was that the top 10 DEGs from EMB and PB samples were selected as representatives, and functions of these genes were annotated. In EMB samples, most of the top 10 DEGs were related to antigen presentation, immune cell activation, and inflammatory process. Genes closely related to immune functions were: IGHG3 (Immunoglobulin Heavy Constant Gamma 3), IGL@ (Immunoglobulin Lambda Locus), CD2, SCYA5 (C-C Motif Chemokine Ligand 5, CCL5), CCL18 (C-C Motif Chemokine Ligand 18), UBD (Ubiquitin D), CSPG2 (Versican), COPG (COPI Coat Complex Subunit Gamma), and NPPB (Natriuretic Peptide B). For example, upregulation of IGHG3 and IGL@ reflects B lymphocytes (B cells) activation, proliferation, differentiation, and immunoglobulins secretion in humoral immunity [12,13]. CD2 is a costimulatory gene highly expressed on activated T and NK cell surfaces [14]. SCYA5 can be released by activated CD8+ T cell, vessel smooth muscle cells, and endothelial cells to induce immune and inflammation reaction by binding to receptors like CCR1, CCR3, and CCR5 [15,16]. CCL18 (C-C Motif Chemokine Ligand 18) not only activates T cells in cellular immune rejection but also promotes immune function of B cells and macrophages in humoral immune response [17]. Additionally, UBD (Ubiquitin D) is a T cell pro-inflammatory cytokine that promotes T cells to excrete TNF-α, inducing cellular apoptosis [18]. Moreover, CSPG2 and COPG are cellular component and function executor proteins, both of which are expressed in endothelial cells, and their increased expression can indicate NK cell response and endothelial activation [19,20]. There are also some genes reflecting heart injury and functional compromise. NPPB mediates cardio-renal homeostasis by regulating the extracellular fluid volume and maintaining the fluid-electrolyte balance, and it increases in patients with congestive heart failure (CHF). SCYA5, mentioned above, is involved in myocardial infarction [21]. In summary, these critical gene changes may reflect the immune landscape and tissue injuries in the EBM samples, and finally reveal tissue pathological conditions, if present.

Next, gene changes in PB groups were evaluated. In PB groups, apart from immune recognition, inflammatory process, and cellular injury, the top 10 DEGs also reflect systematic hypoxia and specific corresponding change of erythrocytes, and these findings are different from those in EBM samples. Firstly, the main genes involved in immune reaction, inflammatory, and organ injury are PBX1 (Pre-B-Cell Leukemia Homeobox 1), IFI27 (interferon Alpha Inducible Protein 27), and RNF182 (ring finger protein 182). Specifically, PBX1 upregulation may indicate B cell activation and expansion [22], and the increased expression of IFI27 and RNF182 not only reflects innate immune response, but also indicates injuries to vessels and myocardium induced by inflammation and ischemia-reperfusion [23,24]. Most importantly, genes indicating hypoxia situation and corresponding changes of erythrocytes were found in PB samples. The expressions of FECH (Ferrochelatase) and BPGM (Bisphosphoglycerate Mutase) usually increase in hypoxia, which then activate and up-regulate hypoxia-responsive genes like hypoxia-induced factor-1 (HIF-1) to compensate for oxygen deficiency [25]. BPGM can decrease the oxygen affinity of hemoglobin in red blood cells, thus helping to compensate for tissue oxygen insufficiency [26]. Meanwhile, XK (X-Linked Kx Blood Group) and HEMGN (Hemogen) are positively associated with erythrocyte and hematopoietic cell generation, proliferation, and differentiation [27,28], which can be stimulated by hypoxia and inflammation [29-31]. These changes in hypoxia- and inflammation-related genes may indicate hypoxia and inflammation in the transplanted heart and can be used as biomarkers.

Next, through enrichment analyses of DEGs from EMB and PB samples, we selected the most centralized GO and KEGG pathway. The top 10 hub genes of PPI networks from EMB and PB were screened out and their involvement in GO and KEGG pathways were analyzed. Firstly, in the EMB group, the top 10 hub genes were enriched in GO of integral component of membrane, plasma membrane, and cell surface as CC (cellular component), immune response, T cell receptor signaling pathway, innate immune response, and inflammatory response as BP (biological process), as well as KEGG pathway of antigen processing and presentation and cell adhesion molecules (CAMs). These categories of GO and KEGG pathways were also reported in recent organ transplantation rejection research [32]. Hub genes, which were enriched in these categories of GO and KEGG pathways mentioned above, also play roles in antigen presentation, immune cell activation, and inflammation. Specifically, HLA class I heavy chain paralogues or component (HLA-B, HLA-C and B2M (Beta-2-Microglobulin) and HLA class II alpha chain paralogues and component (HLA-DQA1, HLA-DRA and CD74) play important roles in the stage of immune recognition, which can recognize exogenous antigen and present them to immune cells like T cells and NK cells [33-36]. TYROBP (Transmembrane Immune Signaling Adaptor) and PSMB9 (Proteasome 20S Subunit Beta 9) participates in the immune activation and immune response processes. They are receptors located on some immune cells, including NK cells, T cells, B cells, dendritic cells, and macrophages, and they are responsible for the activation and regulation of these immune cells [37,38]. In the PB group, categories of mostly enriched GO terms in top 10 hub genes included protein binding (MF), intracellular signal transduction, and negative regulation of apoptotic process (BP). The hub genes enriched in these GO terms above are mainly involved in systematic inflammation, hypoxic stress, cellular injury, and erythrocyte changes. In inflammation, RNF14 (Ring Finger Protein 14) serves as a T cell factor/lymphoid enhancer factor protein in the cellular immune response. Meanwhile, PIKFYVE (Phosphoinositide Kinase and FYVE-Type Zinc Finger Containing) loss in macrophages can damage the lysosomal degradative function, which promotes systemic inflammation [39]. SIAH2 (Siah E3 Ubiquitin Protein Ligase 2) and UBA6 (Ubiquitin Like Modifier Activating Enzyme 6) are involved in hypoxic stress and cellular injury. SIAH2 upregulation can promote stability of HIF-1 responding to hypoxia stress against cell apoptosis [40], while UBA6 downregulation can exacerbate cellular apoptosis in organ ischemia [41]. Finally, ALAS2 (5′-Aminolevulinate Synthase 2), SLC4A1 (Solute Carrier Family 4 Member 1 (Diego Blood Group)), GYPB (GlycophornsB), and AHSP (Alpha Hemoglobin Stabilizing Protein) are genes regulated to bio-activities. Upregulation of these genes can compensate for organ function impairment and may be a marker of the early stages of rejection. In detail, ALAS2 is the specific enzyme located on erythroid mitochondria, SLC4A1 is the integral membrane glycoprotein of the erythrocyte membrane, and GYPB (Glycophorins B) is major sialoglycoproteins of human erythrocytes. Upregulation of these 3 genes can indicate erythrocyte proliferation and activation [42-44]. AHSP (alpha hemoglobin stabilizing protein) is a molecular chaperone which can bind nascent apo-α-globin in the erythroid cell for its development [45], and can also promote NO function of vasodilatation to regulate blood pressure by binding α-globin in arteriolar endothelial cells [46].

In this study, we applied bioinformatics analysis to identify cardiac allograft rejection DEGs and relevant GO, KEGG Pathway, PPI networks, and hub genes from EMB and BP samples. The most significant genes and enrichment were analyzed in-depth as representations. We found that some DEGs from EMB and PB samples are involved in antigen presentation, immune cells activation, inflammatory process, and cellular injuries such as HLA-B, HLA-C, CD74, IGL@, SCYA5, CCL18, IFI27, RNF14, UBD, and UBA6. These genes can be used as biomarkers to identify specific kinds of or the degree of immune activation, so they can guide us to take some action to palliate immune injury before deterioration. In EMB samples, we screened out some genes related to heart tissue injury and cardiac function compromise, including SCYA5 and NPPB. The changed expression of these genes may indicate the deterioration of a cardiac graft, which can remind physicians not only to deal with immune rejection problems but also to address cardiac malfunction [47,48]. In the PB samples, genes containing FECH, BPGM, SIAH2, HEMGN, GYPB ALAS2, SLC4A1, and AHSP regulate hypoxia-induced factor-1 and erythrocyte’s protein synthesis, proliferation, and oxygen affinity as a compensatory mechanism. Changes in expression of these genes sometimes precede the apparent symptoms of immune rejection and cardiac dysfunction. In this way, they can be used as diagnostic early-stage biomarkers and therapeutic targets. Our study suggests how to interpret the nuances of biomarkers of cardiac rejection between EMB and PB samples and to consider the results together. In addition, unlike previous cardiac allograft rejection bioinformatic studies in which only single-center specimens were used or a single kind of sample was analyzed [7,32,49], our present results are more comprehensive. However, a limitation of the present study is that we did not perform in vitro or in vivo experiments on DEGs and their function relating to cardiac graft rejection. Hence, further mechanism studies of DEGs should be conducted to identify biomarkers in a more specific and accurate way to provide an experimental basis for the diagnosis and treatment of cardiac graft rejection at the very early stage.

Conclusions

In this bioinformatics analysis, DEGs from EMB and PB samples of patients with and without cardiac graft rejection were screened and have potential to be candidate biomarkers of rejection diagnosis and treatment. These biomarkers, including CD2, SCYA5, CCL18, UBD, NPPB, HLA-B, HLA-C, PSMB9, B2M, and TYROBP in EMB samples, and PBX1, IFI27, FECH, BPGM, HEMGN, RNF14, PIKFYVE, SIAH2, UBA6, ALAS2, SLC4A1, GYPB, and AHSP in PB samples, can help physicians to better manage cardiac rejection problems. However, more samples should be included and further fundamental investigations should be conducted to identify the translational value of these biomarkers.