## ----set-options, echo=FALSE, cache=FALSE-------------------------------------
options(stringsAsFactors = FALSE, warning = FALSE, message = FALSE)

## ----eval = FALSE-------------------------------------------------------------
#  BiocManager::install("TADCompare")

## -----------------------------------------------------------------------------
library(dplyr)
library(SpectralTAD)
library(TADCompare)

## ----echo = FALSE-------------------------------------------------------------
data("rao_chr22_prim")
row.names(rao_chr22_prim) <- colnames(rao_chr22_prim) <- format(as.numeric(row.names(rao_chr22_prim)), scientific = FALSE)
coords <- 200:275
pheatmap::pheatmap(log10(rao_chr22_prim[coords, coords]), cluster_rows = FALSE, cluster_cols = FALSE)

## ----echo = FALSE-------------------------------------------------------------
coords  <- 50:53
sub_mat <- data.frame(chr = "chr22", start = as.numeric(colnames(rao_chr22_prim[coords, coords])), end   = as.numeric(colnames(rao_chr22_prim[coords, coords])) + 50000, rao_chr22_prim[coords, coords]) 
row.names(sub_mat) = NULL
sub_mat

## ----echo = FALSE-------------------------------------------------------------
data("rao_chr22_prim") 
head(HiCcompare::full2sparse(rao_chr22_prim))

## ----eval=FALSE---------------------------------------------------------------
#  #Read in data
#  primary = read.table('primary.chr22.50kb.txt', header = FALSE)
#  replicate = read.table('replicate.chr22.50kb.txt', header = FALSE)
#  #Run TADCompare
#  tad_diff=TADCompare(primary, replicate, resolution=50000)

## ----eval = FALSE-------------------------------------------------------------
#  # Read in data
#  cool_mat1 <- read.table("Zuin.HEK293.50kb.Control.txt")
#  cool_mat2 <- read.table("Zuin.HEK293.50kb.Depleted.txt")
#  
#  # Convert to sparse 3-column matrix using cooler2sparse from HiCcompare
#  sparse_mat1 <- HiCcompare::cooler2sparse(cool_mat1)
#  sparse_mat2 <- HiCcompare::cooler2sparse(cool_mat2)
#  
#  # Run TADCompare
#  diff_tads = lapply(names(sparse_mat1), function(x) {
#    TADCompare(sparse_mat1[[x]], sparse_mat2[[x]], resolution = 50000)
#  })

## ----eval = FALSE-------------------------------------------------------------
#  # Read in both files
#  mat1 <- read.table("sample1_100000.matrix")
#  bed1 <- read.table("sample1_100000_abs.bed")
#  
#  # Matrix 2
#  
#  mat2 <- read.table("sample2_100000.matrix")
#  bed2 <- read.table("sample2_100000_abs.bed")
#  
#  # Convert to modified bed format
#  sparse_mats1 <- HiCcompare::hicpro2bedpe(mat1,bed1)
#  sparse_mats2 <- HiCcompare::hicpro2bedpe(mat2,bed2)
#  
#  # Remove empty matrices if necessary
#  # sparse_mats$cis = sparse_mats$cis[sapply(sparse_mats, nrow) != 0]
#  
#  
#  # Go through all pairwise chromosomes and run TADCompare
#  sparse_tads = lapply(1:length(sparse_mats1$cis), function(z) {
#    x <- sparse_mats1$cis[[z]]
#    y <- sparse_mats2$cis[[z]]
#  
#    #Pull out chromosome
#    chr <- x[, 1][1]
#    #Subset to make three column matrix
#    x <- x[, c(2, 5, 7)]
#    y <- y[, c(2, 5, 7)]
#    #Run SpectralTAD
#    comp <- TADCompare(x, y, resolution = 100000)
#    return(list(comp, chr))
#  })
#  
#  # Pull out differential TAD results
#  diff_res <- lapply(sparse_tads, function(x) x$comp)
#  # Pull out chromosomes
#  chr      <- lapply(sparse_tads, function(x) x$chr)
#  # Name list by corresponding chr
#  names(diff_res) <- chr

## ----message=FALSE------------------------------------------------------------
library(microbenchmark)
# Reading in the second matrix
data("rao_chr22_rep")
# Converting to sparse
prim_sparse <- HiCcompare::full2sparse(rao_chr22_prim)
rep_sparse  <- HiCcompare::full2sparse(rao_chr22_rep)
# Converting to nxn+3
# Primary
prim_n_n_3 <- data.frame(chr = "chr22",
                         start = as.numeric(colnames(rao_chr22_prim)),
                         end = as.numeric(colnames(rao_chr22_prim))+50000, 
                         rao_chr22_prim)

# Replicate
rep_n_n_3 <- data.frame(chr = "chr22", 
                        start = as.numeric(colnames(rao_chr22_rep)),
                        end = as.numeric(colnames(rao_chr22_rep))+50000,
                        rao_chr22_rep)
# Defining each function
# Sparse
sparse <- TADCompare(cont_mat1 = prim_sparse, cont_mat2 = rep_sparse, resolution = 50000)
# NxN
n_by_n <- TADCompare(cont_mat1 = prim_sparse, cont_mat2 = rep_sparse, resolution = 50000)
# Nx(N+3)
n_by_n_3 <- TADCompare(cont_mat1 = prim_n_n_3, cont_mat2 = rep_n_n_3, resolution = 50000)

# Benchmarking different parameters
bench <- microbenchmark(
# Sparse
sparse <- TADCompare(cont_mat1 = prim_sparse, cont_mat2 = rep_sparse, resolution = 50000),
# NxN
n_by_n <- TADCompare(cont_mat1 = rao_chr22_prim, cont_mat2 = rao_chr22_rep, resolution = 50000),
# Nx(N+3)
n_by_n_3 <- TADCompare(cont_mat1 = prim_n_n_3, cont_mat2 = rep_n_n_3, resolution = 50000), times = 5, unit = "s"
) 

summary_bench <- summary(bench) %>% dplyr::select(mean, median)
rownames(summary_bench) <- c("sparse", "n_by_n", "n_by_n_3")
summary_bench

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

