## ----mock-data----------------------------------------------------------------
suppressPackageStartupMessages(library(GenomicRanges))
library(methylTFRAnnotationHg38)

mock_dir <- file.path(tempdir(), "methylTFRAnnotationHg38_vignette")
dir.create(mock_dir, showWarnings = FALSE)

# Binding sites: a GRangesList with one GRanges per motif. Each range is
# a motif match extended by 200 bp on either side.
mock_tfbs <- GRangesList(
    MA0001.1 = GRanges(
        c("chr1:10001-10411", "chr1:52001-52411", "chr2:3001-3411"),
        strand = c("+", "-", "+")
    ),
    MA0002.1 = GRanges(c("chr1:20001-20415", "chr3:8001-8415"),
        strand = c("-", "+")
    )
)
saveRDS(mock_tfbs, file.path(mock_dir, "altius_tf_bindsites.rds"))

# GC frequency tables: one 5-row matrix per motif. Rows are GC
# quintiles, columns are positions of a 30 nt window along the binding
# site, and every column sums to 1. We use 6 positions to keep it short.
set.seed(1)
mock_freq <- function(n_pos) {
    m <- matrix(runif(5 * n_pos), nrow = 5)
    sweep(m, 2, colSums(m), "/")
}
mock_gcfreq <- list(MA0001.1 = mock_freq(6), MA0002.1 = mock_freq(6))
saveRDS(mock_gcfreq, file.path(mock_dir, "altius_motif_gcfreq.rds"))

# Genome-wide GC: non-overlapping 30 nt windows with their GC fraction
# and GC quintile (1 = lowest, 5 = highest).
mock_gc <- GRanges("chr1", IRanges(start = seq(1, by = 30, length.out = 10),
    width = 30))
mock_gc$GC_bias <- round(runif(10, 0.2, 0.7), 2)
breaks <- quantile(mock_gc$GC_bias, probs = seq(0, 1, 0.2))
mock_gc$GC_bin <- findInterval(mock_gc$GC_bias, breaks,
    rightmost.closed = TRUE)
metadata(mock_gc) <- list(gc_breaks = breaks, bin_scope = "genome",
    step = 30L)
saveRDS(mock_gc, file.path(mock_dir, "genomewide_GC_hg38.rds"))

old_opt <- options(methylTFRAnnotationHg38.datadir = mock_dir)
list.files(mock_dir)

## ----motif-sets---------------------------------------------------------------
availableMotifSets()

## ----distal-------------------------------------------------------------------
file.copy(file.path(mock_dir, "altius_motif_gcfreq.rds"),
    file.path(mock_dir, "jaspar2020_distal_motif_gcfreq.rds"))
file.copy(file.path(mock_dir, "altius_tf_bindsites.rds"),
    file.path(mock_dir, "jaspar2020_tf_bindsites.rds"))
distal_freq <- getGCfreq("jaspar2020_distal")
distal_sites <- getTFbindsites("jaspar2020_distal") # jaspar2020 sites
names(distal_sites)

## ----tfbs---------------------------------------------------------------------
tf_bindsites <- getTFbindsites("altius")
tf_bindsites

## ----tfbs-summary-------------------------------------------------------------
lengths(tf_bindsites)
vapply(width(tf_bindsites), unique, integer(1))

## ----gcfreq-------------------------------------------------------------------
gcfreqs <- getGCfreq("altius")
names(gcfreqs)
dim(gcfreqs$MA0001.1)
round(gcfreqs$MA0001.1, 2)
colSums(gcfreqs$MA0001.1)

## ----genome-gc----------------------------------------------------------------
gc_dist <- getGenomeGC()
gc_dist
table(gc_dist$GC_bin)

## ----genome-gc-breaks---------------------------------------------------------
metadata(gc_dist)$gc_breaks

## ----cleanup------------------------------------------------------------------
options(old_opt)
unlink(mock_dir, recursive = TRUE)
getOption("methylTFRAnnotationHg38.datadir")

## ----session-info-------------------------------------------------------------
sessionInfo()