Title: | Set of tools to make V-plots and compute footprint profiles |
---|---|
Description: | The pattern of digestion and protection from DNA nucleases such as DNAse I, micrococcal nuclease, and Tn5 transposase can be used to infer the location of associated proteins. This package contains useful functions to analyze patterns of paired-end sequencing fragment density. VplotR facilitates the generation of V-plots and footprint profiles over single or aggregated genomic loci of interest. |
Authors: | Jacques Serizay [aut, cre] |
Maintainer: | Jacques Serizay <[email protected]> |
License: | GPL (>= 3) |
Version: | 1.17.0 |
Built: | 2024-10-31 06:42:04 UTC |
Source: | https://github.com/bioc/VplotR |
Genomic loci with a REB1 binding motifs according to http://jaspar.genereg.net/api/v1/matrix/MA0265.1.jaspar. PWM and scanning done with TFBSTools.
data(ABF1_sacCer3)
data(ABF1_sacCer3)
An object of class "GRanges"
.
Rossi, Lai & Pugh 2018 Genome Research
data(ABF1_sacCer3) ABF1_sacCer3
data(ABF1_sacCer3) ABF1_sacCer3
This function re-aligns ranges (typically regulatory elements) to a set of coordinates, either the TSS column or the TSS.fwd and TSS.rev columns. If none are found, the function assumes the ranges are promoters and that the end or the ranges are the TSSs.
alignToTSS(granges, upstream = 0, downstream = 1)
alignToTSS(granges, upstream = 0, downstream = 1)
granges |
A stranded GRanges object with a TSS column or TSS.rev and TSS.fwd columns |
upstream |
How many bases upstream of the TSS should the GRanges object by extended by? [Default: 0] |
downstream |
How many bases downstream of the TSS should the GRanges object by extended by? [Default: 1] |
GRanges aligned to the TSS column or to TSS.rev and TSS.fwd columns, and extended by upstream/downstream bp.
data(ce11_proms) ce11_proms alignToTSS(ce11_proms)
data(ce11_proms) ce11_proms alignToTSS(ce11_proms)
A sample of ATAC-seq fragments from individual worm tissues (Serizay et al. 2020, "Tissue-specific profiling reveals distinctive regulatory architectures for ubiquitous, germline and somatic genes", BiorXiv)
data(ATAC_ce11_Serizay2020)
data(ATAC_ce11_Serizay2020)
An object of class "list"
.
data(ATAC_ce11_Serizay2020) ATAC_ce11_Serizay2020
data(ATAC_ce11_Serizay2020) ATAC_ce11_Serizay2020
A .bam file sample
data(bam_test)
data(bam_test)
An object of class "GRanges"
.
data(bam_test) bam_test
data(bam_test) bam_test
Regulatory elements annotated in C. elegans (ce11) according to Serizay et al. 2020, "Tissue-specific profiling reveals distinctive regulatory architectures for ubiquitous, germline and somatic genes", BiorXiv.
data(ce11_all_REs)
data(ce11_all_REs)
GRanges
Serizay et al. 2020, "Tissue-specific profiling reveals distinctive regulatory architectures for ubiquitous, germline and somatic genes", BiorXiv. (DOI)
data(ce11_all_REs) table(ce11_all_REs$regulatory_class) table(ce11_all_REs$which.tissues)
data(ce11_all_REs) table(ce11_all_REs$regulatory_class) table(ce11_all_REs$which.tissues)
Promoters annotated in C. elegans (ce11) according to Serizay et al. 2020, "Tissue-specific profiling reveals distinctive regulatory architectures for ubiquitous, germline and somatic genes", BiorXiv.
data(ce11_proms)
data(ce11_proms)
An object of class "GRanges"
.
Serizay et al. 2020, "Tissue-specific profiling reveals distinctive regulatory architectures for ubiquitous, germline and somatic genes", BiorXiv. (DOI)
data(ce11_proms) table(ce11_proms$which.tissues)
data(ce11_proms) table(ce11_proms$which.tissues)
This function computes the underlying matrix shown as a heatmap in Vplots. For each pair of coordinates (x: distance from fragment midpoint to center of GRanges of interest; y: fragment size), the function computes how many fragments there are.
computeVmat( bam_granges, granges, cores = 1, xlims = c(-250, 250), ylims = c(50, 300) )
computeVmat( bam_granges, granges, cores = 1, xlims = c(-250, 250), ylims = c(50, 300) )
bam_granges |
GRanges, paired-end fragments |
granges |
GRanges, regions to map the fragments onto |
cores |
Integer, nb of threads to parallelize fragments subsetting |
xlims |
The x limits of the computed Vmat |
ylims |
The y limits of the computed Vmat |
A table object
data(bam_test) data(ce11_all_REs) Vmat <- computeVmat(bam_test, ce11_all_REs) dim(Vmat) Vmat[seq(1,5), seq(1,10)]
data(bam_test) data(ce11_all_REs) Vmat <- computeVmat(bam_test, ce11_all_REs) dim(Vmat) Vmat[seq(1,5), seq(1,10)]
high-score CTCF binding motifs, obtained from JASPAR
data(CTCF_hg38)
data(CTCF_hg38)
An object of class "GRanges"
.
data(CTCF_hg38) CTCF_hg38
data(CTCF_hg38) CTCF_hg38
This function splits bi-directional ranges into + and - stranded ranges. It duplicates the ranges which are '*'.
deconvolveBidirectionalPromoters(granges)
deconvolveBidirectionalPromoters(granges)
granges |
A stranded GRanges object |
GRanges with only '+' and '-' strands. GRanges with '*' strand have been duplicated and split into forward and reverse strands.
data(ce11_all_REs) library(GenomicRanges) proms <- ce11_all_REs[grepl('prom', ce11_all_REs$regulatory_class)] proms table(strand(proms)) proms <- deconvolveBidirectionalPromoters(proms) proms table(strand(proms))
data(ce11_all_REs) library(GenomicRanges) proms <- ce11_all_REs[grepl('prom', ce11_all_REs$regulatory_class)] proms table(strand(proms)) proms <- deconvolveBidirectionalPromoters(proms) proms table(strand(proms))
This function takes fragments and compute the distribution of their sizes over a set or multiple sets of GRanges.
getFragmentsDistribution( fragments, granges_list = NULL, extend_granges = c(-500, 500), limits = c(0, 600), roll = 3, cores = 1 )
getFragmentsDistribution( fragments, granges_list = NULL, extend_granges = c(-500, 500), limits = c(0, 600), roll = 3, cores = 1 )
fragments |
GRanges object containing paired-end fragments. See importPEBamFiles for more details on how to create such object. |
granges_list |
GRanges, can be a list of different sets of GRanges. |
extend_granges |
numeric vector of length 2, how the GRanges should be extended. |
limits |
numeric vector of length 2, only consider fragments within this window of sizes. |
roll |
Integer, apply a moving average of this size |
cores |
Integer, number of threads used to compute fragment size distribution |
A list of tbl, one for each .bam file.
data(bam_test) data(ce11_proms) df <- getFragmentsDistribution( bam_test, ce11_proms, extend_granges = c(-500, 500) ) head(df) which.max(df$y)
data(bam_test) data(ce11_proms) df <- getFragmentsDistribution( bam_test, ce11_proms, extend_granges = c(-500, 500) ) head(df) which.max(df$y)
This function takes bam file paths and read them into GRanges objects. Note: Can be quite lengthy for .bam files with 5+ millions fragments.
importPEBamFiles( files, genome = NULL, where = NULL, max_insert_size = 1000, shift_ATAC_fragments = FALSE, cores = 10, verbose = TRUE )
importPEBamFiles( files, genome = NULL, where = NULL, max_insert_size = 1000, shift_ATAC_fragments = FALSE, cores = 10, verbose = TRUE )
files |
character vector, each element of the vector is the path of an individual .bam file. |
genome |
character, genome ID (e.g. "sacCer3", "ce11", "dm6", "mm10" or "hg38"). |
where |
GRanges, only import the fragments mapping to the input GRanges (can fasten the import process a lot). |
max_insert_size |
Integer, filter out fragments larger than this size. |
shift_ATAC_fragments |
Boolean, if the fragments come from ATAC-seq, one might want to shift the extremities by +5 / -4 bp. |
cores |
Integer, number of cores to use when indexing bam files |
verbose |
Boolean |
A GRanges object containing fragments from the input .bam file.
bamfile <- system.file("extdata", "ex1.bam", package = "Rsamtools") fragments <- importPEBamFiles( bamfile, shift_ATAC_fragments = TRUE ) fragments
bamfile <- system.file("extdata", "ex1.bam", package = "Rsamtools") fragments <- importPEBamFiles( bamfile, shift_ATAC_fragments = TRUE ) fragments
A sample of MNase-seq fragments from yeast (Henikoff et al. 2011, "Epigenome characterization at single base-pair resolution", PNAS)
data(MNase_sacCer3_Henikoff2011)
data(MNase_sacCer3_Henikoff2011)
An object of class "GRanges"
.
data(MNase_sacCer3_Henikoff2011) MNase_sacCer3_Henikoff2011
data(MNase_sacCer3_Henikoff2011) MNase_sacCer3_Henikoff2011
A sample of fragments from multiple MNase-seq experiments performed in yeast (Henikoff et al. 2011, "Epigenome characterization at single base-pair resolution", PNAS), mapping over chrXV:186,400-187,400.
data(MNase_sacCer3_Henikoff2011_subset)
data(MNase_sacCer3_Henikoff2011_subset)
An object of class "GRanges"
.
data(MNase_sacCer3_Henikoff2011_subset) MNase_sacCer3_Henikoff2011_subset
data(MNase_sacCer3_Henikoff2011_subset) MNase_sacCer3_Henikoff2011_subset
This function normalizes a Vmat. Several different approaches have been implemented to normalize the Vmats.
normalizeVmat( Vmat, bam_granges, granges, normFun = c("zscore"), s = 0.99, roll = 1, verbose = TRUE )
normalizeVmat( Vmat, bam_granges, granges, normFun = c("zscore"), s = 0.99, roll = 1, verbose = TRUE )
Vmat |
A Vmat, usually output of computeVmat |
bam_granges |
GRanges, the paired-end fragments |
granges |
GRanges, the regions to map the fragments onto |
normFun |
character. A Vmat should be scaled either by:
|
s |
A float indicating which quantile to use if 'quantile' normalization is chosen |
roll |
integer, to use as the window to smooth the Vmat rows by rolling mean. |
verbose |
Boolean |
A normalized Vmat object
data(bam_test) data(ce11_all_REs) Vmat <- computeVmat(bam_test, ce11_all_REs) Vmat <- normalizeVmat( Vmat, bam_test, ce11_all_REs, normFun = c('libdepth+nloci') )
data(bam_test) data(ce11_all_REs) Vmat <- computeVmat(bam_test, ce11_all_REs) Vmat <- normalizeVmat( Vmat, bam_test, ce11_all_REs, normFun = c('libdepth+nloci') )
A function to compute nucleosome enrichment over a set of GRanges
nucleosomeEnrichment(x, ...)
nucleosomeEnrichment(x, ...)
x |
a GRanges or Vmat |
... |
additional parameters |
list
data(bam_test) data(ce11_proms) n <- nucleosomeEnrichment(bam_test, ce11_proms) n$fisher_test n$plot
data(bam_test) data(ce11_proms) n <- nucleosomeEnrichment(bam_test, ce11_proms) n$fisher_test n$plot
A function to compute nucleosome enrichment over a set of GRanges
## S3 method for class 'GRanges' nucleosomeEnrichment(x, granges, plus1_nuc_only = FALSE, verbose = TRUE, ...)
## S3 method for class 'GRanges' nucleosomeEnrichment(x, granges, plus1_nuc_only = FALSE, verbose = TRUE, ...)
x |
GRanges, paired-end fragments |
granges |
GRanges, loci to map the fragments onto |
plus1_nuc_only |
Boolean, should compute nucleosome enrichment only for +1 nucleosome? |
verbose |
Boolean |
... |
additional parameters |
list
data(bam_test) data(ce11_proms) n <- nucleosomeEnrichment(bam_test, ce11_proms) n$fisher_test n$plot
data(bam_test) data(ce11_proms) n <- nucleosomeEnrichment(bam_test, ce11_proms) n$fisher_test n$plot
A function to compute nucleosome enrichment over a Vmat
## S3 method for class 'Vmat' nucleosomeEnrichment(x, background, plus1_nuc_only = FALSE, ...)
## S3 method for class 'Vmat' nucleosomeEnrichment(x, background, plus1_nuc_only = FALSE, ...)
x |
a computed Vmat. Should be un-normalized. |
background |
a background Vmat. Should be un-normalized. |
plus1_nuc_only |
Boolean, should compute nucleosome enrichment only for +1 nucleosome? |
... |
additional parameters |
list
data(bam_test) data(ce11_proms) V <- plotVmat( bam_test, ce11_proms, normFun = '', return_Vmat = TRUE ) V_bg <- plotVmat( bam_test, sampleGRanges(ce11_proms), normFun = '', return_Vmat = TRUE ) n <- nucleosomeEnrichment(V, V_bg) n$fisher_test n$plot
data(bam_test) data(ce11_proms) V <- plotVmat( bam_test, ce11_proms, normFun = '', return_Vmat = TRUE ) V_bg <- plotVmat( bam_test, sampleGRanges(ce11_proms), normFun = '', return_Vmat = TRUE ) n <- nucleosomeEnrichment(V, V_bg) n$fisher_test n$plot
This function takes paired-end fragments, extract the "cuts" (i.e. extremities) and plot the footprint profile over a set of GRanges.
plotFootprint( frags, targets, split_strand = FALSE, plot_central = TRUE, xlim = c(-75, 75), bin = 1, verbose = 1 )
plotFootprint( frags, targets, split_strand = FALSE, plot_central = TRUE, xlim = c(-75, 75), bin = 1, verbose = 1 )
frags |
GRanges, the paired-end fragments |
targets |
GRanges, the loci to map the fragments onto |
split_strand |
Boolean, should the + and - strand be splitted? |
plot_central |
plot grey rectangle over the loci |
xlim |
numeric vector of length 2, the x limits of the computed Vmat |
bin |
Integer, bin used to smooth the gootprint profile |
verbose |
Integer |
A footprint ggplot
data(bam_test) data(ce11_proms) plotFootprint(bam_test, ce11_proms)
data(bam_test) data(ce11_proms) plotFootprint(bam_test, ce11_proms)
The paired-end fragments overlapping a locus of interest (e.g., binding sites, provided in the 'loci' argument) are shown in red while the remaining fragments mapping to the genomic window are displayed in black. Marginal curves are also plotted on the side of the distribution plot. They highlight the smoothed distribution of the position of paired-end fragment midpoints (top) or of the paired-end fragment length (right)
plotProfile( fragments, window = loc, loci = NULL, annots = NULL, min = 50, max = 200, alpha = 0.5, size = 1, with_densities = TRUE, verbose = TRUE )
plotProfile( fragments, window = loc, loci = NULL, annots = NULL, min = 50, max = 200, alpha = 0.5, size = 1, with_densities = TRUE, verbose = TRUE )
fragments |
GRanges |
window |
character, chromosome location |
loci |
GRanges, optional genomic locus. Fragments overlapping this locus will be in red. |
annots |
GRanges, optional gene annotations |
min |
integer, minimum fragment size |
max |
integer, maximum fragment size |
alpha |
float, transparency value |
size |
float, dot size |
with_densities |
Boolean, should the densities be plotted? |
verbose |
Boolean |
A ggplot
data(bam_test) data(ce11_proms) V <- plotProfile( bam_test, 'chrI:10000-12000', loci = ce11_proms, min = 80, max = 200 )
data(bam_test) data(ce11_proms) V <- plotProfile( bam_test, 'chrI:10000-12000', loci = ce11_proms, min = 80, max = 200 )
See individual methods for further detail
plotVmat(x, ...)
plotVmat(x, ...)
x |
GRanges or list or Vmat |
... |
additional parameters |
A Vmat ggplot
data(bam_test) data(ce11_proms) V <- plotVmat( bam_test, ce11_proms, normFun = 'libdepth+nloci' )
data(bam_test) data(ce11_proms) V <- plotVmat( bam_test, ce11_proms, normFun = 'libdepth+nloci' )
The default plotVmat method generates a ggplot representing a heatmap of fragment density.
## Default S3 method: plotVmat( x, hm = 90, colors = COLORSCALE_VMAT, breaks = NULL, xlim = c(-250, 250), ylim = c(50, 300), main = "", xlab = "Distance from center of elements", ylab = "Fragment length", key = "Score", ... )
## Default S3 method: plotVmat( x, hm = 90, colors = COLORSCALE_VMAT, breaks = NULL, xlim = c(-250, 250), ylim = c(50, 300), main = "", xlab = "Distance from center of elements", ylab = "Fragment length", key = "Score", ... )
x |
A computed Vmat (ideally, should be normalized) |
hm |
Integer, should be between 0 and 100. Used to automatically scale the range of colors (best to keep between 90 and 100) |
colors |
a vector of colors |
breaks |
a vector of breaks. length(breaks) == length(colors) + 1 |
xlim |
vector of two integers, x limits |
ylim |
vector of two integers, y limits |
main |
character, title of the plot |
xlab |
character, x-axis label |
ylab |
character, y-axis label |
key |
character, legend label |
... |
additional parameters |
A Vmat ggplot
data(bam_test) data(ce11_proms) V <- plotVmat( bam_test, ce11_proms, normFun = 'libdepth+nloci', return_Vmat = TRUE ) plotVmat(V)
data(bam_test) data(ce11_proms) V <- plotVmat( bam_test, ce11_proms, normFun = 'libdepth+nloci', return_Vmat = TRUE ) plotVmat(V)
The plotVmat.GRanges() method computes and normalizes a Vmat before passing it to plotVmat.Vmat() method.
## S3 method for class 'GRanges' plotVmat( x, granges, xlims = c(-250, 250), ylims = c(50, 300), normFun = "", s = 0.95, roll = 3, cores = 1, return_Vmat = FALSE, verbose = 1, ... )
## S3 method for class 'GRanges' plotVmat( x, granges, xlims = c(-250, 250), ylims = c(50, 300), normFun = "", s = 0.95, roll = 3, cores = 1, return_Vmat = FALSE, verbose = 1, ... )
x |
GRanges, paired-end fragments |
granges |
GRanges, loci to map the fragments onto |
xlims |
x limits of the computed Vmat |
ylims |
y limits of the computed Vmat |
normFun |
character. A Vmat should be scaled either by:
|
s |
A float indicating which quantile to use if 'quantile' normalization is chosen |
roll |
integer, to use as the window to smooth the Vmat rows by rolling mean. |
cores |
Integer, number of threads to parallelize fragments subsetting |
return_Vmat |
Boolean, should the function return the computed Vmat rather than the plot? |
verbose |
Boolean |
... |
additional parameters |
A Vmat ggplot
data(bam_test) data(ce11_proms) V <- plotVmat( bam_test, ce11_proms, normFun = 'libdepth+nloci', roll = 5 )
data(bam_test) data(ce11_proms) V <- plotVmat( bam_test, ce11_proms, normFun = 'libdepth+nloci', roll = 5 )
The plotVmat.GRanges() method computes and normalizes multiple Vmats before passing them to plotVmat.VmatList() method.
## S3 method for class 'list' plotVmat( x, cores = 1, cores_subsetting = 1, nrow = NULL, ncol = NULL, xlims = c(-250, 250), ylims = c(50, 300), normFun = "libdepth+nloci", s = 0.95, roll = 3, return_Vmat = FALSE, verbose = 1, ... )
## S3 method for class 'list' plotVmat( x, cores = 1, cores_subsetting = 1, nrow = NULL, ncol = NULL, xlims = c(-250, 250), ylims = c(50, 300), normFun = "libdepth+nloci", s = 0.95, roll = 3, return_Vmat = FALSE, verbose = 1, ... )
x |
list Each element of the list should be a list containing paired-end fragments and GRanges of interest. |
cores |
Integer, number of cores to parallelize the plots |
cores_subsetting |
Integer, number of threads to parallelize fragments subsetting |
nrow |
Integer, how many rows in facet? |
ncol |
Integer, how many cols in facet? |
xlims |
x limits of the computed Vmat |
ylims |
y limits of the computed Vmat |
normFun |
character. A Vmat should be scaled either by:
|
s |
A float indicating which quantile to use if 'quantile' normalization is chosen |
roll |
integer, to use as the window to smooth the Vmat rows by rolling mean. |
return_Vmat |
Boolean, should the function return the computed Vmat rather than the plot? |
verbose |
Boolean |
... |
additional parameters |
A list of Vmat ggplots
data(bam_test) data(ce11_proms) list_params <- list( 'germline' = list( bam_test, ce11_proms[ce11_proms$which.tissues == 'Germline'] ), 'muscle' = list( bam_test, ce11_proms[ce11_proms$which.tissues == 'Muscle'] ) ) V <- plotVmat( list_params, normFun = 'libdepth+nloci', roll = 5 )
data(bam_test) data(ce11_proms) list_params <- list( 'germline' = list( bam_test, ce11_proms[ce11_proms$which.tissues == 'Germline'] ), 'muscle' = list( bam_test, ce11_proms[ce11_proms$which.tissues == 'Muscle'] ) ) V <- plotVmat( list_params, normFun = 'libdepth+nloci', roll = 5 )
The plotVmat.Vmat() method forwards the Vmat to plotVmat.default().
## S3 method for class 'Vmat' plotVmat(x, ...)
## S3 method for class 'Vmat' plotVmat(x, ...)
x |
A computed Vmat (ideally, should be normalized) |
... |
additional parameters |
A Vmat ggplot
data(bam_test) data(ce11_proms) V <- plotVmat( bam_test, ce11_proms, normFun = 'libdepth+nloci', return_Vmat = TRUE ) plotVmat(V)
data(bam_test) data(ce11_proms) V <- plotVmat( bam_test, ce11_proms, normFun = 'libdepth+nloci', return_Vmat = TRUE ) plotVmat(V)
The plotVmat.VmatList() method forwards the Vmat to plotVmat.default().
## S3 method for class 'VmatList' plotVmat(x, nrow = NULL, ncol = NULL, dir = "v", ...)
## S3 method for class 'VmatList' plotVmat(x, nrow = NULL, ncol = NULL, dir = "v", ...)
x |
A VmatList (output of plotVmat.list()) |
nrow |
Integer, how many rows in facet? |
ncol |
Integer, how many cols in facet? |
dir |
str, direction of facets? |
... |
additional parameters |
A Vmat ggplot
data(bam_test) data(ce11_proms) list_params <- list( 'germline' = list( bam_test, ce11_proms[ce11_proms$which.tissues == 'Germline'] ), 'muscle' = list( bam_test, ce11_proms[ce11_proms$which.tissues == 'Muscle'] ) ) V <- plotVmat( list_params, normFun = 'libdepth+nloci', roll = 5 )
data(bam_test) data(ce11_proms) list_params <- list( 'germline' = list( bam_test, ce11_proms[ce11_proms$which.tissues == 'Germline'] ), 'muscle' = list( bam_test, ce11_proms[ce11_proms$which.tissues == 'Muscle'] ) ) V <- plotVmat( list_params, normFun = 'libdepth+nloci', roll = 5 )
Genomic loci with a REB1 binding motifs according to http://jaspar.genereg.net/api/v1/matrix/MA0363.1.jaspar. PWM and scanning done with TFBSTools.
data(REB1_sacCer3)
data(REB1_sacCer3)
An object of class "GRanges"
.
Rossi, Lai & Pugh 2018 Genome Research
data(REB1_sacCer3) REB1_sacCer3
data(REB1_sacCer3) REB1_sacCer3
This function takes a given GRanges and returns another GRanges object. The new GRanges has the same number of ranges and the same chromosome, width and strand distributions than the original GRanges.
sampleGRanges( x, n = NULL, width = NULL, exclude = FALSE, avoid_overlap = FALSE )
sampleGRanges( x, n = NULL, width = NULL, exclude = FALSE, avoid_overlap = FALSE )
x |
GRanges object |
n |
Integer, number of sampled GRanges |
width |
Integer, width of sampled GRanges |
exclude |
Boolean, should the original GRanges be excluded? |
avoid_overlap |
Boolean, should the sampled GRanges not be overlapping? |
A GRanges object of length n
data(ce11_proms) sampleGRanges(ce11_proms, 100)
data(ce11_proms) sampleGRanges(ce11_proms, 100)
This function takes a given GRanges and returns another GRanges object. The new GRanges has the same number of ranges and the same chromosome, width and strand distributions than the original GRanges.
## S3 method for class 'GRanges' sampleGRanges( x, n = NULL, width = NULL, exclude = FALSE, avoid_overlap = FALSE )
## S3 method for class 'GRanges' sampleGRanges( x, n = NULL, width = NULL, exclude = FALSE, avoid_overlap = FALSE )
x |
GRanges object |
n |
Integer, number of sampled GRanges |
width |
Integer, width of sampled GRanges |
exclude |
Boolean, should the original GRanges be excluded? |
avoid_overlap |
Boolean, should the sampled GRanges not be overlapping? |
A GRanges object of length n
data(ce11_proms) sampleGRanges(ce11_proms, 100)
data(ce11_proms) sampleGRanges(ce11_proms, 100)
A function to shift GRanges fragments by 5/-4. This is useful when dealing with fragments coming from ATAC-seq.
shiftATACGranges(g, pos_shift = 4, neg_shift = 5)
shiftATACGranges(g, pos_shift = 4, neg_shift = 5)
g |
GRanges of ATAC-seq fragments |
pos_shift |
Integer. How many bases should fragments on direct strand be shifted by? |
neg_shift |
Integer. How many bases should fragments on negative strand be shifted by? |
A GRanges object containing fragments from the input .bam file.
data(bam_test) shiftATACGranges(bam_test)
data(bam_test) shiftATACGranges(bam_test)
This function works on a Vmat (the output of computeVmat()). It shuffles the matrix to randomize the fragment densities.
shuffleVmat(Vmat)
shuffleVmat(Vmat)
Vmat |
A Vmat, usually output of computeVmat |
A shuffled Vmat object
data(bam_test) data(ce11_all_REs) Vmat <- computeVmat(bam_test, ce11_all_REs) Vmat <- shuffleVmat(Vmat)
data(bam_test) data(ce11_all_REs) Vmat <- computeVmat(bam_test, ce11_all_REs) Vmat <- shuffleVmat(Vmat)
Personal ggplot2 theming function, adapted from roboto-condensed at https://github.com/hrbrmstr/hrbrthemes/
theme_ggplot2( grid = TRUE, border = TRUE, base_family = "", base_size = 8, plot_title_family = base_family, plot_title_size = 12, plot_title_face = "plain", plot_title_margin = 5, subtitle_size = 11, subtitle_face = "plain", subtitle_margin = 5, strip_text_family = base_family, strip_text_size = 10, strip_text_face = "bold", caption_size = 9, caption_face = "plain", caption_margin = 3, axis_text_size = base_size, axis_title_family = base_family, axis_title_size = 9, axis_title_face = "plain", axis_title_just = "rt", panel_spacing = grid::unit(2, "lines"), grid_col = "#cccccc", plot_margin = margin(12, 12, 12, 12), axis_col = "#cccccc", axis = FALSE, ticks = FALSE )
theme_ggplot2( grid = TRUE, border = TRUE, base_family = "", base_size = 8, plot_title_family = base_family, plot_title_size = 12, plot_title_face = "plain", plot_title_margin = 5, subtitle_size = 11, subtitle_face = "plain", subtitle_margin = 5, strip_text_family = base_family, strip_text_size = 10, strip_text_face = "bold", caption_size = 9, caption_face = "plain", caption_margin = 3, axis_text_size = base_size, axis_title_family = base_family, axis_title_size = 9, axis_title_face = "plain", axis_title_just = "rt", panel_spacing = grid::unit(2, "lines"), grid_col = "#cccccc", plot_margin = margin(12, 12, 12, 12), axis_col = "#cccccc", axis = FALSE, ticks = FALSE )
grid |
panel grid ('TRUE', 'FALSE', or a combination of 'X', 'x', 'Y', 'y') |
border |
border if 'TRUE' add border |
base_family , base_size
|
base font family and size |
plot_title_family , plot_title_face
|
plot title family, face |
plot_title_size , plot_title_margin
|
plot title size and margin |
subtitle_face , subtitle_size
|
plot subtitle family, face and size |
subtitle_margin |
plot subtitle margin bottom (single numeric value) |
strip_text_family , strip_text_face , strip_text_size
|
facet label font family, face and size |
caption_face , caption_size , caption_margin
|
plot caption family, face, size and margin |
axis_text_size |
font size of axis text |
axis_title_family , axis_title_face , axis_title_size
|
axis title font family, face and size |
axis_title_just |
axis title font justificationk one of '[blmcrt]' |
panel_spacing |
panel spacing (use 'unit()') |
grid_col |
grid color |
plot_margin |
plot margin (specify with [ggplot2::margin]) |
axis_col |
axis color |
axis |
add x or y axes? 'TRUE', 'FALSE', "'xy'" |
ticks |
ticks if 'TRUE' add ticks |
theme A ggplot theme
library(ggplot2) ggplot(mtcars, aes(mpg, wt)) + geom_point() + labs(x="Fuel effiency (mpg)", y="Weight (tons)", title="Seminal ggplot2 scatterplot example") + theme_ggplot2()
library(ggplot2) ggplot(mtcars, aes(mpg, wt)) + geom_point() + labs(x="Fuel effiency (mpg)", y="Weight (tons)", title="Seminal ggplot2 scatterplot example") + theme_ggplot2()