computeLogOdd: Compute the log odd ratio of RIP over background.

In yueli-compbio/RIPSeeker: RIPSeeker: a statistical package for identifying protein-associated transcripts from RIP-seq experiments

Description

The RIPScore is computed as the log odd ratio of the posterior for the RIP state (z_i = 2) over the posterior for the background state (z_i = 1)

Usage

computeLogOdd(nbhGR)

Arguments

nbhGR

GRanges of bins with the value slot saved for the posterior probabilities for the background and RIP state.

Details

To assess the statistical significance of the RIP predictions, we assign each bin a RIPScore defined as the log odd ratio of the posterior for the RIP state (z_i = 2) over the posterior for the background state (z_i = 1). When control is available, the RIPScore is updated as the difference between the RIPScores evaluated separately for RIP and control libraries. The scoring system captures the model confidence for the RIP state of each bin in the RIP library penalized by the false confidence for the RIP state of the same bin in the control library. In addition, RIPScore obviates scaling of read counts. Since sequencing depth usually differs between RIP and control libraries, scaling is necessary if the statistical score were derived from the read count differences. On the other hand, simplistic linear scaling may distort the data.

Value

A vector of log odd scores for each bin in nbhGR.

Author(s)

Yue Li

See Also

seekRIP, scoreMergedBins, logScoreWithoutControl, logScoreWithControl

Examples

# Retrieve system files
extdata.dir <- system.file("extdata", package="RIPSeeker") 

bamFiles <- list.files(extdata.dir, ".bam$", recursive=TRUE, full.names=TRUE)

bamFiles <- grep("PRC2", bamFiles, value=TRUE)

# Retrieve system files
extdata.dir <- system.file("extdata", package="RIPSeeker") 

bamFiles <- list.files(extdata.dir, ".bam$", recursive=TRUE, full.names=TRUE)

bamFiles <- grep("PRC2", bamFiles, value=TRUE)

alignGal <- getAlignGal(bamFiles[1], reverseComplement=TRUE, genomeBuild="mm9")

alignGR <- as(alignGal, "GRanges")

alignGRList <- GRangesList(as.list(split(alignGR, seqnames(alignGR))))

################ run main function for HMM inference on a single chromosome ################
nbhGR <- mainSeekSingleChrom(alignGR=alignGRList$chrX, K = 2, binSize=1e5)

ripscore <- computeLogOdd(nbhGR)

yueli-compbio/RIPSeeker documentation built on May 8, 2019, 2:34 a.m.