## ----include = FALSE----------------------------------------------------------
knitr::opts_chunk$set(
    collapse = TRUE,
    comment = "#>"
)

## -----------------------------------------------------------------------------
library(misha)
gdb.init_examples()

## -----------------------------------------------------------------------------
gtrack.ls() # list tracks in the examples DB
gtrack.info("dense_track") # inspect type/metadata
gtrack.info("sparse_track")

## -----------------------------------------------------------------------------
regions <- gintervals(1, c(0, 250000), c(100000, 260000))
regions

## -----------------------------------------------------------------------------
out <- gextract("dense_track", regions, iterator = 100)
head(out)
log_out <- gextract("log(dense_track + 1)", regions, iterator = 100)
head(log_out)

## -----------------------------------------------------------------------------
gintervals.save("my_intervals_set", regions) # name first, then the intervals
out2 <- gextract("dense_track", gintervals.all(), iterator = "my_intervals_set")
out2

## -----------------------------------------------------------------------------
gvtrack.create("chip.sum", "dense_track", "sum")
out <- gextract("chip.sum", regions, iterator = 200)
head(out)

## -----------------------------------------------------------------------------
gvtrack.create("chip.shifted", "dense_track", "sum")
gvtrack.iterator("chip.shifted", sshift = -100, eshift = 100)
out <- gextract("chip.shifted", regions, iterator = 200)
head(out)

## -----------------------------------------------------------------------------
library(misha)
gdb.init_examples()

# 1) pick scope
regions <- gintervals(1, 0, 50000)

# 2) inspect available tracks
gtrack.ls()

# 3) extract signal with a chosen iterator
chip <- gextract("dense_track", regions, iterator = 100)
head(chip)

# 4) screen high-signal bins (as a simple peak-like filter).
#    Pick the threshold from the data: dense_track only reaches 0.34 in this
#    scope, so a threshold above that would silently screen nothing.
hi_chip <- gscreen("dense_track > 0.2", regions, iterator = 100)
head(hi_chip)

# 5) summarize distribution/coverage
gsummary("dense_track", regions, iterator = 100)

## -----------------------------------------------------------------------------
regions <- gintervals(1, c(1000, 2000), c(1020, 2020))
seqs <- gseq.extract(regions)
seqs

pssm <- matrix(c(
    0.80, 0.05, 0.10, 0.05,
    0.10, 0.10, 0.70, 0.10,
    0.05, 0.80, 0.05, 0.10,
    0.10, 0.10, 0.10, 0.70
), ncol = 4, byrow = TRUE)
colnames(pssm) <- c("A", "C", "G", "T")

scores <- gseq.pwm(seqs, pssm, mode = "lse")
scores

