Rafael Alfaro1
Helios Martínez-Banaclocha1
Santiago Llorente2
Victor Jimenez-Coll1José Antonio Galián1
Carmen Botella1María Rosa Moya-Quiles1
Antonio Parrado1Manuel Muro-Perez1
Alfredo Minguela1
Isabel Legaz3*†
Manuel Muro1*†- 1Immunology Services, University Clinical Hospital Virgen de la Arrixaca-Biomedical Research Institute of Murcia (IMIB), Murcia, Spain
- 2Nephrology Services, University Clinical Hospital Virgen de la Arrixaca-Biomedical Research Institute of Murcia (IMIB), Murcia, Spain
- 3Department of Legal and Forensic Medicine, Biomedical Research Institute (IMIB), University of Murcia, Murcia, Spain
Background: The diagnosis of graft rejection in kidney transplantation (KT) patients is made by evaluating the histological characteristics of biopsy samples. The evolution of omics sciences and bioinformatics techniques has contributed to the advancement in searching and predicting biomarkers, pathways, and new target drugs that allow a more precise and less invasive diagnosis. The aim was to search for differentially expressed genes (DEGs) in patients with/without antibody-mediated rejection (AMR) and find essential cells involved in AMR, new target drugs, protein-protein interactions (PPI), and know their functional and biological analysis.
Material and Methods: Four GEO databases of kidney biopsies of kidney transplantation with/without AMR were analyzed. The infiltrating leukocyte populations in the graft, new target drugs, protein-protein interactions (PPI), functional and biological analysis were studied by different bioinformatics tools.
Results: Our results show DEGs and the infiltrating leukocyte populations in the graft. There is an increase in the expression of genes related to different stages of the activation of the immune system, antigenic presentation such as antibody-mediated cytotoxicity, or leukocyte migration during AMR. The importance of the IRF/STAT1 pathways of response to IFN in controlling the expression of genes related to humoral rejection. The genes of this biological pathway were postulated as potential therapeutic targets and biomarkers of AMR. These biological processes correlated showed the infiltration of NK cells and monocytes towards the allograft. Besides the increase in dendritic cell maturation, it plays a central role in mediating the damage suffered by the graft during AMR. Computational approaches to the search for new therapeutic uses of approved target drugs also showed that imatinib might theoretically be helpful in KT for the prevention and/or treatment of AMR.
Conclusion: Our results suggest the importance of the IRF/STAT1 pathways in humoral kidney rejection. NK cells and monocytes in graft damage have an essential role during rejection, and imatinib improves KT outcomes. Our results will have to be validated for the potential use of overexpressed genes as rejection biomarkers that can be used as diagnostic and prognostic markers and as therapeutic targets to avoid graft rejection in patients undergoing kidney transplantation.
Introduction
The severity and occurrence of rejection in kidney transplant (KT) patients depend on numerous variables that can affect the magnitude and nature of immune responses. Understanding how genetic and molecular factors affect the effector functions of immune cells and donor-specific antibodies (DSA) can better renal stratification receptors based on their immunological risk and thus help the clinician make better decisions to anticipate adverse events (1–5).
In the last two decades, high-throughput technologies such as next-generation sequencing (NGS) and microarrays have been developed. In parallel, international public data repositories have been developed, such as the GEO (Gene Expression Omnibus) database of the National Center for Biotechnology Information (NCBI) (6), the Array Express database of the Institute European Bioinformatics (EBI) (7), in order to store and distribute these data. The large amount of data generated by these technologies makes it necessary to use robust bioinformatics tools that allow us to know in silico how molecules interact and regulate the different biological processes that mediate the biological processes of health and disease.
Therefore, functional analysis tools have been developed to explore and identify critical biological processes. These tools are based on biological knowledge databases such as Gene Ontology (GO) or the Kyoto Encyclopedia of Genes and Genomes (KEGG). GO provides a controlled vocabulary of terms (ontologies) to describe gene products in terms of biological processes, cellular components, and associated molecular functions (8). On the other hand, the KEGG is a reference database for the biological interpretation of metabolism and cellular processes (9). Furthermore, the development of genomic and protein databases has led to many bioinformatics tools to predict the properties of proteins and the genome: splice sites, post-translational modifications, stability, pathogenicity. The integration of these new tools and data offers unprecedented opportunities that promise to accelerate and improve our understanding of biology and the search for biomarkers, and specifically for the field of organ transplantation.
The aim was to search for differentially expressed genes (DEGs) in AMR patients and to analyze the types of populations of leukocytes infiltrated in the graft developing diagnostic models using computational predictions in order to find diagnostic and prognostic markers and as a therapeutic target that allows avoiding graft rejection in patients undergoing kidney transplantation.
Materials and Methods
Gene Expression Data of Gene Expression Omnibus Repository
For this study, the raw gene expression data studies were downloaded from the Gene Expression Omnibus (GEO) database (6). The characteristics of the GSE cohorts are summarized in Table 1. Immunosuppressive treatment in the different GEO studies was similar. Around 40-60% of the transplant recipients in the different studies had triple therapy of corticosteroids + tacrolimus + MMF, while the rest had other combinations in different proportions: corticosteroids + MMF + cyclosporine, tacro + corticosteroids. Regarding induction therapy (thymoglobulin and baxilysimab) there is no information available in the studies. All the bases analyzed in this study have been reviewed and are comparable.
Table 1 Characteristics of the gene expression studies of kidney transplantation outcome obtained from the GEO database.
Venn Diagrams
Venn diagrams (10) were used to determine the intersection and analysis of the four analyzed cohorts’ differentially expressed genes (DEGs).
Functional Analyzes of the KEGG and GO Terms
DEGs in this study were entered into the WebGestalt web tool (6) for the analyzes of Gene Ontology terms (GO) and biological pathways from the English Kyoto Encyclopedia of Genes and Genomes (KEGG) (11). The minimum number of genes to be considered a pathway was adjusted to a value of two.
The Benjamini and Hockberg FDR method adjusted the p-value obtained from each biological pathway with FDR<0.05 considered significant. Only those from the molecular function and biological processes were taken into account for the analyses of GO terms. GO terms with an adjusted p value<0.05 were considered significant (12). Furthermore, the enrichment ratio (ER) is indicated for each GO term and KEGG pathway in this analysis. The ER is defined as the number of genes observed among the number of genes expected for each GO term or KEGG pathway (12).
Distribution of Leukocytes in Renal Graft Biopsy
The xCell algorithm (13) was used to estimate the relative distribution of leukocytes in graft biopsy, which allows the inference of cell subpopulations from transcriptomic data. This algorithm allows the relative proportion of 64 cell types to be inferred, including lymphoid, myeloid, hematopoietic progenitors, stromal cells, and other cell types, such as epithelial cells or melanocytes. In our study, we focus exclusively on lymphoid and myeloid cells.
The gene expression data of GSE36059, GSE44131, GSE50084, GSE93658 was analyzed using the xCell tool (Figure 1), generating a score for each cell population proportional to the relative frequency of that population. The score obtained from each subpopulation and sample will be used later for comparison by traditional statistical methods between transplants with and without AMR.
Figure 1 Design of the bioinformatics study applied to antibody-mediated rejection. NAR, No Acute Rejection; AMR, Antibody-Mediated Rejection; GO, Gene Ontology, KEGG, Kyoto Encyclopedia of Genes and Genomes; PPI, Protein-Protein Interaction, GEO; Gene Expression Omnibus, DEG, Differentially Expressed Genes.
Protein-Protein Interaction Networks
The interaction between the proteins produced by the DEGs obtained in the previous sections was used to generate a protein-protein interaction map (PPI) (14) using the Network Analyst web tool (15). This tool uses the Innate DB database of protein interactions.
The nodes represent the proteins in PPI networks, while the junction lines represent the known interactions between the bound proteins. The constructed network was limited to include only the core proteins (zero-order interactions) to avoid the hairball effect, and better visualize the network. The network’s topological properties were evaluated using the betweenness centrality (BC) and the Degree Centrality (GCy). The GI is the number of connections a node has with other nodes, while BC measures the shortest paths that traverse a node. It indicates the degree to which a node is located between other nodes in the network. Nodes with a high BC and DCy index are essential proteins in cell signaling and control information flow.
The modules of the PPI network are tightly connected areas of the network that are considered relatively independent as they often work together to perform a biological function. Therefore, the modules of the PPI networks usually correspond to metabolic pathways or protein complexes. The module explorer of the Network Analyst tool, based on the random-walk algorithm, was used to identify the modules.
Connections within a module and the rest of the PPI network were analyzed. Values of p<0.05 were considered significant. The significance of each module was calculated using the Wilcoxon test. The significant modules were subsequently functionally analyzed using the KEGG and GO term databases. Those GO terms and KEGG pathways with a value of p<0.05 were considered significant.
Analyzes of Protein Kinases and Transcription Factors
In the search for biomarkers and new therapeutic targets in transplantation, it is crucial to study those genes with altered expression and the entire signaling pathway that regulates them. For this purpose, we employ the Expression2Kinases (X2K) tool that infers the protein kinases and transcription factors that form the regulatory networks upstream of the DEGs obtained previously.
miRNA, B cells were inferred by uploading the list of these genes to the Expression2Kinase (X2K) web tool (15). X2K identifies the transcription factors that regulate DEGs using the ChEA database. Protein kinases associated with transcription factors were obtained through the Kinase Enrichment Analysis module of X2K. Transcription factors and protein kinases with p values <0.05 were considered significant.
Computational Prediction of New Target Drugs in Overexpressed Genes During AMR Kidney Transplantation
The computational prediction of new target drugs overexpressed genes during AMR kidney transplantation used the Gene2Drug web tool (16, 17).
The Gene2drug algorithm searches for drugs or low molecular weight compounds that induce a significant variation in the expression levels of a list of genes of interest. Gene2drug uses gene expression data from the Connectivity Map database (CMap), containing various cell lines treated with 1309 drugs or low molecular weight molecules. CMap identifies drug-associated transcriptomics alterations. CMap show which genes change their expression induced by these compounds (18). The genes introduced in the Gene2drug algorithm had a value of BC≥100 in the topological analysis of the PPI network, presenting a high connection (Table S1). Gene2Drug calculates the English Enrichment Score (ES value) and the p-value corrected by the Benjamini-Hochberg method for the 1309 compounds and drugs of CMap. The ES value ranges from -1 to +1 and represents the degree to which a compound or drug regulates the genes. The ES value uses a generalization of the Kolgomorov-Smirnov method. Compounds with ES values close to +1 will be those that the algorithm predicts will significantly induce the genes of interest, while compounds with values at -1 will be those that will inhibit the genes of interest. Our objective was to search drugs that counteract to focus on those compounds with ES ≤0 values.
Protein-Protein (STRING) and Drug-Protein (STITCH) Interactions Databases
STRING is a database of protein-protein interactions (16), which includes direct (physical) and indirect (functional) associations obtained from experimental data and computational predictions. STRING version 11.0 was used to obtain the interactions between the PDGFRA and PDGFRB receptors and the proteins with the highest BC and GI values from the PPI network previously analyses (GZMB, STAT1, LYN, IRF1, and IRF8).
STITCH is also a database of computational drug-protein interactions (19). Like STRING, it includes direct (physical) and indirect (functional) associations obtained from experimental data and computational predictions (20). STITCH version 5.0 was used to obtain the drug interactions imatinib.
Statistical Analysis
The graphs and statistical analyzes were carried out in the Statistical Package for the Social Sciences (SPSS, version 22, Chicago, IL) and GraphPad Prism (version 6, San Diego, CA) software as well as in the R programming language, used for the latter the environment Integrated Development R Studio version 3.4.
The results have been expressed as the mean ± standard deviation (SD) for quantitative data or percentages for categorical data. For the comparison of categorical variables, the χ2 test or Fisher’s exact test was used. The verification of the normality of the data was carried out using the Kolmogorov-Smirnov test. The Mann-Whitney U test was used to compare two groups with variables that did not adjust to normality. The Kruskal Wallis test and Dunn’s post hoc test with Bonferroni correction for multiple comparisons compared three or more groups. Correlation analyzes were carried out using the Spearman index, as previously published (5, 21, 22). For the longitudinal comparison of two related groups, the Wilcoxon non-parametric test for related samples was used. The Friedman test with Wilcoxon post hoc was used to compare three or more related groups. The evaluation of the sensitivity and specificity of a biomarker was carried out by constructing ROC curves. The discriminatory capacity was evaluated using the AUC (Figure 5). To obtain the optimal cut-off value that maximizes sensitivity and specificity, we use the Youden index. To correct the p-value in multiple comparisons, the Benjamini-Hochberg or Bonferroni method was used. Values of p<0.05, or p-corrected<0.05 in the case of multiple comparisons, were considered significant for all statistical tests.
Results
Description of the Databases Used in the Different Bioinformatic Analyzes
The different databases used in the different bioinformatic analyzes designed to search for biomarkers in KT is shown in Figure 1.
Firstly, a search was conducted in the GEO database for studies performed with kidney graft biopsy samples from transplant recipients with AMR or without acute rejection (NR) (Table 1). Next, the gene expression data in the GEO2R application to obtain the differentially expressed genes (DEGs) in the AMR samples were analyzed. The DEGs list was used in various bioinformatics analyses to obtain essential characteristics at the molecular and cellular level of the AMR, create diagnostic models, and search for new therapeutic options.
Differentially Expressed Genes in Biopsies With AMR
The gene expression data of the selected studies (GSE36059, GSE44131, GSE50084, and GSE93658) were analyzed using the GEO2R web tool establishing two groups for each cohort, one with biopsies from transplants NR and AMR.
The GEO2R tool obtained the Fold-change and adjusted p values using the Benjamini and Hochber FDR method, establishing DEGs between the NR and AMR groups with an FDR value <0.05. The number of DEGs obtained for each study is shown in Table 1. In the study, a total of 1778, 2901, 1387, 2855 DEGs were obtained in GSE36059, GSE44131, GSE50084, and GSE93658, respectively.
We selected those DEGs present in these four studies for further analyses and ruled out genes of little relevance to rejection. In all studies, the deregulation of these genes indicates that they play an essential role in the rejection process and are therefore potential biomarkers. Venn diagram analysis (Figure 2) shows 340 overexpressed and 26 under-expressed DEGs in the AMR group. The complete list of DEGs obtained is shown in Table 2.
Figure 2 Venn diagrams with the number of differentially expressed genes in four different cohorts analyses. (A) The number of genes overexpressed in the antibody-mediated rejection (AMR) group. (B) The number of under-expressed genes in the AMR group in each study and the overlap between studies.
Table 2 Total common DEGs in the AMR group in the four cohorts analyzed.
Biological and Functional Pathway Analyzes of Differentially Expressed Genes
The selected DEGs were subjected to functional analyzes of GO terms and biological pathways of the KEGG. The main biological pathways overrepresented in DEGs are shown in Table S2A, where two biological pathways related to AMR stand out, graft-versus-host disease (GVHD; hsa05332, FDR<0.0001) and allograft rejection (hsa05330, FDR<0.0001).
Regarding the under-expressed DEGs (Table S2B), only two routes were represented, the fatty acid degradation pathway (hsa00071, FDR=0.0022) and the isoleucine, leucine, and valine degradation pathway (hsa00280, FDR=0, 0022).
The functional analyzes of GO terms for biological processes (BP) and molecular functions (MF) are summarized in Table 3.
Table 3 The functional analyzes of GO terms for biological processes and molecular functions in DGEs.
For over-expressed DEGs, the most overrepresented BPs are those related to the immune response activation (GO:0002253, FDR<0.0001) and immune effector processes (GO:0002252, FDR<0.0001). For molecular functions, the most critical terms were pattern recognition receptor activity (GO:0008329, FDR<0.0001) and MHC binding (GO:0042287, FDR<0.0001; Table 3A)
GO term with an FDR <0.05 was obtained for under-represented DEGs due to the low number of genes included in the analysis (Table 3B). According to the unadjusted p values, the main BPs were lipid metabolism (GO:0006629, p=0.0001) and fatty acid oxidation (GO:0019395, p=0.0001). Regarding the molecular functions, they were the oxidoreductase activity of CH-OH groups (GO:0016614, p=0.0002) and the activity of GA-3P DHase (GO: 0043878, p=0.005).
Identification of Interaction Networks Between Proteins
The PPI network obtained was divided into three subnets, of which we analyzed in greater detail the largest one composed of 98 nodes and 154 interactions (Figure 3). The degree of interconnection shows that 25 (25.5%) of the nodes were equal to one, while 73 (74.5%) show a degree ≥2. Of those 73, there are four nodes (STAT1, IRF1, LYN, and IRF8) that show a high GI (≥10): STAT1 (GI=21), IRF1 (GI=19), LYN (GI=10), and IRF8 (GI=10). Regarding the BC, we observe a range from 15 to 2229.9 for a total of 62 (63.2%) of the nodes.
Figure 3 Protein-protein interaction (PPI) network. The size of the nodes represents the degree of interconnection. The color grading represents the centrality of intermediation (BC). Darker colors reflect a high degree of BC. The numbered boxes indicate the functional modules obtained from the analysis of the PPI network.
The nodes with the highest levels of BC were STAT1 (BC=2229.9), IRF1 (BC=997.3), LYN (BC=796.9), and GZMB (BC=781.8). The genes encoding the above proteins were among those over-expressed in the AMR group. Of the proteins encoded by genes under-expressed in AMR, only one was represented in the PPI network as a central gene, ATP5B (GI=2, BC=153.7).
The results showed that proteins with a high degree of GI also show high values of BC, which suggests that, possibly, they are the essential proteins in regulating intracellular signaling. Therefore, they can be essential in the search for therapeutic targets or biomarkers.
The modules of the PPI network are densely connected areas of the network, which are considered relatively independent since they often work together to perform a biological function. Therefore, the modules of the PPI networks usually correspond to metabolic pathways or protein complexes. In our network, we obtained a total of seven modules (Figure 3), although module 6 was not significant (p=0.233) and therefore was not considered for the functional analysis. The proteins of each module are shown in Table 4. The graphic representation of the protein connections of each of the modules is shown in Figure 4.
Table 4 Proteins of different PPI network modules.
Figure 4 Detail of the modules obtained from the PPI network. The figures represent module 0 (A), module 1 (B), module 2 (C), module 3 (D), module 4 (E), module 5 (F), and module 6 (G). Nodes represent the degree of interconnection. The color grading represents the centrality of intermediation (BC). Darker colors reflect a high degree of BC.
Next, a functional analysis was performed with terms from the KEGG database to obtain the main BPs in each module. The results of the functional analysis of the PPI network modules are shown in Table 5. The analysis shows the strong relationship of the modules with effector functions of the immune system such as antigen presentation (Module 0, p=0.0002), NK cell-mediated cytotoxicity (Module 1, p=0.0048), phagocytosis mediated by Fc-γ receptors (Module 3, p=0.0017), and trans-endothelial migration of leukocytes (Module 4, p<0.0001).
Table 5 Biological pathways of the KEGG of the different modules.
Identification of DEG-Associated Transcription Factors and Protein Kinases
Regarding protein kinases, in Table S3, we observed that the most important in the over-expressed genes are MAPK3, CSNK2A1, MAPK14, ERK1, and MAPK1. On the other hand, for the under-expressed genes, the essential protein kinases were DNAPK, CDK1, MAPK14, ERK1, and MAPK1. The MAPK14, MAPK1, and ERK1 protein kinases appear important for over-expressed and under-expressed gene signaling pathways, revealing commonalities upstream of the biological pathways activated during rejection.
For the transcription factors, the most important ones associated with over-expressed genes in the AMR group are SPI, IRF8, RUNX1, RELA, and IRF1, while for under-expressed genes, they were MYC, KLF4, USF1, SOX2, and TCF3 (Table S4), showing a clear divergence in pathways between over- and under-expressed genes. The transcription factors SPI and RUNX1 are the ones that target the highest number of DEGs, with 17.9% (61/340) and 17.3% (59/340), respectively, of the total over-expressed DEGs. Regarding the under-expressed DEGs in the AMR group, the transcription factor USF1 targets 23% (6/26) and KLF4 19.2% (5/26) of the total.
The Cellular Infiltrate in the Graft With Antibody-Mediated Rejection
To evaluate the abundance of immune cells infiltrated in the graft during an AMR, we used the xCell tool (13). This tool is capable of estimating the relative abundance of 64 cell types from transcriptomic data. In our study, we focused on analyzing the 33 immune-type subpopulations that xCell can infer (Figures 5, 6). The transcriptomic data of the cohort GSE36059, composed of 281 biopsies of transplants NR and 65 of transplants with AMR, were entered into xCell. The results obtained of GSE36059 and the rest of the cohorts analyzed were similar are summarized in Table 6.
Figure 5 Cytolytic activity index in AMR. (A) The cytolytic index (CYT) was compared between the NR group without rejection (NR) and AMR. Comparisons were made using the Mann-Whitney U test. Values expressed as mean ± SEM. (B) ROC curve for the diagnosis of AMR using the CYT. As a reference, the non-discrimination diagonal (black line) is represented. The area under the curve (AUC) and statistical significance are indicated on the graph. Values of p<0.05 were considered significant. ***p<0.001. AMR, Antibody-Mediated Rejection. NR, No rejection; AMR, Antibody-Mediated Rejection.
Figure 6 Correlation of the relative abundance of immune subpopulations with the cytolytic index (CYT). The Spearman correlation coefficient is indicated in each cell. The red color indicates a positive correlation, and the blue color indicates a negative correlation.
Table 6 Analysis of leukocyte population results of the different GSEs using xCell tool.
These genes were entered in the Gene2Drug web tool to infer those molecules and drugs with the most remarkable capacity to decrease the expression levels of genes over-expressed in AMR (Figure 7). Gene2Drug searches the CMap database, made up of more than 7,000 drug-induced gene expression patterns and low molecular weight molecules. Gene2Drug assigns each drug an ES (Enrichment Score) score ranging from -1 to +1. Those drugs that induce alterations in gene expression similar to those observed in rejection will have ES values close to 1, while producing an alteration in the opposite direction will have values close to -1. For our purpose, we will choose those drugs with ES values closest to -1 since they will be those that best counteract the expression levels of the altered genes in the AMR. The compounds with statistical significance are shown in Table 7, and only imatinib is currently approved by the Food and Drug Administration (FDA) and the European Medicines Agency (EMA); that is why we decided to deepen the interactions of this drug.
Figure 7 New drug search strategy. During AMR events, an alteration of gene expression profiles occurs. In the search for new drugs, we look for those that induce a transcriptional profile opposite to that observed in AMR to counteract the altered genes and restore normal gene expression levels. AMR, Antibody-Mediated Rejection.
Table 7 Drugs obtained from complete Gene2Drug analysis.
Imatinib is an inhibitor of kinases such as the PDGFRA and PDGFRB receptors, currently used to treat some cancers such as chronic myeloid leukemia or gastrointestinal tumors. As we see in Table 7, imatinib obtains better ES scores than other drugs commonly used in KT patients, such as Tacrolimus (ES= -0.227), Mycophenolic Acid (ES= -0.180), and Methylprednisolone (ES=0.173). Figure 8 shows that imatinib acts directly on CSF1R, a receptor for colony-stimulating factor 1, over-expressed in AMR. Furthermore, indirectly, it can alter the expression of STAT1 and LYN through blocking PDGFRA and PDGFRB, which, as we have seen previously, are central proteins that regulate the expression of a large number of genes involved in AMR. Therefore, imatinib may be a promising drug for the prevention or treatment of AMR.
Figure 8 Interactions of imatinib and PDGFRA and PDGFRB proteins. Imatinib interactions were obtained from the STITCH database (Top). PDGFRA and PDGFRB interactions with important genes in the AMR were obtained from the STRING database (below). AMR, Antibody-Mediated Rejection.
Discussion
The diagnosis of graft rejection in KT patients is made by evaluating the histological characteristics of biopsy samples. The evolution of omics sciences and bioinformatic techniques has contributed to the development of new diagnostic approaches by combining histopathological characteristics with the gene expression profile of the kidney graft, thereby achieving a more precise diagnosis. It is still necessary to improve our knowledge of the molecular mechanisms underlying kidney rejection to achieve non-invasive diagnostic methods or new therapeutic objectives that improve the clinical management of transplanted patients. For this reason, the central genes and biological pathways that regulate the AMR process were studied to understand better the molecular mechanisms that control it and propose possible biomarkers and new therapeutic options.
Transcriptomics studies have proven to be a valuable tool for understanding the molecular processes that govern graft rejection. One of the main limitations of these studies is, on the one hand, the small number of individuals per study due to the high cost of these techniques. On the other, the variability between different analyzes platforms, challenging to interpret and compare results between different research groups.
Although in some biomedical fields, such as oncology, transcriptomics studies are quite widespread, in transplantation, there are still few investigations with this methodology and, for this reason, a meta-analysis on gene expression is scarce in the literature. In transplantation, the meta-analyses of gene expression carried out have focused mainly on the search for expression profiles associated with graft rejection and tolerance (3, 23), in addition to the search for therapeutic objectives, such as for kidney graft fibrosis (24).
In the context of KT, few studies have focused on studying the specific expression patterns of AMR. Kim et al. (25) performed a meta-analysis where they obtained that CXCL10, CXCL9, and GBP1 were the most over-expressed genes with AMR. Unlike the study carried out by Kim et al. (25), our work has focused exclusively on samples derived from biopsies, excluding those obtained from blood and PBMCs, intending to increase the homogeneity of the study, since the expression profiles in blood and graft may differ, leading to discrepancies between studies, making them difficult to compare.
Our study consisted of four cohorts with 466 renal biopsy samples, 137 from transplants with AMR and 329 from transplants with good graft function. The results showed 340 over-expressed genes in the AMR group, while only 26 were under-expressed. DEGs were functionally classified based on the terms GO for annotation of MFs and BPs, as well as based on the biological pathways of the KEGG.
Genes over-expressed in the AMR group were associated with BPs such as immune response activation, immune effector processes, leukocyte-mediated immunity, and molecular functions such as MHC binding or CXCR3 receptor binding. KEGG metabolic pathway analyzes showed biological pathways directly related to transplantation, such as GVHD and allograft rejection. Increased expression of genes related to MHC binding such as HLA, CD74, or TAP1 genes reflects the antigen-presenting solid activity in the allograft. It is known that in rejection processes, there is an increase in the maturation of antigen-presenting cells (APCs) such as DCs (23), which lead to an increase in the expression of genes related to antigen processing pathways such as CD74 or TAP1 (3, 26).
On the other hand, we have obtained an increase in genes related to the CXCR3 receptor. CXCR3 is a chemokine receptor that polarizes CD4 T lymphocytes towards Th1, promoting infiltration of cytotoxic lymphocytes into the graft. Blocking this receptor is effective in inhibiting macrophage infiltration in acute rejection (27). However, some studies indicate that the infiltration of cells that express CXCR3 is only increased in cell-type rejections and not in humoral ones (28). In our study, the over-expression of the CXCR3 gene was only obtained in two of the cohorts (GSE44131 and GSE93658), which indicates that it is not a good biomarker for AMR.
In contrast, the CXCL9 and CXCL10 genes, which code for CXCR3 ligands, were over-expressed in all cohorts in our analyzes, thus being better indicators of rejection. A possible explanation for this discrepancy may lie in the different abundance of receptor and ligand-producing cells into the graft. While the expression of the CXCR3 receptor is mainly limited to cells of the immune system such as lymphocytes and NK cells (29), the CXCL9 and CXCL10 ligands can be produced in large quantities by epithelial cells of the renal tubule in response to certain stimuli such as IFN-γ (30). Previous studies have shown that the CXCL9 and CXCL10 genes’ over-expression was associated with rejection processes, regardless of the tissue’s type of rejection and origin (3, 31, 32). In heart transplant models, CXCL10 inhibition slows rejection and increases graft survival (33) The data from our analyzes confirm those obtained by other researchers, showing that the ligands that regulate the activation of the CXCR3 receptor are a potential biomarker and therapeutic target. In addition, among the over-expressed genes, we also obtained other biological pathways unrelated to transplantation, such as infection by Leishmaniasis and Staphylococcus aureus. These biological pathways are present because the genes involved in the immune response to the allograft are common to those involved in response to infections, such as the HLA or receptor genes of the immunoglobulin Fc fragment (30, 34).
Among the under-expressed genes in the AMR group, the BPs related to the oxidation of fatty acids stand out: catabolic processes of fatty acids, monocarboxylic acids, and lipid oxidation. These results are compatible with those obtained in the analyzes of biological routes of the KEGG where two significantly decreased routes were obtained: degradation of fatty acids and degradation of isoleucine, leucine, and valine, although the differences were not significant after correction for multiple comparisons, probably due to the low number of under-expressed genes obtained in the analyzes. During the activation of the immune system, the proliferation and differentiation processes of lymphocytes lead to metabolic reprogramming. After the lymphocyte activation, there is an inhibition of β-oxidation and an increase in lipid biosynthesis, necessary for cell division (35, 36). These observations are compatible with the results obtained in our study. During allograft rejection, a robust immune response is triggered with increased leukocyte proliferation. This can only be sustained if the metabolic pathways involved in biosynthetic processes are favored, thereby inhibiting the catabolic pathways of essential elements such as lipids and amino acids. Inhibition of lipid synthesis has been shown to interfere with lymphocyte proliferation and slow rejection (37). However, some studies indicate that the microenvironment conditions determine the metabolic status of the effector lymphocytes. In murine models of GVHD, it has been observed that fatty acid oxidation is necessary for effector cell function (38), although most research indicates that carbohydrates remain the primary fuel in immune effector processes (39, 40). This fact becomes evident after immune activation, where the activation of mTORC1 promotes the increase of transporters such as Glut1 to increase the entry of glucose into the cell, necessary to sustain the increase in glycolytic activity. In addition, there is also an increase in other transporters, such as leucine, during lymphocyte activation. Our data suggest that the inhibition of amino acid catabolism, together with the increase in amino acid transporters after immune activation observed in other studies, is a cellular response to maintain the intense biosynthetic activity that occurs during rejection (40, 41). Leucine antagonists or their transporters exert effects similar to those caused after inhibition of mTORC1, altering the differentiation of T cells (41, 42). Our results show that during the processes that mediate AMR, there are changes in immune metabolism compatible with the processes of proliferation and activation of the immune system, and suggest that the pharmacological alteration of the catabolism pathways of some amino acids such as leucine could be a new therapeutic target for the inhibition of lymphocyte proliferation.
Protein-protein interaction networks (PPI) were built to evaluate the relationship between the proteins produced by DEGs and to be able to deduce which are the most important core proteins that control information pathways at the molecular level. The PPI network obtained shows that four proteins control the main biological pathways related to AMR: STAT1, IRF1, LYN, and IRF8. STAT1 is a transcription factor that regulates the expression of genes for cellular response to IFNs and is closely related to the differentiation of T cells towards Th1 (43). The expression of STAT1 has been associated with inflammatory processes in KT related to ischemia-reperfusion damage and BK virus nephropathy (44, 45). Furthermore, STAT1 has been linked to rejection processes in multiple organs in studies conducted in humans (46) and murine models (47). Although STAT1 expression was initially associated with cell-type rejections given its role in Th1 polarization, scientific evidence shows that IFN-related signaling pathways are intimately involved with AMR. That activation of cytotoxic-type effector functions also occurs with intensity in AMR (48). Based on these results, it has been suggested that STAT1 inhibition could be beneficial in slowing down the immune response to graft. Blocking STAT1 with oligonucleotides in murine models of heart transplantation attenuates the recipient’s immune response against the graft, promoting anti-inflammatory effects by inhibiting the expression of MCP-1 (49). Similar results have been obtained in murine models of bone marrow transplantation, where STAT1-deficient mice had decreased IFN-γ production capacity, blocking T cell polarization towards Th1 and reducing rejection events (50).
Furthermore, the STAT1 blockade appears to have more profound effects since it can promote the expansion of Treg (51). However, in other studies, it has been observed that the signaling pathways activated by IFN-γ are also necessary to maintain Treg function, enhancing their ability to prevent rejections. More studies are necessary to understand the true impact of STAT1 blockade for the treatment of rejection.
In this sense, IFN-inducible factors (IRFs, Interferon Regulatory Factor) are transcription factors expressed in response to various stimuli, such as IFN, and have immunoregulatory properties (52). IRF1 is expressed in various organs and various types of immune cells such as NK cells and DCs. IRF1 increases the expression of proinflammatory cytokines such as TNF-α or IL-2 and chemokines such as CXCL10 (30, 53, 54). IRF1 has been shown to act as an early proinflammatory signal during ischemic damage in murine liver transplant models, promoting the secretion of cytokines, such as TNF-α, IFN-γ, or IL-15 (55, 56). Disruption of IRF-1 has been shown to decrease MHC expression in Knock-out (KO) pattern, although all-response rates remained similar to wild-type (57). IRF-1 mediates resistance to necrosis during rejection events, an effect associated with a decrease in IFN-γ (4). In transplant rejection models, the decrease in IRF-1 has been associated with a decrease in the rejection event (58). Other studies carried out with gene expression data obtained from the GEO database have already indicated that the IRF1/STAT1 pathway could be important in rejection processes (46).
On the other hand, the IRF1/STAT1 pathway has also been associated with mechanisms that mediate tolerance in liver transplantation models through the activation of apoptosis (59). On the other hand, IRF8 has a close relationship with B cells differentiation and the secretion of proinflammatory cytokines by the innate immune system (60, 61). IRF8 decreases its expression in liver rejection, although it is a study with a low number of samples (62). In KT, the function of IRF8 has been poorly studied. In a microarray study, IRF8 over-expression has been associated with kidney graft loss in those patients without rejection (63).
Additionally, LYN is a protein kinase of the Src kinase family expressed mainly in hematopoietic cells and played an essential role in regulating the responses of the innate and adaptive immune systems. The study of LYN in the context of KT has not yet been addressed, although the relationship between IFN expression and LYN activation could indicate an important role of the latter in rejection. As we have said previously, IFN pathways appear to be especially active during rejection processes. Plasmacytoid DCs and type I INF-producing APCs increase LYN expression during their activation, acting as a positive regulator in the induction of IFN and the production of cytokines (64). On the other hand, some studies indicate that LYN is involved in the signaling pathways involved in the expression of immunoglobulin Fc receptors such as CD16, or the expression of PIK3CG, a regulator of NK chemotaxis and cytotoxicity, suggesting that the increase in LYN expression could be associated with increased antibody-mediated cytotoxicity (65). Some studies indicate that LYN blockade could prevent the cytotoxic action of NK cells (66). In addition to DCs and NK cells, LYN is also expressed in B cells and is involved in the initiation of BCR signaling and the negative feedback responsible for signal inhibition (67). Inhibitors of this kinase are currently being tested to treat some diseases such as chronic B lymphocytic leukemia (68). Therefore, it could be a possible therapeutic target in transplantation, given its relevant role in activating the B cell and the effector functions of NK cells, although the redundancy of function of the Src family kinases, to which LYN belongs, could limit their activities. In addition, LYN-deficient murine studies produce myeloid hyper-activation, with increased BAFF secretion and a heightened inflammatory state (69). The increase in BAFF has been associated with an increase in the production rate of anti-HLA antibodies. Therefore, the LYN blockade could be counterproductive in the context of KT and should be further studied.
The analysis of transcription factors and protein kinases reveals the importance of IRF1 and IRF8 in regulating genes associated with AMR, like the analysis of PPI networks. However, it also shows the importance of other transcription factors such as RUNX1, RELA, and SPI1. These transcription factors are closely related to the differentiation, proliferation, and activation of cells of the immune system. First, RUNX1 is a hematopoiesis controller that controls the expression of essential genes such as NF-κB and, therefore, inflammatory responses. Second, RELA is the p65 subunit of NF-κB and is essential for its activation (70). SPI1 is a transcription factor of the Ets family expressed in B cells and myeloid cells, crucial in B differentiation and neutrophil activity (71). Bioinformatic analyzes have already reported on the relevance of SPI1 in the regulation of AR processes (70, 72). Thus, PU.1 is a transcription factor encoded by the SPI1 gene and is vital in differentiating and developing macrophages and B cells (73).
Regarding protein kinases, the MAPK family shows a central role in controlling metabolic pathways of both over-and under-expressed genes. The analysis shows that MAPK1 and MAPK14 are important in the expression of genes of both groups. Among the kinases, it is worth mentioning ERK1, a serine-threonine kinase that participates in the Ras-MEK-ERK signaling cascade, involved in multiple processes such as cell adhesion, cell cycle progression, proliferation, and cell survival (72).
During AMR, vasculopathy occurs at the level of the endothelium of the graft vessels. The binding of anti-HLA class I antibodies to the endothelium increases endothelial, smooth muscle proliferation through an ERK1 process, also demonstrated by the generation of anti-non-HLA antibodies such as anti-AT1R (74). The pharmacological blockade of the MEK-ERK pathway effectively reduces chronic allograft nephropathy and the immune response in murine models of KT (75). In cardiac models, ERK inhibition was shown to inhibit the alloresponse, decreasing the leukocyte infiltrate to the graft and IFN production (76). Therefore, ERK1 kinase could be a potential therapeutic target to slow down antibody damage.
The xCell analyzes based on the inference of immune system cells from gene expression data show the important association of NK cell and monocyte infiltration with AMR processes. Therefore, inference data from cell subpopulations infiltrating the graft indicate that NK cells and monocytes are major cellular mediators of graft damage. Therefore, the STAT/IRF molecular pathway acquires great importance since it is directly related to the formation of cytokines, such as IFN and IL-15, or chemistries such as CXCL10, essential for the activation of these cells and their migration towards the graft.
Similar results have been obtained in heart transplantation, where NK cells and monocytes were also associated with AMR (77). Furthermore, several studies carried out in animal transplantation models corroborate the role of these cells in AMR (48, 78). After transplantation, the function of NK cells is exerted through both antibody-dependent mechanisms, such as antibody-dependent cellular cytotoxicity, and non-antibody-dependent mechanisms, through activation by ligands such as osteopontin (79). Although the significance is lost after statistical correction for multiple comparisons, we have observed a clear tendency for patients with AMR to have elevated levels of mast cells. Although mast cells have not been directly related to rejection processes, numerous studies link high levels with fibrosis and chronic graft damage (80, 81). However, its role in transplantation remains controversial since high levels of mast cell-related transcripts have also been found in grafts from tolerant patients (82).
On the other hand, the levels of immature DCs were decreased in patients with AMR. It is known that, during rejection, there is an increase in the antigenic presentation of alloantigens and the maturation and infiltration of APCs (23). Therefore, the maturation of DCs would explain the low relative frequency of immature DCs obtained in our analyzes. In transplantation animal models, inhibition of its maturation has already been shown to suppress graft rejection (83). Therefore, NK cells, monocytes, and DCs are postulated as candidates for their use as biomarkers of AMR. The inhibition of these cells’ maturation and/or the effector functions could be potential therapeutic targets to stop the immune response against the graft effectively.
Despite efforts in searching for new drugs in transplantation, it is still necessary to search for practical strategies for treating or preventing rejection events. In our study, the DEGs derived from the meta-analysis were integrated into bioinformatics applications to suggest new therapeutic options for AMR prevention. Our objective was to find new drugs capable of reversing the expression of the main genes altered during AMR. Of the candidate drugs, only imatinib is currently FDA approved. Imatinib is a tyrosine kinase inhibitor used primarily to treat chronic myeloid leukemia and gastrointestinal tumors (84). With a view to possible use in KT, different studies in murine models have shown the renoprotective effects of imatinib in the long and short term, limiting the progression of glomerulonephritis the prevention of chronic graft nephropathy (85, 86). However, it has not been tested for this purpose in humans. In several clinical cases, imatinib has been used to treat chronic myeloid leukemia after KT has shown to be a safe drug without major side effects (87, 88).
On the other hand, piperlongumine (PL), a natural alkylamide obtained from the Piper longum plant, caught our attention, which has recently attracted the attention of researchers due to its immunomodulatory properties. PL is effective in inhibiting STAT3 and inhibiting proliferation progression in breast cancer cell lines (89). It has also been observed to decrease the proliferation and survival of cell lines derived from hematological neoplasms by decreasing the expression of crucial transcription factors such as NK-κB, STAT1, or MYC (90, 91). In autoimmune diseases, PL has also been observed to suppress DC maturation, decrease proinflammatory cytokines, and inhibit NF-κB (92, 93). Therefore, although PL is still in preclinical studies, its immunomodulatory capacities postulate it as a potential drug for use in transplantation. It should be mentioned that, although our study does not cover the full spectrum of available drugs, limiting themselves to 1309 included in the current version of CMap (18), the strategy used is beneficial for identifying new indications for existing drugs. Future studies in animal models and in vitro will be necessary to confirm the hypotheses arising in this study.
Our study has several limitations; in the first place, the strategy of identifying DEGs by selecting those statistically significant in all the selected GEO studies generates a bias as they are studies with a very uneven number of samples. While the GSE36059 study included 346 samples, the GSE44131 study only 23 samples, so the statistical power of this study will be lower and will act as a limiting factor when identifying DGEs and other genes of interest that may have been lost. Therefore, other statistical methods are necessary to integrate and normalize the data from different studies. Another aspect of assessing is the type of platform used in the studies. In our case, all the studies were carried out using Affymetrix microarrays. Although the GSE36059 was carried out under another version, the same platform allows the data to be more homogeneous and comparable.
The low number of under-expressed genes in the AMR group implies that the functional and topological bioinformatics analyzes have lower power, so the conclusions should be assumed with caution. Finally, inference studies of the immune populations in the graft show the relevant role of NK cells and monocytes, but these algorithms do not discriminate between cell subtypes. It is known that NK cells, in addition to their cytotoxic function, also possess immunoregulatory functions through the secretion of cytokines and chemokines. Therefore, using other tools such as flow cytometry is necessary to identify further the phenotype of infiltrating NKs in the graft. Another of the shortcomings of the study has been not to include samples from cellular rejections. Future studies with cellular rejection samples will be necessary to determine which expression profiles are common to rejection processes and differentiate between subtypes.
Furthermore, the selected GEO studies do not provide information on the acute or chronic nature of the rejection, which means that the analyzes are based on a heterogeneous cohort regarding the post-transplant period in which graft damage occurs. On the other hand, only samples taken at the time of rejection have been studied. In the future, therefore, it would be interesting to carry out longitudinal studies that make it possible to correlate gene expression profiles in the early stages of transplantation with chronic damage, as other previous own studies (94, 95), which allow us to develop prognostic models for the better classification of patients and thereby improve their clinical management. Finally, the results obtained in silico must be confirmed by in vivo and in vitro analysis in the future.
In conclusion, in the present study, we have been able to identify the gene expression profile associated with AMR and the central infiltrating leukocyte populations in the allograft through different bioinformatics tools. The functional analysis shows how the processes related to immune activation and antigenic presentation are fundamental in AMR mechanisms. Furthermore, topological analyses using PPI networks show the importance of the IRF/STAT1 pathways of response to IFN in controlling the expression of genes related to humoral rejection. On the other hand, inference analyzes of leukocyte populations show a crucial role of NK cells and monocytes in graft damage during rejection. Future studies will be necessary to validate the results obtained and assess the potential use of over-expressed genes as rejection biomarkers. In addition, the results suggest that the central proteins obtained in the PPI networks could be potential therapeutic targets to improve the results of kidney transplantation.
Data Availability Statement
The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: https://www.ncbi.nlm.nih.gov/, 93658; https://www.ncbi.nlm.nih.gov/, 44131; https://www.ncbi.nlm.nih.gov/, 50084; https://www.ncbi.nlm.nih.gov/, 36059.
Author Contributions
RA, MM, HM-B, and IL participated in designing and contributed to data processing and supervision, optimization, data analysis, and manuscript writing. IL and RA participated in data analysis for gene expression assays and contributed to writing the manuscript. VJ-C, HM-B, JG, CB, MM-Q, and AP participated in the data processing. RA, SL, and MM-P contributed to gene expression data generation, organization, and manuscript writing. AM provided the study samples, participated in discussions in data analysis strategies. All authors contributed to the article and approved the submitted version.
Funding
Our work was possible thanks to support from Instituto de Salud Carlos III (ISCIII), Spanish Ministry of Economy and Competitiveness. Grant Number PI15/01370, P19/01194, and PI20/00050 and co-funding the European Union with the European Fund of Regional Development (FEDER) with the principle of “A manner to build Europe”.
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.800968/full#supplementary-material
Abbreviations
ACR, Acute Cellular Rejection; ADCC, Antibody-Dependent Cellular Cytotoxicity; AMR, Antibody-Mediated Rejection; APC, Antigen Presenting Cell; AT1R, angiotensin II type 1 receptor; Betweenness Centrality (BC); cPRA, calculated Panel Reactive Antibody; CTLA4, Cytotoxic T-Lymphocyte-Associated Protein 4; CYT, Cytolytic Index; DC, Dendritic Cell; Degree Centrality (DCy); DEGs, Differentially Expressed Genes; DSA, Donor Specific Antibodies; DV, Decision Values; EMA, European Medicines Agency; ES, Enrichment Score; FDA, Food and Drug Administration; GEO, Gene Expression Omnibus; GO, Gene Ontology; GVHD, Graft-versus-host disease; HLA, Human Leukocyte Antigen; KEGG, Kyoto Encyclopedia of Genes and Genomes; PPI, Protein-Protein Interaction; KT, kidney Transplant; MFI, Mean Fluorescence Intensity; MHC, Mayor Histocompatibility Complex; NGS, Next-Generation Sequencing; NR, Not Rejection; PRA, Panel Reactive Antibody; TCR, T-cell Receptor.
References
1. Harris DP, Haynes L, Sayles PC, Duso DK, Eaton SM, Lepak NM, et al. Reciprocal Regulation of Polarized Cytokine Production by Effector B and T Cells. Nat Immunol (2000) 1:475–82. doi: 10.1038/82717
2. Galián JA, Mrowiec A, Muro M. Molecular Targets on B-Cells to Prevent and Treat Antibody-Mediated Rejection in Organ Transplantation. Present Future (2016) 20(7):859–67. doi: 10.1517/14728222.2016.1135904
3. Khatri P, Roedder S, Kimura N, De Vusser K, Morgan AA, Gong Y, et al. A Common Rejection Module (CRM) for Acute Rejection Across Multiple Organs Identifies Novel Therapeutics for Organ Transplantation. J Exp Med (2013) 210:2205–21. doi: 10.1084/JEM.20122709
4. Halloran PF, Pereira AB, Chang J, Matas A, Picton M, de Freitas D, et al. Microarray Diagnosis of Antibody-Mediated Rejection in Kidney Transplant Biopsies: An International Prospective Study (INTERCOM). Am J Transplant (2013) 13:2865–74. doi: 10.1111/AJT.12465
5. Legaz I, Bernardo MV, Alfaro R, Martínez-Banaclocha H, Galián JA, Jimenez-Coll V, et al. PCR Array Technology in Biopsy Samples Identifies Up-Regulated mTOR Pathway Genes as Potential Rejection Biomarkers After Kidney Transplantation. Front Med (2021) 0:547849. doi: 10.3389/FMED.2021.547849
6. Home - GEO - NCBI. Available at: https://www.ncbi.nlm.nih.gov/geo/ (Accessed October 6, 2021).
7. Kolesnikov N, Hastings E, Keays M, Melnichuk O, Tang YA, Williams E, et al. ArrayExpress Update—Simplifying Data Submissions. Nucleic Acids Res (2015) 43:D1113–6. doi: 10.1093/NAR/GKU1057
8. Consortium TGO. Gene Ontology Consortium: Going Forward. Nucleic Acids Res (2015) 43:D1049–56. doi: 10.1093/NAR/GKU1179
9. Kanehisa M, Goto S, Sato Y, Kawashima M, Furumichi M, Tanabe M. Data, Information, Knowledge and Principle: Back to Metabolism in KEGG. Nucleic Acids Res (2014) 42:D199–205. doi: 10.1093/NAR/GKT1076
10. Draw Venn Diagram. Available at: http://bioinformatics.psb.ugent.be/webtools/Venn/ (Accessed October 6, 2021).
11. Kanehisa M, Sato Y, Kawashima M. KEGG Mapping Tools for Uncovering Hidden Features in Biological Data. Protein Sci (2021). doi: 10.1002/PRO.4172
12. Consortium GO. The Gene Ontology (GO) Project in 2006. Nucleic Acids Res (2006) 34:D322–6. doi: 10.1093/NAR/GKJ021
13. Xcell. Available at: https://xcell.ucsf.edu/ (Accessed October 8, 2021).
14. Xia J, Gill EE, Hancock REW. NetworkAnalyst for Statistical, Visual and Network-Based Meta-Analysis of Gene Expression Data. Nat Protoc 2015 106 (2015) 10:823–44. doi: 10.1038/nprot.2015.052
15. Zhou G, Soufan O, Ewald J, Hancock REW, Basu N, Xia J. NetworkAnalyst 3.0: A Visual Analytics Platform for Comprehensive Gene Expression Profiling and Meta-Analysis. Nucleic Acids Res (2019) 47:W234–41. doi: 10.1093/NAR/GKZ240
16. 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:D607–13. doi: 10.1093/NAR/GKY1131
17. Gene2drug. Available at: https://gene2drug.tigem.it/ (Accessed October 8, 2021).
18. Lamb J, Crawford ED, Peck D, Modell JW, Blat IC, Wrobel MJ, et al. The Connectivity Map: Using Gene-Expression Signatures to Connect Small Molecules, Genes, and Disease. Sci (80 ) (2006) 313:1929–35. doi: 10.1126/science.1132939
19. STITCH: Chemical Association Networks. Available at: http://stitch.embl.de/ (Accessed October 8, 2021).
20. Szklarczyk D, Santos A, von Mering C, Jensen LJ, Bork P, Kuhn M. STITCH 5: Augmenting Protein–Chemical Interaction Networks With Tissue and Affinity Data. Nucleic Acids Res (2016) 44:D380–4. doi: 10.1093/NAR/GKV1277
21. Alfaro R, Legaz I, González-Martínez G, Jimenez-Coll V, Martínez-Banaclocha H, Galián JA, et al. Monitoring of B Cell in Kidney Transplantation: Development of a Novel Clusters Analysis and Role of Transitional B Cells in Transplant Outcome. Diagnostics 2021 Vol 11 Page 641 (2021) 11:641. doi: 10.3390/DIAGNOSTICS11040641
22. Alfaro R, Legaz I, Jimenez-Coll V, band JE kaaoui E, Martínez-Banaclocha H, Galián JA, et al. MicroRNA Expression Changes in Kidney Transplant: Diagnostic Efficacy of miR-150-5p as Potential Rejection Biomarker, Pilot Study. J Clin Med 8 (2021) 10:2748. doi: 10.3390/JCM10132748
23. Alegre ML, Lakkis FG, Morelli AE. Antigen Presentation in Transplantation. Trends Immunol (2016) 37:831–43. doi: 10.1016/J.IT.2016.09.003
24. Li L, Greene I, Readhead B, Menon MC, Kidd BA, Uzilov AV, et al. Novel Therapeutics Identification for Fibrosis in Renal Allograft Using Integrative Informatics Approach. Sci Rep 2017 71 (2017) 7:1–14. doi: 10.1038/srep39487
25. Kim I, Kim JH, Han N, Kim S, Kim YS, Oh JM. Gene Expression Profiles for Predicting Antibody−Mediated Kidney Allograft Rejection: Analysis of GEO Datasets. Int J Mol Med (2018) 42:2303–11. doi: 10.3892/IJMM.2018.3798
26. Valiño-Rivas L, Cuarental L, Grana O, Bucala R, Leng L, Sanz A, et al. TWEAK Increases CD74 Expression and Sensitizes to DDT Proinflammatory Actions in Tubular Cells. PloS One (2018) 13:e0199391. doi: 10.1371/JOURNAL.PONE.0199391
27. Kakuta Y, Okumi M, Miyagawa S, Tsutahara K, Abe T, Yazawa K, et al. Blocking of CCR5 and CXCR3 Suppresses the Infiltration of Macrophages in Acute Renal Allograft Rejection. Transplantation (2012) 93:24–31. doi: 10.1097/TP.0B013E31823AA585
28. Segerer S, Böhmig GA, Exner M, Kerjaschki D, Regele H, Schlöndorff D. Role of CXCR3 in Cellular But Not Humoral Renal Allograft Rejection. Transpl Int (2005) 18:676–80. doi: 10.1111/J.1432-2277.2005.00117.X
29. Groom JR, Luster AD. CXCR3 in T Cell Function. Exp Cell Res (2011) 317:620–31. doi: 10.1016/J.YEXCR.2010.12.017
30. Thammavongsa V, Kim HK, Missiakas D, Schneewind O. Staphylococcal Manipulation of Host Immune Responses. Nat Rev Microbiol 2015 139 (2015) 13:529–43. doi: 10.1038/nrmicro3521
31. Chen R, Sigdel TK, Li L, Kambham N, Dudley JT, Hsieh S, et al. Differentially Expressed RNA From Public Microarray Data Identifies Serum Protein Biomarkers for Cross-Organ Transplant Rejection and Other Conditions. PloS Comput Biol (2010) 6:e1000940. doi: 10.1371/JOURNAL.PCBI.1000940
32. Huang H, Xu X, Yao C, Cai M, Qian Y, Wang X, et al. Serum Levels of CXCR3 Ligands Predict T Cell-Mediated Acute Rejection After Kidney Transplantation. Mol Med Rep (2014) 9:45–50. doi: 10.3892/MMR.2013.1753
33. Xu J, Ma T, Deng G, Zhuang J, Li C, Wang S, et al. Inhibition of C−X−C Motif Chemokine 10 Reduces Graft Loss Mediated by Memory CD8+ T Cells in a Rat Cardiac Re−Transplant Model. Exp Ther Med (2018) 15:1560–7. doi: 10.3892/ETM.2017.5585
34. Kautz-Neu K, Schwonberg K, Fischer MR, Schermann AI, von Stebut E. Dendritic Cells in Leishmania Major Infections: Mechanisms of Parasite Uptake, Cell Activation and Evidence for Physiological Relevance. Med Microbiol Immunol 2012 2014 (2012) 201:581–92. doi: 10.1007/S00430-012-0261-2
35. Lochner M, Berod L, Sparwasser T. Fatty Acid Metabolism in the Regulation of T Cell Function. Trends Immunol (2015) 36:81–91. doi: 10.1016/J.IT.2014.12.005
36. Buck MD, O’Sullivan D, Pearce EL. T Cell Metabolism Drives Immunity. J Exp Med (2015) 212:1345–60. doi: 10.1084/JEM.20151159
37. Raha S, Raud B, Oberdörfer L, Castro CN, Schreder A, Freitag J, et al. Disruption of De Novo Fatty Acid Synthesis via Acetyl-CoA Carboxylase 1 Inhibition Prevents Acute Graft-Versus-Host Disease. Eur J Immunol (2016) 46:2233–8. doi: 10.1002/EJI.201546152
38. Byersdorfer CA, Tkachev V, Opipari AW, Goodell S, Swanson J, Sandquist S, et al. Effector T Cells Require Fatty Acid Metabolism During Murine Graft-Versus-Host Disease. Blood (2013) 122:3230–7. doi: 10.1182/BLOOD-2013-04-495515
39. Amorim LM, van Tong H, Hoan NX, Vargas L de B, Ribeiro EM de SF, Petzl-Erler ML, et al. KIR-HLA Distribution in a Vietnamese Population From Hanoi. Hum Immunol (2018) 79:93–100. doi: 10.1016/j.humimm.2017.11.011
40. Macintyre AN, Gerriets VA, Nichols AG, Michalek RD, Rudolph MC, Deoliveira D, et al. The Glucose Transporter Glut1 is Selectively Essential for CD4 T Cell Activation and Effector Function. Cell Metab (2014) 20:61. doi: 10.1016/J.CMET.2014.05.004
41. Sinclair LV, Rolf J, Emslie E, Shi Y-B, Taylor PM, Cantrell DA. Control of Amino-Acid Transport by Antigen Receptors Coordinates the Metabolic Reprogramming Essential for T Cell Differentiation. Nat Immunol 2013 145 (2013) 14:500–8. doi: 10.1038/ni.2556
42. Hidayat S, Yoshino KI, Tokunaga C, Hara K, Matsuo M, Yonezawa K. Inhibition of Amino acid-mTOR Signaling by a Leucine Derivative Induces G1 Arrest in Jurkat Cells. Biochem Biophys Res Commun (2003) 301:417–23. doi: 10.1016/S0006-291X(02)03052-8
43. Oestreich KJ, Weinmann AS. Transcriptional Mechanisms That Regulate T Helper 1 Cell Differentiation. Curr Opin Immunol (2012) 24:191–5. doi: 10.1016/J.COI.2011.12.004
44. Kemmner S, Bachmann Q, Steiger S, Lorenz G, Honarpisheh M, Foresto-Neto O, et al. STAT1 Regulates Macrophage Number and Phenotype and Prevents Renal Fibrosis After Ischemia-Reperfusion Injury. Am J Physiol Renal Physiol (2019) 316(2):F277–91. doi: 10.1152/AJPRENAL.00004.2018
45. Jia L, Fu W, Jia R, Wu L, Li X, Jia Q, et al. Identification of Potential Key Protein Interaction Networks of BK Virus Nephropathy in Patients Receiving Kidney Transplantation. Sci Rep (2018) 8:1–8. doi: 10.1038/s41598-018-23492-2
46. Spivey TL, Uccellini L, Ascierto ML, Zoppoli G, De Giorgi V, Delogu LG, et al. Gene Expression Profiling in Acute Allograft Rejection: Challenging the Immunologic Constant of Rejection Hypothesis. J Transl Med (2011) 9:1–22. doi: 10.1186/1479-5876-9-174
47. Famulski KS, Einecke G, Reeve J, Ramassar V, Allanach K, Mueller T, et al. Changes in the Transcriptome in Allograft Rejection: IFN-γ-Induced Transcripts in Mouse Kidney Allografts. Am J Transplant (2006) 6:1342–54. doi: 10.1111/J.1600-6143.2006.01337.X
48. Hidalgo LG, Sis B, Sellares J, Campbell PM, Mengel M, Einecke G, et al. NK Cell Transcripts and NK Cells in Kidney Biopsies From Patients With Donor-Specific Antibodies: Evidence for NK Cell Involvement in Antibody-Mediated Rejection. Am J Transplant (2010) 10:1812–22. doi: 10.1111/J.1600-6143.2010.03201.X
49. Stojanovic T, Wagner AH, Wang S, Kiss E, Rockstroh N, Bedke J, et al. STAT-1 Decoy Oligodeoxynucleotide Inhibition of Acute Rejection in Mouse Heart Transplants. Basic Res Cardiol (2009) 104:719–29. doi: 10.1007/S00395-009-0028-0
50. Mariotti J, Foley J, Ryan K, Buxhoeveden N, Kapoor V, Amarnath S, et al. Graft Rejection as a Th1-Type Process Amenable to Regulation by Donor Th2-Type Cells Through an Interleukin-4/STAT6 Pathway. Blood (2008) 112:4765–75. doi: 10.1182/BLOOD-2008-05-154278
51. Ma H, Lu C, Ziegler J, Liu A, Sepulveda A, Okada H, et al. Absence of Stat1 in Donor CD4+ T Cells Promotes the Expansion of Tregs and Reduces Graft-Versus-Host Disease in Mice. J Clin Invest (2011) 121:2554–69. doi: 10.1172/JCI43706
52. Mogensen TH. IRF and STAT Transcription Factors - From Basic Biology to Roles in Infection, Protective Immunity, and Primary Immunodeficiencies. Front Immunol (2019) 0:3047. doi: 10.3389/FIMMU.2018.03047
53. Ogasawara K, Hida S, Azimi N, Tagaya Y, Sato T, Yokochi-Fukuda T, et al. Requirement for IRF-1 in the Microenvironment Supporting Development of Natural Killer Cells. Nat 1998 3916668 (1998) 391:700–3. doi: 10.1038/35636
54. Yeruva S, Ramadori G, Raddatz D. NF-κb-Dependent Synergistic Regulation of CXCL10 Gene Expression by IL-1β and IFN-γ in Human Intestinal Epithelial Cell Lines. Int J Color Dis (2007) 23:305–17. doi: 10.1007/S00384-007-0396-6
55. Wang Y, John R, Chen J, Richardson JA, Shelton JM, Bennett M, et al. IRF-1 Promotes Inflammation Early After Ischemic Acute Kidney Injury. J Am Soc Nephrol (2009) 20:1544–55. doi: 10.1681/ASN.2008080843
56. Yokota S, Yoshida O, Dou L, Spadaro AV, Isse K, Ross MA, et al. IRF-1 Promotes Liver Transplant Ischemia/Reperfusion Injury via Hepatocyte IL-15/IL-15rα Production. J Immunol (2015) 194:6045–56. doi: 10.4049/JIMMUNOL.1402505
57. Afrouzian M, Ramassar V, Urmson J, Zhu L-F, Halloran PF. Transcription Factor IRF-1 in Kidney Transplants Mediates Resistance to Graft Necrosis During Rejection. J Am Soc Nephrol (2002) 13:1199–209. doi: 10.1097/01.ASN.0000013302.11876.A5
58. Solomon M, Flodström-Tullberg M, Sarvetnick N. Beta-Cell Specific Expression of Suppressor of Cytokine Signaling-1 (SOCS-1) Delays Islet Allograft Rejection by Down-Regulating Interferon Regulatory Factor-1 (IRF-1) Signaling. Transpl Immunol (2011) 24:181–8. doi: 10.1016/j.trim.2010.11.007
59. Cordoba SP, Wang C, Williams R, Li J, Smit L, Sharland A, et al. Gene Array Analysis of a Rat Model of Liver Transplant Tolerance Identifies Increased Complement C3 and the STAT-1 / IRF-1 Pathway During Tolerance Induction. Publ Am Assoc Study Liver Dis Int Liver Transplant Soc (2006) 12:636–43. doi: 10.1002/lt
60. Zhao J, Hee JK, Li H, Huang B, Yang M, Zhu C, et al. IRF-8/Interferon (IFN) Consensus Sequence-Binding Protein is Involved in Toll-Like Receptor (TLR) Signaling and Contributes to the Cross-Talk Between TLR and IFN-γ Signaling Pathways. J Biol Chem (2006) 281:10073–80. doi: 10.1074/jbc.M507788200
61. Shin DM, Lee CH, Morse HC. IRF8 Governs Expression of Genes Involved in Innate and Adaptive Immunity in Human and Mouse Germinal Center B Cells. PloS One (2011) 6:e27384. doi: 10.1371/journal.pone.0027384
62. Nasiri M, Geramizadeh B, Nabavizadeh SH, Male-Hosseini SA, Karimi MH, Saadat I. mRNA Expression of Interferon Regulatory Factors During Acute Rejection of Liver Transplants in Patients With Autoimmune Hepatitis. Int J Organ Transplant Med (2018) 9:34–40.
63. Fu Q, Liao M, Feng C, Tang J, Liao R, Wei L, et al. Profiling of mRNA of Interstitial Fibrosis and Tubular Atrophy With Subclinical Inflammation in Recipients After Kidney Transplantation. Aging (Albany NY) (2019) 11:5215–31. doi: 10.18632/aging.102115
64. Dallari S, Macal M, Loureiro ME, Jo Y, Swanson L, Hesser C, et al. Src Family Kinases Fyn and Lyn Are Constitutively Activated and Mediate Plasmacytoid Dendritic Cell Responses. Nat Commun (2017) 8:1–16. doi: 10.1038/ncomms14830
65. Mota G, Moldovan I, Calugaru A, Hirt M, Kozma E, Galatiuc C, et al. Interaction of Human Immunoglobulin G With CD16 on Natural Killer Cells: Ligand Clearance, Fcγriiia Turnover and Effects of Metalloproteinases on Fcγriiia-Mediated Binding, Signal Transduction and Killing. Scand J Immunol (2004) 59:278–84. doi: 10.1111/j.0300-9475.2004.01398.x
66. Oykhman P, Timm-McCann M, Xiang RF, Islam A, Li SS, Stack D, et al. Requirement and Redundancy of the Src Family Kinases Fyn and Lyn in Perforin-Dependent Killing of Cryptococcus Neoformans by NK Cells. Infect Immun (2013) 81:3912–22. doi: 10.1128/IAI.00533-13
67. Manuscript A. Lyn Deficiency Affects B Cell Maturation as Well as Survival. Eur J Inmunol (2012) 42:511–21. doi: 10.1002/eji.201141940.Lyn
68. Kadia T, Delioukina ML, Kantarjian HM, Keating MJ. A Pilot Phase II Study of the Lyn Kinase Inhibitor Bafetinib in Patients With Relapsed or Refractory B Cell Chronic Lymphocytic Leukemia. Blood (2011) 118:2858–8. doi: 10.1182/blood.V118.21.2858.2858
69. Okuda T, Nishimura M, Nakao M, Fujita Y. RUNX1/AML1: A Central Player in Hematopoiesis. Int J Hematol (2001) 74:252–7. doi: 10.1007/BF02982057
70. Molinero LL, Alegre ML. Role of T Cell-Nuclear Factor κb in Transplantation. Transplant Rev (2012) 26:189–200. doi: 10.1016/j.trre.2011.07.005
71. Lloberas J, Soler C, Celada A. The Key Role of PU.1/SPI-1 in B Cells, Myeloid Cells and Macrophages. Immunol Today (1999) 20:184–9. doi: 10.1016/S0167-5699(99)01442-5
72. Roskoski R. ERK1/2 MAP Kinases: Structure, Function, and Regulation. Pharmacol Res (2012) 66:105–43. doi: 10.1016/j.phrs.2012.04.005
73. Oikawa T, Yamada T, Kihara-Negishi F, Yamamoto H, Kondoh N, Hitomi Y, et al. The Role of Ets Family Transcription Factor PU.1 in Hematopoietic Cell Differentiation, Proliferation and Apoptosis. Cell Death Differ (1999) 6:599–608. doi: 10.1038/sj.cdd.4400534
74. Hiemann NE, Meyer R, Wellnhofer E, Schoenemann C, Heidecke H, Lachmann N, et al. Non-HLA Antibodies Targeting Vascular Receptors Enhance Alloimmune Response and Microvasculopathy After Heart Transplantation. Transplantation (2012) 94:919–24. doi: 10.1097/TP.0b013e3182692ad2
75. Wang S, Jiang J, Guan Q, Wang H, Nguan CYC, Jevnikar AM, et al. Reduction of Chronic Allograft Nephropathy by Inhibition of Extracellular Signal-Regulated Kinase 1 and 2 Signaling. Am J Physiol Ren Physiol (2008) 295:F672–679. doi: 10.1152/ajprenal.90285.2008
76. Wang S, Guan Q, Diao H, Lian D, Zhong R, Jevnikar AM, et al. Prolongation of Cardiac Allograft Survival by Inhibition of ERK1/2 Signaling in a Mouse Model. Transplantation (2007) 83:323–32. doi: 10.1097/01.tp.0000251374.49225.19
77. Loupy A, Duong Van Huyen JP, Hidalgo L, Reeve J, Racapé M, Aubert O, et al. Gene Expression Profiling for the Identification and Classification of Antibody-Mediated Heart Rejection. Circulation (2017) 135:917–35. doi: 10.1161/CIRCULATIONAHA.116.022907
78. Bergler T, Jung B, Bourier F, Kühne L, Banas MC. Infiltration of Macrophages Correlates With Severity of Allograft Rejection and Outcome in Human Kidney Transplantation. PloS One (2016) 11:e0156900. doi: 10.1371/journal.pone.0156900
79. Zhang ZX, Huang X, Jiang J, Lau A, Yin Z, Liu W, et al. Natural Killer Cells Mediate Long-Term Kidney Allograft Injury. Transplantation (2015) 99:916–24. doi: 10.1097/TP.0000000000000665
80. Papadimitriou JC, Drachenberg CB, Ramos E, Ugarte R, Haririan A. Mast Cell Quantitation in Renal Transplant Biopsy Specimens as a Potential Marker for the Cumulative Burden of Tissue Injury. Transplant Proc (2013) 45:1469–71. doi: 10.1016/j.transproceed.2013.01.078
81. Pardo J, Diaz L, Errasti P, Idoate M, De Alava E, Sola I, et al. Mast Cells in Chronic Rejection of Human Renal Allografts. Virchows Arch (2000) 437:167–72. doi: 10.1007/s004280000211
82. Elieh D, Komi A, Ribatti D. Mast Cell-Mediated Mechanistic Pathways in Organ Transplantation. Eur J Pharmacol (2019) 857:172458. doi: 10.1016/j.ejphar.2019.172458
83. Yang C, Zhang Y, Wang J, Li L, Wang L, Hu M, et al. Novel Cyclic Helix B Peptide Inhibits Dendritic Cell Maturation During Amelioration of Acute Kidney Graft Rejection Through Jak-2/Stat3/Socs1. Cell Death Dis (2015) 6:e1993. doi: 10.1038/cddis.2015.338
84. Miura M. Therapeutic Drug Monitoring of Imatinib. Ann Oncol (2015) 26:vii38. doi: 10.1093/annonc/mdv432.01
85. Iyoda M, Shibata T, Wada Y, Kuno Y, Shindo-Hirai Y, Matsumoto K, et al. Long- and Short-Term Treatment With Imatinib Attenuates the Development of Chronic Kidney Disease in Experimental Anti-Glomerular Basement Membrane Nephritis. Nephrol Dial Transplant (2013) 28:576–84. doi: 10.1093/ndt/gfs414
86. Savikko J, Rintala JM, Rintala SE, Koskinen PK, Von Willebrand E. Early Short-Term Imatinib Treatment is Sufficient to Prevent the Development of Chronic Allograft Nephropathy. Nephrol Dial Transplant (2011) 26:3026–32. doi: 10.1093/ndt/gfq790
87. Thierry A, Dreyfus B, Bridoux F, Ayache RA, Milin S, Guilhot F, et al. Long-Term Molecular Efficacy and Safety of Imatinib in a Patient With Chronic Myeloid Leukaemia After Renal Transplantation. Nephrol Dial Transplant (2007) 22:1791–2. doi: 10.1093/ndt/gfm052
88. Mulder KE, Egorin MJ, Sawyer MB. Renal Dysfunction in a Renal Transplant Patient Treated Concurrently With Cyclosporine and Imatinib. Invest New Drugs (2012) 30:2400–2. doi: 10.1007/s10637-011-9769-3
89. Bharadwaj U, Eckols TK, Kolosov M, Kasembeli MM, Adam A, Torres D, et al. Drug-Repositioning Screening Identi Fi Ed Piperlongumine as a Direct STAT3 Inhibitor With Potent Activity Against Breast Cancer. Oncogene (2015) 34:1341–53. doi: 10.1038/onc.2014.72
90. Han SS, Han S, Kamberos NL. Piperlongumine Inhibits the Proliferation and Survival of B-Cell Acute Lymphoblastic Leukemia Cell Lines Irrespective of Glucocorticoid Resistance. Biochem Biophys Res Commun (2014) 452:669–75. doi: 10.1016/j.bbrc.2014.08.131
91. Han SS, Son DJ, Yun H, Kamberos NL, Janz S. Piperlongumine Inhibits Proliferation and Survival of Burkitt Lymphoma. Vitro Leuk Res (2013) 37:146–54. doi: 10.1016/j.leukres.2012.11.009
92. Xiao Y, Shi M, Qiu Q, Huang M, Zeng S, Zou Y, et al. Piperlongumine Suppresses Dendritic Cell Maturation by Reducing Production of Reactive Oxygen Species and Has Therapeutic Potential for Rheumatoid Arthritis. J Inmunol (2016) 196(12):4925–34. doi: 10.4049/jimmunol.1501281
93. Gu SM, Yun J, Son DJ, Kim HY, Nam KT, Kim HD, et al. Piperlongumine Attenuates Experimental Autoimmune Encephalomyelitis Through Inhibition of NF-KappaB Activity. Free Radic Biol Med (2017) 103:133–45. doi: 10.1016/j.freeradbiomed.2016.12.027
94. Boix F, Legaz I, Minhas A, Alfaro R, Mrowiec A, Mart H, et al. Identification of Peripheral CD154+ T Cells and HLA-DRB1 as Biomarkers of Acute Cellular Rejection in Adult Liver Transplant Recipients. Clin Exp Immunol (2020) 203(2):315–28. doi: 10.1111/cei.13533
Keywords: bioinformatics tool, biomarkers, acute rejection, kidney transplant, new target drugs
Citation: Alfaro R, Martínez-Banaclocha H, Llorente S, Jimenez-Coll V, Galián JA, Botella C, Moya-Quiles MR, Parrado A, Muro-Perez M, Minguela A, Legaz I and Muro M (2021) Computational Prediction of Biomarkers, Pathways, and New Target Drugs in the Pathogenesis of Immune-Based Diseases Regarding Kidney Transplantation Rejection. Front. Immunol. 12:800968. doi: 10.3389/fimmu.2021.800968
Received: 24 October 2021; Accepted: 29 November 2021;
Published: 15 December 2021.
Edited by:
Antonij Slavcev, Institute for Clinical and Experimental Medicine (IKEM), CzechiaReviewed by:
Rusan Ali Catar, Charité University Medicine Berlin, GermanyPeixiang Lan, Huazhong University of Science and Technology, China
Copyright © 2021 Alfaro, Martínez-Banaclocha, Llorente, Jimenez-Coll, Galián, Botella, Moya-Quiles, Parrado, Muro-Perez, Minguela, Legaz and Muro. 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: Manuel Muro, manuel.muro@carm.es; Isabel Legaz, isalegaz@um.es
†These authors have contributed equally to this work and share last authorship