normR-diffR: Difference Calling between conditional ChIP-seq data in normR...
In your-highness/normR: Normalization and difference calling in ChIP-seq data

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

Difference calling between treatment (ChIP-seq 1) and control (ChIP-seq 2) in normR is done by fitting background and two conditional enrichments simultaenously. Therefore, a mixture of three binomials is fit to the data with Expectation Maximization (EM). After convergence of the EM, the fitted background component is used to calculate significance for treatment and control count pair. Based on this statistic, user can extract significantly enriched/depleted regions in a condition with a desired significance level. These regions can be further analyzed within R or exported (see NormRFit-class). Furthermore, diffR calculates a standardized conditional-specific enrichment given the fitted background component. See also Details

diffR(treatment, control, genome, ...)

## S4 method for signature 'integer,integer,GenomicRanges'
diffR(treatment, control, genome,
  procs = 1L, verbose = TRUE, eps = 1e-05, iterations = 10,
  minP = 0.05)

## S4 method for signature 'character,character,GenomicRanges'
diffR(treatment, control, genome,
  countConfig = countConfigSingleEnd(), procs = 1L, verbose = TRUE,
  eps = 1e-05, iterations = 10, minP = 0.05)

## S4 method for signature 'character,character,data.frame'
diffR(treatment, control, genome,
  countConfig = countConfigSingleEnd(), procs = 1L, verbose = TRUE,
  eps = 1e-05, iterations = 10, minP = 0.05)

## S4 method for signature 'character,character,character'
diffR(treatment, control,
  genome = "", countConfig = countConfigSingleEnd(), procs = 1L,
  verbose = TRUE, eps = 1e-05, iterations = 10, minP = 0.05)

`treatment`	An `integer` vector of treatment counts or a `character` pointing to the treatment bam file. In the latter case an "`treatment`.bai" index file should exist in the same folder.
`control`	An `integer` vector of control counts or a character pointing to the control bam file. In the latter case an "`control`.bai" index file should exist in the same folder.
`genome`	Either `NULL` (default), a `character` specifying a USCS genome identifier, a data.frame consisting of two columns or a GenomicRanges specifying the genomic regions (see Details).
`...`	Optional arguments for the respective implementations of `diffR`.
`procs`	An `integer` giving the number of parallel threads to use.
`verbose`	A `logical` indicating whether verbose output is desired.
`eps`	A `numeric` specifying the T Filter threshold and the threshold for EM convergence, i.e. the minimal difference in log-likelihood in two consecutive steps.
`iterations`	An `integer` specifying how many times the EM is initialized with random model parameters.
`minP`	An `integer` controlling the threshold for the T method when filtering low power regions, i.e. regions with low counts.
`countConfig`	A `NormRCountConfig` object specifying bam counting parameters for read count retrieval. See Details.

Supplied count vectors for treatment and control should be of same length and of type integer.

For convenience, read count vectors can be obtained directly from bam files. In this case, please specify a bam file for treatment and control each and a genome. Bam files should be indexed using samtools (i.e. samtools index file file.bai). Furthermore, bam files should contain a valid header with given chromosome names. If genome == NULL(default), chromosome names will be read from treatment bamheader. Please be aware that bamheader might contain irregular contigs and chrM which influence the fit. Also be sure that treatment and control contain the same chromosomes. Otherwise an error will be thrown. If genome is a character, fetchExtendedChromInfoFromUCSC is used to resolve this to a valid UCSC genome identifier (see https://genome.ucsc.edu/cgi-bin/hgGateway for available genomes). In this case, only assembled molecules will be considered (no circular). Please check if your bam files obey this annotation. If genome is a data.frame, it represents the chromosome specification. The first column will be used as chromosome ID and the second column will be used as the chromosome lengths. If genome is a GenomicRanges, it should contain the equally sized genomic loci to count in, e.g. promoters. The binsize in the supplied NormRCountConfig is ignore in this case.

bamCountConfig is an instance of class NormRCountConfig specifying settings for read counting on bam files. You can specify the binsize, minimum mapping quality, shifting of read ends etc.. Please refer to NormRFit-class for details.

A NormRFit container holding results of the fit with type diffR.

Johannes Helmuth helmuth@molgen.mpg.de

NormRFit-class for functions on accessing and exporting the diffR fit. NormRCountConfig-class for configuration of the read counting procedure (binsize, mapping quality,...).

require(GenomicRanges)

### enrichR(): Calling Enrichment over Input
#load some example bamfiles
input <- system.file("extdata", "K562_Input.bam", package="normr")
chipK4 <- system.file("extdata", "K562_H3K4me3.bam", package="normr")
#region to count in (example files contain information only in this region)
gr <- GRanges("chr1", IRanges(seq(22500001, 25000000, 1000), width = 1000))
#configure your counting strategy (see BamCountConfig-class)
countConfiguration <- countConfigSingleEnd(binsize = 1000,
                                           mapq = 30, shift = 100)
#invoke enrichR to call enrichment
enrich <- enrichR(treatment = chipK4, control = input,
                  genome = gr,  countConfig = countConfiguration,
                  iterations = 10, procs = 1, verbose = TRUE)
#inspect the fit
enrich
summary(enrich)

## Not run:
#write significant regions to bed
#exportR(enrich, filename = "enrich.bed", fdr = 0.01)
#write normalized enrichment to bigWig
#exportR(enrich, filename = "enrich.bw")
## End(**Not run**)

### diffR(): Calling differences between two conditions
chipK36 <- system.file("extdata", "K562_H3K36me3.bam", package="normr")
diff <- diffR(treatment = chipK36, control = chipK4,
              genome = gr,  countConfig = countConfiguration,
              iterations = 10, procs = 1, verbose = TRUE)
summary(diff)

### regimeR(): Identification of broad and peak enrichment
regime <- regimeR(treatment = chipK36, control = input, models = 3,
                  genome = gr,  countConfig = countConfiguration,
                  iterations = 10, procs = 1, verbose = TRUE)
summary(regime)