Deduct noise by comparing two signal tracks in bedGraph. Note: All regions on the same chromosome in the bedGraph file should be continuous so only bedGraph files from MACS3 are accpetable.

bdgcmp(
  tfile,
  cfile,
  sfactor = 1,
  pseudocount = 0,
  method = c("ppois", "qpois", "subtract", "logFE", "FE", "logLR", "slogLR", "max"),
  oprefix = character(),
  outputfile = list(),
  outdir = ".",
  log = TRUE,
  verbose = 2L
)

Arguments

tfile

Treatment bedGraph file, e.g. *_treat_pileup.bdg from MACSv2. REQUIRED

cfile

Control bedGraph file, e.g. *_control_lambda.bdg from MACSv2. REQUIRED

sfactor

Scaling factor for treatment and control track. Keep it as 1.0 or default in most cases. Set it ONLY while you have SPMR output from MACS3 callpeak, and plan to calculate scores as MACS3 callpeak module. If you want to simulate 'callpeak' w/o '--to-large', calculate effective smaller sample size after filtering redudant reads in million (e.g., put 31.415926 if effective reads are 31,415,926) and input it for '-S'; for 'callpeak --to-large', calculate effective reads in larger sample. DEFAULT: 1.0

pseudocount

The pseudocount used for calculating logLR, logFE or FE. The count will be applied after normalization of sequencing depth. DEFAULT: 0.0, no pseudocount is applied.

method

Method to use while calculating a score in any bin by comparing treatment value and control value. Available choices are: ppois, qpois, subtract, logFE, logLR, and slogLR. They represent Poisson Pvalue (-log10(pvalue) form) using control as lambda and treatment as observation, q-value through a BH process for poisson pvalues, subtraction from treatment, linear scale fold enrichment, log10 fold enrichment(need to set pseudocount), log10 likelihood between ChIP-enriched model and open chromatin model(need to set pseudocount), symmetric log10 likelihood between two ChIP-enrichment models, or maximum value between the two tracks. Default option is ppois.",default="ppois".

oprefix

The PREFIX of output bedGraph file to write scores. If it is given as A, and method is 'ppois', output file will be A_ppois.bdg. Mutually exclusive with -o/--ofile.

outputfile

Output filename. Mutually exclusive with --o-prefix. The number and the order of arguments for --ofile must be the same as for -m.

outdir

The output directory.

log

Whether to capture logs.

verbose

Set verbose level. 0: only show critical message, 1: show additional warning message, 2: show process information, 3: show debug messages. If you want to know where are the duplicate reads, use 3. DEFAULT:2

Value

macsList object.

Examples

eh <- ExperimentHub::ExperimentHub()
CHIP <- eh[["EH4558"]]
#> see ?MACSdata and browseVignettes('MACSdata') for documentation
#> loading from cache
CTRL <- eh[["EH4563"]]
#> see ?MACSdata and browseVignettes('MACSdata') for documentation
#> loading from cache
p1 <- pileup(CHIP, outdir = tempdir(),
             outputfile = "pileup_ChIP_bed.bdg", format = "BED")
#> INFO  @ 17 Nov 2023 13:42:08: [706 MB] # Existing file b'/tmp/RtmpEPj3Vo/pileup_ChIP_bed.bdg' will be replaced! 
#> INFO  @ 17 Nov 2023 13:42:08: [706 MB] # read alignment files... 
#> INFO  @ 17 Nov 2023 13:42:08: [706 MB] # read treatment tags... 
#> INFO  @ 17 Nov 2023 13:42:08: [708 MB] tag size is determined as 101 bps 
#> INFO  @ 17 Nov 2023 13:42:08: [708 MB] # tag size = 101 
#> INFO  @ 17 Nov 2023 13:42:08: [708 MB] # total tags in alignment file: 49622 
#> INFO  @ 17 Nov 2023 13:42:08: [708 MB] # Pileup alignment file, extend each read towards downstream direction with 200 bps 
#> INFO  @ 17 Nov 2023 13:42:08: [708 MB] # Done! Check pileup_ChIP_bed.bdg 
#> 
p2 <- pileup(CTRL, outdir = tempdir(),
             outputfile = "pileup_CTRL_bed.bdg", format = "BED")
#> INFO  @ 17 Nov 2023 13:42:08: [708 MB] # Existing file b'/tmp/RtmpEPj3Vo/pileup_CTRL_bed.bdg' will be replaced! 
#> INFO  @ 17 Nov 2023 13:42:08: [708 MB] # read alignment files... 
#> INFO  @ 17 Nov 2023 13:42:08: [708 MB] # read treatment tags... 
#> INFO  @ 17 Nov 2023 13:42:08: [708 MB] tag size is determined as 101 bps 
#> INFO  @ 17 Nov 2023 13:42:08: [708 MB] # tag size = 101 
#> INFO  @ 17 Nov 2023 13:42:08: [708 MB] # total tags in alignment file: 50837 
#> INFO  @ 17 Nov 2023 13:42:08: [708 MB] # Pileup alignment file, extend each read towards downstream direction with 200 bps 
#> INFO  @ 17 Nov 2023 13:42:08: [708 MB] # Done! Check pileup_CTRL_bed.bdg 
#> 
c1 <- bdgcmp(p1$outputs, p2$outputs, outdir = tempdir(),
             oprefix = "bdgcmp", pseudocount = 1, method = "FE")
#> INFO  @ 17 Nov 2023 13:42:08: [708 MB] Read and build treatment bedGraph... 
#> INFO  @ 17 Nov 2023 13:42:08: [708 MB] Read and build control bedGraph... 
#> INFO  @ 17 Nov 2023 13:42:08: [708 MB] Build ScoreTrackII... 
#> INFO  @ 17 Nov 2023 13:42:08: [708 MB] Calculate scores comparing treatment and control by 'FE'... 
#> INFO  @ 17 Nov 2023 13:42:09: [708 MB] Write bedGraph of scores... 
#> INFO  @ 17 Nov 2023 13:42:09: [708 MB] Finished 'FE'! Please check '/tmp/RtmpEPj3Vo/bdgcmp_FE.bdg'! 
#>