Introduction.R
In chromVAR: Chromatin Variation Across Regions

## ---- message = FALSE, warning = FALSE----------------------------------------
library(chromVAR)
library(motifmatchr)
library(Matrix)
library(SummarizedExperiment)
library(BiocParallel)
set.seed(2017)

## -----------------------------------------------------------------------------
register(MulticoreParam(8))

## -----------------------------------------------------------------------------
register(MulticoreParam(8, progressbar = TRUE))

## -----------------------------------------------------------------------------
register(SnowParam(workers = 1, type = "SOCK"))

## -----------------------------------------------------------------------------
register(SerialParam())

## -----------------------------------------------------------------------------
peakfile <- system.file("extdata/test_bed.txt", package = "chromVAR")
peaks <- getPeaks(peakfile, sort_peaks = TRUE)

## -----------------------------------------------------------------------------
bamfile <- system.file("extdata/test_RG.bam", package = "chromVAR")
fragment_counts <- getCounts(bamfile, 
                             peaks, 
                             paired =  TRUE, 
                             by_rg = TRUE, 
                             format = "bam", 
                             colData = DataFrame(celltype = "GM"))

## -----------------------------------------------------------------------------
data(example_counts, package = "chromVAR")
head(example_counts)

## ---- message = FALSE---------------------------------------------------------
library(BSgenome.Hsapiens.UCSC.hg19)
example_counts <- addGCBias(example_counts, 
                            genome = BSgenome.Hsapiens.UCSC.hg19)
head(rowData(example_counts))

## -----------------------------------------------------------------------------
#find indices of samples to keep
counts_filtered <- filterSamples(example_counts, min_depth = 1500, 
                                 min_in_peaks = 0.15, shiny = FALSE)

## ----filter_plot, fig.width = 6, fig.height = 6-------------------------------
#find indices of samples to keep
filtering_plot <- filterSamplesPlot(example_counts, min_depth = 1500, 
                                    min_in_peaks = 0.15, use_plotly = FALSE)
filtering_plot

## -----------------------------------------------------------------------------
ix <- filterSamples(example_counts, ix_return = TRUE, shiny = FALSE)

## -----------------------------------------------------------------------------
counts_filtered <- filterPeaks(counts_filtered, non_overlapping = TRUE)

## -----------------------------------------------------------------------------
motifs <- getJasparMotifs()

## -----------------------------------------------------------------------------
yeast_motifs <- getJasparMotifs(species = "Saccharomyces cerevisiae")

## -----------------------------------------------------------------------------
library(motifmatchr)
motif_ix <- matchMotifs(motifs, counts_filtered, 
                        genome = BSgenome.Hsapiens.UCSC.hg19)

## -----------------------------------------------------------------------------
kmer_ix <- matchKmers(6, counts_filtered, 
                      genome = BSgenome.Hsapiens.UCSC.hg19)

## -----------------------------------------------------------------------------
dev <- computeDeviations(object = counts_filtered, annotations = motif_ix)

## -----------------------------------------------------------------------------
bg <- getBackgroundPeaks(object = counts_filtered)

## -----------------------------------------------------------------------------
dev <- computeDeviations(object = counts_filtered, annotations = motif_ix,
                         background_peaks = bg)

## ----variability, fig.width = 6, fig.height = 6-------------------------------
variability <- computeVariability(dev)

plotVariability(variability, use_plotly = FALSE) 

## -----------------------------------------------------------------------------
tsne_results <- deviationsTsne(dev, threshold = 1.5, perplexity = 10)

## -----------------------------------------------------------------------------
tsne_plots <- plotDeviationsTsne(dev, tsne_results, 
                                 annotation_name = "TEAD3", 
                                   sample_column = "Cell_Type", 
                                 shiny = FALSE)
tsne_plots[[1]]
tsne_plots[[2]]

## -----------------------------------------------------------------------------
Sys.Date()

## -----------------------------------------------------------------------------
sessionInfo()