## ----style, echo=FALSE, results="asis", message=FALSE---------------------- BiocStyle::markdown() knitr::opts_chunk$set(tidy = FALSE, warning = FALSE, message = FALSE) ## ----echo=FALSE, results='hide', message=FALSE----------------------------- library(BaalChIP) ## ---- eval=TRUE------------------------------------------------------------ wd <- system.file("test",package="BaalChIP") #system directory samplesheet <- file.path(wd, "exampleChIP1.tsv") hets <- c("MCF7"="MCF7_hetSNP.txt", "GM12891"="GM12891_hetSNP.txt") res <- BaalChIP(samplesheet=samplesheet, hets=hets) ## ----quick1, eval=FALSE---------------------------------------------------- # res <- BaalChIP(samplesheet=samplesheet, hets=hets) # res <- BaalChIP.run(res, cores=4, verbose=TRUE) #cores for parallel computing ## ----quick2, eval=FALSE---------------------------------------------------- # # #create BaalChIP object # samplesheet <- file.path(wd, "exampleChIP1.tsv") # hets <- c("MCF7"="MCF7_hetSNP.txt", "GM12891"="GM12891_hetSNP.txt") # res <- BaalChIP(samplesheet=samplesheet, hets=hets) # # #Now, load some data # data(blacklist_hg19) # data(pickrell2011cov1_hg19) # data(UniqueMappability50bp_hg19) # # #run one at the time (instead of BaalChIP.run) # res <- alleleCounts(res, min_base_quality=10, min_mapq=15, verbose=FALSE) # res <- QCfilter(res, # RegionsToFilter=c("blacklist_hg19", "pickrell2011cov1_hg19"), # RegionsToKeep="UniqueMappability50bp_hg19", # verbose=FALSE) # res <- mergePerGroup(res) # res <- filter1allele(res) # res <- getASB(res, # Iter=5000, # conf_level=0.95, # cores=4, RMcorrection = TRUE, # RAFcorrection=TRUE) ## -------------------------------------------------------------------------- gDNA <- list("TaT1"= c(file.path(wd, "bamFiles/TaT1_1_gDNA.test.bam"), file.path(wd, "bamFiles/TaT1_2_gDNA.test.bam")), "AMOVC"= c(file.path(wd, "bamFiles/AMOVC_1_gDNA.test.bam"), file.path(wd, "bamFiles/AMOVC_2_gDNA.test.bam"))) ## ----quick3, eval=FALSE---------------------------------------------------- # wd <- system.file("test",package="BaalChIP") #system directory # # samplesheet <- file.path(wd, "exampleChIP2.tsv") # hets <- c("TaT1"="TaT1_hetSNP.txt", "AMOVC"="AMOVC_hetSNP.txt") # res2 <- BaalChIP(samplesheet=samplesheet, hets=hets, CorrectWithgDNA=gDNA) # res2 <- BaalChIP.run(res2, cores=4, verbose=TRUE) #cores for parallel computing ## -------------------------------------------------------------------------- samplesheet <- read.delim(file.path(wd,"exampleChIP1.tsv")) samplesheet ## -------------------------------------------------------------------------- samplesheet <- read.delim(file.path(wd,"exampleChIP2.tsv")) samplesheet ## -------------------------------------------------------------------------- head(read.delim(file.path(system.file("test",package="BaalChIP"),"MCF7_hetSNP.txt"))) ## ---- eval=TRUE------------------------------------------------------------ samplesheet <- file.path(wd,"exampleChIP1.tsv") hets <- c("MCF7"="MCF7_hetSNP.txt", "GM12891"="GM12891_hetSNP.txt") res <- BaalChIP(samplesheet=samplesheet, hets=hets) res ## -------------------------------------------------------------------------- BaalChIP.get(res, what="samples") ## ---- eval=TRUE------------------------------------------------------------ #run alleleCounts res <- alleleCounts(res, min_base_quality=10, min_mapq=15, verbose=FALSE) ## -------------------------------------------------------------------------- data(blacklist_hg19) data(pickrell2011cov1_hg19) data(UniqueMappability50bp_hg19) ## ---- eval=TRUE------------------------------------------------------------ #run QC filter res <- QCfilter(res, RegionsToFilter=list("blacklist"=blacklist_hg19, "highcoverage"=pickrell2011cov1_hg19), RegionsToKeep=list("UniqueMappability"=UniqueMappability50bp_hg19), verbose=FALSE) res <- mergePerGroup(res) res <- filter1allele(res) ## ---- eval=FALSE, message=FALSE, error=FALSE, warning=FALSE---------------- # # res <- filterIntbias(res, # simul_output="directory_name", # tmpfile_prefix="prefix_name", # simulation_script = "local", # alignmentSimulArgs=c("picard-tools-1.119", # "bowtie-1.1.1", # "genomes_test/male.hg19", # "genomes_test/maleByChrom") # verbose=FALSE) ## -------------------------------------------------------------------------- preComputed_output <- system.file("test/simuloutput",package="BaalChIP") list.files(preComputed_output) ## ---- eval=TRUE, message=FALSE--------------------------------------------- res <- filterIntbias(res, simul_output=preComputed_output, tmpfile_prefix="c67c6ec6c433", skipScriptRun=TRUE, verbose=FALSE) ## -------------------------------------------------------------------------- res <- mergePerGroup(res) ## -------------------------------------------------------------------------- res <- filter1allele(res) ## ---- eval=FALSE, message=FALSE-------------------------------------------- # res <- getASB(res, Iter=5000, conf_level=0.95, cores = 4, # RMcorrection = TRUE, # RAFcorrection=TRUE) ## ---- eval=TRUE, echo=FALSE, include=FALSE--------------------------------- data(baalObject) res <- BaalObject ## -------------------------------------------------------------------------- res ## -------------------------------------------------------------------------- result <- BaalChIP.report(res) head(result[["MCF7"]]) ## -------------------------------------------------------------------------- data(ENCODEexample) ENCODEexample ## -------------------------------------------------------------------------- a <- summaryQC(ENCODEexample) ## -------------------------------------------------------------------------- summaryQC(ENCODEexample)[["filtering_stats"]] ## -------------------------------------------------------------------------- summaryQC(ENCODEexample)[["average_stats"]] ## ----plotENCODE1----------------------------------------------------------- plotQC(ENCODEexample, what="barplot_per_group") ## ----plotENCODE2----------------------------------------------------------- plotQC(ENCODEexample, what="boxplot_per_filter") ## ----plotENCODE3----------------------------------------------------------- plotQC(ENCODEexample, what="overall_pie") ## ----plotSimul------------------------------------------------------------- plotSimul(ENCODEexample) ## -------------------------------------------------------------------------- #ENCODE example a <- BaalChIP.get(ENCODEexample, "assayedVar")[["MCF7"]] a[1:5,1:5] #FAIRE exmaple a <- BaalChIP.get(FAIREexample, "assayedVar")[["MDA134"]] a[1:5,] ## -------------------------------------------------------------------------- summaryASB(ENCODEexample) ## -------------------------------------------------------------------------- summaryASB(FAIREexample) ## ----adjENCODE, fig.width=12----------------------------------------------- adjustmentBaalPlot(ENCODEexample) ## ----adjFAIRE1, fig.width=7, fig.height=2.5-------------------------------- adjustmentBaalPlot(FAIREexample) ## ----adjFAIRE2, fig.width=7, fig.height=2.5-------------------------------- adjustmentBaalPlot(FAIREexample, col=c("cyan4","chocolate3")) ## -------------------------------------------------------------------------- biastable <- BaalChIP.get(ENCODEexample, "biasTable") head(biastable[["K562"]]) ## -------------------------------------------------------------------------- biastable <- BaalChIP.get(FAIREexample, "biasTable") head(biastable[["T47D"]]) ## -------------------------------------------------------------------------- result <- BaalChIP.report(FAIREexample)[["T47D"]] #show ASB SNPs result[result$isASB==TRUE,] ## ----echo=FALSE------------------------------------------------------------ sessionInfo()