htSeqTools (version 1.18.0)

alignPeaks: Align peaks in a ChIP-Seq experiment by removing the strand specific bias.

Description

Align peaks in a ChIP-Seq experiment by removing the shift between reads aligned to the plus and the minus strands.

Usage

alignPeaks(x, strand, npeaks = 1000, bandwidth = 150, mc.cores=1)

Arguments

x

A RangedDataList, RangedData or an IRangesList object containing the aligned reads in each chromosome.

strand

Strand that each read was aligned to. If x is of class RangedDataList, strand can be a character vector of length 1 indicating the name of the field in x indicating the strand, i.e. x[[1]][[strand]] contains the strand information.

npeaks

Number of peaks to be used to estimate the shift size.

bandwidth

Only reads with distance less than bandwidth between them and their closest gene are used to estimate the shift size.

mc.cores

Number of cores to be used for parallel computing (passed on to mclapply). Only used if x is of class RangedDataList.

Value

A CompressedIRangesList object with all reads shifted so that the strand specific bias is no longer present.

Methods

signature(x = "IRangesList", strand = "list")

Each element in x corresponds to a chromosome, and each range gives the start/end of a sequence. strand indicates the strand for the ranges in x.

signature(x = "RangedData", strand = "character")

x gives read start and end positions, and strand gives the name of the variable in values(x) containing the strand information.

signature(x = "RangedDataList", strand = "character")

The method for RangedData is applied to each element in x separately, as each element may have a different strand-specific bias.

Details

The procedure detects the npeaks highest peaks (using reads from both strands simultaneously). Then it selects reads which are less than bandwidth base pairs away from any of the peaks. Then it computes (a) the average distance between reads on the plus strand and the closest peak, (b) the same distance for reads on the minus strand. The mean difference between (a) and (b) is the estimated shift size. Reads on the plus strand are shifted to the right, whereas reads on the minus strands are shifted to the left.

Examples

#Generate 1000 reads containing strand-specific bias
st <- runif(1000,1,250)
strand <- rep(c('+','-'),each=500)
st[strand=='-'] <- st[strand=='-'] + runif(500,50,100)
x <- RangedData(IRanges(st,st+38),strand=strand)
#Estimate and remove the bias
xalign <- alignPeaks(x, strand='strand', npeaks=1)