DIANA-miRPath v4.0: expanding target-based miRNA functional analysis in cell-type and tissue contexts

1microRNAs (miRNAs) are central post-transcriptional modulators of gene expression and protein production that exhibit abundant yet distinct expression profiles across tissues and cell-types. Robust high throughput methods directly capture functional miRNA Recognition Elements (MREs) for one or more miRNAs in the 3’ untranslated region (UTR) and/or the coding (CDS) sequence of most transcripts, providing valuable insights into miRNA regulation of most biological processes. The relative abundance of miRNAs and their targets creates unique RNA landscapes in different cell types, tissues and conditions, in which distinct miRNA regulatory roles are exerted.

2De novo miRNA target prediction methods and collections of miRNA interactions supported by robust high-throughput experimental methods paved the way for the development of applications conducting target-based analysis of miRNA functions. Previous versions of DIANA-miRPath provided a statistical framework for evaluation of combined miRNA effects on Gene Ontology (GO) terms and KEGG pathways, based on predicted and experimentally supported targets.

3DIANA-miRPath v4.0 is the first relevant application to integrate robust experimental and predicted interactions with phenotypic and gene expression contexts, towards tailored miRNA functional analyses. In the current version, analyses have been dissected into separate modules that each answer different biological questions (Figure 1). Upon providing one or multiple miRNAs and selecting resources of miRNA interactions and terms/pathways/ontologies (hereafter, ‘terms’), enrichment analysis can be performed using the miRNA targets union, intersections, or more complex designs (Figure 1A). miRPath v4.0 extends the functionality of previous versions towards context-specific interrogations by applying nonparametric analysis of targeted and non-targeted term components in cell types, tissues and diseases. Nonparametric testing of log-scaled gene expression fold changes has been successfully utilized in the past to assess the functional efficacy of miRNA binding events. In the current scenario, it is repurposed to enable the identification of contexts and states in which the distribution of the targets is found repressed/de-regulated, transforming miRPath analyses from context-agnostic to context-aware endeavours (Figure 1B). Distribution testing is performed against internal libraries of expression sets, including tissues (healthy or cancerous) and cell types (adult or fetal), outputting results per term, sorted by level of significance. The option to also upload and use custom datasets expands the applicability of context-specific testing even more. In the absence of specific miRNAs of interest, explorations focused instead on terms can be initiated via Term-centric analysis, uncovering potential common miRNA regulators across one or multiple terms (Figure 1C). A plethora of interaction, term and context options are available by default in DIANA-miRPath v4.0 (Figure 1D).

4Potential biases in miRNA functional analysis, relevant to the multiplicity of miRNA targeting, as well as the 3'-UTR length and MRE differences among transcripts within terms, have been discussed in the past. In DIANA-miRPath v3.0, permutation-based testing, which is computationally intensive and yields non-granular p-values, was provided as an alternative testing method in order to close this gap. In the current version, it is replaced by Wallenius non-central hypergeometric testing, parameterized to account for bias in targeting events among tested terms.

5DIANA-miRPath v4.0 provides a novel module that caters miRNA-specific functional analysis of positive/negative selection results obtained using CRISPR knock-out screen libraries in stressors of interest. Via the relevant module, users can provide tab-separated results from MAGeCK analysis (RRA or MLE) of libraries targeting both genes and miRNAs, such as GeCKO v2.0 human, and obtain significantly enriched terms, miRNAs and targets that are positively/negatively selected in conditions under study.

6MIENTURNET, miRPathDB 2.0 and miR2GO are three webservers that are also devoted to miRNA functional analysis. In MIENTURNET users can also produce miRNA-mRNA interaction networks however a term-based analysis is not supported. miRPathDB 2.0 provides a set of downstream information for the queried miRNA (e.g. stem-loop, mature and seed sequence) for both miRBase and miRCarta annotations but supports only human and mouse species and single miRNA queries. On the contrary, the miRNA-centric module of miRPath-v4.0 performs functional enrichment analysis allowing as input both single and multiple miRNAs, to uncover the potential combinatorial effects of miRNAs in terms. miR2GO provides an interface for miRNA functional analysis that also evaluates the functional effects of mutations in miRNA seeds. Even though this feature is not catered as a built-in miRPath capacity, users may still conduct such investigations by making use of the Upload Dashboard module.

7miRPath-v4.0 extends the functionality of its previous version (DIANA-miRPath-v3.0) and offers numerous improvements over all the aforementioned implementations. miRPath-v4.0 provides the largest set of resources for both experimentally supported and predicted miRNA-gene interactions, biological terms and pathways and gene expression profiles. Importantly, it is the only method to our knowledge to incorporate gene expression profiles with single-cell resolution, as demonstrated in Supplementary Note S1, Supplementary Figures S1–S4, and Supplementary Tables S1–S3. Two newly introduced key advantages of miRPath-v4.0 are the dedicated workflow for miRNA-centric functional analysis of CRISPR-KO Screen data, for which a specific use-case scenario is provided in Supplementary Note S2 and Supplementary Tables S4–S5, and the Upload Dashboard, where users can unlock the full capacity of miRPath with custom annotations (i.e. miRNAs, genes, terms) for the species/phenotypes of their choice. The case that specific terms are of interest, yet no miRNAs are provided to initiate a miRNA-centric analysis is shown in Supplementary Note S3 and Supplementary Tables S6–S7, which provide a Term-centric miRPath application. The latest version of miRPath is the only relevant implementation to provide support for both miRBase and MirGeneDB miRNA annotation schemes. It hosts the largest number of supported species, offers an option to account for bias in targeting events and provides a completely redesigned, intuitive interface to enhance user experience and enable its ingenious use.

8Enrichment analysis is performed by parameterizing function kegga from limma package to user-selected resources. Depending on the scenario, the union or intersection(s) of miRNA targets is calculated and subjected to testing. In the case of Pathways union, the targets of each miRNA are tested separately and significantly enriched instances are merged per-term using Fisher's method, as implemented in metap package. The option to account for targeting differences among tested terms through Wallenius non-central hypergeometric testing is integrated using packages GOseq and BiasedUrn.

9Post-enrichment analysis compares miRNA–targeted versus non-targeted pathway components across a selected context set (Section ‘Data collections’ below) using either a Mann Whitney U test (default) or a two-sample Kolmogorov-Smirnov test, as implemented in base R stats package (functions wilcox.test and ks.test), to assess difference in the mean ranks of miRNA–targeted and non-targeted genes, or if they are drawn from the same distribution, respectively. By default, one-tailed tests are used, in order to assess the repressive effect of miRNAs, however two-sided testing can also be selected to accommodate, for example, scenarios where the abundance of miRNAs under question is limited.

10Primary term analysis of the top positively or negatively selected miRNAs (user-defined) is performed taking the union of their targets and employing hypergeometric testing, as described in, ‘Main analysis engine’. For significantly enriched pathways, two downstream options are provided: the default option performs competitive, rank-based CAMERA testing (proposed and described by Wu and Smyth) to identify significant left/right trends of targets within the opposite ranking (e.g. for positively selected miRNAs the negative ranking of genes is assessed); the second option also requires setting the number of top opposite-ranked genes and performs one-sided Fisher's exact test to evaluate the enrichment on the rank's extreme in miRNA targets.

11Gene identifiers were retrieved via Biomart from Ensembl v102. miRNA annotations were retrieved from miRBase v22.1 and also matched to instances from MirGeneDB v2.1.

12Experimentally verified interactions were retrieved from reference resources DIANA-TarBase v8, DIANA-LncBase v3, and miRTarBase 2022. Context-aware interaction sets were obtained by analysing publicly available CLIP-Seq datasets in 35 human and 6 mouse tissues with microCLIP. Predicted interactions were obtained directly from TargetScan 8 and by performing de novo DIANA-microT-CDS executions on Ensembl v102 and miRBase v22.1 annotation. microT interactions include for the first time the capacity to query host targets for virally encoded miRNAs regarding 20 viruses infecting human (EBV, KSHV, HSV1-2, HCMV, SFV, HBV, HIV1, HHV6B, MCPV, JCV, BKV, SV40 and TTV), mouse (MGHV, MCMV) and chicken (MDV1-2, HVT, ILTV). For TarBase, LncBase, miRTarBase and TargetScan, gene IDs were converted to Ensembl v102 using Ensembl ID History Converter tool or Biomart (miRTarBase) and miRNA names were manually converted to miRBase v22.1. Supplemental confidence information, including annotation on direct/indirect proof of interaction by TarBase, high/low interaction confidence annotation by miRTarBase, and interaction scores by target prediction tools, was retained for collections.

13Annotation of GO terms against Gene IDs was retrieved through Ensembl. Pathview package was used to retrieve KEGG-to-Entrez Gene ID relationships, which were annotated to Ensembl Gene IDs (Biomart). Reactome pathway annotations were retrieved and Ensembl Gene IDs were retained and converted to Ensembl v102. Expertly curated Molecular Signatures Database (MSigDB) gene sets regarding Hallmark (H), Oncogenic, (C7), Immunologic (C8) and Cell type (C9) signatures were retrieved using the msigdbr package and annotated to Ensembl Gene IDs (Biomart). Protein family annotations were retrieved from PFam.

14Adult and fetal single-cell RNA-Seq data was collected from the Han et al. human cell atlas. Log2(Fold Change) (log2FC) gene expression values were obtained for each cell type by estimating the percent of individual cells exhibiting at least one read within each annotated cell type and contrasting it to the median of percentages of the remaining cell types. Transcripts-Per-Million (TPM-) level gene expression data were obtained for each cancer project in TCGA via FireBrowse (https://firebrowse.org). Adjacent healthy and blood-derived normal samples were filtered out. Healthy tissue expression data (TPM-level) were retrieved from the GTEx Consortium. TPM-level data were used to denote a per-gene per-tissue log2FC, similarly to the single-cell data, depicting the percentage of samples presenting at least 10 TPM in a tissue relative to the median of the percentages of the remaining tissues. This metric was adopted as a balanced estimator to depict consistent and cell-/tissue-specific non-zero expression. For all three resources (single-cell atlas, TCGA, GTEx), gene names were converted to Ensembl v102 Gene IDs using Ensembl ID History Converter tool and R.

15The user interface of miRPath has been upgraded to provide further flexibility without sacrificing user-friendliness. Figure 2 briefly showcases the input menu and results from the tab enabling miRNA-centric functional analysis. All interaction, term and context resources included in DIANA-miRPath v4.0 were converted to Mongo collections and stored locally in a MongoDB schema (v4.4.2). DIANA-miRPath v4 was set up as an R Shiny application (https://shiny.rstudio.com/), and deployed using ShinyProxy (https://www.shinyproxy.io). The data access layer is implemented using mongolite package. The presentation layer relies primarily on packages shiny, shinyalert, shinyjs, markdown, tippy, while visualizations are created making use of packages ggplot2, limma, Pathview and pheatmap (v1.0.12, Kolde, 2019). R Markdown was used to create a Help Section, while all functionalities requiring user input are accompanied with example data to further aid user experience.

16DIANA-miRPath v4.0 constitutes an upgrade both in terms of interface, which now holds a modular application-like character, and of content. Available visualization options were also expanded. Importantly, a User Upload Dashboard is introduced to further augment tailored use and flexibility. By it, the full arsenal of miRPath modules is unlocked for application with custom, user-provided annotations e.g. non-model organisms, modified/edited miRNAs, or even other non-coding RNAs. DIANA-miRPath v4.0 constitutes a valuable resource, delivering numerous computational methods to initiate or consolidate miRNA-related investigations.

17DIANA-miRPath is freely available, without registration or login, at http://www.microrna.gr/miRPathv4.