nbea: Network-based enrichment analysis (NBEA)

Description Usage Arguments Details Value Author(s) References See Also Examples

View source: R/nbea.R

Description

This is the main function for network-based enrichment analysis. It implements and wraps existing implementations of several frequently used methods and allows a flexible inspection of resulting gene set rankings.

Usage

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
nbeaMethods()

nbea(
  method = EnrichmentBrowser::nbeaMethods(),
  se,
  gs,
  grn,
  prune.grn = TRUE,
  alpha = 0.05,
  perm = 1000,
  padj.method = "none",
  out.file = NULL,
  browse = FALSE,
  assay = "auto",
  ...
)

Arguments

method

Network-based enrichment analysis method. Currently, the following network-based enrichment analysis methods are supported: ‘ggea’, ‘spia’, ‘pathnet’, ‘degraph’, ‘topologygsa’, ‘ganpa’, ‘cepa’, ‘netgsa’, and ‘neat’. Default is 'ggea'. This can also be the name of a user-defined function implementing a method for network-based enrichment analysis. See Details.

se

Expression dataset. An object of class SummarizedExperiment. Mandatory minimal annotations:

  • colData column storing binary group assignment (named "GROUP")

  • rowData column storing (log2) fold changes of differential expression between sample groups (named "FC")

  • rowData column storing adjusted (corrected for multiple testing) p-values of differential expression between sample groups (named "ADJ.PVAL")

Additional optional annotations:

  • colData column defining paired samples or sample blocks (named "BLOCK")

  • metadata slot named "annotation" giving the organism under investigation in KEGG three letter code (e.g. "hsa" for Homo sapiens)

  • metadata slot named "dataType" indicating the expression data type ("ma" for microarray, "rseq" for RNA-seq)

gs

Gene sets. Either a list of gene sets (character vectors of gene IDs) or a text file in GMT format storing all gene sets under investigation.

grn

Gene regulatory network. Either an absolute file path to a tabular file or a character matrix with exactly *THREE* cols; 1st col = IDs of regulating genes; 2nd col = corresponding regulated genes; 3rd col = regulation effect; Use '+' and '-' for activation/inhibition.

prune.grn

Logical. Should the GRN be pruned? This removes duplicated, self, and reversed edges. Defaults to TRUE.

alpha

Statistical significance level. Defaults to 0.05.

perm

Number of permutations of the expression matrix to estimate the null distribution. Defaults to 1000. If using method=‘ggea’, it is possible to set 'perm=0' to use a fast approximation of gene set significance to avoid permutation testing. See Details.

padj.method

Method for adjusting nominal gene set p-values to multiple testing. For available methods see the man page of the stats function p.adjust. Defaults to 'none', i.e. leaves the nominal gene set p-values unadjusted.

out.file

Optional output file the gene set ranking will be written to.

browse

Logical. Should results be displayed in the browser for interactive exploration? Defaults to FALSE.

assay

Character. The name of the assay for enrichment analysis if se is a SummarizedExperiment with *multiple assays*. Defaults to "auto", which automatically determines the appropriate assay based on data type provided and enrichment method selected. See details.

...

Additional arguments passed to individual nbea methods. This includes currently:

  • beta: Log2 fold change significance level. Defaults to 1 (2-fold).

For SPIA and NEAT:

  • sig.stat: decides which statistic is used for determining significant DE genes. Options are:

    • 'p' (Default): genes with adjusted p-value below alpha.

    • 'fc': genes with abs(log2(fold change)) above beta

    • '&': p & fc (logical AND)

    • '|': p | fc (logical OR)

    • 'xxp': top xx % of genes sorted by adjusted p-value

    • 'xxfc' top xx % of genes sorted by absolute log2 fold change.

For GGEA:

  • cons.thresh: edge consistency threshold between -1 and 1. Defaults to 0.2, i.e. only edges of the GRN with consistency >= 0.2 are included in the analysis. Evaluation on real datasets has shown that this works best to distinguish relevant gene sets. Use consistency of -1 to include all edges.

  • gs.edges: decides which edges of the grn are considered for a gene set under investigation. Should be one out of c('&', '|'), denoting logical AND and OR. respectively. Accordingly, this either includes edges for which regulator AND / OR target gene are members of the investigated gene set.

Details

'ggea': gene graph enrichment analysis, scores gene sets according to consistency within the given gene regulatory network, i.e. checks activating regulations for positive correlation and repressing regulations for negative correlation of regulator and target gene expression (Geistlinger et al., 2011). When using 'ggea' it is possible to estimate the statistical significance of the consistency score of each gene set in two different ways: (1) based on sample permutation as described in the original publication (Geistlinger et al., 2011) or (2) using an approximation in the spirit of Bioconductor's npGSEA package that is much faster.

'spia': signaling pathway impact analysis, combines ORA with the probability that expression changes are propagated across the pathway topology; implemented in Bioconductor's SPIA package (Tarca et al., 2009).

'pathnet': pathway analysis using network information, applies ORA on combined evidence for the observed signal for gene nodes and the signal implied by connected neighbors in the network; implemented in Bioconductor's PathNet package.

'degraph': differential expression testing for gene graphs, multivariate testing of differences in mean incorporating underlying graph structure; implemented in Bioconductor's DEGraph package.

'topologygsa': topology-based gene set analysis, uses Gaussian graphical models to incorporate the dependence structure among genes as implied by pathway topology; implemented in CRAN's topologyGSA package.

'ganpa': gene association network-based pathway analysis, incorporates network-derived gene weights in the enrichment analysis; implemented in CRAN's GANPA package.

'cepa': centrality-based pathway enrichment, incorporates network centralities as node weights mapped from differentially expressed genes in pathways; implemented in CRAN's CePa package.

'netgsa': network-based gene set analysis, incorporates external information about interactions among genes as well as novel interactions learned from data; implemented in CRAN's NetGSA package.

'neat': network enrichment analysis test, compares the number of links between differentially expressed genes and a gene set of interest to the number of links expected under a hypergeometric null model; proposed by Signorelli et al. (2016) and implemented in CRAN's neat package.

It is also possible to use additional network-based enrichment methods. This requires to implement a function that takes 'se', 'gs', and 'grn' and as arguments and returns a numeric matrix 'res.tbl' with a mandatory column named 'PVAL' storing the resulting p-value for each gene set in 'gs'. The rows of this matrix must be named accordingly (i.e. rownames(res.tbl) == names(gs)). See examples.

Using a SummarizedExperiment with *multiple assays*:

For the typical use case within the EnrichmentBrowser workflow this will be a SummarizedExperiment with two assays: (i) an assay storing the *raw* expression values, and (ii) an assay storing the *norm*alized expression values as obtained with the normalize function.

In this case, assay = "auto" will *auto*matically determine the assay based on the data type provided and the enrichment method selected. For usage outside of the typical workflow, the assay argument can be used to provide the name of the assay for the enrichment analysis.

Value

nbeaMethods: a character vector of currently supported methods;

nbea: if(is.null(out.file)): an enrichment analysis result object that can be detailedly explored by calling eaBrowse and from which a flat gene set ranking can be extracted by calling gsRanking. If 'out.file' is given, the ranking is written to the specified file.

Author(s)

Ludwig Geistlinger <Ludwig.Geistlinger@sph.cuny.edu>

References

Geistlinger et al. (2011) From sets to graphs: towards a realistic enrichment analysis of transcriptomic systems. Bioinformatics, 27(13), i366-73.

Tarca et al. (2009) A novel signaling pathway impact analysis. Bioinformatics, 25(1):75-82.

Signorelli et al. (2016) NEAT: an efficient network enrichment analysis test. BMC Bioinformatics, 17:352.

See Also

Input: readSE, probe2gene, getGenesets to retrieve gene set definitions from databases such as GO and KEGG. compileGRN to construct a GRN from pathway databases.

Output: gsRanking to rank the list of gene sets. eaBrowse for exploration of resulting gene sets.

Other: sbea to perform set-based enrichment analysis. combResults to combine results from different methods.

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
    # currently supported methods
    nbeaMethods()

    # (1) expression data: 
    # simulated expression values of 100 genes
    # in two sample groups of 6 samples each
    se <- makeExampleData(what="SE")
    se <- deAna(se)

    # (2) gene sets:
    # draw 10 gene sets with 15-25 genes
    gs <- makeExampleData(what="gs", gnames=names(se))

    # (3) make 2 artificially enriched sets:
    sig.genes <- names(se)[rowData(se)$ADJ.PVAL < 0.1]
    gs[[1]] <- sample(sig.genes, length(gs[[1]])) 
    gs[[2]] <- sample(sig.genes, length(gs[[2]]))   
   
    # (4) gene regulatory network 
    grn <- makeExampleData(what="grn", nodes=names(se))
    
    # (5) performing the enrichment analysis
    ea.res <- nbea(method="ggea", se=se, gs=gs, grn=grn)

    # (6) result visualization and exploration
    gsRanking(ea.res, signif.only=FALSE)

    # using your own function as enrichment method
    dummyNBEA <- function(se, gs, grn)
    {
        sig.ps <- sample(seq(0,0.05, length=1000),5)
        insig.ps <- sample(seq(0.1,1, length=1000), length(gs)-5)
        ps <- sample(c(sig.ps, insig.ps), length(gs))
        score <- sample(1:100, length(gs), replace=TRUE)
        res.tbl <- cbind(score, ps)
        colnames(res.tbl) <- c("SCORE", "PVAL")
        rownames(res.tbl) <- names(gs)
        return(res.tbl[order(ps),])
    }

    ea.res2 <- nbea(method=dummyNBEA, se=se, gs=gs, grn=grn)
    gsRanking(ea.res2) 

EnrichmentBrowser documentation built on Dec. 12, 2020, 2 a.m.