Skip to content

Commit

Permalink
Expand dba rules (#44)
Browse files Browse the repository at this point in the history
* feat: add cftoolkit v0.6.0 w rmarkdown

* fix: inherit from docker v0.5.0

* feat: move to static images for diff bind and uropa rules

* fix: rework diff bind and uropa scripts for modularity and new rules

* fix: fix hooks, split diff bind and uropa into multiple rule hierarchical rules

* chore: remove planned new docker version

* fix: full path for config

* fix: set cwd at rule instead of script

* fix: cd to output dir in diff bind rule

* fix: make sif cache aware of env

* fix: review changes

* fix: syntax error with rmarkdown::render

* fix: images not image

* fix: remove references to container sif files

* fix: pathing in rmarkdown and syntax in prep_diffbind
  • Loading branch information
rroutsong authored Aug 27, 2024
1 parent 27a0e29 commit 3325348
Show file tree
Hide file tree
Showing 8 changed files with 493 additions and 290 deletions.
107 changes: 72 additions & 35 deletions bin/DiffBind_v2_ChIPseq.Rmd
Original file line number Diff line number Diff line change
Expand Up @@ -16,34 +16,44 @@ params:
---

<style type="text/css">
body{
font-size: 12pt;
}
body {
font-size: 12pt;
}
</style>

```{r, include=FALSE, warning=FALSE, message=FALSE}
## grab args
dateandtime<-format(Sys.time(), "%a %b %d %Y - %X")
# inputs
dateandtime <- format(Sys.time(), "%a %b %d %Y - %X")
csvfile <- params$csvfile
outbase <- dirname(csvfile)
contrasts <- params$contrasts
peakcaller <- params$peakcaller
```
**Groups being compared:**
*`r contrasts`*
**Peak sources:**
*`r peakcaller`*
**Report generated:**
*`r dateandtime`*
# file output suffixes
cp_bed <- "_Diffbind_consensusPeaks.bed"
edger_txt <- "_Diffbind_EdgeR.txt"
deseq2_txt <- "_Diffbind_Deseq2.txt"
edger_bed <- "_Diffbind_EdgeR.bed"
deseq2_bed <- "_Diffbind_Deseq2.bed"
deseq2_bed_fullist <- "_Diffbind_Deseq2_fullList.txt"
edger_bed_fullist <- "_Diffbind_EdgeR_fullList.txt"
```{r setup, echo=FALSE, warning=FALSE,message=FALSE}
# knitr options
knitr::opts_chunk$set(echo = FALSE, include=TRUE, message=FALSE, warning=FALSE, error=FALSE)
# libraries
suppressMessages(library(DT))
suppressMessages(library(DiffBind))
suppressMessages(library(parallel))
```

**Groups being compared:**
*`r contrasts`*
**Peak sources:**
*`r peakcaller`*
**Report generated:**
*`r dateandtime`*

# Peak Data
Read in sample sheet information and peak information
```{r samples}
Expand Down Expand Up @@ -88,13 +98,13 @@ if ( grepl("narrow",samples$samples$Peaks[1]) ) {
print ("Narrow peak calling tool.")
print ("Differential peaks are 250bp upstream and downstream of the summits.")
} else if ( grepl("broad",samples$samples$Peaks[1]) ) {
summits <- FALSE
print ("Broad peak calling tool.")
print ("Differential peaks are consensus peaks.")
summits <- FALSE
print ("Broad peak calling tool.")
print ("Differential peaks are consensus peaks.")
} else {
summits <- FALSE
print ("Indeterminate peak calling tool.")
print ("Differential peaks are consensus peaks.")
summits <- FALSE
print ("Indeterminate peak calling tool.")
print ("Differential peaks are consensus peaks.")
}
```

Expand All @@ -105,9 +115,9 @@ if (summits == TRUE) {
DBdataCounts <- dba.count(samples)
}
print(DBdataCounts)
outfile2 <- paste0(contrasts, "-", peakcaller,"_Diffbind_consensusPeaks.bed")
consensus2 <- dba.peakset(DBdataCounts,bRetrieve=T)
consensus2$name <- paste0("Peak",1:length(consensus2))
outfile2 <- paste0(contrasts, "-", peakcaller, cp_bed)
consensus2 <- dba.peakset(DBdataCounts, bRetrieve=T)
consensus2$name <- paste0("Peak", 1:length(consensus2))
#rtracklayer::export(consensus2,outfile2)
```

Expand Down Expand Up @@ -210,11 +220,24 @@ try(dba.plotHeatmap(DBAnalysisEdgeR,contrast=1,method = DBA_EDGER,correlations=F
## Top 500 differentially bound peaks {.tabset .tabset-fade}
### DeSeq2 {-}
```{r Deseq2Report}
outfile <- paste0(contrasts, "-", peakcaller, "_Diffbind_Deseq2.txt")
outfile2 <- paste0(contrasts, "-", peakcaller, "_Diffbind_Deseq2.bed")
outfile <- paste0(contrasts, "-", peakcaller, deseq2_txt)
outfile2 <- paste0(contrasts, "-", peakcaller, deseq2_bed)
DBReportDeseq2$name <- paste0("Peak",1:length(DBReportDeseq2))
try(rtracklayer::export(DBReportDeseq2, outfile2),silent=TRUE)
write.table(DBReportDeseq2, outfile, quote=F, sep="\t", row.names=F)
tryDeseqExport <- function(DBReportDeseq2, outfile2) {
tryCatch(
{
rtracklayer::export(DBReportDeseq2, outfile2)
},
error = function(cond) {
print("ERROR: Failed to export DeSeq bed file `rtracklayer::export(DBReportDeseq2, outfile2)`, output blank file")
write.table(outfile2, file='empty', col.names=FALSE)
}
)
}
tryDeseqExport(DBReportDeseq2, file.path(outbase, outfile2))
write.table(DBReportDeseq2, file.path(outbase, outfile), quote=F, sep="\t", row.names=F)
D2i <- length(DBReportDeseq2)
if (D2i == 0) {
i=1
Expand All @@ -226,17 +249,31 @@ if (D2i == 0) {
try(DT::datatable(data.frame(DBReportDeseq2)[1:i,], rownames=F),silent=TRUE)
report2 <- dba.report(DBAnalysisDeseq2,method = DBA_DESEQ2,th=100,bNormalized=T,bFlip=FALSE,precision=0,bCalled=T)
outfile3 <- paste0(contrasts, "-", peakcaller, "_Diffbind_Deseq2_fullList.txt")
write.table(report2, outfile3, quote=F, sep="\t", row.names=F)
outfile3 <- paste0(contrasts, "-", peakcaller, deseq2_bed_fullist)
write.table(report2, file.path(outbase, outfile3), quote=F, sep="\t", row.names=F)
```

### EdgeR {-}
```{r EdgeRReport}
outfile <- paste0(contrasts, "-", peakcaller,"_Diffbind_EdgeR.txt")
outfile2 <- paste0(contrasts, "-", peakcaller,"_Diffbind_EdgeR.bed")
outfile <- paste0(contrasts, "-", peakcaller, edger_txt)
outfile2 <- paste0(contrasts, "-", peakcaller, edger_bed)
DBReportEdgeR$name <- paste0("Peak",1:length(DBReportEdgeR))
try(rtracklayer::export(DBReportEdgeR, outfile2),silent=TRUE)
write.table(DBReportEdgeR, outfile, quote=F, sep="\t", row.names=F)
tryEdgeRExport <- function(edger_report, fout) {
tryCatch(
{
rtracklayer::export(edger_report, fout)
},
error = function(cond) {
print("ERROR: Failed to export EdgeR bed file `rtracklayer::export(edger_report, fout))`, output blank file")
write.table(fout, file='empty', col.names=FALSE)
}
)
}
tryEdgeRExport(DBReportEdgeR, file.path(outbase, outfile2))
write.table(DBReportEdgeR, file.path(outbase, outfile), quote=F, sep="\t", row.names=F)
Ei <- length(DBReportEdgeR)
if (Ei == 0) {
i=1
Expand All @@ -248,8 +285,8 @@ if (Ei == 0) {
try(DT::datatable(data.frame(DBReportEdgeR)[1:i,], rownames=F),silent=TRUE)
report2 <- dba.report(DBAnalysisEdgeR,method = DBA_EDGER,th=100,bNormalized=T,bFlip=FALSE,precision=0,bCalled=T)
outfile3 <- paste0(contrasts, "-", peakcaller, "_Diffbind_EdgeR_fullList.txt")
write.table(report2, outfile3, quote=F, sep="\t", row.names=F)
outfile3 <- paste0(contrasts, "-", peakcaller, edger_bed_fullist)
write.table(report2, file.path(outbase, outfile3), quote=F, sep="\t", row.names=F)
```

## R tool version information
Expand Down
Loading

0 comments on commit 3325348

Please sign in to comment.