Skip to main content

ORIGINAL RESEARCH article

Front. Immunol. , 15 December 2021

Sec. Autoimmune and Autoinflammatory Disorders

Volume 12 - 2021 | https://doi.org/10.3389/fimmu.2021.794800

Novel Insights Into Gene Signatures and Their Correlation With Immune Infiltration of Peripheral Blood Mononuclear Cells in Behcet’s Disease

Haoting Zhan,&#x;Haoting Zhan1,2†Haolong Li,&#x;Haolong Li1,2†Linlin Cheng,Linlin Cheng1,2Songxin Yan,Songxin Yan1,2Wenjie ZhengWenjie Zheng3Yongzhe Li,*Yongzhe Li1,2*
  • 1Department of Clinical Laboratory, Peking Union Medical College Hospital, Peking Union Medical College and Chinese Academy of Medical Sciences, Beijing, China
  • 2Department, State Key Laboratory of Complex, Severe and Rare Diseases, Peking Union Medical College Hospital, Chinese Academy of Medical Science and Peking Union Medical College, Beijing, China
  • 3Department of Rheumatology and Clinical Immunology, Peking Union Medical College Hospital, Chinese Academy of Medical Sciences & Peking Union Medical College, National Clinical Research Center for Dermatologic and Immunologic Diseases (NCRC-DID), Ministry of Science & Technology, State Key Laboratory of Complex Severe and Rare Diseases, Key Laboratory of Rheumatology and Clinical Immunology, Ministry of Education, Beijing, China

Background: Behcet’s disease (BD) is a chronic inflammatory disease that involves systemic vasculitis and mainly manifests as oral and genital ulcers, uveitis, and skin damage as the first clinical symptoms, leading to gastrointestinal, aortic, or even neural deterioration. There is an urgent need for effective gene signatures for BD’s early diagnosis and elucidation of its underlying etiology.

Methods: We identified 82 differentially expressed genes (DEGs) in BD cases compared with healthy controls (HC) after combining two Gene Expression Omnibus datasets. We performed pathway analyses on these DEGs and constructed a gene co-expression network and its correlation with clinical traits. Hub genes were identified using a protein–protein interaction network. We manually selected CCL4 as a central hub gene, and gene-set enrichment and immune cell subset analyses were applied on patients in high- and low-CCL4 expression groups. Meanwhile, we validated the diagnostic value of hub genes in differentiating BD patients from HC in peripheral blood mononuclear cells using real-time PCR.

Results: Twelve hub genes were identified, and we validated the upregulation of CCL4 and the downregulation of NPY2R mRNA expression. Higher expression of CCL4 was accompanied by larger fractions of CD8 + T cells, natural killer cells, M1 macrophages, and activated mast cells. Receiver operator characteristic curves showed good discrimination between cases and controls based on the expression of these genes.

Conclusion: CCL4 and NPY2R could be diagnostic biomarkers for BD that reveal inflammatory status and predict vascular involvement in BD, respectively.

Introduction

Behcet’s disease (BD), considered a variable vessel vasculitis, is characterized by a wide spectrum of clinical features including mucocutaneous lesions, non-granulomatous uveitis, lower extremity vein/cerebral venous sinus thrombosis, aortic aneurysms, and neurological involvement (1). The diagnostic criteria of the international study group for BD requires oral aphthous as the entry symptom, with the concurrence of genital ulcers, ocular or vascular involvement, as well as pathergy test positivity, which renders an extreme specificity (97.5%) (2). Nonetheless, its dynamic manifestation with disease relapse and remission periods, distinct population and sexual clinical expression, and disease type subsets pose grand challenges in identifying novel biomarkers for the diagnosis of BD.

Genetic risk factors impact BD pathogenesis, as the occurrence of BD is 5.78 times higher in HLA B51/B5 carriers than in noncarriers (3). BD susceptibility loci, for instance IL-10, STAT, and IL-23R polymorphisms (4, 5), indicate the predominant pathogenetic role of inflammation in the disease etiology. MicroRNA (miRNA), influenced by inflammatory signaling pathways, regulates mRNA expression through pre- and post-translational processes, and miRNA variants have shed light on disease status (6). A recent study in Turkey identified the FAS rs1800682 polymorphism and mir146a as protective factors against BD (7). Ibrahim and Jadideslam confirmed the association between mir146a expression and disease activity, with involvement of the eyes and vasculature (8, 9). Downregulation of miR-155 and that of miR146a likely exert their effects by disturbing Th1/Th17 cell and cytokine homeostasis (10). Profiling analysis of miRNA expression has also identified epigenetic pathways that are enriched in vascular biology, the immune response, and apoptosis (11). Long non-coding RNAs (lncRNAs) are also involved in BD etiology, with the homozygous T allele of rs9517723 in lncRNA LOC107984558 correlated with ocular and central nervous system (CNS) lesions, the induction of UBAC2 expression, and overactivation of the ubiquitination-related pathway (12). However, we could not only concentrate on individual gene and ignore the potential interactions between them which might be the backbone of intracellular communications and etiology pathways.

Hub genes are highly connected nodes between co-expression network, exploring of which is vital for seeking diagnostic biomarkers or therapeutic drug targets and clarify the biological regulation process in disease conditions. In this study, we utilized weighted gene co-expression (WGCNA) network and protein-protein interaction (PPI) network to identify hub genes. Since WGCNA analysis could provide a systematic insight into co-expression gene sets and their relationships with clinical phenotypes, we invoked this algorithm to establish a scale-free network using adjacency matrix from which topological overlap matrix and corresponding dissimilarity were converted to differentiate gene co-expression modules (13). For subsequent analysis, we rigorously chose the turquoise module as a candidate module which gains highest correlation coefficient and significant corresponding p value with clinical features among patients suffering Behcet’s disease and healthy controls. Recently, a few of studies have comprehensively utilize WGCNA to filtrate synergistically altered gene sets in autoimmune diseases (1418). Up till now, amassed structure and function of single protein fragmented our comprehension of their intracellular interactions and pathogenesis pathways. By conducting PPI program, we could easily attain a tighter protein-based co-expressed DEG channel and a more integrated interaction network, the connectivity of which is ranked by combined score calculated from online sources including text-mining, experiments, databases and genomic context information, etc. (19). To predict hub genes bridged through networks, we employed MCODE plugin to visualize intensively related modules and CytoHubba package to screen out hub nodes avail of functional enrichment.

Despite the body of research elucidating the genetic variation contributing to BD, a deep exploration of its molecular mechanisms, interaction networks, and key pathways is lacking. In the present study, we aimed to uncover hub genes acting as biomarkers for the diagnosis, classification, and biological functions and networks of BD using bioinformatics methods and a validation cohort.

Materials and Methods

Dataset Acquisition and Normalization

Two microarray profiling datasets from peripheral blood mononuclear cells (PBMCs) (GSE70403 and GSE17114) were downloaded from the Gene Expression Omnibus (GEO) (http://www.ncbi.nlm.nih.gov/geo/) with the screening terms “Behcet Syndrome,” “Triple Symptom Complex,” “Behçet Disease,” “Behcet’s Disease,” “Behcet Disease,” “Adamantiades Behcet Disease,” and the organism “Homo Sapiens.” We merged the two datasets and eliminated batch effects utilizing the R package “sva” (20). Consequently, 14 healthy controls (HCs) from GSE17114 and 56 BD samples (15 from GSE17114 and 41 from GSE70403) were combined.

Differentially Expressed Gene and Enrichment Analysis

Differentially expressed gene (DEG) analysis was conducted using the “limma” package in R with the thresholds of an adjusted P < 0.05 and |log Fold Change (FC)| > 0.5 (21). To identify the biological functions and corresponding pathways involving these DEGs, we performed Gene Ontology (GO) enrichment and visualized the biological process, cellular component (CC), and molecular function (MF) categories and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis of signaling pathways using the R packages “clusterProfiler” and “enrichplot” (22). The significance cutoff criteria were P < 0.05 and FDR < 0.05.

Identification of Hub Genes and Construction of a Protein–Protein Interaction Network

We clustered PBMC-specific gene modules and established a co-expression network using the R package “WGCNA.” The obvious superiority of WGCNA analysis is to cluster DEGs into co-expression modules and ultimately to identify the specific gene module which is most relevant to clinical phenotype of Behcet’s disease. After shearing samples below the abline (h=20000), we constructed the co-expression network with soft thresholding power to obtain a higher level of scale free R2 and mean connectivity (13). In dynamic tree cut and module identification section, we altered 10 as the minimum number of gene modules. Utilizing the clinical traits data containing BD patients and healthy controls from GSE70403 and GSE17114, the gene significance (GS) and module membership (MM) were calculated. The relationship between the gene modules and clinical traits was represented in the form of a heatmap, from which we selected a module that was most clinically relevant.

With a threshold score of 0.150, we constructed a protein–protein interaction (PPI) network consisting of the DEGs using the STRING database (19). We have realized from previous reference (23) that the combined score from STRING database is meant to express an approximate confidence of the association between proteins being true based on every channel of evidence (depended on genomic context information, co-expression, text-mining, experiments and databases). Thus, we have ranked the first 10 interactions by the combined score in the PPI network. To attain a higher degree of connectivity, we applied CytoHubba and then extracted the key Hubba nodes ranked by maximum clique centrality (MCC) in Cytoscape software (24). Using the MCODE plugin (25), we sought dense regions divided by the critical value of a degree cut-off = 2, node score cut-off = 0.2, k-core = 2, and max. Depth = 100.

Furthermore, we united the first 10 genes ranked by the combined score in the PPI network and node genes filtered with an MCC > 5 in CytoHubba analysis, together with genes exported using the MCODE program, and 12 hub genes were ultimately identified for principal component analysis (PCA) in the GSE17114 dataset where Oğuz et al. have divided BD patients into three subtypes containing isolated mucocutaneous manifestations (MB), ocular involvement (OB), and large vein thrombosis (VB) subtypes in accordance with the major clinical symptoms provided by Xavier et al. (26, 27) (detailed methodology is provided in online Supplementary Materials).

Correlation of Hub Genes and Vasculitis

Genecards (https://www.genecards.org/) was screened using the keyword “vasculitis.” We ranked the relevance score given by website and selected top 15 genes, then revisited the normalized gene expression data merged from Method 2.1 to see which genes from Genecards have the intersections with GEO datasets. Thus, we deleted 3 genes which were not documented in normalized gene expression data merged from GEO datasets (ADA2, HLA-DRB, C4A). Finally, the top 12 genes (PRTN3, PTPN22, CTLA4, DNASE1L3, MPO, MEFV, HLA-B, HLA-DPA1, HLA-DPB1, IL-10, TNF and CRP) were selected. According to the relevance score, the top 12 genes were selected. Moreover, we visualized the correlation heatmap using the R package “corrplot” (28).

Gene-Set Enrichment Analysis

To explore the biological functions of hub genes, we downloaded Gene-Set Enrichment Analysis (GSEA) software (http://software.broadinstitute.org/gsea/index.jsp). The patients were divided into high/low-expression groups according to the median expression of the central hub genes. Afterward, we selected “c2.cp.kegg.v6.1.symbols.gmt” as the reference gene-set in the Molecular Signatures Database (29) and diagrammed the KEGG enrichment plot of distinct expression groups by normalized P < 0.05 and FDR < 0.2.

Correlation Analysis of Gene Signatures and Immune Cells

CYBERSORT is a deconvolution algorithm used to determine the relative proportions of 22 human immune cells from bulk expression data (30). Using the R packages “e1071,” “preprocessCore,” and “limma,” we performed immune cell profiling among the BD samples. Samples with a P < 0.05 were included. To estimate the diverse immune cell patterns in different central hub gene expression groups, we used “vioplot” in R studio.

Diagnostic Prediction and Validation of Hub Genes

To validate the expression pattern and diagnostic value of the selected hub genes, we collected blood samples from 16 BD patients and 16 age- and sex-matched HCs under a protocol approved by the Medical Ethics Committee of Peking Union Medical College Hospital with informed consent acquired from the enrolled subjects (demographic and clinical data of the validation population and controls in online Supplementary Materials, Table S3). We then isolated PBMCs by Ficoll density gradient centrifugation. Further, total RNA was extracted from PBMCs using the Trizol Reagent, and cDNA was reverse-transcribed following the manufacturer’s instructions (Mei5 Biotechnology, Co., Ltd). For amplification, SYBER green (Mei5 Biotechnology, Co., Ltd) was utilized for real-time PCR on a Roche LightCycler 480 system. Finally, the relative expression levels were determined using the 2−∆∆Ct formula after normalization with endogenous GAPDH mRNA expression. The primer sequences were as follows: CCL4 forward CTTTTCTTACACCGCGAGGA and reverse GCTTGCTTCTTTTGGTTTGG, AGTR2 forward TTCCCTTCCATGTTCTGACC and reverse AAACACACTGCGGAGCTTCT, NPY2R forward CAAAACTTCTCCTCCAGTCCCC and reverse AAGGTGGGAGGATCAGAGATGG, ZIC1 forward GCATCCCAGTTCGCTGCGCAAA and reverse GGAGACACGATGGTGGGAGGCG, and GADPH forward CCTCAAGATCATCAGCAAT and reverse CCATCCACAGTCTTCTGGGT. Receiver operator characteristic (ROC) curves were rendered using the R package “pROC,” and DEGs between BD and HC were visualized using “ggpubr.”

Results

Screening for DEGs Between BD and HC

The workflow illustrated in Figure 1 demonstrates the procedure of exploring hub genes and diagnostic markers in BD. Two eligible datasets GSE70403 and GSE17114 with 56 BD and 14 HC combined were included. In total, 82 genes were differentially expressed in BD compared with HC (details provided in Table S1). Meanwhile, 11 DEGs were upregulated, and 71 were downregulated (Figures 2A, B).

FIGURE 1
www.frontiersin.org

Figure 1 Flow diagram of exploring hub genes and diagnostic markers in BD.

FIGURE 2
www.frontiersin.org

Figure 2 Differentially expressed genes displayed by volcano plot (A) and heatmap (B).

GO and KEGG Annotation of DEGs

We investigated the pathogenesis of BD by applying functional and pathway analysis. Several DEGs were prominently enriched in “clathrin-coated endocytic vesicle membrane” and “clathrin-coated endocytic vesicle” terms in the context of CC. For MF, DEGs were preferentially included in “IgG binding,” “immunoglobulin binding,” and “carbohydrate binding” demonstrating the immunopathogenic background of BD (Figure 3A). Various disease-associated KEGG pathways were hyperactivated, including “Fc gamma R-mediated phagocytosis,” “Staphylococcus aureus infection,” and “Leishmaniasis” (Figure 3B), indicating that the pathogenesis of BD might be triggered by innate immune activation and immunological response after infection (Figures 3C, D).

FIGURE 3
www.frontiersin.org

Figure 3 Gene Ontology and KEGG annotation of DEGs (A) Histogram of Gene Ontology enrichment analysis: several DEGs were prominently enriched in “clathrin-coated endocytic vesicle membrane” and “clathrin-coated endocytic vesicle” terms in the context of cellular component (CC). For molecular function (MF), DEGs were preferentially included in “IgG binding”, “immunoglobulin binding”, and “carbohydrate binding” which demonstrates the immunopathogenic background of BD. (B) Histogram of KEGG pathway analysis: various disease-associated KEGG pathways were hyperactivated, including “Fc gamma R-mediated phagocytosis,” “Staphylococcus aureus infection,” and “Leishmaniasis” which indicates the pathogenesis of BD might be triggered by innate immune activation and immunological response after infection. (C, D) Circle plot of Gene Ontology and KEGG annotations.

Identification of Hub Genes

Five clusters of modules were constructed with WGCNA co-expression network analysis (Figures 4A, B and Table S2). A sample dendrogram and trait heatmap calculated using Pearson’s correlation coefficients indicated that there were no outlier samples (Figure 4C). When incorporating the expression profile with clinical traits, we found that the turquoise module was significantly associated with clinical phenotypes and contained 72 genes (r = 0.60, P < 0.05) (Figures 4D, E). A Venn plot intersecting the DEGs and the turquoise module narrowed the hub gene candidates down to 71 (Figure 4F).

FIGURE 4
www.frontiersin.org

Figure 4 Identification of hub genes (A) Gene dendrogram and the modules of DEGs correlated to BD in PBMC of 56 patients and 14 HCs. (B) Clustering of module eigengenes (C) Sample dendrogram and trait heatmap (D) Module-trait relationships: every module has its correlation coefficient and corresponding p value. (E) Module membership in turquoise module: the gene significance (GS) and module membership (MM) were calculated and their correlation was 0.93 (P <0.001) which indicates the strong relationship between module genes and clinical traits. One dot represents one gene in turquoise module, and dots at upper left corner were highly related to both the turquoise module and clinical traits. (F) Venn plot intersecting the DEGs and the turquoise module: there are 71 intersected DEGs for subsequent PPI network analysis.

Next, we input these 71 genes into PPI network analysis to investigate the interrelationship between the proteins encoded by these genes and discovered 70 interactions (Figure 5A). We have extracted the first 10 interactions by the combined score in the PPI network. CCL4, NPY2R, AGTR2, TAS2R1, ASB14, ASB17, C1orf110, SOX14, MAGEA1, NPAS4, EYA1 and HOXA11 was elected (Set 1). We filtered the top 10 genes ordered by a combined score and an MCC > 5 filtered by CytoHubba, including ZIC1, CTXN3, NPY2R, AGTR2, LRRC3B, EYA1, HOXA11, SOX6, CCL4, CAMKV, ANGPTL3, TAS2R1, SLC6A3 and SOX14 (Set 2). (Figure 5B and Table 1), along with subnetwork gene nodes exported using the MCODE program where NPY2R, LRRC3B, AGTR2, CTXN3, CAMKV, CCL4 and TAS2R1 was clustered as subset 1 (shown in Figure 5C-ii) and HOXA11, EYA1, SOX6 as subset 2 (shown in Figure 5C-i), they are all named as Set 3. (Figure 5C and Table 2). Eventually, 12 hub genes (AGTR2, CAMKV, CTXN3, EYA1, HOXA11, LRRC3B, NPY2R, SOX14, SOX6, TAS2R1, ZIC1, and CCL4) were identified (shown in Figure S1: a union of the intersection of set 2 and 3 plus the intersection of set 1 and 2). ZIC1 is the top 1 significant node gene ranked by MCC, therefore we could not neglect it for further analysis and also add it into hub genes. These hub genes were identified for PCA analysis in GSE17114 to divide BD patients based on clinical manifestations (26). These hub genes could well discriminate isolated mucocutaneous manifestations (MB), ocular involvement (OB), and large vein thrombosis (VB) patients diagnosed with BD (Figure 5D).

FIGURE 5
www.frontiersin.org

Figure 5 PPI network analysis (A) Interrelationship between the proteins encoded by DEGs from STRING database. The thickness of connecting lines suggests the strength of the association between proteins. (B) Key Hubba nodes ranked by maximum clique centrality (MCC): the top 10 candidate hub genes ranked by MCC linkage degrees were ZIC1, CTXN3, NPY2R, AGTR2, LRRC3B, EYA1, HOXA11, SOX6, CCL4 and CAMKV. (C) Subnetwork gene nodes exported using the MCODE program: the hub nodes linking highly connected modules were NPY2R, LRRC3B, AGTR2, CTXN3, CAMKV, CCL4, TAS2R1, HOXA11, EYA1 and SOX6. (D) principal component analysis using the GSE17114 dataset: these hub genes could discriminate isolated mucocutaneous manifestations (MB), ocular involvement (OB), and large vein thrombosis (VB) patients diagnosed with BD from in GSE17114.

TABLE 1
www.frontiersin.org

Table 1 Hubba nodes ranked by maximum clique centrality (MCC) in CytoHubba.

TABLE 2
www.frontiersin.org

Table 2 DEGs ranked by combined score in PPI network.

Correlation Between Hub Genes and Vasculitis

Genecards was searched, and the top 12 genes associated with vasculitis were exported (PRTN3, PTPN22, CTLA4, DNASE1L3, MPO, MEFV, HLA-B, HLA-DPA1, HLA-DPB1, IL10, TNF, and CRP). The correlation heatmap revealed tight relationships among these 12 hub genes with Pearson’s r ranging from 0.99 to 1 (P < 0.05); however, except for CCL4 and TNF (r = 0.56, P < 0.05), there was no relationship between the hub genes and the output genes from the Genecards website (Figure 6).

FIGURE 6
www.frontiersin.org

Figure 6 Correlation heatmap between hub Genes and vasculitis. The top 12 genes associated with vasculitis were exported (PRTN3, PTPN22, CTLA4, DNASE1L3, MPO, MEFV, HLA-B, HLA-DPA1, HLA-DPB1, IL10, TNF, and CRP) from Genecards together with 12 hub genes from 12 hub genes (AGTR2, CAMKV, CTXN3, EYA1, HOXA11, LRRC3B, NPY2R, SOX14, SOX6, TAS2R1, ZIC1, and CCL4) from above WGCNA and PPI analysis.

GSEA Analysis of Central Hub Genes

Through an intensive literature search in PubMed, we selected CCL4 as the central hub gene because of its known pathogenic role in autoimmune diseases including rheumatoid arthritis (RA) (31), multiple sclerosis (MS) (32) and Behcet’s disease (33). Accordingly, our patients were categorized into high- and low-expression groups of CCL4 for GSEA function and pathway analysis. Twelve KEGG pathways, “natural killer cell-mediated cytotoxicity,” “glycosaminoglycan biosynthesis chondroitin surface,” “antigen processing and presentation,” “regulation of autophagy,” “graft versus host disease,” “type I diabetes,” “T cell receptor signaling pathway,” “cytosolic DNA sensing pathway,” “TGF-β signaling pathway,” “chemokine signaling pathway,” “cytokine-to-cytokine receptor interaction,” and “toll-like receptor signaling pathway,” were enriched in the CCL4 high-expression BD patients (Figure 7).

FIGURE 7
www.frontiersin.org

Figure 7 GSEA analysis of central hub genes. BD patients were categorized into high- and low-expression groups of CCL4 for GSEA function and pathway analysis from which twelve KEGG pathways were enriched.

Landscape of Immune Cell Proportions

We calculated the relative percentage of immune cells among BD patients and HC (Figure 8A). In comparison with the low-expression group, the high-expression group had a higher fraction of CD8 + T cells (P = 0.019), Natural killer (NK) cells (both resting and activated, P = 0.041 and 0.010, respectively), M1 macrophages (P = 0.044), and activated mast cells (P = 0.041). Inversely, monocytes (P = 0.007), M0 macrophages (P = 0.006), and resting mast cells (P = 0.008) were decreased in BD patients with high CCL4 expression (Figure 8B).

FIGURE 8
www.frontiersin.org

Figure 8 Immune cell proportions among BD patients and HC according to the high/low expression of central hub genes. (A) The relative percentage of immune cells are presented. (B) Violin plot shows the CCL4 high-expression group had a higher fraction of CD8 + T cells (P = 0.019), Natural killer (NK) cells (both resting and activated, P = 0.041 and 0.010, respectively), M1 macrophages (P = 0.044), and activated mast cells (P = 0.041) in contrast with low-expression group.

Diagnostic Performance of Hub Genes and Validation

We then performed ROC analysis to detect the performance of hub genes in differentiating patients with BD from HC. The area under the curve (AUC) ranged from 0.81 to 0.88 (Figure 9A). Additionally, CCL4 was highly expressed in our validation cohort (P = 0.0234), whereas NPY2R was downregulated (P = 0.024) (Figures 9B, C and Table 3); however, no significant difference was found in AGTR2 and ZIC1 expression among BD and HC individuals (data not presented).

FIGURE 9
www.frontiersin.org

Figure 9 Diagnostic performance of hub genes and validation (A) ROC analysis to detect the performance of hub genes in differentiating patients with BD from HC. (B) ROC validation and differential expression among BD and HC of CCL4. (C) ROC validation and differential expression among BD and HC of NPY2R.

TABLE 3
www.frontiersin.org

Table 3 mRNA expression in PBMC of BD and HC from validation cohort.

Discussion

In this current study, we knit in-silico results and the in vitro data closely then disclosed the biological functions as well as pathogenic pathways involved in Behçet’s disease (BD). To begin with, we downloaded and merged 2 GEO datasets GSE70403 and GSE17114 by microarray profiling of peripheral blood mononuclear cells from BD patients. WGCNA algorithm identified the most relevant gene set to clinical phenotype of Behçet’s disease which named as turquoise module. Subsequently, PPI network (CytoHubba and MOCDE plugin) was employed to explore the co-expression interactions between DEGs from which hub nodes were extracted (AGTR2, CAMKV, CTXN3, EYA1, HOXA11, LRRC3B, NPY2R, SOX14, SOX6, TAS2R1, ZIC1, and CCL4). It turns out that these above hub genes could be used to discriminate BD patients into MB, OB, and VB subtypes in accordance with the major clinical symptoms provided by Xavier et al. in GSE 17114 dataset which additionally demonstrated a reliable antidiastole performance (26, 27). Wondering which hub gene is previously considered as the potential pathogenic manipulator for autoimmune disease in references, we retrieved PubMed and only to find CCL4, as a central hub gene, stands out in RA, MS and BD etiology. Therefore, we classified the BD patients from microarray datasets following the expression level of CCL4 for GSEA function, twelve KEGG pathways intensively enriched in intracellular and intercellular biological communications among immune cells. Hence, CYBERSORT was utilized to visualize the immune cell profiling based on the expression level of CCL4 in BD patients. The upregulated CCL4 could recruit a higher proportion of CD8 + T cells (P = 0.019), NK cells (both resting and activated, P = 0.041 and 0.010, respectively), M1 macrophages (P = 0.044), and activated mast cells (P = 0.041) into inflammatory sites (skin, eyes and even vascular), which likely exacerbates proinflammatory situation, immune cell function and intensify disease activity in BD patients. Last but not least, the performance of hub genes exerts a great ability in discerning BD and HC (AUC ranging from 0.81-0.88) where CCL4 and NPY2R are proved as a moderate diagnostic biomarker in disease recognition from healthy subjects in our validation cohort (clinical characteristics are provided in the online Supplementary Material, Table S1).

Taken together, we have identified the immune-activated, inflammation-induced, and carbohydrate-mediated pathogenesis of BD and identified DEGs and pathways involved in its etiology. Our results indicate that the identified hub genes could be of diagnostic value and reveal new insights into the etiological analysis of BD.

Subsequent validation disclosed that CCL4 expression is significantly upregulated in BD patients compared to HC. Additionally, we found that there was a higher proportion of CD8 + T cells, NK cells, M1 macrophages, and activated mast cells infiltrating in inflammation sites in the subset of BD patients with high CCL4 expression. CCL4 is a proinflammatory chemokine that extensively participates in autoimmune pathogenesis, with immunoregulatory effects on trafficking immune cells into sites of inflammation. For instance, CCL4/MIP-1β secreted by IL-2-activated NK cells amplified IFN-γ release and invoked other inflammatory cells like neutrophils, macrophages, and T cells, which are responsive to chemokines and aggravate the inflammatory process (34). In RA, CCL4 is vital for guiding CD56dimCD16bright NK cells into the synovium (35); in addition, it is positively correlated with IL-36 agonists, IL-36RA, and IL-38 cytokines, which intensify joint inflammation (36).

Additionally, CCL4 could indicate the deterioration of BD, as extremely high expression of CCL4 mRNA was found in patients who suffer from thrombosis and CNS involvement in our validation cohort. A research published by Schecter et al. demonstrated that in human smooth muscle cells enriched in plaque area CCL4/MIP-1β is abundant and introduce the elevation of tissue factor that involves in arterial coagulation together with thrombotic cascade (37). From snap shots of CCL4/MIP-1β participation in clotting formation, there is supposed to exist a hidden mechanism in causing thrombosis of BD patients. Previous evidence supports that upregulated CCL4 could be considered to be a circulating inflammatory cytokine marker in primary and systemic lupus erythematosus (SLE)-related autoimmune hemolytic anemia, indicating erythroid compensation and disease severity due to its positive relationship with reticulocyte count (38). Being a major macrophage attractant, universal involvement of CCL4 has been found in diabetes mellitus, since its circulating concentration increased regardless of disease stage, which may contribute to destruction of islet cells and the progression of diabetes. By providing an inflammatory microenvironment, CCL4 induces cell adhesion to endothelial cells through oxidative stress and tends to augment atherosclerosis and vascular injury (39). Additionally, elevated CCL4/MIP-1β in the cerebrospinal fluid is a characteristic of IgG4 anti‐neurofascin 155-positive chronic inflammatory demyelinating polyneuropathy, which implies exacerbation of protein leakage, damage of encephalic and CCL4/MIP-1β decreases with combined immune-treatment response (40). Della-Torre et al. demonstrated that B cells and fibroblasts might orchestrate inflammatory infiltration through CCL4, which in turn amplified the recruitment of monocytes, CD4 + cytotoxic T cells through CCR5, and IgG4-related disease lesions (41). In discriminating between relapsing–remitting and progressive clinical phenotypes of multiple sclerosis (MS), the plasma CCL4/MIP-1β level showed a good diagnostic value (AUC = 0.702, P < 0.0002) and was a protective factor preventing disease progression (odds ratio = 0.9873, P = 0.0171) (42). However, the pathogenic impact of CCL4 in BD has not been reported before. As the similar pathogenesis existing among autoimmune diseases, we could deduce that the expression of CCL4 might be a circulating inflammatory chemokine which aggravates inflammation and disease progression via inducing immune cell infiltration and damaging surrounding tissue in BD patients.

Neuropeptide Y receptor subtype 2 (NPY2R) is a G-protein binding receptor for neurotransmitter NPY, which triggers vascular smooth cell proliferation and vasoconstriction, and results in ischemic angiogenesis (43). In our study, NPY2R mRNA expression in PBMCs of BD patients, especially in those with severe symptoms including thrombosis and gastrointestinal involvement, was distinctly lower than that in HC. A genome-wide association study demonstrated that the NPY2R rs1902491 polymorphism was associated with a severe diabetic retinopathy subgroup without end-stage renal disease (44); however, no significant association was replicated in a small cohort with rapidly proliferative diabetic retinopathy in Lithuania (45). This may imply NPY2R is potentially correlated with microvascular leakage and obstruction which could stands for vascular destruction within BD patients. The gene encoding the type 2 angiotensin II (Ang II) receptor (AGT2R) plays a pivotal role in atherosclerosis, vascular inflammation, and remodeling and has anti-inflammatory and antiproliferative influences, and it exhibited a female advantage by counteracting the AT1 receptor in mice (46). In SLE, AGT2R polymorphisms containing the A allele (CA  +  AA) in exon 3 or G allele (AG  +  GG) in intron 1 are of apparent significance compared with HC (47). AGT2R expression is also strongly enhanced in RA synovial tissue, an exogenous agonist that could ameliorate inflammation (48). On the contrary, AGT2R blockers are postulated to dampen the clinical situation and laboratory metrics in RA patients (49). In idiopathic pulmonary fibrosis, AGT2R is highly expressed in interstitial fibroblasts and serves as a determinant of fibroblast proliferation and migration (50). In 883 Austrian MS patients and 972 controls. The ZIC1 rs1841770 polymorphism was associated with an earlier disease onset of MS (51), but a Spanish replication study did not support it as an MS predisposition factor (52). Nevertheless, we did not find any statistical difference in AGT2R and ZIC1 mRNA expression between BD and HC; this may be attributed to the small population of our validation cohort.

NK cells are innate immune lymphocytes that have a cytolytic immune-regulatory effect in autoinflammatory disease, particularly in BD. We confirmed that in patients with upregulated CCL4 expression, there is a higher fraction of NK cells infiltrating the inflammatory zone. However, conflicting studies have reported decreased, normal, and increased NK cells in BD compared with those in HC (53). It has been hypothesized that impaired NK cytotoxicity and degranulation facilitate and active disease progress (54, 55). The NK1/NK2 ratio is significantly increased, resulting in CD16 + IFN-γ + NK1-induced IFN-γ secretion and suppressed NK regulation in mucocutaneous BD. Moreover, increasing IFN-γ produced by NK cells reflects active disease and relapse of disease, which demonstrated the dominant interaction between NK cells and interferons (56). Accordingly, NKMHC-I modules encoding HLA-B*51 could also contribute to NK cell inhibition through binding killer immunoglobin-like receptors, further supporting hereditary factors as key drivers of BD etiology.

Hyperactivation of monocytes and neutrophils is a feature of BD (57). Our results identified a high level of monocytes in CCL4 low expression group, instead of neutrophils in BD. Monocytes in BD patients appear to have intensified oxidative burst activity compared with those in RA patients and HC (58). The enriched pathways JAK/STAT, IL-6, and interferon signaling are upregulated in CD14 + monocytes isolated from BD patients (59), whereas toll-like receptors 2, 4, and 5 are highly expressed in monocytes, causing a cascade of activation of nuclear factor-κB and production of TNF-α, unveiling a proinflammatory mechanism of BD (60). Bacterial infection might be a trigger of BD pathogenesis. Lipopolysaccharide (LPS) stimulation increases TNF-α secretion in vivo in monocytes of active BD rather than in quiescent ones (61). A recent study demonstrated a perturbed FcγRs balance in monocytes from BD patients that represented downregulated inhibitory FcγRIIb and upregulated activating FcγRIII, the level of which fluctuated along with BD treatment and disease severity. Moreover, LPS and IgG-stimulated monocytes produced more IL-6 than those from HC (57). IgG binding activates FcγRs and FcγRIIb, mainly leading to endocytosis and internalization of immune complexes; therefore, these studies corroborate our finding of “IgG binding,” “clathrin-coated endocytic vesicle membrane,” and “clathrin-coated endocytic vesicle” in our GO analysis and the enriched pathways “Fc gamma R-mediated phagocytosis” and “Staphylococcus aureus infection” in the KEGG analysis.

M1 macrophages stimulate an inflammatory reaction by triggering cytokines such as TNF-α and IL-6, resulting in Th1 reactivity, whereas BD patients with skin lesions conferred M1 macrophage dominance in systemic sclerosis (SSc). It confers that M1 macrophages is more like a predominant pathogenesis factor in BD than in SSc. A higher expression of CCR1 ligands for CCL3/MIP-1α than that in HC again demonstrated that a skewed M1 and M2 macrophage balance led to MIP-1α-induced inflammation (62). Macrophage incubation with active BD serum facilitated CD64 (FcγRI) positivity and IL-8 mRNA expression compared to inactive BD, suggesting an overwhelming proinflammatory response and aberrant endothelial alteration (63). Mast cells are prominently elevated in the skin lesions of active Behcet’s disease compared with those uninvolved skin (64). It seems that mast cells are attracted by chemotaxis influenced by complement cascades in BD patients with high-expressed CCL4, which is analogical to basophil hyperactivation process.

The profound role of IL-17A and CD8 + T cell-mediated immune response is pronounced in BD with papulopustular skin lesions. Investigation into Tc17 cell activation confirmed that CD8 + T cells are infiltrated in peripheral vascular walls and the subcutis in BD, and, concurrently, IL-17A + CD8 + T cells were elevated in BD compared to psoriasis vulgaris lesions and were a major source of IL-17A rather than CD4 + T cells (65). IFN-γ-producing CD8 + T cells are enhanced in the peripheral blood of active BD patients and are even higher in patients with pulmonary manifestations and neurological involvement. Furthermore, the Tc1/Tc2 ratio is dramatically increased in the active and remission phases of BD (66). Increased IFN-γ and CD8 + Leu7 + T cells together with an increased titer of IgG antibodies to HSV-1 demonstrated an antigen-driven T cell-dominant etiology characteristic of BD (67), which is in accordance with our GSEA pathway analysis results.

Our GO analysis highlighted the “carbohydrate binding” pathway, and accumulating evidence suggests that galectins regulate immune cell homeostasis, in which galectin-1 inhibits T cell proliferation, enhances apoptotic activity, and thus suppresses inflammation, whereas galectin-3 amplifies the inflammatory response and survival of T lymphocytes, which subsequentially activates NADPH oxidase, the superoxidative reaction of neutrophils, and recruits monocytes (68). This coincides with our finding of “alcohol dehydrogenase [NAD(P)+] activity” among the MF GO terms. The combination of lectins and epithelial cells and lectins binding to specific tissues triggers antigen humoral responses leading to autoimmune rheumatic responses, for instance, in RA or type I diabetes (69). S-nitrosated mannose binding lectin (SNO-MBL) is upregulated in RA patients and was dysfunctional in depositing C4, phagocytosing bacteria, and binding to apoptotic cells, which consequently caused enhancement of anti-SNO-MBL antibodies (70).

There are several limitations in the present study. First, although we have verified significant difference between BD and HC group in mRNA expressions of CCL4 and NPY2R, we could not associate our results with disease activity index (Behçet’s Disease Current Activity Form) attributed to lack of integrated clinical information of each patient. Second, we aimed to identify a diagnostic marker in discerning BD patients and healthy subjects, thus we did not add disease controls such as takayasu arteritis, inflammatory bowel disease, ANCA-associated systemic vasculitis and even patients with either coronary or venous thrombosis which shares semblable clinical feature with BD. Third, via MCODE program we have discovered the protein interaction network between CCL4 and NPY2R, future analysis should more concentrate on underlying up/down regulation between them which might be the mechanism leading to vasculitis moreover vascular deteriorations in BD. Further studies using more comprehensive data and larger validation cohort should be performed to elucidate the association of CCL4 and NPY2R with BD activity as well as clinical phenotypes.

In conclusion, as a proinflammatory factor synthesized in the inflammation zone, higher mRNA expression of CCL4 might be an upstream switch to attract more immune cells into inflammatory sites, these cells in turn positively amplify cytokine secretion and exacerbates BD disease stage. In vasculitis related to BD, NPY2R expression is potentially impaired and may be a marker of vascular complications. Future studies are required to elucidate the potential mechanisms underlying the dysregulated expression of target genes and their correlation with the immune microenvironment and clinical phenotypes in BD.

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 including human participants were approved by Medical Ethics Committee of Peking Union Medical College Hospital with informed consensus acquired from enrolled subjects. The patients/participants provided their written informed consent to participate in this study.

Author Contributions

YL conceived and designed the research. HZ extracted data, performed software analysis and visualized graphs and tables. LC, SY, HL, and HZ categorized and supervised graphs and tables. WZ provided the clinical data of each patient. HZ and HL went on validation experiments. HZ wrote the paper. All authors are accountable for all aspects of the study, and attest to the accuracy and integrity of the results. All authors contributed to the article and approved the submitted version.

Funding

This research was supported by grants from the National Key Research and Development Program of China (2018YFE0207300), the National Natural Science Foundation of China Grants (81871302,81871299) and Beijing Key Clinical Specialty for Laboratory Medicine - Excellent Project (No. ZK201000).

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/fimmu.2021.794800/full#supplementary-material

References

1. Yazici H, Seyahi E, Hatemi G, Yazici Y. Behçet Syndrome: A Contemporary View. Nat Rev Rheumatol (2018) 14(2):107–19. doi: 10.1038/nrrheum.2017.208

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Davatchi F. Diagnosis/Classification Criteria for Behcet’s Disease. Patholog Res Int (2012) 2012:607921. doi: 10.1155/2012/607921

PubMed Abstract | CrossRef Full Text | Google Scholar

3. de Menthon M, Lavalley MP, Maldini C, Guillevin L, Mahr A. HLA-B51/B5 and the Risk of Behçet’s Disease: A Systematic Review and Meta-Analysis of Case-Control Genetic Association Studies. Arthritis Rheum (2009) 61(10):1287–96. doi: 10.1002/art.24642

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Hou S, Yang Z, Du L, Jiang Z, Shu Q, Chen Y, et al. Identification of a Susceptibility Locus in STAT4 for Behçet’s Disease in Han Chinese in a Genome-Wide Association Study. Arthritis Rheum (2012) 64(12):4104–13. doi: 10.1002/art.37708

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Mizuki N, Meguro A, Ota M, Ohno S, Shiota T, Kawagoe T, et al. Genome-Wide Association Studies Identify IL23R-IL12RB2 and IL10 as Behçet’s Disease Susceptibility Loci. Nat Genet (2010) 42(8):703–6. doi: 10.1038/ng.624

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Mendell JT, Olson EN. MicroRNAs in Stress Signaling and Human Disease. Cell (2012) 148(6):1172–87. doi: 10.1016/j.cell.2012.02.005

PubMed Abstract | CrossRef Full Text | Google Scholar

7. Yenmis G, Soydas T, Ekmekci CG, Yazici Guvercin AC, Kucuk OS, Sultuybek GK. Fas and MicroRNAs Variations as a Possible Risk for Behçet Disease. J Clin Rheumatol (2020) 27(8):306–10. doi: 10.1097/rhu.0000000000001254

CrossRef Full Text | Google Scholar

8. Ibrahim W, Sakr BR, Obaya E, Ghonem H. MicroRNA-146a Expression and MicroRNA-146a Rs2910164 Polymorphism in Behcet’s Disease Patients. Clin Rheumatol (2019) 38(2):397–402. doi: 10.1007/s10067-018-4191-2

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Jadideslam G, Ansarin K, Sakhinia E, Babaloo Z, Abhari A, Alipour S, et al. Expression Levels of MiR-21, MiR-146b and MiR-326 as Potential Biomarkers in Behcet’s Disease. biomark Med (2019) 13(16):1339–48. doi: 10.2217/bmm-2019-0098

PubMed Abstract | CrossRef Full Text | Google Scholar

10. Ahmadi M, Yousefi M, Abbaspour-Aghdam S, Dolati S, Aghebati-Maleki L, Eghbal-Fard S, et al. Disturbed Th17/Treg Balance, Cytokines, and MiRNAs in Peripheral Blood of Patients With Behcet’s Disease. J Cell Physiol (2019) 234(4):3985–94. doi: 10.1002/jcp.27207

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Puccetti A, Pelosi A, Fiore PF, Patuzzo G, Lunardi C, Dolcino M. MicroRNA Expression Profiling in Behçet’s Disease. J Immunol Res (2018) 2018:2405150. doi: 10.1155/2018/2405150

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Yamazoe K, Meguro A, Takeuchi M, Shibuya E, Ohno S, Mizuki N. Comprehensive Analysis of the Association Between UBAC2 Polymorphisms and Behçet’s Disease in a Japanese Population. Sci Rep (2017) 7(1):742. doi: 10.1038/s41598-017-00877-3

PubMed Abstract | CrossRef Full Text | Google Scholar

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

CrossRef Full Text | Google Scholar

14. Shao X, Wang B, Mu K, Li L, Li Q, He W, et al. Key Gene Co-Expression Modules and Functional Pathways Involved in the Pathogenesis of Graves’ Disease. Mol Cell Endocrinol (2018) 474:252–9. doi: 10.1016/j.mce.2018.03.015

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Li X, LeBlanc J, Elashoff D, McHardy I, Tong M, Roth B, et al. Microgeographic Proteomic Networks of the Human Colonic Mucosa and Their Association With Inflammatory Bowel Disease. Cell Mol Gastroenterol Hepatol (2016) 2(5):567–83. doi: 10.1016/j.jcmgh.2016.05.003

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Ma C, Lv Q, Teng S, Yu Y, Niu K, Yi C, et al. Identifying Key Genes in Rheumatoid Arthritis by Weighted Gene Co-Expression Network Analysis. Int J Rheum Dis (2017) 20(8):971–9. doi: 10.1111/1756-185x.13063

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Yao Q, Song Z, Wang B, Qin Q, Zhang JA. Identifying Key Genes and Functionally Enriched Pathways in Sjögren’s Syndrome by Weighted Gene Co-Expression Network Analysis. Front Genet (2019) 10:1142. doi: 10.3389/fgene.2019.01142

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Zhou S, Lu H, Xiong M. Identifying Immune Cell Infiltration and Effective Diagnostic Biomarkers in Rheumatoid Arthritis by Bioinformatics Analysis. Front Immunol (2021) 12:726747. doi: 10.3389/fimmu.2021.726747

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Szklarczyk D, Franceschini A, Wyder S, Forslund K, Heller D, Huerta-Cepas J, et al. STRING V10: Protein-Protein Interaction Networks, Integrated Over the Tree of Life. Nucleic Acids Res (2015) 43(Database issue):D447–52. doi: 10.1093/nar/gku1003

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Leek JT, Scharpf RB, Bravo HC, Simcha D, Langmead B, Johnson WE, et al. Tackling the Widespread and Critical Impact of Batch Effects in High-Throughput Data. Nat Rev Genet (2010) 11(10):733–9. doi: 10.1038/nrg2825

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Diboun I, Wernisch L, Orengo CA, Koltzenburg M. Microarray Analysis After RNA Amplification can Detect Pronounced Differences in Gene Expression Using Limma. BMC Genomics (2006) 7:252. doi: 10.1186/1471-2164-7-252

PubMed Abstract | CrossRef Full Text | Google Scholar

22. Yu G, Wang LG, Han Y, He QY. Clusterprofiler: An R Package for Comparing Biological Themes Among Gene Clusters. Omics (2012) 16(5):284–7. doi: 10.1089/omi.2011.0118

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Szklarczyk D, Gable AL, Lyon D, Junge A, Wyder S, Huerta-Cepas J, et al. STRING V11: Protein-Protein Association Networks With Increased Coverage, Supporting Functional Discovery in Genome-Wide Experimental Datasets. Nucleic Acids Res (2019) 47(D1):D607–d13. doi: 10.1093/nar/gky1131

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Chin CH, Chen SH, Wu HH, Ho CW, Ko MT, Lin CY. Cytohubba: Identifying Hub Objects and Sub-Networks From Complex Interactome. BMC Syst Biol (2014) 8 Suppl 4(Suppl 4):S11. doi: 10.1186/1752-0509-8-s4-s11

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Bader GD, Hogue CW. An Automated Method for Finding Molecular Complexes in Large Protein Interaction Networks. BMC Bioinf (2003) 4:2. doi: 10.1186/1471-2105-4-2

CrossRef Full Text | Google Scholar

26. Oğuz AK, Yılmaz ST, Oygür Ç, Çandar T, Sayın I, Kılıçoğlu SS, et al. Behçet’s: A Disease or a Syndrome? Answer From an Expression Profiling Study. PloS One (2016) 11(2):e0149052. doi: 10.1371/journal.pone.0149052

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Xavier JM, Krug T, Davatchi F, Shahram F, Fonseca BV, Jesus G, et al. Gene Expression Profiling and Association Studies Implicate the Neuregulin Signaling Pathway in Behçet’s Disease Susceptibility. J Mol Med (Berl) (2013) 91(8):1013–23. doi: 10.1007/s00109-013-1022-4

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Wei T SV. Corrplot: Visualization of a Correlation Matrix. R Package Version 0.77 (2016). Available at: http://cran.r-project.org/package=corrplot.

Google Scholar

29. Liberzon A, Birger C, Thorvaldsdóttir H, Ghandi M, Mesirov JP, Tamayo P, et al. The Molecular Signatures Database (Msigdb) Hallmark Gene Set Collection. Cell Syst (2015) 1(6):417–25. doi: 10.1016/j.cels.2015.12.004

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Newman AM, Liu CL, Green MR, Gentles AJ, Feng W, Xu Y, et al. Robust Enumeration of Cell Subsets From Tissue Expression Profiles. Nat Methods (2015) 12(5):453–7. doi: 10.1038/nmeth.3337

PubMed Abstract | CrossRef Full Text | Google Scholar

31. Zhang Y, Aldridge J, Vasileiadis GK, Edebo H, Ekwall AH, Lundell AC, et al. Recombinant Adiponectin Induces the Production of Pro-Inflammatory Chemokines and Cytokines in Circulating Mononuclear Cells and Fibroblast-Like Synoviocytes From non-Inflamed Subjects. Front Immunol (2020) 11:569883. doi: 10.3389/fimmu.2020.569883

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Dhaiban S, Al-Ani M, Elemam NM, Maghazachi AA. Targeting Chemokines and Chemokine Receptors in Multiple Sclerosis and Experimental Autoimmune Encephalomyelitis. J Inflammation Res (2020) 13:619–33. doi: 10.2147/jir.S270872

CrossRef Full Text | Google Scholar

33. Verrou KM, Vlachogiannis NI, Ampatziadis-Michailidis G, Moulos P, Pavlopoulos GA, Hatzis P, et al. Distinct Transcriptional Profile of Blood Mononuclear Cells in Behcet’s Disease: Insights Into the Central Role of Neutrophil Chemotaxis. Rheumatol (Oxford) (2021) 60(10):4910–9. doi: 10.1093/rheumatology/keab052

CrossRef Full Text | Google Scholar

34. Rolin J, Sand KL, Knudsen E, Maghazachi AA. FTY720 and SEW2871 Reverse the Inhibitory Effect of S1P on Natural Killer Cell Mediated Lysis of K562 Tumor Cells and Dendritic Cells But Not on Cytokine Release. Cancer Immunol Immunother (2010) 59(4):575–86. doi: 10.1007/s00262-009-0775-7

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Dalbeth N, Gundle R, Davies RJ, Lee YC, McMichael AJ, Callan MF. CD56bright NK Cells are Enriched at Inflammatory Sites and can Engage With Monocytes in a Reciprocal Program of Activation. J Immunol (2004) 173(10):6418–26. doi: 10.4049/jimmunol.173.10.6418

PubMed Abstract | CrossRef Full Text | Google Scholar

36. Boutet MA, Bart G, Penhoat M, Amiaud J, Brulin B, Charrier C, et al. Distinct Expression of Interleukin (IL)-36α, β and γ, Their Antagonist IL-36Ra and IL-38 in Psoriasis, Rheumatoid Arthritis and Crohn’s Disease. Clin Exp Immunol (2016) 184(2):159–73. doi: 10.1111/cei.12761

PubMed Abstract | CrossRef Full Text | Google Scholar

37. Schecter AD, Calderon TM, Berman AB, McManus CM, Fallon JT, Rossikhina M, et al. Human Vascular Smooth Muscle Cells Possess Functional CCR5. J Biol Chem (2000) 275(8):5466–71. doi: 10.1074/jbc.275.8.5466

PubMed Abstract | CrossRef Full Text | Google Scholar

38. Wu B, Wang W, Zhan Y, Li F, Zou S, Sun L, et al. CXCL13, CCL4, and sTNFR as Circulating Inflammatory Cytokine Markers in Primary and SLE-Related Autoimmune Hemolytic Anemia. J Transl Med (2015) 13:112. doi: 10.1186/s12967-015-0474-4

PubMed Abstract | CrossRef Full Text | Google Scholar

39. Chang TT, Chen JW. Emerging Role of Chemokine CC Motif Ligand 4 Related Mechanisms in Diabetes Mellitus and Cardiovascular Disease: Friends or Foes? Cardiovasc Diabetol (2016) 15(1):117. doi: 10.1186/s12933-016-0439-9

PubMed Abstract | CrossRef Full Text | Google Scholar

40. Ogata H, Zhang X, Yamasaki R, Fujii T, Machida A, Morimoto N, et al. Intrathecal Cytokine Profile in Neuropathy With Anti-Neurofascin 155 Antibody. Ann Clin Transl Neurol (2019) 6(11):2304–16. doi: 10.1002/acn3.50931

PubMed Abstract | CrossRef Full Text | Google Scholar

41. Della-Torre E, Rigamonti E, Perugino C, Baghai-Sain S, Sun N, Kaneko N, et al. B Lymphocytes Directly Contribute to Tissue Fibrosis in Patients With Igg(4)-Related Disease. J Allergy Clin Immunol (2020) 145(3):968–81.e14. doi: 10.1016/j.jaci.2019.07.004

PubMed Abstract | CrossRef Full Text | Google Scholar

42. Tejera-Alhambra M, Casrouge A, de Andrés C, Seyfferth A, Ramos-Medina R, Alonso B, et al. Plasma Biomarkers Discriminate Clinical Forms of Multiple Sclerosis. PloS One (2015) 10(6):e0128952. doi: 10.1371/journal.pone.0128952

PubMed Abstract | CrossRef Full Text | Google Scholar

43. Zukowska Z, Grant DS, Lee EW. Neuropeptide Y: A Novel Mechanism for Ischemic Angiogenesis. Trends Cardiovasc Med (2003) 13(2):86–92. doi: 10.1016/s1050-1738(02)00232-3

PubMed Abstract | CrossRef Full Text | Google Scholar

44. Grassi MA, Tikhomirov A, Ramalingam S, Lee KE, Hosseini SM, Klein BE, et al. Replication Analysis for Severe Diabetic Retinopathy. Invest Ophthalmol Vis Sci (2012) 53(4):2377–81. doi: 10.1167/iovs.11-8068

PubMed Abstract | CrossRef Full Text | Google Scholar

45. Derkac I, Januleviciene I, Sepetiene R, Valiauga R, Velickiene D, et al. The Association of CEP135 Rs4865047 and NPY2R Rs1902491 Single Nucleotide Polymorphisms (SNPs) With Rapid Progression of Proliferative Diabetic Retinopathy in Patients With Type 1 Diabetes Mellitus. Med Sci Monit (2018) 24:8891–8. doi: 10.12659/msm.909803

PubMed Abstract | CrossRef Full Text | Google Scholar

46. Okumura M, Iwai M, Ide A, Mogi M, Ito M, Horiuchi M. Sex Difference in Vascular Injury and the Vasoprotective Effect of Valsartan are Related to Differential AT2 Receptor Expression. Hypertension (2005) 46(3):577–83. doi: 10.1161/01.Hyp.0000178564.14464.80

PubMed Abstract | CrossRef Full Text | Google Scholar

47. Shoaib RMS, Yahia S, Elsaid A, Abdel-Malak C, Hammad A. Angiotensin II Type 2 Receptor Gene Polymorphisms and Serum Angiotensin-Converting Enzyme Level in Egyptian Children With Systemic Lupus Erythematosus. Lupus (2019) 28(2):223–33. doi: 10.1177/0961203318820707

PubMed Abstract | CrossRef Full Text | Google Scholar

48. Terenzi R, Manetti M, Rosa I, Romano E, Galluccio F, Guiducci S, et al. Angiotensin II Type 2 Receptor (AT2R) as a Novel Modulator of Inflammation in Rheumatoid Arthritis Synovium. Sci Rep (2017) 7(1):13293. doi: 10.1038/s41598-017-13746-w

PubMed Abstract | CrossRef Full Text | Google Scholar

49. Perry ME, Chee MM, Ferrell WR, Lockhart JC, Sturrock RD. Angiotensin Receptor Blockers Reduce Erythrocyte Sedimentation Rate Levels in Patients With Rheumatoid Arthritis. Ann Rheum Dis (2008) 67(11):1646–7. doi: 10.1136/ard.2007.082917

PubMed Abstract | CrossRef Full Text | Google Scholar

50. Königshoff M, Wilhelm A, Jahn A, Sedding D, Amarie OV, Eul B, et al. The Angiotensin II Receptor 2 is Expressed and Mediates Angiotensin II Signaling in Lung Fibrosis. Am J Respir Cell Mol Biol (2007) 37(6):640–50. doi: 10.1165/rcmb.2006-0379TR

PubMed Abstract | CrossRef Full Text | Google Scholar

51. Schmied MC, Zehetmayer S, Reindl M, Ehling R, Bajer-Kornek B, Leutmezer F, et al. Replication Study of Multiple Sclerosis (MS) Susceptibility Alleles and Correlation of DNA-Variants With Disease Features in a Cohort of Austrian MS Patients. Neurogenetics (2012) 13(2):181–7. doi: 10.1007/s10048-012-0316-y

PubMed Abstract | CrossRef Full Text | Google Scholar

52. Cavanillas ML, Fernández O, Comabella M, Alcina A, Fedetz M, Izquierdo G, et al. Replication of Top Markers of a Genome-Wide Association Study in Multiple Sclerosis in Spain. Genes Immun (2011) 12(2):110–5. doi: 10.1038/gene.2010.52

PubMed Abstract | CrossRef Full Text | Google Scholar

53. Vandenhaute J, Wouters CH, Matthys P. Natural Killer Cells in Systemic Autoinflammatory Diseases: A Focus on Systemic Juvenile Idiopathic Arthritis and Macrophage Activation Syndrome. Front Immunol (2019) 10:3089. doi: 10.3389/fimmu.2019.03089

PubMed Abstract | CrossRef Full Text | Google Scholar

54. Onder M, Bozkurt M, Gürer MA, Gülekon A, Sezgin P, Imir T, et al. Natural Cellular Cytotoxicity in Behçet’s Disease. J Dermatol (1994) 21(4):239–43. doi: 10.1111/j.1346-8138.1994.tb01729.x

PubMed Abstract | CrossRef Full Text | Google Scholar

55. Hamzaoui K, Berraies A, Kaabachi W, Ammar J, Hamzaoui A. Pulmonary Manifestations in Behçet Disease: Impaired Natural Killer Cells Activity. Multidiscip Respir Med (2013) 8(1):29. doi: 10.1186/2049-6958-8-29

PubMed Abstract | CrossRef Full Text | Google Scholar

56. Cosan F, Aktas Cetin E, Akdeniz N, Emrence Z, Cefle A, Deniz G. Natural Killer Cell Subsets and Their Functional Activity in Behçet’s Disease. Immunol Invest (2017) 46(4):419–32. doi: 10.1080/08820139.2017.1288240

PubMed Abstract | CrossRef Full Text | Google Scholar

57. Huang L, Yu X, Li L, Liu J, Wu X, Zeng Y, et al. Aberrant Fcγriib and Fcγriii Expression on Monocytes From Patients With Behçet’s Disease. Clin Immunol (2020) 219:108549. doi: 10.1016/j.clim.2020.108549

PubMed Abstract | CrossRef Full Text | Google Scholar

58. Gogus F, Fresko I, Elbir Y, Eksioglu-Demiralp E, Direskeneli H. Oxidative Burst Response to Monosodium Urate Crystals in Patients With Behçet’s Syndrome. Clin Exp Rheumatol (2005) 23(4 Suppl 38):S81–5.

PubMed Abstract | Google Scholar

59. Tulunay A, Dozmorov MG, Ture-Ozdemir F, Yilmaz V, Eksioglu-Demiralp E, Alibaz-Oner F, et al. Activation of the JAK/STAT Pathway in Behcet’s Disease. Genes Immun (2015) 16(2):170–5. doi: 10.1038/gene.2014.64

PubMed Abstract | CrossRef Full Text | Google Scholar

60. van der Houwen TB, Dik WA, Goeijenbier M, Hayat M, Nagtzaam NMA, van Hagen M, et al. Leukocyte Toll-Like Receptor Expression in Pathergy Positive and Negative Behçet’s Disease Patients. Rheumatol (Oxford) (2020) 59(12):3971–9. doi: 10.1093/rheumatology/keaa251

CrossRef Full Text | Google Scholar

61. Slobodin G, Toukan Y, Rosner I, Rozenbaum M, Boulman N, Pavlotzky E, et al. LPS-Stimulated Production of TNF-Alpha by Peripheral Blood Monocytes in Patients With Behcet’s Disease. Clin Rheumatol (2007) 26(5):764–7. doi: 10.1007/s10067-006-0371-6

PubMed Abstract | CrossRef Full Text | Google Scholar

62. Nakano H, Kirino Y, Takeno M, Higashitani K, Nagai H, Yoshimi R, et al. GWAS-Identified CCR1 and IL10 Loci Contribute to M1 Macrophage-Predominant Inflammation in Behçet’s Disease. Arthritis Res Ther (2018) 20(1):124. doi: 10.1186/s13075-018-1613-0

PubMed Abstract | CrossRef Full Text | Google Scholar

63. Yamashita N. Hyperreactivity of Neutrophils and Abnormal T Cell Homeostasis: A New Insight for Pathogenesis of Behçet’s Disease. Int Rev Immunol (1997) 14(1):11–9. doi: 10.3109/08830189709116841

PubMed Abstract | CrossRef Full Text | Google Scholar

64. Lichtig C, Haim S, Hammel I, Friedman-Birnbaum R. The Quantification and Significance of Mast Cells in Lesions of Behçet’s Disease. Br J Dermatol (1980) 102(3):255–9. doi: 10.1111/j.1365-2133.1980.tb08137.x

PubMed Abstract | CrossRef Full Text | Google Scholar

65. Vural S, Kerl K, Doğan PE, Vollmer S, Puchta U, He M, et al. Lesional Activation of T(C) 17 Cells in Behçet’s Disease and Psoriasis Supports HLA-Class I-Mediated Autoimmune Responses. Br J Dermatol (2021). doi: 10.1111/bjd.20643

PubMed Abstract | CrossRef Full Text | Google Scholar

66. Houman H, Hamzaoui A, Ben Ghorbal I, Khanfir MS, Feki M, Hamzaoui K. Tc1/Tc2 Ratio in the Inflammatory Process in Patients With Behçet’s Disease. Mediators Inflammation (2004) 13(4):247–53. doi: 10.1080/09629350400003167

CrossRef Full Text | Google Scholar

67. Hamzaoui K, Ayed K, Slim A, Hamza M, Touraine J. Natural Killer Cell Activity, Interferon-Gamma and Antibodies to Herpes Viruses in Patients With Behçet’s Disease. Clin Exp Immunol (1990) 79(1):28–34. doi: 10.1111/j.1365-2249.1990.tb05122.x

PubMed Abstract | CrossRef Full Text | Google Scholar

68. Rabinovich GA, Baum LG, Tinari N, Paganelli R, Natoli C, Liu FT, et al. Galectins and Their Ligands: Amplifiers, Silencers or Tuners of the Inflammatory Response? Trends Immunol (2002) 23(6):313–20. doi: 10.1016/s1471-4906(02)02232-9

PubMed Abstract | CrossRef Full Text | Google Scholar

69. Vojdani A, Afar D, Vojdani E. Reaction of Lectin-Specific Antibody With Human Tissue: Possible Contributions to Autoimmunity. J Immunol Res (2020) 2020:1438957. doi: 10.1155/2020/1438957

PubMed Abstract | CrossRef Full Text | Google Scholar

70. Gupta B, Raghav SK, Das HR. S-Nitrosylation of Mannose Binding Lectin Regulates its Functional Activities and the Formation of Autoantibody in Rheumatoid Arthritis. Nitric Oxide (2008) 18(4):266–73. doi: 10.1016/j.niox.2008.01.009

PubMed Abstract | CrossRef Full Text | Google Scholar

Keywords: Behcet’s disease, gene expression omnibus, biomarkers, diagnosis, immune cell infiltration

Citation: Zhan H, Li H, Cheng L, Yan S, Zheng W and Li Y (2021) Novel Insights Into Gene Signatures and Their Correlation With Immune Infiltration of Peripheral Blood Mononuclear Cells in Behcet’s Disease. Front. Immunol. 12:794800. doi: 10.3389/fimmu.2021.794800

Received: 14 October 2021; Accepted: 23 November 2021;
Published: 15 December 2021.

Edited by:

Michele Maria Luchetti Gentiloni, Marche Polytechnic University, Italy

Reviewed by:

Kamel Hamzaoui, Tunis El Manar University, Tunisia
Esın Cetın Aktas, Istanbul University, Turkey
Lesley Ann Bergmeier, Queen Mary University of London, United Kingdom

Copyright © 2021 Zhan, Li, Cheng, Yan, Zheng and Li. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Yongzhe Li, eW9uZ3poZWxpcHVtY2hAMTI2LmNvbQ==

These authors have contributed equally to this work and share first authorship

Disclaimer: 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.

Research integrity at Frontiers

Man ultramarathon runner in the mountains he trains at sunset

94% of researchers rate our articles as excellent or good

Learn more about the work of our research integrity team to safeguard the quality of each article we publish.


Find out more