## ---- echo=FALSE--------------------------------------------------------------
# These settings make the vignette prettier
knitr::opts_chunk$set(results="hold", collapse=FALSE, message=FALSE)
#refreshPackage("GenomicDistributions")
#devtools::build_vignettes("code/GenomicDistributions")
#devtools::test("code/GenomicDistributions")

## ---- eval=FALSE--------------------------------------------------------------
#  if (!requireNamespace("BiocManager", quietly = TRUE))
#  install.packages("BiocManager")
#  BiocManager::install("GenomicDistributions")

## ---- echo=TRUE, results="hide", message=FALSE, warning=FALSE-----------------
library("GenomicDistributions")
queryFile = system.file("extdata", "vistaEnhancers.bed.gz", package="GenomicDistributions")
query = rtracklayer::import(queryFile)

## ----chrom-plots-single-------------------------------------------------------
# First, calculate the distribution:
x = calcChromBinsRef(query, "hg19")

# Then, plot the result:
plotChromBins(x)

## ----chrom-plots-multi--------------------------------------------------------
# Let's fudge a second region set by shifting the first one over 
query2 = GenomicRanges::shift(query, 1e6)
queryList = GRangesList(vistaEnhancers=query, shifted=query2)
x2 = calcChromBinsRef(queryList, "hg19")
plotChromBins(x2)

## ----tss-distribution, fig.cap="TSS plot. Distribution of query regions relative to TSSs", fig.small=TRUE----
# Calculate the distances:
TSSdist = calcFeatureDistRefTSS(query, "hg19")

# Then plot the result:
plotFeatureDist(TSSdist, featureName="TSS")

## ----tss-distribution-multi, fig.cap="TSS plots with multiple region sets"----

TSSdist2 = calcFeatureDistRefTSS(queryList, "hg19")
plotFeatureDist(TSSdist2, featureName="TSS")


## ----tiled-distance-plot, fig.cap="Tiled feature distance plot. Distribution of query regions relative to arbitrary features", fig.small=TRUE----
plotFeatureDist(TSSdist2, featureName="TSS", tile=TRUE)

## ----Build features-----------------------------------------------------------
featureExample = GenomicRanges::shift(query, round(rnorm(length(query), 0,1000)))

## ----distance-to-features-plot, fig.cap="Feature distance plot. Distribution of query regions relative to arbitrary features", fig.small=TRUE----
fdd = calcFeatureDist(query, featureExample)
plotFeatureDist(fdd)

## ----percentage-partition-plot, fig.cap="Partition distribution plot. Percentage distribution of query regions across genomic features", fig.small=TRUE----
gp = calcPartitionsRef(query, "hg19")
plotPartitions(gp)

## ----multiple-percentage-partition-plot, fig.cap="Partition distribution plot for multiple query region sets.", fig.small=TRUE----
gp2 = calcPartitionsRef(queryList, "hg19")
plotPartitions(gp2)

## ----multiple-raw-partition-plot, fig.cap="Raw partition distribution plot for multiple regionsets", fig.small=TRUE----
# Plot the results:
plotPartitions(gp2, numbers=TRUE)

## ----corrected-partition-plot, fig.cap="Expected partition distribution plot. Distribution of query regions across genomic features relative to the expected distribution of those features.", fig.small=TRUE----
ep = calcExpectedPartitionsRef(query, "hg19")
plotExpectedPartitions(ep)

## ----multiple-corrected-partition-plot, fig.cap="Expected partition distribution plots for multiple query region sets.", fig.small=TRUE----
ep2 = calcExpectedPartitionsRef(queryList, "hg19")
plotExpectedPartitions(ep2)

## ----cumulative-partition-plot, fig.cap="Cumulative partition distribution plot. Cumulative distribution of query regions across genomic features.", fig.small=TRUE----
cp = calcCumulativePartitionsRef(query, "hg19")
plotCumulativePartitions(cp)

## ----multiple-cumulative-partition-plot, fig.cap="Cumulative partition distribution plots for multiple query region sets.", fig.small=TRUE----
cp2 = calcCumulativePartitionsRef(queryList, "hg19")
plotCumulativePartitions(cp2)

## ----specificity-plot-single, fig.height = 6, fig.cap="Specificity of chromatin accessibility across cell types."----
exampleCellMatrixFile = system.file("extdata", "example_cell_matrix.txt", package="GenomicDistributions")
cellMatrix = data.table::fread(exampleCellMatrixFile)
op = calcOpenSignal(query, cellMatrix)
plotOpenSignal(op)

## ----specificity-plot-multi, fig.height = 7, fig.cap="Specificity of chromatin accessibility across cell types."----
op2 = calcOpenSignal(queryList, cellMatrix)
plotOpenSignal(op2)

## ----neighbor-distance-plots, fig.cap="Distances between neighboring intervals of a regionset", fig.small=TRUE----

# Calculate the distances 
nd = calcNeighborDist(query)

# Plot the distribution
plotNeighborDist(nd)

## ----multiple-neighbor-distance-plots, fig.cap="Neighboring regions distance for multiple regionsets", fig.small=TRUE----

# Create a list of GRanges objects
s = system.file("extdata", "setB_100.bed.gz", package="GenomicDistributions")
setB = rtracklayer::import(s)
queryList2 = GRangesList(vistaEnhancers=query, setB=setB)

# Calculate the distances
dist2 = calcNeighborDist(queryList2)

# Plot the distribution
plotNeighborDist(dist2)

## ----gc-content-plots, fig.cap="GC content plot. Probability density function of GC percentage", eval=FALSE----
#  # Calculate the GC content
#  gc1 = calcGCContentRef(query, "hg19")
#  
#  # Plot the GC distribution
#  plotGCContent(gc1)

## ----gc-content-plots-multi, fig.cap="Multiple GC content plots.", eval=FALSE----
#  gc2 = calcGCContentRef(queryList, "hg19")
#  plotGCContent(gc2)

## ----width-distribution-single, fig.cap="Width distribution plot. Frequency of widths in this query"----
# Calculate the widths
qt1 = calcWidth(query)

# Plot the width distribution
plotQTHist(qt1)

## ----width-distribution-multi, fig.cap="Multiple width distribution plots.", fig.small=TRUE----
qt2 = calcWidth(queryList)
plotQTHist(qt2)

## ----width-distribution-colors, fig.cap="Width distribution plot with color options.", fig.small=TRUE----
plotQTHist(qt1, bins=7, quantile = .015, EndBarColor = 'black', MiddleBarColor='darkblue')