RMotifScanPair {esATAC}R Documentation

Search Motif Position in Given Regions

Description

Search motif position in given genome regions according PWM matrix.

Usage

atacMotifScanPair(atacProc, peak1 = NULL, peak2 = NULL, background = NULL,
  genome = NULL, motifPWM = NULL, min.score = "85%", scanO.dir = NULL,
  n.cores = NULL, prefix = NULL, use.SavedPWM = NULL, ...)

## S4 method for signature 'ATACProc'
atacMotifScanPair(atacProc, peak1 = NULL, peak2 = NULL,
  background = NULL, genome = NULL, motifPWM = NULL, min.score = "85%",
  scanO.dir = NULL, n.cores = NULL, prefix = NULL, use.SavedPWM = NULL,
  ...)

motifscanpair(peak1 = NULL, peak2 = NULL, background = NULL,
  genome = NULL, motifPWM = NULL, min.score = "85%", scanO.dir = NULL,
  n.cores = NULL, prefix = NULL, use.SavedPWM = NULL, ...)

Arguments

atacProc

ATACProc-class object scalar. It has to be the return value of upstream process: atacpeakComp.

peak1

peak file path.

peak2

peak file path.

background

background peak file path.

genome

A DNAString object.

motifPWM

list scalar. Default: from setConfigure. Every element in the list contains a motif PWM matrix. e.g. pwm <- list("CTCF" = CTCF_PWMmatrix)

min.score

The minimum score for counting a match. Can be given as a character string containing a percentage (e.g. "85 possible score or as a single number.

scanO.dir

Character scalar. the output file directory. This function will use the index in motifPWM as the file name to save the motif position information in separate files.

n.cores

How many core to run this function. Default: from setConfigure.

prefix

prefix for Output file. Order: peak1, peak2, backgroud.

use.SavedPWM

use local motif position information. This data is download or generate by users. It must be a rds file (generated by R) and the information saved as GRanges. Once this parameter is used, parameters "genome", "motifPWM", "min.score", "n.cores" will be set to NULL.

...

Additional arguments, currently unused.

Details

This function scan motif position in a given genome regions.

Value

An invisible ATACProc-class object scalar for downstream analysis.

Author(s)

Wei Zhang

See Also

atacpeakComp

Examples



library(R.utils)
library(BSgenome.Hsapiens.UCSC.hg19)
p1bz <- system.file("extdata", "Example_peak1.bed.bz2", package="esATAC")
p2bz <- system.file("extdata", "Example_peak2.bed.bz2", package="esATAC")
peak1_path <- as.vector(bunzip2(filename = p1bz,
destname = file.path(getwd(), "Example_peak1.bed"),
ext="bz2", FUN=bzfile, overwrite=TRUE , remove = FALSE))
peak2_path <- as.vector(bunzip2(filename = p2bz,
destname = file.path(getwd(), "Example_peak2.bed"),
ext="bz2", FUN=bzfile, overwrite=TRUE, remove = FALSE))
peakcom.output <- peakcomp(bedInput1 = peak1_path, bedInput2 = peak2_path,
olap.rate = 0.1)

pwm <- readRDS(system.file("extdata", "motifPWM.rds", package="esATAC"))
## Not run: 
output <- atacMotifScanPair(atacProc = peakcom.output,
genome = BSgenome.Hsapiens.UCSC.hg19,
motifPWM = pwm)

## End(Not run)


[Package esATAC version 1.0.23 Index]