peaky: Obtains peak locations and heights based on adjusted...

In cqgd/pky: Capture Hi-C calling algorithm

Description

Fits additive models of peaks of varying strengths in various locations to the adjusted readcounts via RJMCMC, and stores these models on disk.

Usage

peaky(
  bait,
  omega_power,
  iterations = 1e+06,
  min_interactions = 20,
  log_file = NA
)

Arguments

bait	Data table containing the putative interactions of a bait, having the columns 'baitID', 'dist', and 'residual'. These report the bait ID, its distance to putative preys, and the adjusted readcounts for its interactions with them, respectively.
omega_power	Expected decay of adjusted read counts around a truly interacting prey. See details.
iterations	Number of models to parametrize. Greated numbers should lead to increased reproducibility.
log_file	Path to a log file.

Value

The output directory.

Details

The steepness of the function to be fitted to putative peaks is determined by ω according to β exp(-|ω * d|), where β represents peak height and d the distance from the center of the peak in bp.

Examples

base = system.file("extdata",package="peaky")
interactions_file = paste0(base,"/counts.tsv")
fragments_file = paste0(base,"/fragments.bed")
interactions = data.table::fread(interactions_file)
fragments = data.table::fread(fragments_file)
## Not run: 
BI = bin_interactions(interactions, fragments, bins=5)
models = by(BI$interactions, BI$interactions$dist.bin, model_bin, subsample_size=1000)
residuals = lapply(models, "[[", "residuals")
bins = split(BI$interactions, BI$interactions$dist.bin)

BTS = split_baits(bins, residuals)

peaky(BTS[baitID==618421], omega_power=-4.7, iterations=10e3)

## End(Not run)

cqgd/pky documentation built on Dec. 13, 2020, 3:32 a.m.