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

## -----------------------------------------------------------------------------
library(rigvf)
gene_variants(gene_name = "GCK")

## -----------------------------------------------------------------------------
gene_variants(gene_name = "GCK", page=1L)
gene_variants(gene_name = "GCK", limit=50L)

## -----------------------------------------------------------------------------
gene_variants(gene_name = "GCK", log10pvalue="gt:5.0")
gene_variants(gene_name = "GCK", effect_size="gt:0.5")

## -----------------------------------------------------------------------------
res <- gene_elements(gene_id = "ENSG00000187961", verbose = TRUE)
res

res |>
    dplyr::select(`regions`) |>
    tidyr::unnest_wider(`regions`)

## -----------------------------------------------------------------------------
db_gene_variants("ENSG00000106633", threshold = 0.85)

## -----------------------------------------------------------------------------
aql <- system.file(package = "rigvf", "aql", "gene_variants.aql")
readLines(aql) |> noquote()

## ----message=FALSE------------------------------------------------------------
library(plyranges)
rng <- data.frame(seqnames="chr1", start=10e6+1, end=10.1e6) |>
  as_granges()

e <- elements(rng, limit=200L) |>
  plyranges::filter(source != "ENCODE_EpiRaction")
e

## -----------------------------------------------------------------------------
tiles <- tile_ranges(rng, width=10e3) %>%
  plyranges::select(-partition) %>%
  plyranges::mutate(id = letters[seq_along(.)])

# count overlaps of central basepair of elements in tiles
e |>
  plyranges::anchor_center() |>
  plyranges::mutate(width = 1) |>
  plyranges::join_overlap_left(tiles) |>
  tibble::as_tibble() |>
  dplyr::count(id)

## -----------------------------------------------------------------------------
# up to 200 variants:
v <- gene_variants(gene_name = "GCK", limit=200L, verbose=TRUE) |>
  dplyr::select(-c(gene, chr, source, source_url)) |>
  tidyr::unnest_wider(`sequence variant`) |>
  dplyr::rename(seqnames = chr) |>
  dplyr::mutate(start = pos + 1, end = pos + 1) |>
  as_granges()

## ----message=FALSE------------------------------------------------------------
library(plotgardener)
par <- pgParams(
  chrom = "chr7",
  chromstart = 44.1e6,
  chromend = 44.25e6,
  assembly = "hg38",
  just = c("left", "bottom")
)

## -----------------------------------------------------------------------------
v_for_plot <- v |>
  plyranges::select(snp = rsid, p = log10pvalue, effect_size)

## -----------------------------------------------------------------------------
library(TxDb.Hsapiens.UCSC.hg38.knownGene)
library(org.Hs.eg.db)

## ----plotgardener, fig.alt="Manhattan of IGVF annotated variants", fig.dim = c(5, 4)----
pageCreate(width = 5, height = 4, showGuides = FALSE)
plotGenes(
  params = par, x = 0.5, y = 3.5, width = 4, height = 1
)
plotGenomeLabel(
  params = par,
  x = 0.5, y = 3.5, length = 4,
  just = c("left", "top")
)
mplot <- plotManhattan(
  params = par, x = 0.5, y = 2.5, width = 4, height = 2,
  v_for_plot, trans = "",
  sigVal = -log10(5e-8), sigLine = TRUE, col = "grey", lty = 2
)
annoYaxis(
    plot = mplot, at=0:4 * 4, axisLine = TRUE, fontsize = 8
)
annoXaxis(
    plot = mplot, axisLine = TRUE, label = FALSE
)
plotText(
    params = par,
    label = "-log10(p-value)", x = 0.2, y = 2, rot = 90,
    fontsize = 8, fontface = "bold",
    default.units = "inches"
)

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