Commit 15b7a26f authored by Jonathan Gesell's avatar Jonathan Gesell

Create branch and initial attempts to add UMI-Tools as an option for deduplication.

parents 59b330c3 712d035d
Pipeline #2847 failed with stage
in 0 seconds
## RNASeq Analysis Worklow
```
module load nextflow/0.24.1-SNAPSHOT
module load nextflow
nextflow run workflow/main.nf
```
This diff is collapsed.
......@@ -3,6 +3,7 @@ update.packages()
install.packages("sqldf",dep=TRUE)
install.packages("gmp",dep=TRUE)
install.packages("shinythemes",dep=TRUE)
install.packages(c('gplots','lattice','latticeExtra','vegan','labdsv','cluster','ggplot2'))
install.packages("Vennerable", repos="http://R-Forge.R-project.org",type='source')
source("http://bioconductor.org/biocLite.R")
......
library(shiny)
library(Vennerable)
library(qusage)
library(DT)
library(ggplot2)
library(ballgown)
library(sqldf)
library(reshape2)
library("gplots")
shinyServer(function(input, output, session) {
data.dir <- Sys.getenv('outputDir')
rda.dir <- data.dir
symsyn <- read.table(file='symbol2synonym.txt',header=FALSE)
names(symsyn) <- c('symbol','synonyms')
source('functions.R', local = TRUE)
source('tools/intro.R', local = TRUE)
source('tools/gc.R', local = TRUE)
source('tools/qc.R', local = TRUE)
source('tools/altsplice.R', local = TRUE)
source('tools/dea.R', local = TRUE)
source('tools/gsea.R', local = TRUE)
#source('tools/fastqc.R', local = TRUE)
source('tools/gc_ui.R', local = TRUE)
source('tools/qc_ui.R', local = TRUE)
source('tools/dea_ui.R', local = TRUE)
source('tools/gsea_ui.R', local = TRUE)
#source('tools/fastqc_ui.R', local = TRUE)
source('tools/altsplice_ui.R', local = TRUE)
})
library(shiny)
library(Vennerable)
library(qusage)
library(DT)
library(ggplot2)
library(ballgown)
library(sqldf)
library(reshape2)
library("gplots")
shinyServer(function(input, output, session) {
data.dir <- Sys.getenv('outputDir')
#data.dir <-"../../vizapp/workflow_322_output"
rda.dir <- data.dir
symsyn <- read.table(file='symbol2synonym.txt',header=FALSE)
names(symsyn) <- c('symbol','synonyms')
source('functions.R', local = TRUE)
source('tools/intro.R', local = TRUE)
source('tools/intro_ui.R', local = TRUE)
source('tools/gc.R', local = TRUE)
source('tools/qc.R', local = TRUE)
source('tools/altsplice.R', local = TRUE)
source('tools/dea.R', local = TRUE)
source('tools/gsea.R', local = TRUE)
#source('tools/fastqc.R', local = TRUE)
source('tools/gc_ui.R', local = TRUE)
source('tools/qc_ui.R', local = TRUE)
source('tools/dea_ui.R', local = TRUE)
source('tools/gsea_ui.R', local = TRUE)
#source('tools/fastqc_ui.R', local = TRUE)
source('tools/altsplice_ui.R', local = TRUE)
})
get.bgdata <- function(var) {
rdafile <- paste(data.dir,'bg.rda',sep='/')
load(rdafile)
genetrx <- indexes(bg)$t2g
genetrx$SYMBOL <- geneNames(bg)
geneid <- as.character(unique(genetrx[genetrx$SYMBOL %in% var$symsearch,]$g_id))
if (exists("var$enssearch") && var$enssearch != '') {
geneid <- as.character(unique(genetrx$g_id[grep(enssearch,genetrx$g_id)]))
}
genebg <- subset(bg,paste0("gene_id=='",geneid,"'"))
rownames(genetrx) <- transcriptNames(bg)
fpkm <- texpr(bg,meas='FPKM')
keep <- genetrx[genetrx$g_id == geneid,]$t_id
cttbl <- fpkm[keep,]
grps <- pData(bg)$group
trxnames <- genetrx[genetrx$g_id == geneid,]
dtm <- melt(t(cttbl))
names(dtm) <- c('sample','transcript','value')
dtm$grps <- rep(grps,length(keep))
if (length(keep) < 2) {
dtm$transcript <- rep(keep,length(keep))
}
dtm$grptrx <- paste(dtm$grps,dtm$transcript,sep='.')
test <- stattest(gown=genebg, pData=pData(bg), feature='transcript',covariate='group', libadjust=FALSE,getFC=TRUE)
if (length(keep) > 1) {
agg = collapseTranscripts(gene=geneid, gown=bg, k=var$kct, method='kmeans')
test <- stattest(gowntable=agg$tab, pData=pData(bg), feature='transcript_cluster',covariate='group', libadjust=FALSE,getFC=TRUE)
}
return(list(stattbl=test,gid=geneid,obj=genebg,cttbl=dtm,tname=trxnames))
}
get.bgdata <- function(var) {
rdafile <- paste(data.dir, 'bg.rda', sep = '/')
load(rdafile)
genetrx <- indexes(bg)$t2g
genetrx$SYMBOL <- geneNames(bg)
geneid <- as.character(unique(genetrx[genetrx$SYMBOL %in% var$symsearch, ]$g_id))
if (exists("var$enssearch") && var$enssearch != '') {
geneid <- as.character(unique(genetrx$g_id[grep(enssearch, genetrx$g_id)]))
}
genebg <- subset(bg, paste0("gene_id=='", geneid, "'"))
rownames(genetrx) <- transcriptNames(bg)
fpkm <- texpr(bg, meas = 'FPKM')
keep <- genetrx[genetrx$g_id == geneid, ]$t_id
cttbl <- fpkm[keep, ]
grps <- pData(bg)$group
trxnames <- genetrx[genetrx$g_id == geneid, ]
dtm <- melt(t(cttbl))
names(dtm) <- c('sample', 'transcript', 'value')
dtm$grps <- rep(grps, length(keep))
if (length(keep) < 2) {
dtm$transcript <- rep(keep, length(keep))
}
dtm$grptrx <- paste(dtm$grps, dtm$transcript, sep = '.')
test <-
stattest(
gown = genebg,
pData = pData(bg),
feature = 'transcript',
covariate = 'group',
libadjust = FALSE,
getFC = TRUE
)
if (length(keep) > 1) {
agg = collapseTranscripts(
gene = geneid,
gown = bg,
k = var$kct,
method = 'kmeans'
)
test <-
stattest(
gowntable = agg$tab,
pData = pData(bg),
feature = 'transcript_cluster',
covariate = 'group',
libadjust = FALSE,
getFC = TRUE
)
}
return(list(
stattbl = test,
gid = geneid,
obj = genebg,
cttbl = dtm,
tname = trxnames
))
}
find_sym <- function (sym) {
if (!(toupper(sym) %in% toupper(symsyn$symbol))) {
syns <- symsyn[grep(input$symsearch,symsyn$synonym,ignore.case=TRUE),]$symbol
synlist <- paste(as.character(syns),collapse=',')
if (length(syns) > 1) {
paste("Please Use Official Gene Symbols",synlist,sep=':')
}else {"Please Use Official Gene Symbols"}
}else {NULL}
}
getgeneid <- eventReactive(input$altButton,{
if (input$symsearch == '' & input$enssearch != '') {
if (!(toupper(sym) %in% toupper(symsyn$symbol))) {
syns <-
symsyn[grep(input$symsearch, symsyn$synonym, ignore.case = TRUE), ]$symbol
synlist <- paste(as.character(syns), collapse = ',')
if (length(syns) > 1) {
paste("Please Use Official Gene Symbols", synlist, sep = ':')
} else {
"Please Use Official Gene Symbols"
}
} else {
NULL
}
}
getgeneid <- eventReactive(input$altButton, {
if (input$symsearch == '' & input$enssearch != '') {
validate (find_sym(input$symsearch))
}
get.bgdata(input)
}
get.bgdata(input)
})
output$plot.cluster <- renderPlot({
par(oma=c(4,4,1,1))
par(oma = c(0, 0, 0, 0))
gid <- getgeneid()$gid
bg <- getgeneid()$obj
tname <- getgeneid()$tname
if (nrow(tname) > 1) {
plotLatentTranscripts(gene=gid, gown=bg, k=input$kct, method='kmeans', returncluster=FALSE)
plotLatentTranscripts(
gene = gid,
gown = bg,
k = input$kct,
method = 'kmeans',
returncluster = FALSE
)
}
}, height = "auto", width = 'auto')
output$dlplotcluster <- renderUI({
if (!is.null(getgeneid()$gid) &
!is.null(getgeneid()$obj) & !is.null(getgeneid()$tname)) {
downloadButton('Downloadpcluster', 'Download PNG')
}
})
output$Downloadpcluster <- downloadHandler(
file = function() {
paste('plotcluster', 'png', sep = ".")
},
content <- function(file) {
png(
file,
width = 8 * 120,
height = 4 * 120,
res = 120,
pointsize = 5
)
if (is.null(getgeneid()$gid) &
is.null(getgeneid()$obj) & is.null(getgeneid()$tname)) {
return()
}
par(oma = c(4, 4, 1, 1))
gid <- getgeneid()$gid
bg <- getgeneid()$obj
tname <- getgeneid()$tname
if (nrow(tname) > 1) {
plot <-
plotLatentTranscripts(
gene = gid,
gown = bg,
k = input$kct,
method = 'kmeans',
returncluster = FALSE
)
}
print(plot)
dev.off()
},
contentType = "image/png"
)
output$gene.stat <- DT::renderDataTable({
getgeneid()$stattbl
},escape=FALSE)
}, escape = FALSE, options = list(dom = 'ft'))
output$trx.name <- DT::renderDataTable({
getgeneid()$tname
},escape=FALSE)
}, escape = FALSE , options = list(dom = 'ft'))
output$downloadTrx <- renderUI({
if (!is.null(getgeneid()$tname)) {
downloadButton('DownloadTrx', 'Download CSV')
}
})
output$DownloadTrx <- downloadHandler(
file <- paste('Gene.trxname.txt'),
content = function(file) {
write.table(
getgeneid()$tname,
file,
quote = FALSE,
row.names = FALSE,
sep = '\t'
)
}
)
output$plot.means <- renderPlot({
par(oma=c(4,4,1,1))
par(cex.main=0.75)
par(oma = c(1, 4, 1, 8))
par(cex.main = 0.75)
gid <- getgeneid()$gid
bg <- getgeneid()$obj
tname <- getgeneid()$tname
if (nrow(tname) > 1) {
plotMeans(gid, bg, groupvar='group', meas='cov', colorby='transcript')
plotMeans(gid,
bg,
groupvar = 'group',
meas = 'cov',
colorby = 'transcript')
}
}, height = 900, width = "auto")
output$dlplotmean <- renderUI({
if (!is.null(getgeneid()$gid) & !is.null(getgeneid()$obj) & !is.null(getgeneid()$tname)) {
downloadButton('Downloadmean', 'Download PNG')
}
}, height = 900, width = 900)
})
output$Downloadmean <- downloadHandler(
file = function() {
paste('mean', 'png', sep = ".")
},
content <- function(file) {
png(
file,
width = 8 * 120,
height = 4 * 180,
res = 140,
pointsize = 5
)
if (is.null(getgeneid()$gid) & is.null(getgeneid()$obj) & is.null(getgeneid()$tname)) {
return()
}
countTable <- getgeneid()$cttbl
par(oma = c(1, 4, 1, 8))
par(cex.main = 0.75)
gid <- getgeneid()$gid
bg <- getgeneid()$obj
tname <- getgeneid()$tname
if (nrow(tname) > 1) {
plot <- plotMeans(gid,
bg,
groupvar = 'group',
meas = 'cov',
colorby = 'transcript')
}
print(plot)
dev.off()
},
contentType = "image/png"
)
output$trx.gene <- renderPlot({
countTable <- getgeneid()$cttbl
par(oma=c(4,4,1,1))
p <- ggplot(countTable,aes(x=grptrx,y=log2(value+1))) + geom_boxplot(aes(fill = factor(grptrx))) + geom_jitter(height = 0) + theme(legend.position="left",axis.text.x=element_text(angle=45,hjust=1, vjust=1),legend.key.height=unit(0.5,"line"),legend.text=element_text(size=8),legend.title=element_blank()) + ylab("Relative Abundance (FPKM)") + xlab("")
par(oma = c(4, 1, 4, 1))
p <-
ggplot(countTable, aes(x = grptrx, y = log2(value + 1))) + geom_boxplot(aes(fill = factor(grptrx))) +
geom_jitter(height = 0) + theme_bw() +
theme(
legend.position = "left",
axis.text.x = element_text(
angle = 45,
hjust = 1,
vjust = 1
),
legend.key.height = unit(0.5, "line"),
legend.text = element_text(size = 8),
legend.title = element_blank()
) + ylab("Relative Abundance (FPKM)") + xlab("")
print(p)
}, height = 'auto', width = 'auto')
output$dltrxgene <- renderUI({
if (!is.null(getgeneid()$cttbl)) {
downloadButton('Downloadtregene', 'Download PNG')
}
})
output$Downloadtregene <- downloadHandler(
file = function() {
paste('tregene', 'png', sep = ".")
},
content <- function(file) {
png(
file,
width = 8 * 160,
height = 4 * 120,
res = 120,
pointsize = 5
)
if (is.null(getgeneid()$cttbl)) {
return()
}
countTable <- getgeneid()$cttbl
par(oma = c(4, 1, 4, 1))
p <-
ggplot(countTable, aes(x = grptrx, y = log2(value + 1))) + geom_boxplot(aes(fill = factor(grptrx))) +
geom_jitter(height = 0) + theme_bw() +
theme(
legend.position = "left",
axis.text.x = element_text(
angle = 45,
hjust = 1,
vjust = 1
),
legend.key.height = unit(0.5, "line"),
legend.text = element_text(size = 8),
legend.title = element_blank()
) + ylab("Relative Abundance (FPKM)") + xlab("")
print(p)
dev.off()
},
contentType = "image/png"
)
output$ui_altsplice <- renderUI({
fluidPage(
includeCSS("www/style.css"),
sidebarLayout(
sidebarPanel(
uiOutput("dir.splice"),
textInput("symsearch", "Search By Gene Symbol", 'IL1B'),
textInput("enssearch", "Search By ENS ID",''),
numericInput("kct",label = "Number of Clusters",value = 2),
actionButton("altButton", "GO")
),
mainPanel(
tabsetPanel(
tabPanel("Gene Compare",
dataTableOutput('gene.stat'),
dataTableOutput('trx.name'),
plotOutput("plot.cluster"),
plotOutput("trx.gene"),
plotOutput("plot.means")
)
)
)
)
)
list(fluidPage(
includeCSS("www/style.css"),
fluidRow(
sidebarPanel(
uiOutput("dir.splice"),
textInput("symsearch", "Search By Gene Symbol", 'IL1B'),
textInput("enssearch", "Search By ENS ID", ''),
numericInput("kct", label = "Number of Clusters", value = 2),
actionButton("altButton", "GO", class = "btn btn-primary btn-bg centerbtn")
),
column(7,
tabsetPanel(
tabPanel(
"Gene Compare",
br(),br(),
dataTableOutput('gene.stat'),
br(),br(),
dataTableOutput('trx.name'),
br(),
uiOutput('downloadTrx'),
plotOutput("plot.cluster"),
uiOutput('dlplotcluster'),
br(),br(),
plotOutput("trx.gene"),
uiOutput('dltrxgene'),
plotOutput("plot.means", height = 900, width = 'auto'),
uiOutput('dlplotmean')
)
))
)
))
})
This diff is collapsed.
#######################################
# Shiny interface for data functions
#######################################
# data ui and tabs
output$ui_dea <- renderUI({
list(
includeCSS("www/style.css"),
sidebarLayout(
sidebarPanel(
numericInput("fc.thresh",
"LogFold Change Threshold:", 1),
numericInput("gene.thresh",
"Gene Mean Threshold:", 5),
numericInput("pval.thresh",
"P-Value Threshold:", 0.05),
selectInput("adjust", "Choose P-Value Correction", choices=c("raw","FDR",'BONF'),selected='FDR'),
uiOutput("pick.dea"),
selectInput(
"heatmap", "HeatMap",
c(Top = "top",
HallmarkGeneSet = "hgeneset",
CustomGeneSet= "cgeneset")
),
conditionalPanel(
condition = "input.heatmap == 'cgeneset'",
uiOutput("selectgenes")
),
conditionalPanel(
condition = "input.heatmap == 'top'",
numericInput("numgenes","Number Top Genes:", 50)
),
conditionalPanel(
condition = "input.heatmap == 'hgeneset'",
uiOutput("pick.pathway")
),
actionButton("deButton", "Go")
),
mainPanel(
tabsetPanel(
tabPanel("Differential Gene Set Comparison",
downloadButton('downloadC', 'Download CSV'),
dataTableOutput('dge.c')
),
tabPanel("Heatmap",
#downloadButton('downloadpdf', 'Download Heatmap'),
h1("HeatMap Comparison"),
h3("Top 50 User Defined Genes"),
plotOutput("plot.heatmap"),
h3("All Differentially Expressed Genes"),
imageOutput("hm.comp",width="1200px",height="1200px"),
textOutput("hmcomp.desc")
)
)
)
)
)
})
#######################################
# Shiny interface for data functions
#######################################
# data ui and tabs
output$ui_dea <- renderUI({
list(fluidPage(
includeCSS("www/style.css"),
fluidRow(
sidebarPanel(
numericInput("fc.thresh",
"LogFold Change Threshold:", 1),
numericInput("gene.thresh",
"Gene Mean Threshold:", 5),
numericInput("pval.thresh",
"P-Value Threshold:", 0.05),
selectInput(
"adjust",
"Choose P-Value Correction",
choices = c("raw", "FDR", 'BONF'),
selected = 'FDR'
),
uiOutput("pick.dea"),
selectInput(
"heatmap",
"HeatMap",
c(
Top = "top",
HallmarkGeneSet = "hgeneset",
CustomGeneSet = "cgeneset"
)
),
conditionalPanel(condition = "input.heatmap == 'cgeneset'",
uiOutput("selectgenes")),
conditionalPanel(condition = "input.heatmap == 'top'",
numericInput("numgenes", "Number Top Genes:", 50)),
conditionalPanel(condition = "input.heatmap == 'hgeneset'",
uiOutput("pick.pathway")),
radioButtons("cluster", label = "Cluster",
choices = list("Display" = 1, "Hide" = 2),
selected = 1),
actionButton("deButton", "Go", class = "btn btn-primary btn-bg centerbtn")
),
column(7,
tabsetPanel(
tabPanel(
"Differential Gene Set Comparison",
br(),br(),
dataTableOutput('dge.c', width = '100%'),
uiOutput('downloadC')
),
tabPanel(
"Heatmap Comparison",
#downloadButton('downloadpdf', 'Download Heatmap'),
h3(textOutput("hmcomp.top50")),
plotOutput("plot.heatmap"),