1. Laboratory of Systems Tumor Immunology, Department of Dermatology, Friedrich-Alexander-Universität Erlangen-Nürnberg (FAU) and Universitätsklinikum Erlangen, Erlangen, Germany.
2. RNA Group, Department of Dermatology, Friedrich-Alexander-Universität Erlangen-Nürnberg (FAU) and Universitätsklinikum Erlangen, Erlangen, Germany.
3. Department of Human Genetics, Universitätsklinikum Erlangen, Erlangen, Germany.
4. Division of Molecular Immunology, Department of Medicine 3, Universitätsklinikum Erlangen, Erlangen, Germany.
5. Deutsches Zentrum Immuntherapie (DZI), Erlangen, Germany.
6. Comprehensive Cancer Center (CCC) Erlangen, Erlangen, Germany.
Dendritic cells (DCs) are professional antigen-presenting cells that induce and regulate adaptive immunity by presenting antigens to T cells. Due to their coordinative role in adaptive immune responses, DCs have been used as cell-based therapeutic vaccination against cancer. The capacity of DCs to induce a therapeutic immune response can be enhanced by re-wiring of cellular signalling pathways with microRNAs (miRNAs).
Methods: Since the activation and maturation of DCs is controlled by an interconnected signalling network, we deploy an approach that combines RNA sequencing data and systems biology methods to delineate miRNA-based strategies that enhance DC-elicited immune responses.
Results: Through RNA sequencing of IKKβ-matured DCs that are currently being tested in a clinical trial on therapeutic anti-cancer vaccination, we identified 44 differentially expressed miRNAs. According to a network analysis, most of these miRNAs regulate targets that are linked to immune pathways, such as cytokine and interleukin signalling. We employed a network topology-oriented scoring model to rank the miRNAs, analysed their impact on immunogenic potency of DCs, and identified dozens of promising miRNA candidates, with miR-15a and miR-16 as the top ones. The results of our analysis are presented in a database that constitutes a tool to identify DC-relevant miRNA-gene interactions with therapeutic potential (
Conclusions: Our approach enables the systematic analysis and identification of functional miRNA-gene interactions that can be experimentally tested for improving DC immunogenic potency.
Keywords: dendritic cell, therapeutic vaccination in cancer, systems biology, network biology, microRNAs
Dendritic cells (DCs) play an important role in regulating adaptive immunity by presenting antigens to T cells . Due to the unique function of DCs in the coordination of adaptive immune responses, they have been tested as cell-based vaccination against tumours [2-5]. To obtain immunogenic potency ex vivo, monocyte-derived DCs need to go through a complex maturation process, in which DCs are exposed to a monocyte-conditioned medium or a cocktail of cytokines [6,7]. These treatments result in various phenotypic changes in DCs, such as upregulation of co-stimulatory surface markers (e.g., CD80 and CD40) and secretion of pro-inflammatory cytokines (e.g., IL-12 and TNFα). The matured DCs loaded with cancer antigens are infused into patients and trigger a selective immune response by migrating into the peripheral lymphatic tissue, where they encounter and activate tumour-specific T cells .
The capacity of DCs to induce an immune response can be improved by molecular engineering. Pfeiffer and co-workers enhanced DCs through electroporation with mRNA encoding a constitutively active variant of IKKβ (caIKK), a kinase upstream of NF-κB that is a key regulator of the immune response [9,10]. Specifically, the kinase phosphorylates IκB resulting in desequestration of the transcription factor (TF) NF-κB and its translocation into the nucleus, where it regulates expression of immune-related genes such as cytokines [11,12]. The engineered caIKK promotes constant activation of NF-κB signalling, and the cells expressing it (hereafter labelled caIKK-DCs) displayed increased expression of co-stimulatory molecules and pro-inflammatory cytokines [10,13]. When those cells were in addition loaded with the Melan-A melanoma antigen, they were able to induce repeated expansion of Melan-A-specific cytotoxic T cells with a memory phenotype . Besides, caIKK-DCs displayed an increased ability to activate autologous NK cells . They are currently under evaluation as vaccine in a phase-I clinical trial for the treatment of uveal melanoma patients (NCT04335890).
Our hypothesis is that DCs can be further improved using non-coding RNAs, in particular microRNAs (miRNAs), interacting with key regulators of DC activation and maturation. miRNAs are a class of small endogenous non-coding RNAs with a length of 22-25 nucleotides. Through the inhibition and modulation of the transcription and translation of specific protein-coding genes, miRNAs can alter the basal state of cells and the outcome of stimulatory events [15,16]. Increasing evidence shows that miRNAs play a crucial role in the development and function of DCs [17-19]. They serve as important regulators of complex networks by targeting key signalling genes to regulate proliferation and cell death as well as homeostasis . It has also been found that miRNAs are pivotal in both adaptive and innate immunity, e.g., by controlling the differentiation of immune cell subsets and their immunological functions . In particular, miRNAs can modulate the immune response by inducing apoptosis, affecting homeostasis, and changing the cytokine profile of DCs . Further, one can use miRNAs, alone or in combination, in therapeutic setups to inhibit expression of selected genes in cancer and other targeted cells [23-26].
To facilitate the re-wiring of DCs, it is crucial to understand the intracellular regulatory processes involved in DC maturation and activation. However, the regulatory networks eliciting the activation and maturation of DCs involve multiple interconnected signalling and transcriptional circuits, and their understanding and proper manipulation requires the combined use of high-throughput data and systems biology methods [27,28]. We here present a systems biology approach to understanding the role that miRNAs play in regulating the function of DCs in immunotherapy (Figure 1), and exploit this knowledge to enhance their potential to stimulate an immune response using miRNAs.
In this study, we chose the caIKK-DCs as an ideal model system to identify miRNAs that are involved in DC activation via NF-κB signalling and can boost pro-inflammatory signals. We think the identified miRNAs can enable the DCs to repetitively stimulate T cell expansion. To this end, we performed RNA sequencing (RNA-seq) to obtain the transcriptomic profile (i.e., protein-coding genes and miRNAs) of caIKK-DCs in relation to standard DCs. Next, we analysed miRNA-gene interactions at the pathway level and reconstructed regulatory networks underlying immunological functions of DCs. We then performed network-based prioritization of miRNAs by integrating their expression profiles and their strength of association with other protein-coding genes.
Our analysis identified dozens of miRNA candidates in the regulation of caIKK-DCs, with miR-15a-5p and miR-16-5p as prominent examples. We showed that both miRNAs may exert a strong regulatory effect on genes involved in NF-κB signalling and also target chemokines and cytokines regulating T-cell responses. Moreover, we delineated molecular mechanisms through which the miRNAs alter the immunogenic potency of caIKK-DCs. The results of our analysis are available in a web database that facilitates their exploration and visualization (https://www.synmirapy.net/dc-optimization), thereby providing researchers with a tool to select functional miRNA-gene interactions with therapeutic potential in DCs for experimental investigation.
Monocyte-derived DCs were generated as previously described . In brief, blood samples from seven healthy donors were collected after approval was granted by the responsible institutional review board (Ethikkommission der Friedrich-Alexander-Universität Erlangen-Nürnberg, Ref. no. 4158) and written informed consent was obtained. Peripheral blood mononuclear cells were isolated from whole blood using density gradient centrifugation. Monocytes were extracted from the non-adherent fraction by plastic adherence, were cultured in DC medium (RPMI; Lonza, Verviers, Belgium) containing 1% non-autologous human plasma (Sigma-Aldrich, St Louis, USA), 2 mM L‑glutamine (Lonza), and 20 mg/L gentamicin (Lonza), and were differentiated into DCs by application of 800 U/mL GM‑CSF (Miltenyi Biotec, Bergisch Gladbach, Germany) and 250 IU/mL IL‑4 (Miltenyi Biotec) on days 1, 3, and 5. After six days in culture, DCs were matured for 24 hours using a cytokine cocktail consisting of 200 IU/mL IL‑1β (CellGenix, Freiburg, Germany), 1 000 U/mL IL‑6 (CellGenix), 10 ng/mL TNF (Beromun, Boehringer Ingelheim Pharma, Germany), and 1 μg/mL PGE2 (Pfizer, Zurich, Switzerland).
Generation of in vitro transcribed RNA with mMESSAGE mMACHINE™ T7 ULTRA transcription-kits (Thermofisher scientific, Waltham, USA) and the electroporation of cocktail-matured DCs (cmDCs) was carried out as previously described . For transcriptome analyses, cmDCs were electroporated using 30 µg RNA encoding constitutively active IKKβ (caIKK) since its introduction into mature DCs has been shown to improve activation of T cells  and natural killer cells . DCs electroporated with RNA encoding melanoma antigen recognized by T cells 1 (Melan-A) were used as control, and such DCs were shown to have no influence on the DCs' transcriptome profile . After electroporation, DCs were cultured in DC medium containing GM‑CSF and IL‑4 at the concentrations indicated above.
Total RNA including small RNAs was extracted four hours after electroporation using the RNeasy Plus Mini Kit (QIAGEN GmbH, Hilden, Germany) and the generated samples were sequenced using an Illumina HiSeq-2500. Demultiplexed reads were filtered for ribosomal RNAs, transfer RNAs, mitochondrial rRNAs, and mitochondrial tRNAs. The reads were aligned to the human reference genome (hg19) using STAR (v2.5.2b) and assigned to genes using Subread (v1.5.2). Only uniquely mapping reads that could unambiguously be assigned to a single gene were considered for analysis (Supplementary Table S1).
For miRNA expression quantification, we performed a quality check of the RNA-seq reads with FastQC , mapped the short sequences to the human reference genome (hg19) using BWA , and calculated raw read counts of mature miRNAs that are known and annotated in miRBase v21  (Supplementary Table S2; See Supplementary Materials for details).
Before differential expression analysis, we aggregated read counts of Ensembl identifiers that represent the same gene and discarded genes with less than 5 read counts in any sample to increase power for detecting differentially expressed genes [34,35]. Next, we used DESeq2  in R version 3.6.3  to assess differential expression for protein-coding genes and miRNAs. Then, we performed independent filtering on the results to remove genes that have no or little chance of showing significant evidence (Supplementary Table S3 and S4). Specifically, the independent filtering uses the mean of normalized counts as a threshold to optimize the number of adjusted p-values ≤ 0.05 . If the normalized expression of a gene was lower than the threshold, it was discarded. The Benjamini-Hochberg method was then used on the set of remaining genes to correct for multiple comparisons . Genes with adjusted p-values ≤ 0.05 were regarded as significantly differentially expressed.
We extracted all curated pathways from the Reactome pathway knowledgebase (release 68)  together with their hierarchical and biological classification according to the database developers. We retraced Reactome's pathway hierarchy by assigning every pathway from Homo sapiens to its corresponding root categories, such as signalling transduction and immune system (see Supplementary Materials for details). As a result, we obtained a table of Reactome pathways matched to the 26 root categories (Supplementary Figure S3 and Supplementary Table S5).
We applied a competitive gene set test to perform gene set enrichment analyses for Reactome pathways. The algorithm CAMERA  tests whether the genes in the set are lowly or highly ranked in terms of differential expression relative to genes not in the set, with a positive gene set score indicating a shared tendency for upregulation of the corresponding genes, and vice versa (see Supplementary Materials for details). All genes identified as differentially expressed from our RNA-seq data were used as the background gene list for the enrichment analysis. All obtained p-values were corrected using the Benjamini-Hochberg method. Pathways with false discovery rate (FDR) ≤ 0.05 were regarded as significantly up- (positive score) or down-regulated (negative score) in our comparison of caIKK-DCs with controls (Supplementary Table S5). The gene set enrichment analysis was performed using the CAMERA implementation in the package limma  in R.
A systems biology approach to study miRNA regulation in DCs. The study starts with preparing of donor DCs followed by cocktail maturation and subsequent electroporation of mRNA encoding Melan-A (control DCs) or a constitutively active variant of IKKβ (caIKK-DCs). The obtained RNA-seq data are processed and analysed for annotating and quantifying protein-coding genes and miRNAs. The identified differentially expressed genes in DCs are used for pathway enrichment analyses and reconstruction of gene regulatory networks. A network topology-oriented scoring model is employed to prioritize miRNAs in different pathway categories of DCs. Finally, a literature review of the top ranking miRNAs in immune signalling pathways elaborates their potential function for improving the immunogenic potency of caIKK-DCs.
We downloaded functional interactions from the Reactome database (release 68). The collection includes protein-protein interactions, transcriptional regulation, gene co-expression, protein domain interaction, gene ontology annotations and text-mined protein interactions, which cover almost half of the human proteome . There are different types of directional molecular interactions including: activation, inhibition or repression, and co-expression or complex formation. The biochemical reactions covered are phosphorylation and ubiquitination. We processed the list to transform bidirectional interactions into their two unidirectional constituents. The result list contained 435,043 unidirectional interactions among 13,852 protein-coding genes.
To derive miRNA-gene interactions, we first obtained conserved and non-conserved miRNA binding sites as predicted by Targetscan version 7.2 . Then, we filtered the predicted interactions with experimental evidence from miRTarBase version 8.0  and starBase version 2.0 . By doing so, we obtained a list of miRNA-gene interactions (Supplementary Table S6) that contain putative miRNA binding sites with experimental support, such as high-throughput experiments (e.g., RNA-seq, microarray, and AGO CLIP-sequencing) and/or low-throughput experiments (e.g., q-PCR, reporter assay, and Western blot).
The obtained lists of protein-coding genes' functional interactions and miRNA-gene interactions were used to reconstruct gene regulatory networks from Reactome pathways. Specifically, for a pathway of interest, we built a network around its participating genes and calculated the pairwise Pearson correlation coefficients between interaction partners from their normalized count values in caIKK-DCs. The normalized counts were obtained using the regularized logarithm method . In the reconstructed networks, we used the Pearson correlation coefficients to filter out interactions that disagree with their regulation type. We assumed that positive interactions (i.e., activation) require positive Pearson correlations and negative interactions (i.e., inhibition) negative Pearson correlations between interacting molecules. Interactions annotated as gene co-expression or formations of protein complexes were kept, assuming that the involved genes can affect each other's expression or activity in both directions.
Furthermore, we added annotation to the reconstructed networks' components in the form of differential expression profiles (i.e., fold-change and FDR), types of genes (e.g., protein-coding gene or miRNA), gene interaction types (e.g., functional interaction or post-transcriptional regulation), gene interaction strengths as denoted by the Pearson correlation coefficients introduced above, and immune categories of genes. Immune categories of protein-coding genes were annotated using curated data from the Immport database . Data from the TcoF database were used to identify TFs in our networks .
We prioritized genes in a network using SANTA in R . The algorithm determines a score of relative importance for each node in a network through a clustering model that accounts for network topology (distances between nodes) and node weights (in our case, a measure of differential expression called perturbation). Briefly, a gene is assigned a high score when itself and its close neighbours in the network have a higher-than-average node weight. The closeness, or distance, between genes is calculated by finding the shortest path through the network.
The node weight is given by the gene's perturbation (i.e., -log10(adj-p) • |log2(fold-change)|). Both adjusted p-value and log2 fold-change of the gene were taken from the differential gene expression analysis. The distances between neighbouring nodes were calculated as 1 - |p|, where p represents the Pearson correlation coefficient between the two interacting genes. Higher correlation coefficients (i.e., higher interaction strengths) correspond to lower edge length and thus shorter distance between the nodes. The calculated score was used to prioritize genes (see Supplementary Materials for details). As our networks contain both miRNAs and protein-coding genes that have different types of interactions, miRNAs and protein-coding genes were ranked separately.
To generate the curated DC network, we made use of a previously published network of macrophage pathways , since macrophages and DCs are generally considered to be quite similar. We manually added pathways for antigen processing and presentation that were not present in the macrophage map through a comprehensive database search in Reactome. The enriched DC network is a restructuring (see  for details on the algorithm) of the curated version that also incorporates miRNA-gene interactions identified in this study. The reconstructed network is accessible at https://vcells.net/dendritic-cell.
We used the TriplexRNA database  to identify miRNAs that can cooperate with significantly differentially expressed miRNAs in caIKK-DCs to repress protein-coding genes of interest. The obtained RNA triplexes were further filtered using pre-computed equilibrium concentrations and minimum free energies. We kept the RNA triplexes with equilibrium concentrations ≥ 50 nM and minimum free energies ≤ -25 kcal/mol and regarded the participating miRNAs as cooperative partners to repress protein-coding genes.
The gene regulatory networks for significantly enriched pathways were drawn using ggraph  and igraph  in R. Heat maps were plotted using Complexheatmap  in R. Scatter and bar plots were drawn using ggplot2  in R. Sankey diagram was drawn using networkD3  in R. The clustered Reactome pathways were visualized using Cytoscape version 3.72 .
To characterize the gene expression profiles induced by a constitutively active IKKβ in DCs, we collected blood samples from seven healthy donors and generated monocyte-derived DCs (see Materials and Methods). The DCs were matured with a cytokine cocktail of IL-1β, IL-6, TNF, and PGE2. The cmDCs were split for the subsequent experiments. One half was electroporated with mRNA encoding constitutively active IKKβ (caIKK, encoded by an engineered IKBKB). The other half was electroporated with mRNA encoding the melanoma antigen Melan-A (encoded by MLANA). The Melan-A protein was used because it is non-functional in DCs as it is naturally involved in melanocyte-specific pathways [57,58]. Additionally, we have previously shown that maturated DCs transfected with Melan-A mRNA did not show any significant changes in their transcriptomic profiles compared to mock-electroporated DCs . To confirm that caIKK electroporation of cmDCs induced activation in DCs, while electroporation of cmDCs with the control RNA did not, we analysed their phenotypes and cytokine secretion profiles. In caIKK-DCs, we verified the secretion of pro-inflammatory cytokines such as IL-6, IL-8, IL-12, and TNFα (Supplementary Figure S1A), upregulation of co-stimulatory surface markers such as CD25, CD40, CD70, CD80, CD86, and OX40L (Supplementary Figure S1B), and showed that the DC preparations used for RNA sequencing were successfully transfected as indicated by the expression of CD25 and CD70 in the seven donors (Supplementary Figure S1C).
Four hours after electroporation, RNA was isolated and assessed via bulk RNA sequencing (RNA-seq). We chose this early time point because we were interested in mRNA levels, which are expected to quickly respond to the activation of NF-κB as a result of continuous IKKβ activation. From the RNA-seq data, we identified 63 protein-coding genes and 44 miRNAs that were significantly differentially expressed (DE) between caIKK-DCs and controls (Supplementary Table S3 and S4; see Materials and Methods). Among the protein-coding genes, MLANA (encoding Melan-A) and IKBKB (encoding IKKβ) were the most down- and upregulated in caIKK-DCs, respectively (Supplementary Figure S2A). This is in line with the fact that the mRNA content of the two genes was artificially altered in the respective populations and can be considered a quality control for the experimental results. For the miRNAs, miR-146a/b and miR-155 were upregulated in caIKK-DCs, in consistence with them being transcriptional targets of NF-κB  and being upregulated in mature DCs . By performing principal component analysis, we assessed the clustering tendency in the RNA-seq data. Controls and caIKK-DCs showed better separation when restricting the input to the measured miRNAs rather than the whole transcriptome (Supplementary Figure S2B). In addition, the DE miRNAs unequivocally separated the caIKK-DCs from the controls in hierarchical clustering (Supplementary Figure S2C). These results suggested that caIKK-DCs harbour a distinct miRNA expression profile.
To understand the molecular function of the identified DE genes in the caIKK-DCs, we performed gene set enrichment analysis using the Reactome pathway database. The database contains more than 2,000 cellular pathways curated from 30,721 peer-reviewed publications and classified into 26 root categories , thereby enabling a systematic and comprehensive analysis of DE genes. The 26 categories consist of a set of pathways that are annotated to be hierarchically and functionally linked. We calculated enrichment scores for each pathway which reflect the degree to which its corresponding gene set tends to be up- or downregulated in caIKK-DCs (Figure 2A; see Materials and Methods).
Gene set enrichment analysis on the differential expression analysis between caIKK-enhanced and control DCs. (A) Enrichment plots of two exemplary gene sets. The genes' weighted log2 fold-change (i.e., log2 fold-change divided by the standard error of the log2 fold-change) obtained from the differentially gene expression analysis are sorted from the smallest to the largest (the barcode plot). Genes from TNF signalling (red bars) accumulate at the upregulated end while genes from ligand-receptor interaction (blue bars) in Hedgehog signalling accumulate at the downregulated end. The curves show the local enrichment score of the vertical bars in the barcode plot. For the red curve, parts above the dashed line signify enrichment while parts below the line signify depletion. For the blue one, parts below the dashed line signify enrichment while parts above the line signify depletion. (B) Network visualization of the gene set enrichment results for Reactome's 26 root categories (nodes annotated with texts). Categories are connected when they have shared genes (edges). The size of a node denotes the number of pathways that belongs to the respective category. The node colour represents the number of significantly enriched pathways in a category. A grey node means that no significant pathways were identified in the category. The width and colour of an edge represent the number of shared genes between the two connected categories. A detailed map containing all pathways and their corresponding enrichment scores can be found in Supplementary Figure S3.
We found that the DE genes in caIKK-DCs are significantly enriched in 195 Reactome pathways, most of which belonged to the Reactome categories signal transduction and immune system (Figure 2B; Supplementary Figure S3 and Table S5). In the category of immune system, 65 out of 182 pathways were identified as significantly enriched, including cytokine signalling pathways and pathways associated with innate and adaptive immune response. This suggested that the continuous activation of IKKβ in DCs has a generalized effect on DC-mediated immune responses. In addition, we identified 12 enriched pathways that are directly associated with NF-κB activation and signalling (see Supplementary Table S5), such as NF-κB activation by IκB kinase complexes. This was consistent with the current model of the canonical NF-κB activation pathway, in which IKKs phosphorylate IκB resulting in desequestration of NF-κB and its translocation into the nucleus, where it regulates expression of immune-related cytokine genes and others [11,12]. All these NF-κB pathways had positive enrichment scores, indicating that the genes involved are more likely to be up-regulated in the caIKK-DCs. The results were in line with our expectation that the caIKK-DCs can trigger a stronger immune response as a result of NF-κB desequestration by constitutive activation of IKKβ.
To identify potential miRNA-gene interactions regulating the immunogenic potency of DCs, we first obtained putative miRNA-gene interactions for the significantly DE miRNAs (see Materials and Methods). We kept the putative interactions that are validated by experiments. For each identified miRNA-gene interaction, we then computed the Pearson correlation coefficient between miRNA and target gene expression. The interactions with negative correlation were regarded as reliable and functional, as miRNAs canonically repress translation initiation or stimulate mRNA degradation  and miRNA-mediated gene activation usually results from indirect regulation mechanisms . The data showed that 36 out of the 44 miRNAs are involved in the regulation of protein-coding genes belonging to 195 enriched pathways of the 26 Reactome root categories (Figure 3A).
In individual pathway categories, the number of molecules (i.e., protein-coding genes, DNA/RNA, drugs, and chemical compounds) ranged from 27 to 2727, and genes identified by our RNA-seq data covered between 51% and 95% of the molecules in the respective category. For all categories except digestion and absorption, the DE miRNAs were identified to regulate 1 to 103 protein-coding genes (denoted by the numbers on the heat map grid cells in Figure 3A). Some miRNAs were found to regulate the expression of dozens of protein-coding genes in more than ten categories, suggesting that they act as regulatory hubs in caIKK-DCs. For example, miR-15a-5p, miR-16-5p, miR-20a-5p, and miR-424-5p can potentially regulate more than 80 protein-coding genes in the category signal transduction (Figure 3B), and they also have more than 60 targets in immune system. In contrast, miR-15a-3p and miR-9-3p target only BCL2 and MAPK1 in immune system, suggesting a specific role for them in regulating cytokine signalling in caIKK-DCs (Supplementary Figure S4). Furthermore, we found that some miRNAs target a high fraction of enriched pathways belonging to specific Reactome root categories (denoted by the colour of heat map grid cells in Figure 3A). For instance, miR-16-5p regulates 61 out of 64 enriched immune system pathways and 34 out of 36 enriched pathways in signal transduction. This suggested that it plays a vital role in regulating the immunogenic potency of caIKK-DCs. On the other hand, some categories included abundant enriched pathways that are regulated by multiple DE miRNAs. Interesting cases were pathways associated with protein metabolism, RNA metabolism, programmed cell death, and cell cycle. This result suggests that the DE miRNAs in caIKK-DCs are involved in regulating synthesis, processing, and modification of mRNAs and proteins and can also participate in other biological processes, such as cell cycle and cell apoptosis. Taken together, the DE miRNAs in caIKK-DCs target and potentially coordinate the activity of immune-relevant pathways in a pleiotropic fashion.
The ubiquitous, pleiotropic, and concerted gene regulation by miRNAs makes it challenging to quantify the relative impact of each individual miRNA. To prioritize the DE miRNAs according to their potential to act synergistically with NF-κB in DC activation, we applied a network-based method that integrates their expression and interaction profiles.
First, we reconstructed one gene regulatory network for each of the 26 pathway categories. The reconstructed networks were composed of miRNA-gene interactions and functional interactions among protein-coding genes. Interactions were discarded when the sign of their Pearson correlation coefficient of expression disagreed with their regulation type, such as inhibition or activation (see Materials and Methods). Depending on the category, the size of the corresponding networks varied from 1,915 genes and 57,520 interactions (for signal transduction) to 30 genes and 153 interactions (for mitophagy). To prioritize the network components involved in regulating the immunogenic potency of DCs, we used a clustering model  to calculate a node score (Figure 4A; see Materials and Methods).
The score ranked IKBKB, whose expression was greatly increased by mRNA electroporation, as the top protein-coding gene in 51 out of 59 networks in which it is involved (Supplementary Table S7; Supplementary Figure S5). This result is consistent with our expectation that the intentionally modulated gene in experiments is prioritized, and thus demonstrating the ability of the model to identify crucial regulatory genes in the experiments. In the two prominent categories signal transduction and immune system, the NF-κB family and genes related with immune signalling or antigen processing and presentation tended to rank higher than other genes (Supplementary Figure S6). This result again justified the ability of the model to prioritize important genes in networks, as members of the NF-κB family are downstream targets of IKBKB while signalling and antigen presentation genes are supposed to be crucial regulating immune function of DCs.
Furthermore, we analysed the data to identify crucial miRNAs for each Reactome root category. As shown in Figure 4B, miRNAs with higher node weights (i.e., stronger perturbation) generally ranked higher in a category, as miRNA scores and node weights showed a positive correlation, ranging from 0.19 to 0.98. Specifically, miR-503-5p, miR-503-3p, and miR-146-5p had the highest perturbation in the DE miRNAs, and they ranked top in 22 out of 26 categories. However, the interaction profile also plays a role, as for example in signal transduction, the three top-ranking miRNAs miR-101-3p, miR-16-5p, and miR-15a-5p had lower perturbation than miR-146-5p but interacted with more protein-coding genes. In addition, the three above miRNAs and miR-144-3p ranked top in immune system, most probably due to the reason that they regulate a large number of protein-coding genes associated with immune signalling pathways. To facilitate the visualization of our results, we integrated the data and the identified miRNA-gene interactions into a comprehensive, manually curated regulatory network including key pathways in DC priming and activation according to the literature (see Materials and Methods; https://vcells.net/dendritic-cell).
miRNA targeting profiles in Reactome pathway categories. (A) Overview of the 36 DE miRNAs (in columns) targeting the 26 Reactome pathway categories (in rows). On the heat map grid, the number of protein-coding genes targeted by the respective miRNA is given, and the colour represents the category's number of significantly enriched pathways that are regulated by the miRNA. For example, if an entry shows n with a colour corresponding to m on the figure legend, it means a miRNA regulates n targets in m significantly enriched pathways of a category. The bar plot on the left indicates the per-category percentage of the protein-coding genes (black bars) that were found in our RNA-seq data, with the total number of molecules in the category given next to it. The bar plot on the right indicates the percentage of enriched pathways (black bars) per category, with the total number of pathways in the category given next to it. The figures at the bottom tabulate how many categories a miRNA regulates. The top annotation shows statistics from the differential expression analysis for the miRNAs (i.e., fold-change in log2 scale and FDR). (B) The network shows miRNA-gene interactions in the category signal transduction. The four miRNAs (miR-16, miR-15a, miR-20a, and miR-424) that have the largest number of targets were selected. The node size is proportional to the node degree. The node colour represents a gene's fold-change in log2 scale. The node shape denotes the type of a gene, including protein-coding (square) and miRNA (circle), with their names shown in blue and black labels, respectively. TFs are drawn as diamonds with their names shown in purple font. The colour of node borders represents different categories of annotated immune genes, with the gene names given as labels. The edge colour shows Pearson correlation coefficients between the expression of miRNAs and that of their targets.
Ranking of miRNA relevance based on expression and regulatory network neighbourhood. (A) Network nodes were scored with an algorithm that uses the guilt-by-association principle to rank genes. In other words, a gene inside of or close to a cluster of important genes is potentially more important than a gene that is further away. In our case, the importance of a gene for a phenotype is quantified by their perturbation in expression (denoted by node colours). In a network, the distance of the gene in question (red or blue border) to other genes is calculated via the weighted shortest-path method. The range of observed distances is then subdivided into discrete bins (denoted by circles in the figure) and an estimate of neighbourhood importance calculated for each bin, i.e., for all nodes up to the respective distance. The area under the curve (AUC) in the plot of bins (d) vs neighbourhood importance is used as the gene's score. In the example, gene k is much closer to genes with high node weights (i.e., their perturbation in caIKK-DCs is large) than gene j. As a result, the blue area is bigger than the red area, and thus gene k ranks higher than j. (B) Heat map of miRNA ranking in pathway categories. The columns of the matrix indicate the 36 DE miRNAs sorted by perturbation (i.e., node weights), and the rows of the matrix indicate 25 categories of Reactome pathways. The category digestion and absorption is not shown, as we did not identify functional and reliable interactions among its genes. On the heat map grid, the rank of a miRNA in a category is given as a number, and the colour represents the number of protein-coding genes targeted by it. For example, if an entry shows k with a colour corresponding to n on the figure legend, it means that the miRNA ranks kth (1st is the highest ranking) and regulates n targets in the category. A white grid cell means that the miRNA has no targets in the category and thus no ranking. The top annotation shows node weights of the miRNAs. The numbers in parentheses on the left side list how many genes and edges the reconstructed regulatory network of the category possessed. The box plots on the left show the distribution of edge weights (denoted by Pearson correlation coefficients between genes) in the networks. The bar plots on the right show the Pearson correlation coefficients between a miRNA's perturbation and its score. The numbers at the bottom show the number of times miRNAs ranked 1st in the pathway categories.
Taken together, the reconstructed regulatory networks underlying different cell functions allowed us to identify important miRNA regulators based on their expression and interaction profiles. The miRNAs with the highest scores possibly exert regulatory functions, and manipulation of their expression levels may enhance the immunogenic potency of DCs.
To characterize the functional role that miRNAs play in caIKK-DCs, we delineated landscapes of miRNA-gene interactions in the significantly enriched pathways that were found in corresponding categories (see Figure 5 and https://www.synmirapy.net/DC-optimization). The interaction landscapes are a way of systematically mapping relevant gene interactions, and in our case they served as a tool for identifying functional miRNA-gene interactions in DC priming and activation. As we were particularly interested in identifying miRNAs that can enhance the caIKK-DCs' immunogenic potency, we focused on analysing miRNA-mediated gene regulation in the category immune system. In this category, we identified hundreds of miRNA-gene interactions in significantly enriched pathways, including toll-like receptor, cytokine, and interleukin signalling as well as MHC processing and presentation. All of these pathways had positive enrichment scores, indicating that the involved genes tended to be upregulated in caIKK-DCs according to our analysis. Most protein-coding genes used as indicators of DC activation and maturation [5,60,63,64] were found to be upregulated in the enriched pathways (Supplementary Table S8). The activation of NF-κB signalling led to upregulation of surface proteins that can prime T cells (e.g., CD40, CD70, CD80, and CD86), chemokines (e.g., CCL3 and CXCL10) that are necessary for T-cell migration, TNF superfamily members that can induce crosstalk between T cells and DCs (e.g., TNF, TNFRSF4, and TNFSF9), and cytokines that are responsible for stimulating proliferation and activation of T cells (e.g., CXCL8, IL6, IL12A, and IL12B).
Landscape of miRNA-mediated DC gene regulation in immune signalling pathways. The heat map has two components that share a set of columns corresponding to 98 DE protein-coding genes that are targeted by the DE miRNAs. In the upper component, rows represent pathways from the category immune system, and grid cell colours indicate whether a protein-coding gene is involved in the pathway (grey), involved in the pathway and an immune gene (grey grid cells with red borders), or not involved in the pathway (white). The figures next to a pathway name indicate how many DE immune genes (left) and how many protein-coding genes found in our RNA-seq data (right) belong to it. The top annotation highlights genes with different characteristic immune function using a colour code. The annotation on the right side shows the statistics of the gene set enrichment analysis including the enrichment score and the FDR. The bar plot between the heat map components shows the log2 fold-change of the genes in caIKK-DCs (blue: downregulated; red: upregulated). In the lower component, the rows represent the ranking miRNAs in immune system (from high to low) and the grid cells show the regulative influence of a protein-coding gene by a miRNA, which is estimated by the Pearson correlation coefficients between their expression profiles. If a gene is a known immune gene, the corresponding grid cell has a red border. The numbers in the parentheses next to the miRNA names show the number of DE immune genes and the number of DE protein-coding genes that are regulated by a miRNA. The right annotation shows the results of the differential expression analysis including the log2 fold-change of miRNA expressions and their FDRs. For lack of space, we show only enriched pathways with more than 30 protein-coding genes picked up in the RNA-seq data, and in each pathway, only a subset of protein-coding genes that are estimated to be strongly influenced by the miRNAs (Pearson correlation ≤ -0.3) are shown. The complete landscape of miRNA-gene interactions in immune system is shown in Supplementary Figure S7.
Potential miRNAs for improving immunogenic potency of caIKK-DCs. The Sankey diagram contains three columns made up of nodes representing the DE miRNAs from our RNA-seq data and their cooperating miRNAs, protein-coding genes targeted by the miRNAs, and DC phenotypes associated with the protein-coding genes, respectively. miRNA pairs that were identified to cooperatively repress a protein-coding gene are connected by a brace. Colours of miRNAs and protein-coding genes indicate whether or not they were significantly DE in caIKK-DCs. Arrows in miRNA nodes indicate how the expression of miRNAs should be manipulated to obtain DCs with higher immunogenic potency: upregulation (↑), downregulation (↓), and no suggestion due to potentially conflicting effects of the miRNA (-). Connections between miRNAs and protein-coding genes show regulative influence of protein-coding genes by miRNAs (strong: Pearson correlation ≤ -0.5; weak: -0.3 ≤ Pearson correlation < -0.5). The connections between protein-coding genes and phenotypes denote how a gene regulates a phenotype according to literature. For instance, miR-424-3p and miR-224-5p target IRF4 that is known to positively regulate differentiation of DCs. The two miRNAs cooperatively repress the protein-coding gene, but the observed downregulation of miR-424-3p results in a decreased inhibitory effect on the expression of IRF4. A detailed discussion of the results can be found in the main text. The corresponding miRNA-gene interactions in immune system as well as annotated gene-phenotype associations can be found in Supplementary Table S9 and Table S10, respectively.
Furthermore, our data showed that the identified DE miRNAs have a regulative influence (represented by Pearson correlation ≤ -0.3) on protein-coding genes associated with NF-κB activation, cytokines, chemokines, and TFs that are associated with immunophenotypes of DCs (Figure 6). Some of the DE miRNAs were found to cooperate with other miRNAs to regulate the expression of a protein-coding gene (see Materials and Methods). This mechanism, known as miRNA cooperativity, is characterized by more efficient inhibitory effects on the target's expression compared to the regulation by individual miRNAs [23,25,26]. Moreover, for most of the identified miRNAs, our analysis proposed specific modulation of their expression levels to improve immunogenic potency of DCs. However, in some cases, such as miR- 34a-5p and miR-20a-5p, up- or downregulating their expression levels may result in contradictory effects in DC-mediated immune response, thereby requiring further analysis and experimental investigation. In the following paragraphs, we illustrate and discuss specific functions of the miRNAs in caIKK-DCs.
Chemokines direct cell migration via induction of chemotaxis. For example, CCL5 and CXCL10 improve CD8+ T-cell infiltration  and CCL20 plays a role in recruiting regulatory T cells and T helper (Th) 17 cells [66,67]. CXCL10 is repressed by miR-16-5p and CCL20 is targeted by miR-21a-5p with its cooperating miR-25-3p. This suggested that the miRNAs exert an inhibitory function in the recruitment of T cells. miR-21a-5p also targets IL12A, a subunit of the inflammatory cytokine IL-12 that is necessary for CD8+ T-cell clonal expansion, function and memory [63,68].
Control of DC survival is necessary for maintaining their homeostasis [69,70]. We showed that miR-15a-5p with its cooperating miRNAs (i.e., miR-156-5p and miR-876-5p) and miR-20a-5p target BCL2, and miR-101-3p targets MCL1. The repression of the BCL2 family of anti-apoptosis genes by these miRNAs suggested their ability to undermine the survival mechanism of DCs.
miR-424-5p and miR-224-5p can co-repress IRF4, which is a member of the interferon-regulatory family and can regulate differentiation of specific DCs that can induce Th 2 cell responses . miR-20a-5p and miR-144-3p regulate the MAPK signalling pathway by targeting MAPK1 and MAP3K8 respectively, and these MAP kinases have been found to activate the IKK complex that triggers NF-κB activation [72,73] and also to regulate release of TNFα by DCs . miR-34a-5p has a strong regulative influence on CD44 whose presence is important for the immune synapse between DCs and T cells that subsequently regulates T-cell activation  and apoptosis . miR-34a-5p also targets TNFAIP8 whose knockdown in DCs has been found to promote DC maturation and activation followed by increased proliferation and differentiation of T cells . miR-9-5p can cooperate with miR-139-5p to repress CXCR4 that is required for DC migration into the skin's draining lymph nodes . miR-142-3p with its cooperating miR-429 and miR-142-5p target the small GTP-binding protein RAC1 that controls the formation of dendrites in mature DCs and their migration toward T cells .
Some identified DE miRNAs target protein-coding genes involved in regulating the DC-mediated secretion of cytokines that are important for the T-cell response. The repression of NFATC3 by miR-424-5p and its cooperating miR-370-3p suggested a regulating influence on the production of IL-2 that is involved in T-cell priming [80,81]. STX3 has been shown to play a role in trafficking of IL-6 or MIP-1α in DCs and thus regulating their secretion  and is targeted by let-7e-5p, miR-146a-5p, and miR-146b-5p with its cooperative partner miR-519d-3p. The deficiency of IRAK1 in plasmacytoid DCs abrogates IFNα production, leading to a remodulation of T cell function [83-85], and IRAK1 is a target of miR-142-3p.
Finally, miR-16-5p and miR-15a-5p can cooperate with miR-203a-3p to repress IL-15, an interleukin which can induce T-cell proliferation, enhance cytolytic effector cells including natural killer and cytotoxic T cells, and reinforce B-cell stimulation . A recent in vivo study has shown that an IL15-enhanced DC vaccine is a potent delayer of tumour growth, improves mouse survival, and induces a stronger Th1-skewed T-cell response . The two miRNAs also target the receptor TNFSF9 (also known as CD137), whose stimulation in DCs by its ligand CD137L can lead to secretion of IL-6 and IL-12 and induce T-cell proliferation [63,88,89]. In addition, miR-16-5p and miR-15a-5p were identified to strongly repress IKBKB itself. Since both miRNAs were found to be downregulated in caIKK-DCs, this implied a positive feedback loop in NF-κB signalling as the miR-15/16 cluster is a transcriptional target of NFKB1 . The results suggested both miRNAs as promising candidates for improving the immunogenic potency of caIKK-DCs, as they not only have the ability to strengthen NF-κB activation but also to improve DC-induced immune responses through regulating cytokines and chemokines.
We applied a systems biology approach to investigate the regulatory functions of miRNAs in caIKK-DCs. Due to the promiscuous binding of miRNAs, it is challenging to identify relevant miRNA-gene interactions for experimental validation and cell re-engineering [91,92]. Our approach, which integrates transcriptomic profiling, networks of curated signalling pathways, and a prioritisation score, allows the systematic identification of condition-specific miRNA-gene interactions.
Through RNA sequencing of monocyte-derived DCs matured with a cytokine cocktail and electroporated with caIKK mRNA, we identified DE protein-coding genes and miRNAs in the caIKK-DCs. The identified DE miRNAs correctly separated the caIKK-DCs from the control, suggesting a well-defined transcriptional response to caIKK that is consistent with our understanding that miRNAs act as post-transcriptional regulators of expression in DC differentiation and function . Among the identified miRNAs there were several, such as miR-15a-5p, miR-16-5p, miR-20a-5p, and miR-424-5p, which target a considerable number of genes. Such hubs have been shown to be important regulators, as they represent sites of signalling convergence in gene regulatory networks and coordinate cell development and function [94-96]. In contrast, other DE miRNAs, such as miR-15a-3p and miR-9-3p, exert a narrow function by regulating the expression of specific protein-coding genes in the caIKK-DCs.
Integration of the transcriptomic response into the curated pathways from Reactome provides an understanding of the functional changes at the pathway level. The gene set enrichment analysis highlighted cytokine, interleukin, and toll-like receptor signalling pathways that are involved in regulating various aspects of innate and adaptive immune responses . Such results may be compromised when other pathway databases such as KEGG  and WikiPathways  are employed, as the relevant pathways and molecular interactions in the pathways are different from Reactome [100,101]. To circumvent this issue, one possibility could be to extract the overlapped networks between the different databases of pathways; however, this is not always possible due to the differences in annotation of genes and interactions. An alternative option is to integrate the data and the detected miRNA-gene interactions into comprehensive, manually curated regulatory networks based on the current literature on DC regulation. This way, one can put the newly discovered relevant interactions into the context of the existing knowledge and facilitate the mining and interpretation of the omics data [102,103]. However, when used inappropriately, knowledge-based networks mainly rediscover existing knowledge but may overlook insights gained from the evaluation of an all-encompassing network.
We reconstructed regulatory networks from Reactome pathways and used them to rank genes and miRNAs according to their predicted impact on DC function. Systematic computation of such a ranking supports and facilitates experimental efforts, allowing them to focus on the most promising candidates. Gene prioritization algorithms have been widely used in recent times to rank genes in networks [104,105]. For instance, the PageRank algorithm designed to analyse the relative importance of websites was adapted to identify crucial genes in biological networks , and diffusion-based methods were used on dense networks to prioritize genes . We used a gene prioritization algorithm that utilizes the guilt-by-association principle to rank genes based on their own perturbation, i.e., differential expression profile, and their weighted distances to other perturbed genes in a network. The algorithm prioritized dozens of miRNAs, of which miR-16-5p and miR-15a-5p are the top candidates to regulate the immunogenic potency of DCs.
Finally, an in-depth analysis of the identified miRNA-gene interactions in immune signalling pathways showed diverse roles of the DE miRNAs in regulating DC-mediated immune response. For instance, miR-16-5p and miR-15a-5p may have strong regulatory influence on IKBKB that activates NF-κB and on TNFSF9 that controls cytokine secretion of DCs; both miRNAs could have weak inhibitory effects on BCL2 that maintains DC homeostasis and on cytokines (such as CXCL10 and IL15) that regulates T-cell response but may cooperate with other miRNAs to more efficiently repress the protein-coding genes. While these predictions were made based on validated miRNA binding sites and negative correlation between the expression levels of miRNAs and their targets, they cannot quantify the strength of individual repression effects . The results suggested both miRNAs as potential candidates for improving immunogenic potency of caIKK-DCs through strengthening NF-κB stimulation and also synergistically regulating other genes related with immunogenic potency.
For most identified crucial miRNAs, our analysis suggested up- or downregulation of their expression levels to improve immunogenic potency of DCs, but in some cases the pleiotropic nature of miRNAs in regulating gene expression makes it difficult to decide how to experimentally modulate their expression. In addition, it is worth noting that the results reflect the early transcriptional response that may differ from that in the long-term. From an experimental perspective, the next step would be to analyse the kinetics of the expressions of miRNA and mRNA after the activation of NF-κB. Further, it remains to be tested how co-electroporating the selected miRNAs, or artificial antagonists thereof, with caIKK or introducing them into the cells after a delay will alter the DCs' phenotype and immunogenic potency.
Taken together, our approach enables the systematic analysis and identification of functional miRNA-gene interactions that can be experimentally tested for improving DC immunogenic potency. Since the results were produced in a computationally reproducible manner and were stored in a public database, experimental tests of the predictions can be performed in the future not only by our group but also by other researchers working on DC-based cancer immunotherapy. Additionally, since the approach is not specific for DCs, it can be adapted to study miRNAs in other immune cells and relevant immunotherapies.
DC: Dendritic cell; cmDCs: cocktail-matured dendritic cells; miRNA: microRNA; caIKK: constitutively active IKKβ; TF: transcription factor; RNA-seq: RNA sequencing; Melan-A: melanoma antigen recognized by T cells 1; DE: differentially expressed; FDR: false discovery rate; AGO-CLIP: argonaute-crosslinking and immunoprecipitation; q-PCR: quantitative polymerase chain reaction; Th: T helper.
We would like to thank Mr. Fitsumbirhan Mehari for his assistance in the project. We would like to thank the following funding agencies for supporting our research: ELAN program of Universitätsklinikum Erlangen (16-08-16-1-Lai to X.L.); German Federal Ministry of Education and Research (e:Bio-MelEVIR 031L0073A and e:Med-MelAutim 01ZX1905A to J.V.); STAEDTLER Stiftung (ww/eh 30/16 to J.V.); Manfred-Roth Stiftung (to J.V.). We also acknowledge support by Deutsche Forschungsgemeinschaft and Friedrich-Alexander-Universität Erlangen-Nürnberg within the funding program Open Access Publishing.
RNA-seq data have been deposited in the ArrayExpress database under accession number E-MTAB-9666.
The database for selecting functional miRNA-gene interactions with therapeutic potential in DCs for experimental investigation is available at https://www.synmirapy.net/DC-optimization.
The curated DC network containing identified miRNA-gene interactions and the RNA-seq data are available at https://vcells.net/dendritic-cell.
Conceptualization: XL, JV; Data curation: XL, MC, ME, TJ, SU, CL; Formal analysis: XL; Funding acquisition: XL, JV; Investigation: XL, FD, KFG, SU, AE, JD, NS; Methodology: XL; Project administration: XL; Resources: XL, FD, AE, JW, HMJ, JD, NS, JV; Software: XL, MC, TJ; Supervision: XL, JV; Validation: XL, ME; Visualization: XL; Writing-original draft: XL, FD, MC, ME, JV; Writing-review & editing: XL, ME, JD, NS, JV. All authors have read and agreed to the published version of the manuscript.
The authors have declared that no competing interest exists.
1. Banchereau J, Steinman RM. Dendritic cells and the control of immunity. Nature. 1998;392:245-52
2. Janikashvili N, Larmonier N, Katsanis E. Personalized dendritic cell-based tumor immunotherapy. Immunotherapy. 2010;2:57-68
3. Kalinski P, Muthuswamy R, Urban J. Dendritic cells in cancer immunetherapy: vaccines and combination immunotherapies. Expert Rev Vaccines. 2013;12:285-95
4. Perez CR, De Palma M. Engineering dendritic cell vaccines to improve cancer immunotherapy. Nat Commun. 2019;10:5408
5. Dörrie J, Schaft N, Schuler G, Schuler-Thurner B. Therapeutic Cancer Vaccination with Ex vivo RNA-Transfected Dendritic Cells-An Update. Pharmaceutics. 2020 12
6. Steinman RM, Hemmi H. Dendritic cells: translating innate to adaptive immunity. Curr Top Microbiol Immunol. 2006;311:17-58
7. Constantino J, Gomes C, Falcão A, Neves BM, Cruz MT. Dendritic cell-based immunotherapy: a basic review and recent advances. Immunol Res. 2017;65:798-810
8. Tacken PJ, de Vries IJM, Torensma R, Figdor CG. Dendritic-cell immunotherapy: from ex vivo loading to in vivo targeting. Nat Rev Immunol. 2007;7:790-802
9. Hayden MS, Ghosh S. NF-κB in immunobiology. Cell Res. 2011;21:223-44
10. Pfeiffer IA, Hoyer S, Gerer KF. et al. Triggering of NF-κB in cytokine-matured human DCs generates superior DCs for T-cell priming in cancer immunetherapy. Eur J Immunol. 2014;44:3413-28
11. Liu T, Zhang L, Joo D, Sun S-C. NF-κB signaling in inflammation. Signal Transduct Target Ther. 2017 2
12. Taniguchi K, Karin M. NF-κB, inflammation, immunity and cancer: coming of age. Nat Rev Immunol. 2018;18:309-24
13. Gerer KF, Erdmann M, Hadrup SR. et al. Preclinical evaluation of NF-κB-triggered dendritic cells expressing the viral oncogenic driver of Merkel cell carcinoma for therapeutic vaccination. Ther Adv Med Oncol. 2017;9:451-64
14. Bosch NC, Voll RE, Voskens CJ. et al. NF-κB activation triggers NK-cell stimulation by monocyte-derived dendritic cells. Ther Adv Med Oncol. 2019;11:1758835919891622
15. Bartel DP. MicroRNAs: target recognition and regulatory functions. Cell. 2009;136:215-33
16. Bartel DP. Metazoan MicroRNAs. Cell. 2018;173:20-51
17. Wang J, Ye H, Zhang D. et al. Cancer-derived Circulating MicroRNAs Promote Tumor Angiogenesis by Entering Dendritic Cells to Degrade Highly Complementary MicroRNAs. Theranostics. 2017;7:1407-21
18. Turner ML, Schnorfeil FM, Brocker T. MicroRNAs regulate dendritic cell differentiation and function. J Immunol. 2011;187:3911-7
19. Zhou H, Wu L. The development and function of dendritic cell populations and their regulation by miRNAs. Protein Cell. 2017;8:501-13
20. Inui M, Martello G, Piccolo S. MicroRNA control of signal transduction. Nat Rev Mol Cell Biol. 2010;11:252-63
21. Smyth LA, Boardman DA, Tung SL, Lechler R, Lombardi G. MicroRNAs affect dendritic cell function and phenotype. Immunology. 2015;144:197-205
22. Lu C, Huang X, Zhang X. et al. miR-221 and miR-155 regulate human dendritic cell development, apoptosis, and IL-12 production through targeting of p27kip1, KPC1, and SOCS-1. Blood. 2011;117:4293-303
23. Lai X, Schmitz U, Gupta SK. et al. Computational analysis of target hub gene repression regulated by multiple and cooperative miRNAs. Nucleic Acids Res. 2012;40:8818-34
24. Bertoli G, Cava C, Castiglioni I. MicroRNAs: New Biomarkers for Diagnosis, Prognosis, Therapy Prediction and Therapeutic Tools for Breast Cancer. Theranostics. 2015;5:1122-43
25. Lai X, Gupta SK, Schmitz U. et al. MiR-205-5p and miR-342-3p cooperate in the repression of the E2F1 transcription factor in the context of anticancer chemotherapy resistance. Theranostics. 2018;8:1106-20
26. Lai X, Eberhardt M, Schmitz U, Vera J. Systems biology-based investigation of cooperating microRNAs as monotherapy or adjuvant therapy in cancer. Nucleic Acids Res. 2019;47:7753-66
27. Hao Shi null, Yan K-K, Ding L, Qian C, Chi H, Yu J. Network Approaches for Dissecting the Immune System. iScience. 2020;23:101354
28. Szeto GL, Finley SD. Integrative Approaches to Cancer Immunotherapy. Trends Cancer. 2019;5:400-10
29. Schaft N, Dörrie J, Thumann P. et al. Generation of an optimized polyvalent monocyte-derived dendritic cell vaccine by transfecting defined RNAs after rather than before maturation. Journal of Immunology. 2005;174:3087-97
30. Hoyer S, Gerer KF, Pfeiffer IA. et al. Electroporated Antigen-Encoding mRNA is not a Danger Signal to Human Mature Monocyte-Derived Dendritic Cells. J Immunol Res. 2015;2015:952184
31. Andrews S. Babraham Bioinformatics - FastQC A Quality Control tool for High Throughput Sequence Data. 2017. http://www.bioinformatics.babraham.ac.uk/projects/fastqc/.
32. Li H. Aligning sequence reads, clone sequences and assembly contigs with BWA-MEM. arXiv:13033997. 2013 http://arxiv.org/abs/1303.3997
33. Kozomara A, Griffiths-Jones S. miRBase: annotating high confidence microRNAs using deep sequencing data. Nucleic Acids Research. 2014;42:D68-73
34. Bourgon R, Gentleman R, Huber W. Independent filtering increases detection power for high-throughput experiments. Proc Natl Acad Sci USA. 2010;107:9546-51
35. Tam S, Tsao M-S, McPherson JD. Optimization of miRNA-seq data preprocessing. Brief Bioinformatics. 2015;16:950-63
36. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15:550
37. R Core Team. R: The R Project for Statistical Computing. 2017. https://www.r-project.org.
38. Benjamini Y, Hochberg Y. Controlling the False Discovery Rate: A Practical and Powerful Approach to Multiple Testing. Journal of the Royal Statistical Society: Series B (Methodological). 1995;57:289-300
39. Fabregat A, Jupe S, Matthews L. et al. The Reactome Pathway Knowledgebase. Nucleic Acids Res. 2018;46:D649-55
40. Wu D, Smyth GK. Camera: a competitive gene set test accounting for inter-gene correlation. Nucleic Acids Res. 2012;40:e133-e133
41. Ritchie ME, Phipson B, Wu D. et al. limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res. 2015;43:e47
42. Wu G, Feng X, Stein L. A human functional protein interaction network and its application to cancer data analysis. Genome Biol. 2010;11:R53
43. Agarwal V, Bell GW, Nam J-W, Bartel DP. Predicting effective microRNA target sites in mammalian mRNAs. eLife. 2015;4:e05005
44. Chou C-H, Shrestha S, Yang C-D. et al. miRTarBase update 2018: a resource for experimentally validated microRNA-target interactions. Nucleic Acids Res. 2018;46:D296-302
45. Li J-H, Liu S, Zhou H, Qu L-H, Yang J-H. starBase v2.0: decoding miRNA-ceRNA, miRNA-ncRNA and protein-RNA interaction networks from large-scale CLIP-Seq data. Nucleic Acids Res. 2014;42:D92-97
46. Bhattacharya S, Dunn P, Thomas CG. et al. ImmPort, toward repurposing of open access immunological assay data for translational and clinical research. Sci Data. 2018;5:180015
47. Schmeier S, Alam T, Essack M, Bajic VB. TcoF-DB v2: update of the database of human and mouse transcription co-factors and transcription factor interactions. Nucleic Acids Res. 2017;45:D145-50
48. Cornish AJ, Markowetz F. SANTA: quantifying the functional content of molecular networks. PLoS Comput Biol. 2014;10:e1003808
49. Wentker P, Eberhardt M, Dreyer FS. et al. An Interactive Macrophage Signal Transduction Map Facilitates Comparative Analyses of High-Throughput Data. J Immunol. 2017;198:2191-201
50. Schmitz U, Lai X, Winter F, Wolkenhauer O, Vera J, Gupta SK. Cooperative gene regulation by microRNA pairs and their identification using a computational workflow. Nucleic Acids Res. 2014;42:7539-52
51. Pedersen TL. ggraph: An Implementation of Grammar of Graphics for Graphs and Networks. 2020. https://CRAN.R-project.org/package=ggraph.
52. Csardi G, Nepusz T. igraph - Network analysis software. 2015. https://igraph.org/r/.
53. Gu Z, Eils R, Schlesner M. Complex heatmaps reveal patterns and correlations in multidimensional genomic data. Bioinformatics. 2016;32:2847-9
54. Wickham H. ggplot2: Elegant Graphics for Data Analysis. 2nd ed. Springer International Publishing. 2016 https://www.springer.com/gp/book/9783319242750
56. Shannon P, Markiel A, Ozier O. et al. Cytoscape: a software environment for integrated models of biomolecular interaction networks. Genome Res. 2003;13:2498-504
57. Ceppi M, Ruggli N, Tache V, Gerber H, McCullough KC, Summerfield A. Double-stranded secondary structures on mRNA induce type I interferon (IFN alpha/beta) production and maturation of mRNA-transfected monocyte-derived dendritic cells. J Gene Med. 2005;7:452-65
58. Hoashi T, Watabe H, Muller J, Yamaguchi Y, Vieira WD, Hearing VJ. MART-1 is required for the function of the melanosomal matrix protein PMEL17/GP100 and the maturation of melanosomes. J Biol Chem. 2005;280:14006-16
59. Ma X, Becker Buscaglia LE, Barker JR, Li Y. MicroRNAs in NF-kappaB signaling. J Mol Cell Biol. 2011;3:159-66
60. Jin P, Han TH, Ren J. et al. Molecular signatures of maturing dendritic cells: implications for testing the quality of dendritic cell therapies. J Transl Med. 2010;8:4
61. Gebert LFR, MacRae IJ. Regulation of microRNA function in animals. Nat Rev Mol Cell Biol. 2019;20:21-37
62. Tan H, Huang S, Zhang Z, Qian X, Sun P, Zhou X. Pan-cancer analysis on microRNA-associated gene activation. EBioMedicine. 2019;43:82-97
63. Summers deLuca L, Gommerman JL. Fine-tuning of dendritic cell biology by the TNF superfamily. Nat Rev Immunol. 2012;12:339-51
64. Locy H, Melhaoui S, Maenhout SK, KrisThielemans. Dendritic Cells: The Tools for Cancer Treatment. In. 2018
65. Ji R-R, Chasalow SD, Wang L. et al. An immune-active tumor microenvironment favors clinical response to ipilimumab. Cancer Immunol Immunother. 2012;61:1019-31
66. Comerford I, Bunting M, Fenix K. et al. An immune paradox: how can the same chemokine axis regulate both immune tolerance and activation?: CCR6/CCL20: a chemokine axis balancing immunological tolerance and inflammation in autoimmune disease. Bioessays. 2010;32:1067-76
67. Lee AYS, Eri R, Lyons AB, Grimm MC, Korner H. CC Chemokine Ligand 20 and Its Cognate Receptor CCR6 in Mucosal T Cell Immunology and Inflammatory Bowel Disease: Odd Couple or Axis of Evil?. Front Immunol. 2013;4:194
68. Xiao Z, Casey KA, Jameson SC, Curtsinger JM, Mescher MF. Programming for CD8 T cell memory development requires IL-12 or type I IFN. J Immunol. 2009;182:2786-94
69. Olsson Åkefeldt S, Maisse C, Belot A. et al. Chemoresistance of human monocyte-derived dendritic cells is regulated by IL-17A. PLoS ONE. 2013;8:e56865
70. Carrington EM, Zhang J-G, Sutherland RM. et al. Prosurvival Bcl-2 family members reveal a distinct apoptotic identity between conventional and plasmacytoid dendritic cells. Proc Natl Acad Sci USA. 2015;112:4044-9
71. Gao Y, Nish SA, Jiang R. et al. Control of T helper 2 responses by transcription factor IRF4-dependent dendritic cells. Immunity. 2013;39:722-32
72. Moynagh PN. The NF-kappaB pathway. J Cell Sci. 2005;118:4589-92
73. Hoesel B, Schmid JA. The complexity of NF-κB signaling in inflammation and cancer. Mol Cancer. 2013;12:86
74. Paardekooper LM, Bendix MB, Ottria A. et al. Hypoxia potentiates monocyte-derived dendritic cells for release of tumor necrosis factor α via MAP3K8. Biosci Rep. 2018 38
75. Hegde VL, Singh NP, Nagarkatti PS, Nagarkatti M. CD44 mobilization in allogeneic dendritic cell-T cell immunological synapse plays a key role in T cell activation. J Leukoc Biol. 2008;84:134-42
76. Termeer C, Averbeck M, Hara H. et al. Targeting dendritic cells with CD44 monoclonal antibodies selectively inhibits the proliferation of naive CD4+ T-helper cells by induction of FAS-independent T-cell apoptosis. Immunology. 2003;109:32-40
77. Luan Y-Y, Yao R-Q, Tong S, Dong N, Sheng Z-Y, Yao Y-M. Effect of tumor necrosis factor-α induced protein 8 like-2 on immune function of dendritic cells in mice following acute insults. Oncotarget. 2016;7:30178-92
78. Kabashima K, Shiraishi N, Sugita K. et al. CXCL12-CXCR4 engagement is required for migration of cutaneous dendritic cells. Am J Pathol. 2007;171:1249-57
79. Benvenuti F, Hugues S, Walmsley M. et al. Requirement of Rac1 and Rac2 expression by mature dendritic cells for T cell priming. Science. 2004;305:1150-3
80. Zanoni I, Granucci F. Regulation and dysregulation of innate immunity by NFAT signaling downstream of pattern recognition receptors (PRRs). Eur J Immunol. 2012;42:1924-31
81. Bao M, Wang Y, Liu Y. et al. NFATC3 promotes IRF7 transcriptional activity in plasmacy-toid dendritic cells. J Exp Med. 2016;213:2383-98
82. Collins LE, DeCourcey J, Rochfort KD, Kristek M, Loscher CE. A role for syntaxin 3 in the secretion of IL-6 from dendritic cells following activation of toll-like receptors. Front Immunol. 2014;5:691
83. Gottipati S, Rao NL, Fung-Leung W-P. IRAK1: a critical signaling mediator of innate immunity. Cell Signal. 2008;20:269-76
84. McNab F, Mayer-Barber K, Sher A, Wack A, O'Garra A. Type I interferons in infectious disease. Nat Rev Immunol. 2015;15:87-103
85. Finotti G, Tamassia N, Cassatella MA. Interferon-λs and Plasmacytoid Dendritic Cells: A Close Relationship. Front Immunol. 2017;8:1015
86. Anguille S, Lion E, Van den Bergh J. et al. Interleukin-15 dendritic cells as vaccine candidates for cancer immunotherapy. Hum Vaccin Immunother. 2013;9:1956-61
87. Mookerjee A, Graciotti M, Kandalaft LE. IL-15 and a Two-Step Maturation Process Improve Bone Marrow-Derived Dendritic Cell Cancer Vaccine. Cancers (Basel). 2019 11
88. Wilcox RA, Chapoval AI, Gorski KS. et al. Cutting edge: Expression of functional CD137 receptor by dendritic cells. J Immunol. 2002;168:4262-7
89. Lagali NS, Badian RA, Liu X. et al. Dendritic cell maturation in the corneal epithelium with onset of type 2 diabetes is associated with tumor necrosis factor receptor superfamily member 9. Sci Rep. 2018;8:14248
90. Shin VY, Jin H, Ng EKO. et al. NF-κB targets miR-16 and miR-21 in gastric cancer: involvement of prostaglandin E receptors. Carcinogenesis. 2011;32:240-5
91. Pla A, Zhong X, Rayner S. miRAW: A deep learning-based approach to predict microRNA targets by analyzing whole microRNA transcripts. PLoS Comput Biol. 2018;14:e1006185
92. Pinzón N, Li B, Martinez L. et al. microRNA target prediction programs predict many false positives. Genome Res. 2017;27:234-45
93. Amon L, Lehmann CHK, Baranska A, Schoen J, Heger L, Dudziak D. Transcriptional control of dendritic cell development and functions. Int Rev Cell Mol Biol. 2019;349:55-151
94. Ye H, Liu X, Lv M. et al. MicroRNA and transcription factor co-regulatory network analysis reveals miR-19 inhibits CYLD in T-cell acute lymphoblastic leukemia. Nucleic Acids Res. 2012;40:5201-14
95. Bracken CP, Scott HS, Goodall GJ. A network-biology perspective of microRNA function and dysfunction in cancer. Nat Rev Genet. 2016;17:719-32
96. Winkler I, Bitter C, Winkler S. et al. Identification of Pparγ-modulated miRNA hubs that target the fibrotic tumor microenvironment. Proc Natl Acad Sci USA. 2020;117:454-63
97. Hayden MS, West AP, Ghosh S. NF-kappaB and the immune response. Oncogene. 2006;25:6758-80
98. Kanehisa M, Furumichi M, Tanabe M, Sato Y, Morishima K. KEGG: new perspectives on genomes, pathways, diseases and drugs. Nucleic Acids Res. 2017;45:D353-61
99. Slenter DN, Kutmon M, Hanspers K. et al. WikiPathways: a multifaceted pathway database bridging metabolomics to other omics research. Nucleic Acids Res. 2018;46:D661-7
100. Mubeen S, Hoyt CT, Gemünd A, Hofmann-Apitius M, Fröhlich H, Domingo-Fernández D. The Impact of Pathway Database Choice on Statistical Enrichment Analysis and Predictive Modeling. Front Genet. 2019;10:1203
101. Geistlinger L, Csaba G, Santarelli M. et al. Toward a gold standard for benchmarking gene set enrichment analysis. Brief Bioinformatics. 2020
102. Khan FM, Marquardt S, Gupta SK. et al. Unraveling a tumor type-specific regulatory core underlying E2F1-mediated epithelial-mesenchymal transition to predict receptor protein signatures. Nat Commun. 2017;8:198
103. Dreyer FS, Cantone M, Eberhardt M. et al. A web platform for the network analysis of high-throughput data in melanoma and its use to investigate mechanisms of resistance to anti-PD1 immunotherapy. Biochim Biophys Acta. 2018;1864:2315-28
104. Moreau Y, Tranchevent L-C. Computational tools for prioritizing candidate genes: boosting disease gene discovery. Nat Rev Genet. 2012;13:523-36
105. Guala D, Sonnhammer ELL. A large-scale benchmark of gene prioritization methods. Sci Rep. 2017;7:46598
106. Iván G, Grolmusz V. When the Web meets the cell: using personalized PageRank for analyzing protein interaction networks. Bioinformatics. 2011;27:405-7
107. Hsu C-L, Huang Y-H, Hsu C-T, Yang U-C. Prioritizing disease candidate genes by a gene interconnectedness-based approach. BMC Genomics. 2011;12(Suppl 3):S25
108. McGeary SE, Lin KS, Shi CY. et al. The biochemical basis of microRNA targeting efficacy. Science. 2019 366
Corresponding authors: Xin Lai, Tel: +49(0)91318545888, E-mail: xin.laide; Julio Vera, Tel: +49(0)91318545876, E-mail: julio.vera-gonzalezde.