motif_enrichment()
calculates motif enrichment relative to a set of
background sequences using Analysis of Motif Enrichment (AME) from
memes.
Usage
motif_enrichment(
peak_input,
motif,
genome_build,
out_dir = tempdir(),
verbose = FALSE,
meme_path = NULL,
...
)
Arguments
- peak_input
Either a path to the narrowPeak file or a GRanges peak object generated by
read_peak_file()
.- motif
An object of class
universalmotif
.- genome_build
The genome build that the peak sequences should be derived from.
- out_dir
Location to save the 0-order background file along with the AME output files.
- verbose
A logical indicating whether to print verbose messages while running the function. (default = FALSE)
- meme_path
path to "meme/bin/" (default:
NULL
). Will use default search behavior as described incheck_meme_install()
if unset.- ...
Arguments passed on to
memes::runAme
method
default: fisher (allowed values: fisher, ranksum, pearson, spearman, 3dmhg, 4dmhg)
sequences
logical(1)
add results fromsequences.tsv
tosequences
list column to returned data.frame. Valid only if method = "fisher". See AME outputs webpage for more information (Default: FALSE).silent
whether to suppress stdout (default: TRUE), useful for debugging.
Value
A list containing a AME results data frame and a numeric referring to the proportion of peaks with a motif.
Examples
if (memes::meme_is_installed()) {
data("CTCF_TIP_peaks", package = "MotifPeeker")
data("motif_MA1102.3", package = "MotifPeeker")
res <- motif_enrichment(
peak_input = CTCF_TIP_peaks,
motif = motif_MA1102.3,
genome_build =
BSgenome.Hsapiens.UCSC.hg38::BSgenome.Hsapiens.UCSC.hg38,
)
print(res)
}
#> $tp
#> [1] 73.00 40.11
#>
#> $fp
#> [1] 286.00 26.19
#>
#> $positive_peaks
#> GRanges object with 73 ranges and 7 metadata columns:
#> seqnames ranges strand | name score
#> <Rle> <IRanges> <Rle> | <character> <numeric>
#> TIPseq_peak_63 chr10 69801154-69801320 * | TIPseq_peak_63 109
#> TIPseq_peak_90 chr10 70559156-70559588 * | TIPseq_peak_90 122
#> TIPseq_peak_157 chr10 73529363-73529638 * | TIPseq_peak_157 106
#> TIPseq_peak_167 chr10 73957451-73957745 * | TIPseq_peak_167 73
#> TIPseq_peak_39 chr10 68988538-68988838 * | TIPseq_peak_39 111
#> ... ... ... ... . ... ...
#> TIPseq_peak_105 chr10 71270297-71270663 * | TIPseq_peak_105 138
#> TIPseq_peak_51 chr10 69507836-69508269 * | TIPseq_peak_51 206
#> TIPseq_peak_169 chr10 73966072-73966376 * | TIPseq_peak_169 70
#> TIPseq_peak_44 chr10 69232402-69232647 * | TIPseq_peak_44 45
#> TIPseq_peak_46 chr10 69392521-69392898 * | TIPseq_peak_46 132
#> signalValue pValue qValue peak summit
#> <numeric> <numeric> <numeric> <integer> <integer>
#> TIPseq_peak_63 7.19424 15.0082 10.96220 54 69801208
#> TIPseq_peak_90 8.60956 16.5013 12.29430 312 70559468
#> TIPseq_peak_157 8.79456 14.6750 10.66540 105 73529468
#> TIPseq_peak_167 6.58212 11.0334 7.35935 125 73957576
#> TIPseq_peak_39 9.10857 15.2180 11.14710 168 68988706
#> ... ... ... ... ... ...
#> TIPseq_peak_105 9.39681 18.29610 13.87650 113 71270410
#> TIPseq_peak_51 11.20730 25.79100 20.65960 119 69507955
#> TIPseq_peak_169 6.33687 10.63470 7.00211 168 73966240
#> TIPseq_peak_44 5.62082 7.95049 4.55375 110 69232512
#> TIPseq_peak_46 9.25369 17.56690 13.23280 225 69392746
#> -------
#> seqinfo: 1 sequence from an unspecified genome; no seqlengths
#>