PeakSEA — Peak-Set Enrichment Analysis

The goal of PeakSEA is to test for enrichment in public peak set data using preranked peak data as input.

NB: I renamed the package from PSEA to PeakSEA to avoid confusion with other packages.

For example ATACseq peak can be tested with DiffBind for differential opening, and ranking metrics, such as log Fold change can be used. We can then test for enrichment in opening (up) or closing (down) with publicly available databases, eg LOLAcore for epigenetic annotations or causalDB for finemapped credible intervals of disease loci derived from GWAS dataset.

PSEA returns a gseResult object that is compatible with DOSE package routines.

Installation

You can install the development version of PeakSEA with:

# install.packages("pak")
pak::pak("ImmuneAxisa/PeakSEA")
# install.packages("devtools")
devtools::install_github("ImmuneAxisa/PeakSEA")

Example

library(PeakSEA)
#> 
library(plyranges)
#> Loading required package: BiocGenerics
#> 
#> Attaching package: 'BiocGenerics'
#> The following objects are masked from 'package:stats':
#> 
#>     IQR, mad, sd, var, xtabs
#> The following objects are masked from 'package:base':
#> 
#>     anyDuplicated, aperm, append, as.data.frame, basename, cbind,
#>     colnames, dirname, do.call, duplicated, eval, evalq, Filter, Find,
#>     get, grep, grepl, intersect, is.unsorted, lapply, Map, mapply,
#>     match, mget, order, paste, pmax, pmax.int, pmin, pmin.int,
#>     Position, rank, rbind, Reduce, rownames, sapply, saveRDS, setdiff,
#>     table, tapply, union, unique, unsplit, which.max, which.min
#> Loading required package: IRanges
#> Warning: package 'IRanges' was built under R version 4.4.2
#> Loading required package: S4Vectors
#> Loading required package: stats4
#> 
#> Attaching package: 'S4Vectors'
#> The following object is masked from 'package:utils':
#> 
#>     findMatches
#> The following objects are masked from 'package:base':
#> 
#>     expand.grid, I, unname
#> Loading required package: GenomicRanges
#> Loading required package: GenomeInfoDb
#> Warning: package 'GenomeInfoDb' was built under R version 4.4.2
#> 
#> Attaching package: 'plyranges'
#> The following object is masked from 'package:IRanges':
#> 
#>     slice
#> The following object is masked from 'package:stats':
#> 
#>     filter
library(dplyr)
#> 
#> Attaching package: 'dplyr'
#> The following objects are masked from 'package:plyranges':
#> 
#>     between, n, n_distinct
#> The following objects are masked from 'package:GenomicRanges':
#> 
#>     intersect, setdiff, union
#> The following object is masked from 'package:GenomeInfoDb':
#> 
#>     intersect
#> The following objects are masked from 'package:IRanges':
#> 
#>     collapse, desc, intersect, setdiff, slice, union
#> The following objects are masked from 'package:S4Vectors':
#> 
#>     first, intersect, rename, setdiff, setequal, union
#> The following objects are masked from 'package:BiocGenerics':
#> 
#>     combine, intersect, setdiff, union
#> The following objects are masked from 'package:stats':
#> 
#>     filter, lag
#> The following objects are masked from 'package:base':
#> 
#>     intersect, setdiff, setequal, union

gr <- data.frame(
seqnames = "chr1",
start = seq(0,9999,100),
width = 30,
ranking_metric = -50:49
  ) %>%
  as_granges()

gr_set <- shift_right(gr, rep_len(c(20,50),100)) %>%
  mutate(term = "peakSet1")

psea(gr, ranking_metric, gr_set)
#> using 'fgsea' for GSEA analysis, please cite Korotkevich et al (2019).
#> preparing geneSet collections...
#> GSEA analysis...
#> no term enriched under specific pvalueCutoff...
#> #
#> # Gene Set Enrichment Analysis
#> #
#> #...@organism     UNKNOWN 
#> #...@setType      UNKNOWN 
#> #...@geneList     Named int [1:100] 49 48 47 46 45 44 43 42 41 40 ...
#>  - attr(*, "names")= chr [1:100] "chr1_9900_9929" "chr1_9800_9829" "chr1_9700_9729" "chr1_9600_9629" ...
#> #...nPerm     
#> #...pvalues adjusted by 'BH' with cutoff <0.05 
#> #...0 enriched terms found
#> 'data.frame':    0 obs. of  8 variables:
#>  $ ID             : chr 
#>  $ Description    : chr 
#>  $ setSize        : int 
#>  $ enrichmentScore: num 
#>  $ NES            : num 
#>  $ pvalue         : num 
#>  $ p.adjust       : num 
#>  $ qvalue         : logi 
#> #...Citation
#> G Yu. Thirteen years of clusterProfiler. The Innovation. 2024, 5(6):100722