circularRandomizeRegions: Circular Randomize Regions

View source: R/circularRandomizeRegions.R

circularRandomizeRegionsR Documentation

Circular Randomize Regions

Description

Given a set of regions A and a genome, this function returns a new set of regions created by applying a random spin to each chromosome.

Usage

circularRandomizeRegions(A, genome="hg19", mask=NULL, max.mask.overlap=NULL, max.retries=10, verbose=TRUE, ...)

Arguments

A

The set of regions to randomize. A region set in any of the accepted formats by toGRanges (GenomicRanges, data.frame, etc...)

genome

The reference genome to use. A valid genome object. Either a GenomicRanges or data.frame containing one region per whole chromosome or a character uniquely identifying a genome in BSgenome (e.g. "hg19", "mm10" but not "hg"). Internally it uses getGenomeAndMask.

mask

The set of regions specifying where a random region can not be (centromeres, repetitive regions, unmappable regions...). A region set in any of the accepted formats by toGRanges (GenomicRanges,data.frame, ...). If NULL it will try to derive a mask from the genome (currently only works is the genome is a character string) and if NA it will explicitly give an empty mask.

max.mask.overlap

numeric value

max.retries

numeric value

verbose

a boolean.

...

further arguments to be passed to or from methods.

Details

This randomization strategy is useful when the spatial relation between the regions in the RS is important and has to be conserved.

Value

It returns a GenomicRanges object with the regions resulting from the randomization process.

See Also

randomizeRegions, toDataframe, toGRanges, getGenome, getMask, getGenomeAndMask, characterToBSGenome, maskFromBSGenome, resampleRegions, createRandomRegions

Examples

A <- data.frame("chr1", c(1, 10, 20, 30), c(12, 13, 28, 40))

mask <- data.frame("chr1", c(20000000, 100000000), c(22000000, 130000000))

genome <- data.frame(c("chr1", "chr2"), c(1, 1), c(180000000, 20000000))

circularRandomizeRegions(A)

circularRandomizeRegions(A, genome=genome, mask=mask, per.chromosome=TRUE, non.overlapping=TRUE)


bernatgel/regioneR documentation built on Sept. 10, 2023, 12:03 a.m.