Skip to content

Instantly share code, notes, and snippets.

Embed
What would you like to do?
comparisons for GSEA
# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #
# GLOBAL #
# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #
if (!require(shiny)) install.packages("shiny")
if (!require(shinydashboard)) install.packages("shinydashboard")
require(shiny)
require(shinydashboard)
require(stringr)
require(data.table)
require(dplyr)
# primary samplesheet
primary_ss_csv <- "/Volumes/MOLECPATH/shared/TPU/Runs/NovaSeq/815/run_815_SAMPLESHEET.csv"
primary_ss <- fread(primary_ss_csv)
sample_info_dt <- unique(primary_ss[,.(Sample_ID, Sample_Name)])
# specify files and dir for analysis
analysis_template <- "GSEA_report_RNASeq_template.Rmd"
projdir <- "../"
# specify all the factors to be used in comparisons
factors <- list(
level1 = c("786-O","A498","U2OS"),
level2 = c("WT","PBRM1"),
level3 = c("IFN")
)
# parse sample names (i.e. remove replicate numbers etc.)
# helper functions
replace_NAs_DT <- function(DT, replacement='', cols=colnames(DT)) {
for (j in cols)
set(DT,which(is.na(DT[[j]])),j,replacement)
}
paste_cols_DT <- function(DT, sep='', cols=colnames(DT), colname='newcolumn') {
DT[,I(colname):=Reduce(function(...) paste(..., sep=sep), .SD[,..cols])]
}
# extract factors from sample names
fv <- unlist(unname(factors))
factors_dt <- sample_info_dt[,tstrsplit(Sample_Name,split=" ")][
,lapply(.SD,str_extract,pattern=paste(fv,collapse="|"))]
# remove NAs for empty factors
replace_NAs_DT(factors_dt)
# make unique groups
paste_cols_DT(factors_dt, sep=' ', colname='group')
# replace spaces in group labels
sample_info_dt[,sampleConditions:=factors_dt$group][
,sampleConditions:=factors_dt$group][,sampleConditions:=str_trim(sampleConditions)][
,sampleConditions:=str_replace_all(sampleConditions,"\\s+","_")]
# separate comparison levels
sample_info_dt[,names(factors):=lapply(names(factors), function(x)
str_extract(sampleConditions, paste(factors[[x]],collapse="|")))]
# get all the groups that will be compared
comparisons <- function(DT, factors, group_level, compare_level, inter=NULL)
{
groups <- unlist(factors[group_level])
compare <- unlist(factors[compare_level])
if (is.numeric(group_level)) group_level <- names(factors)[group_level]
if (is.numeric(compare_level)) compare_level <- names(factors)[compare_level]
if (!any((compare_level %in% group_level))) {
DT_list <- split(DT, by=group_level)
not_used <- !(names(factors) %in% c(group_level, compare_level))
out_list <- lapply(DT_list, function(x) {
DT <- x[,names(factors)[not_used]:=NULL]
group <- unique(unlist(DT[,..group_level]))
if (length(group)>1)
group <- paste(group,collapse="_")
treatment <- unique(unlist(DT[,..compare_level]))
if (any(is.na(treatment)))
treatment <- treatment[!is.na(treatment)]
if (length(treatment)==1) {
control_gp <- group
treated_gp <- sprintf("%s_%s",group,treatment)
} else {
control_gp <- sprintf("%s_%s",group,treatment[1])
treated_gp <- sprintf("%s_%s",group,treatment[2])
}
Exp_Name <- sprintf("%s_vs_%s",control_gp,treated_gp)
Exp_Name
#list(Exp_Name=Exp_Name, control_gp=control_gp, treated_gp=treated_gp, sample_info=DT)
})
} else if (any(compare_level %in% group_level)) {
out_list <- lapply(compare_level, function(cl) {
constant_level <- group_level[group_level != cl]
DT_list <- split(DT, by=constant_level)
not_used <- !(names(factors) %in% c(group_level, compare_level))
lapply(DT_list, function(x) {
ns <- x[,..group_level]
replace_NAs_DT(DT=ns, replacement='')
names <- unique(paste_cols_DT(ns, sep='_', colname='name')[,name])
names <- str_replace(str_trim(names), pattern='_$', replacement='')
dm <- matrix(nrow=length(names),ncol=length(names))
for (i in seq_along(names)) {
for (j in seq_along(names)) {
g1 <- names[i]
g2 <- names[j]
dm[i,j] <- sprintf("%s_vs_%s", g1, g2)
}
}
dm[upper.tri(dm) & !is.na(dm)]
})
})
}
names(out_list) <- compare_level
out_list
}
# function to generate reports
render_report <- function(Exp_Name, sample_info, analysis_template) {
rmarkdown::render(
analysis_template,
params = list(
Exp_Name = Exp_Name,
sample_info = sample_info
),
output_file = paste0("GSEA_report_",Exp_Name,".docx")
)
}
# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #
# SERVER #
# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #
server <- function(input, output, session) {
# table of level options
output$levels <- renderTable({
cbind.data.frame(
levels = rep(names(factors),sapply(factors, length),each=TRUE),
values = c(
c("786-O","A498","U2OS"),
c("WT","PBRM1"),
c("IFN")
)
)
})
# table of comparisons i.e. experiments
clist <- reactive({
comparisons(sample_info_dt,factors=factors,group_level=input$group_level,compare_level=input$compare_level)
})
output$comparisons <- renderTable({
arrange(cbind.data.frame(Exp_Name = unname(unlist(clist()))), Exp_Name)
})
# analyses
eventReactive(
eventExpr = input$analysis,
valueExpr = cat("Mock analysis for comparisons starting with ", output$comparisons[1,]),
ignoreNULL = FALSE
)
# # report download
# output$downloa_report <- downloadHandler(
#
# filename = function() {
# paste('report', sep = '.', 'docx') #switch(input$format, PDF='pdf', HTML='html', Word='docx'))
# },
#
# content = function(file) {
# src <- normalizePath('report.Rmd')
# # temporarily switch to the temp dir, in case you do not have write
# # permission to the current working directory
# owd <- setwd(tempdir())
# on.exit(setwd(owd))
# file.copy(src, 'report.Rmd', overwrite=TRUE)
# out <- render_report(
# Exp_Name = input$Exp_Name, sample_info = input$sample_info, analysis_template = switch(input$analysis, GSEA='')
# )
# file.rename(out, file)
# }
# )
}
# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #
# UI #
# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #
ui <- pageWithSidebar(
headerPanel(' '),
sidebarPanel(
tableOutput('levels'),
selectInput('group_level', 'grouping levels', names(factors), selected=names(factors)[1], multiple=TRUE, selectize=TRUE),
selectInput('compare_level', 'comparison levels', names(factors), selected=names(factors)[1], multiple=TRUE, selectize=TRUE),
submitButton("analysis")
# downloadButton('download_report')
),
mainPanel(
tableOutput('comparisons')
)
)
shinyApp(ui = ui, server = server)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
You can’t perform that action at this time.