## ----, echo = FALSE, message = FALSE-------------------------------------
library(markdown)
options(markdown.HTML.options = c(options('markdown.HTML.options')[[1]], "toc"))

library(knitr)
knitr::opts_chunk$set(
    error = FALSE,
    tidy  = FALSE,
    message = FALSE,
    fig.align = "center",
    fig.width = 6,
    fig.height = 6)
options(markdown.HTML.stylesheet = "custom.css")

options(width = 100)

## ------------------------------------------------------------------------
suppressPackageStartupMessages(library(HilbertCurve))
library(circlize)
set.seed(12345)

## ----, fig.width = 12, fig.height = 3, echo = FALSE----------------------
grid.newpage()
pushViewport(viewport(layout = grid.layout(nr = 1, nc = 4)))
pushViewport(viewport(layout.pos.row = 1, layout.pos.col = 1))
HilbertCurve(1, 100, level = 2, reference = TRUE, newpage = FALSE, title = "level = 2")
upViewport()
pushViewport(viewport(layout.pos.row = 1, layout.pos.col = 2))
HilbertCurve(1, 100, level = 3, reference = TRUE, newpage = FALSE, title = "level = 3")
upViewport()
pushViewport(viewport(layout.pos.row = 1, layout.pos.col = 3))
HilbertCurve(1, 100, level = 4, reference = TRUE, newpage = FALSE, title = "level = 4")
upViewport()
pushViewport(viewport(layout.pos.row = 1, layout.pos.col = 4))
HilbertCurve(1, 100, level = 5, reference = TRUE, arrow = FALSE, newpage = FALSE, title = "level = 5")
upViewport()
upViewport()

## ----, fig.width = 7-----------------------------------------------------
library(HilbertVis)
pos = HilbertVis::hilbertCurve(5)
mat = as.matrix(dist(pos))
library(ComplexHeatmap)

ht = Heatmap(mat, name = "dist", cluster_rows = FALSE, cluster_columns = FALSE, 
    show_row_names = FALSE, show_column_names = FALSE, 
    heatmap_legend_param = list(title = "euclidean_dist"))
draw(ht, padding = unit(c(5, 5, 5, 2), "mm"))
decorate_heatmap_body("dist", {
    grid.segments(c(0.25, 0.5, 0.75, 0, 0, 0), c(0, 0, 0, 0.25, 0.5, 0.75), 
	      c(0.25, 0.5, 0.75, 1, 1, 1), c(1, 1, 1, 0.25, 0.5, 0.75), gp = gpar(lty = 2))
    grid.text(rev(c(256, 512, 768, 1024)), 0, c(0, 256, 512, 768)/1024, just = "bottom", 
        rot = 90, gp = gpar(fontsize = 10))
    grid.text(c(1, 256, 512, 768, 1024), c(1, 256, 512, 768, 1024)/1024, 1, just = "bottom",
        gp = gpar(fontsize = 10))
})

## ----, eval = FALSE------------------------------------------------------
## hc = HilbertCurve(...)
## hc_points(hc, ...)
## hc_rect(hc, ...)
## hc_segments(hc, ...)
## hc_text(hc, ...)

## ----, eval = FALSE------------------------------------------------------
## hc = HilbertCurve(1, 100, level = 4, reference = TRUE)

## ------------------------------------------------------------------------
x = sort(sample(100, 20))
s = x[1:10*2 - 1]
e = x[1:10*2]
ir = IRanges(s, e)
ir

## ------------------------------------------------------------------------
hc = HilbertCurve(1, 100, level = 4, reference = TRUE)
hc_points(hc, ir)

## ------------------------------------------------------------------------
hc = HilbertCurve(1, 100, level = 4, reference = TRUE)
hc_points(hc, ir, np = 3, gp = gpar(fill = rand_color(length(ir))))

## ------------------------------------------------------------------------
hc = HilbertCurve(1, 100, level = 4, reference = TRUE)
hc_points(hc, ir, np = NULL, size = unit(runif(length(ir)), "cm"), pch = 16)

## ------------------------------------------------------------------------
hc = HilbertCurve(1, 100, level = 4, reference = TRUE)
hc_segments(hc, ir)

## ------------------------------------------------------------------------
hc = HilbertCurve(1, 100, level = 4, reference = TRUE)
hc_rect(hc, ir)

## ------------------------------------------------------------------------
hc = HilbertCurve(1, 100, level = 4, reference = TRUE)
labels = sample(letters, length(ir), replace = TRUE)
hc_text(hc, ir, labels = labels)

## ------------------------------------------------------------------------
hc = HilbertCurve(1, 100, level = 4)
hc_segments(hc, IRanges(1, 100))   # background line
hc_rect(hc, ir)
hc_points(hc, ir, np = 3)
hc_text(hc, ir, labels = labels, gp = gpar(fontsize = 16, col = "blue"))

## ------------------------------------------------------------------------
hc = HilbertCurve(0.1, 0.8, level = 4, reference = TRUE)
hc_points(hc, x1 = c(0.15, 0.55), x2 = c(0.25, 0.65))

## ------------------------------------------------------------------------
value = runif(length(ir))
col_fun = colorRamp2(c(0, 1), c("white", "red"))
cm = ColorMapping(col_fun = col_fun)
legend = color_mapping_legend(cm, plot = FALSE, title = "value")
hc = HilbertCurve(1, 100, reference = TRUE, title = "points", legend = legend)
hc_points(hc, ir, np = 3, gp = gpar(fill = col_fun(value)))

## ------------------------------------------------------------------------
hc = HilbertCurve(1, 100, level = 9, mode = "pixel")
hc_layer(hc, ir)

## ------------------------------------------------------------------------
hc = HilbertCurve(1, 100, level = 9, mode = "pixel", title = "pixel mode", legend = legend)
hc_layer(hc, ir, col = col_fun(value), grid_line = 2)
hc_layer(hc, x1 = 75, x2 = 85, col = "#0000FF10")

## ----, eval = FALSE------------------------------------------------------
## hc_png(hc, file = "test.png")

## ------------------------------------------------------------------------
col = rainbow(100)
hc = HilbertCurve(1, 100, level = 5)
ir = IRanges(start = 1:99, end = 2:100)
hc_points(hc, ir, np = 4, gp = gpar(col = col, fill = col))

## ------------------------------------------------------------------------
chr1_len = 249250621
load(paste0(system.file("extdata", package = "HilbertCurve"), "/refseq_chr1.RData"))
hc = HilbertCurve(1, chr1_len, level = 5)
hc_segments(hc, IRanges(start = 1, end = chr1_len), gp = gpar(col = "grey"))
hc_segments(hc, ranges(g), gp = gpar(lwd = unit(6, "mm"), col = rand_color(length(g))))

## ------------------------------------------------------------------------
load(paste0(system.file("extdata", package = "HilbertCurve"), "/mouse_net.RData"))
hc = HilbertCurve(1, chr1_len, level = 6)
ir1 = ranges(mouse)
ir2 = setdiff(IRanges(1, chr1_len), ir1)
ir = c(ir1, ir2)
col = c(rep("red", length(ir1)), rep("yellow", length(ir2)))
hc_points(hc, ir, np = 3, gp = gpar(col = col, fill = col))

load(paste0(system.file("extdata", package = "HilbertCurve"), "/zebrafish_net.RData"))
hc = HilbertCurve(1, chr1_len, level = 6)
ir1 = ranges(zebrafish)
ir2 = setdiff(IRanges(1, chr1_len), ir1)
ir = c(ir1, ir2)
col = c(rep("red", length(ir1)), rep("yellow", length(ir2)))
hc_points(hc, ir, np = 3, gp = gpar(col = col, fill = col))

## ----, eval = FALSE------------------------------------------------------
## df = read.table(pipe("awk '$1==\"chr1\"' /path-to/hg19_gc_percent_window1000.bed"))
## col_fun = colorRamp2(c(0, 500, 1000), c("green", "#FFFFCC", "red"))
## png("gc_percent_chr1_points.png", width = 500, height = 500)
## hc = HilbertCurve(1, chr1_len, level = 6)
## hc_points(hc, IRanges(df[[2]], df[[3]]), np = 3, gp = gpar(fill = col_fun(df[[5]]), col = col_fun(df[[5]])))
## hc_rect(hc, reduce(ranges(g)), gp = gpar(fill = "#00000020", col = "#00000020"))
## dev.off()

## ----, eval = FALSE------------------------------------------------------
## hc = HilbertCurve(1, chr1_len, level = 9, mode = "pixel")
## hc_layer(hc, IRanges(df[[2]], df[[3]]), col = col_fun(df[[5]]))
## hc_layer(hc, reduce(ranges(g)), col = "#00000020")
## hc_png(hc, file = "gc_percent_chr1.png")

## ----, eval = FALSE------------------------------------------------------
## library(GetoptLong)
## 
## for(mark in c("H3K27ac", "H3K36me3", "H3K4me3", "H3K9me3")) {
## 	df = read.table(pipe(qq("awk '$5>0 && $1==\"chr1\"' /path-to/UCSD.Lung.@{mark}.STL002.bed")), sep = "\t")
## 	col_fun = colorRamp2(c(0, quantile(df[[5]], 0.99)), c("white", "red"))
## 	hc = HilbertCurve(1, chr1_len, level = 9, mode = "pixel")
## 	hc_layer(hc, IRanges(df[[2]], df[[3]]), col = col_fun(df[[5]]))
##   hc_layer(hc, reduce(ranges(g)), col = "#00000010")
## 	hc_png(hc, file = qq("@{mark}_chr1.png"))
## }

## ----, eval = FALSE------------------------------------------------------
## df = read.table(pipe(qq("awk '$5>0 && $1==\"chr1\"' /path-to/UCSD.Lung.mRNA-Seq.STL002.bed")), sep = "\t")
## df[[5]] = log(df[[5]] + 1)
## col_fun = colorRamp2(c(0, quantile(df[[5]], 0.99)), c("white", "red"))
## hc = HilbertCurve(1, chr1_len, level = 9, mode = "pixel")
## hc_layer(hc, IRanges(df[[2]], df[[3]]), col = col_fun(df[[5]]))
## hc_layer(hc, reduce(ranges(g)), col = "#00000010")
## hc_png(hc, file = qq("RNASeq_coverage_chr1.png"))

## ----, eval = FALSE------------------------------------------------------
## df = read.table(pipe("awk '$1==\"chr1\"' /path-to/UCSD.Lung.Bisulfite-Seq.STL002.bed"), sep = "\t")
## col_fun = colorRamp2(c(0, 0.5, 1), c("blue", "white", "red"))
## hc = HilbertCurve(1, chr1_len, level = 9, mode = "pixel")
## hc_layer(hc, IRanges(df[[2]], df[[3]]), col = col_fun(df[[5]]), mean_mode = "absolute")
## hc_png(hc, file = "methylation_chr1.png")

## ----, eval = FALSE------------------------------------------------------
## chr.len = read.chromInfo()$chr.len
## sum_chr = 0
## s = numeric(length(ratio))
## e = numeric(length(ratio))
## for(chr in paste0("chr", c(1:22, "X", "Y"))) {
##     l = as.character(seqnames(ratio)) == chr
##     s[l] = start(ratio[l]) + sum_chr
##     e[l] = end(ratio[l]) + sum_chr
##     sum_chr = sum_chr + chr.len[chr]
## }
## 
## ir = IRanges(start = round(s/1000),
##              end = round(e/1000))
## hc = HilbertCurve(0, round(max(chr_cumsum)/1000), level = 9, mode = "pixel")
## hc_layer(hc, ir, col = ..., grid_line = 3)

## ----, eval = FALSE------------------------------------------------------
## chr_cumsum = cumsum(chr.len[paste0("chr", c(1:22, "X", "Y"))])
## chr_ir = IRanges(start = round(c(1, chr_cumsum[-length(chr_cumsum)]+1)/1000),
##                  end = round(chr_cumsum/1000),
##                  names = names(chr_cumsum))
## hc = HilbertCurve(0, round(max(chr_cumsum)/1000), level = 8)
## hc_rect(hc, chr_ir, gp = gpar(fill = rand_color(length(chr_ir), transparency = 0.5), col = NA))
## hc_text(hc, chr_ir, names(chr_ir), gp = gpar(fontsize = 16))

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