fixed some gene set issues mainly

R/main.r

@@ -38,88 +38,118 @@ cGSEAcore = function(ElistObject, contrast, geneset = "H", specie = "Mus musculu
   res = list()
   time = list()
 
-  if(gsva == FALSE && ssgsea == FALSE && zscore == FALSE && plage == FALSE){
-    cat("Please set at least one of the methods listed below to TRUE to run cGSEA : \n(mandatory for logFC calculations)\n - gsva\n - ssgsea\n - zscore\n - plage")
-    return()
-  }
 
   if(camera == TRUE){
     print("RUNNING CAMERA")
-    camerastart = Sys.time()
-    res$camera = runcamera(voom.results = ElistObject, contrast = contrast, genesetIdx = genesetIdx, num.workers = num.workers, verbose = verbose)
-    camerastop = Sys.time()
-    time$camera = difftime(camerastop, camerastart, units = c("secs"))
+    try({
+        camerastart = Sys.time()
+        res$camera = runcamera(voom.results = ElistObject, contrast = contrast, genesetIdx = genesetIdx, num.workers = num.workers, verbose = verbose)
+        camerastop = Sys.time()
+        time$camera = difftime(camerastop, camerastart, units = c("secs"))
+
+        })
+
   }
 
   if (gage == TRUE){
     print("RUNNING GAGE")
-    gagestart = Sys.time()
-    res$gage = rungage(voom.results = ElistObject, contrast = contrast, genesetIdx = genesetIdx, num.workers = num.workers, verbose = verbose)
-    gagestop = Sys.time()
-    time$gage = difftime(gagestop,gagestart, units = c("secs"))
+    try({
+        gagestart = Sys.time()
+        res$gage = rungage(voom.results = ElistObject, contrast = contrast, genesetIdx = genesetIdx, num.workers = num.workers, verbose = verbose)
+        gagestop = Sys.time()
+        time$gage = difftime(gagestop,gagestart, units = c("secs"))
+        })
+
   }
   if (globaltest == TRUE){
     print("RUNNING GLOBALTEST")
-    globalteststart = Sys.time()
-    res$globaltest = runglobaltest(voom.results = ElistObject, contrast = contrast, genesetIdx = genesetIdx, num.workers = num.workers, verbose = verbose)
-    globalteststop = Sys.time()
-    time$globaltest = difftime(globalteststop,globalteststart, units = c("secs"))
+    try({
+        globalteststart = Sys.time()
+        res$globaltest = runglobaltest(voom.results = ElistObject, contrast = contrast, genesetIdx = genesetIdx, num.workers = num.workers, verbose = verbose)
+        globalteststop = Sys.time()
+        time$globaltest = difftime(globalteststop,globalteststart, units = c("secs"))
+        })
+
   }
   if (gsva == TRUE){
     print("RUNNING GSVA")
-    gsvastart = Sys.time()
-    res$gsva = rungsva(method = "gsva", voom.results = ElistObject, contrast = contrast,  genesetIdx = genesetIdx, num.workers = num.workers, verbose = verbose)
-    gsvastop = Sys.time()
-    time$gsva = difftime(gsvastop,gsvastart, units = c("secs"))
+    try({
+        gsvastart = Sys.time()
+        res$gsva = rungsva(method = "gsva", voom.results = ElistObject, contrast = contrast,  genesetIdx = genesetIdx, num.workers = num.workers, verbose = verbose)
+        gsvastop = Sys.time()
+        time$gsva = difftime(gsvastop,gsvastart, units = c("secs"))
+        })
+
   }
   if (ssgsea == TRUE){
     print("RUNNING SSGSEA")
-    ssgseastart = Sys.time()
-    res$ssgsea = rungsva(method = "ssgsea", voom.results = ElistObject, contrast = contrast,  genesetIdx = genesetIdx, num.workers = num.workers, verbose = verbose)
-    ssgseastop = Sys.time()
-    time$ssgsea = difftime(ssgseastop,ssgseastart, units = c("secs"))
+    try({
+        ssgseastart = Sys.time()
+        res$ssgsea = rungsva(method = "ssgsea", voom.results = ElistObject, contrast = contrast,  genesetIdx = genesetIdx, num.workers = num.workers, verbose = verbose)
+        ssgseastop = Sys.time()
+        time$ssgsea = difftime(ssgseastop,ssgseastart, units = c("secs"))
+        })
+
   }
   if (zscore == TRUE){
     print("RUNNING ZSCORE")
-    zscorestart = Sys.time()
-    res$zscore = rungsva(method = "zscore", voom.results = ElistObject, contrast = contrast,  genesetIdx = genesetIdx, num.workers = num.workers, verbose = verbose)
-    zscorestop = Sys.time()
-    time$zscore = difftime(zscorestop,zscorestart, units = c("secs"))
+    try({
+        zscorestart = Sys.time()
+        res$zscore = rungsva(method = "zscore", voom.results = ElistObject, contrast = contrast,  genesetIdx = genesetIdx, num.workers = num.workers, verbose = verbose)
+        zscorestop = Sys.time()
+        time$zscore = difftime(zscorestop,zscorestart, units = c("secs"))
+        })
+
   }
   if (plage == TRUE){
-    print("RUNNING PLAGE")
-    plagestart = Sys.time()
-    res$plage = rungsva(method = "plage", voom.results = ElistObject, contrast = contrast,  genesetIdx = genesetIdx, num.workers = num.workers, verbose = verbose)
-    plagestop = Sys.time()
-    time$plage = difftime(plagestop,plagestart, units = c("secs"))
+    print("RUNNING PLAGE ")
+    try({
+        plagestart = Sys.time()
+        res$plage = rungsva(method = "plage", voom.results = ElistObject, contrast = contrast,  genesetIdx = genesetIdx, num.workers = num.workers, verbose = verbose)
+        plagestop = Sys.time()
+        time$plage = difftime(plagestop,plagestart, units = c("secs"))
+        })
+
   }
   if (ora == TRUE){
     print("RUNNING ORA")
-    orastart = Sys.time()
-    res$ora = runora(voom.results = ElistObject, contrast = contrast,  genesetIdx = genesetIdx, num.workers = num.workers, verbose = verbose)
-    orastop = Sys.time()
-    time$ora = difftime(orastop, orastart, units = c("secs"))
+    try({
+        orastart = Sys.time()
+        res$ora = runora(voom.results = ElistObject, contrast = contrast,  genesetIdx = genesetIdx, num.workers = num.workers, verbose = verbose)
+        orastop = Sys.time()
+        time$ora = difftime(orastop, orastart, units = c("secs"))
+        })
+
   }
   if (padog == TRUE){
     print("RUNNING PADOG")
-    padogstart = Sys.time()
-    res$padog = runpadog(voom.results = ElistObject, contrast = contrast,  genesetIdx = genesetIdx, num.workers = num.workers, verbose = verbose)
-    padogstop = Sys.time()
-    time$padog = difftime(padogstop, padogstart, units = c("secs"))
+    try({
+        padogstart = Sys.time()
+        res$padog = runpadog(voom.results = ElistObject, contrast = contrast,  genesetIdx = genesetIdx, num.workers = num.workers, verbose = verbose)
+        padogstop = Sys.time()
+        time$padog = difftime(padogstop, padogstart, units = c("secs"))
+        })
+
   }
   if (roast == TRUE){
     print("RUNNING ROAST")
-    roaststart = Sys.time()
-    res$roast = runroast(voom.results = ElistObject, contrast = contrast,  genesetIdx = genesetIdx, num.workers = num.workers, verbose = verbose)
-    roaststop = Sys.time()
-    time$roast = difftime(roaststop,roaststart, units = c("secs"))
+    try({
+        roaststart = Sys.time()
+        res$roast = runroast(voom.results = ElistObject, contrast = contrast,  genesetIdx = genesetIdx, num.workers = num.workers, verbose = verbose)
+        roaststop = Sys.time()
+        time$roast = difftime(roaststop,roaststart, units = c("secs"))
+        })
+
   }
   if (safe == TRUE){
     print("RUNNING SAFE")
-    safestart = Sys.time()
-    res$safe = runsafe(voom.results = ElistObject, contrast = contrast,  genesetIdx = genesetIdx, num.workers = num.workers, verbose = verbose)
-    safestop = Sys.time()
-    time$safe = difftime(safestop, safestart, units = c("secs"))
+    try({
+        safestart = Sys.time()
+        res$safe = runsafe(voom.results = ElistObject, contrast = contrast,  genesetIdx = genesetIdx, num.workers = num.workers, verbose = verbose)
+        safestop = Sys.time()
+        time$safe = difftime(safestop, safestart, units = c("secs"))
+        })
+
   }
   # if (setrank == TRUE){
   #   print("RUNNING SETRANK")

R/runners.r

@@ -177,7 +177,11 @@ rungage <- function(voom.results, contrast, genesetIdx,
                 mc.cores=num.workers)
     #stop("here")
     names(gage.results) = colnames(contrast)
+    for (i in length(names(gage.results))){
+      gage.results[[i]][which(gage.results[[i]]$set.size == 1),"p.value"] = 1
+    }
     return(gage.results)
+    # return(tmp)
 }
 
 rungage.contrast <- function(args){
@@ -190,7 +194,7 @@ rungage.contrast <- function(args){
     gage.results = gage(exprs=args$logCPM,
             gsets=args$gsets,
             ref=args$group1,
-            samp=args$group2, same.dir=FALSE,
+            samp=args$group2, set.size = c(0,2000), same.dir=FALSE,
             compare = "unpaired")$greater[, 1:5]
     # returns PropDown/PropUp ==> proportion of genes that are
 # down/up-regulated
@@ -539,7 +543,7 @@ runpadog.contrast <- function(args){
 
     padog.results = padog(esetm=args$logCPM,
             group=group, paired=FALSE,
-            gslist=args$gsets, NI=100,
+            gslist=args$gsets, NI=100, Nmin = 1,
             verbose=FALSE)
 
     padog.results = padog.results[order(
@@ -644,7 +648,7 @@ runsafe <- function(voom.results, contrast, genesetIdx,
 
     groupData = prepareTwoGroupsData(voom.results, contrast, genesetIdx,
             min.samples = 3, verbose)
-    capture.output(C.mat <- getCmatrix(keyword.list=groupData$gsets,
+    capture.output(C.mat <- getCmatrix(keyword.list=groupData$gsets, min.size = 1,
                     present.genes=rownames(groupData$data[[1]]$logCPM)))
 #    safe.results = vector("list", ncol(contrast))
     args.all = list()
@@ -682,6 +686,8 @@ runsafe.contrast <- function(args){
     safe.results = safe.toptable(temp,
             number=length(args$gsets),
             description = FALSE)
+    # safe.results = na.omit(safe.results)
+    # rownames(safe.results) = names(args$genesetIdx)
     safe.results[, "P.value"] = as.numeric(
             safe.results[, "P.value"])
     safe.results[, "Adj.p.value"] = as.numeric(

R/utils.r

@@ -151,6 +151,8 @@ genesetToIdx <- function(geneset = "H", specie = "Mus musculus", entrezGenesIds)
             load(url("http://bioinf.wehi.edu.au/software/MSigDB/mouse_c2_v5p2.rdata"))
             Mm.c2.subset = Mm.c2[grep("REACTOME",attributes(Mm.c2)$names)]
             genesetToIdx = ids2indices(Mm.c2.subset, entrezGenesIds)
+        } else if (class(geneset) == "list"){
+            genesetToIdx = ids2indices(geneset, entrezGenesIds, remove.empty=TRUE)
         } else {
             cat(paste(geneset , "is not a valid MSigDB geneset. Available genesets : \n - H \n - C2"))
             return(NULL)
@@ -167,6 +169,8 @@ genesetToIdx <- function(geneset = "H", specie = "Mus musculus", entrezGenesIds)
             load(url("http://bioinf.wehi.edu.au/software/MSigDB/human_c2_v5p2.rdata"))
             Hs.c2.subset = Hs.c2[grep("REACTOME",attributes(Hs.c2)$names)]
             genesetToIdx = ids2indices(Hs.c2.subset, entrezGenesIds)
+        } else if (class(geneset) == "list"){
+            genesetToIdx = ids2indices(geneset, entrezGenesIds, remove.empty=TRUE)
         } else {
             cat(paste(geneset , "is not a valid MSigDB geneset. Available genesets : \n - H \n - C2"))
             return(NULL)
@@ -343,10 +347,10 @@ cGSEAOutputTable = function(cGSEAcoreOutput, contrastLevel = contrastLevel, gene
   for (method in names(cGSEAcoreOutput$res)){
     method_stats = rank_getter(geneSetname = geneSetLevels,method_result = cGSEAcoreOutput[["res"]][[method]][[contrastLevel]])
 
-    dt[[paste(method,"Rank",sep="_")]] = method_stats[["ranks"]]
+    dt[[paste(method,"Rank",sep="_")]] = as.numeric(method_stats[["ranks"]])
     # ranks[[method]] = rank_getter(geneSetname = geneSetLevels,method_result = cGSEAOutput[["res"]][[method]][[contrastLevel]])
     #p.values
-    dt[[paste(method,"p.value",sep="_")]] = method_stats[["pvals"]]
+    dt[[paste(method,"p.value",sep="_")]] = as.numeric(method_stats[["pvals"]])
   }
 
   # df = data.frame(matrix(unlist(dt), nrow = length(geneSetLevels), byrow = FALSE),stringsAsFactors=FALSE)
@@ -582,8 +586,10 @@ prepareData = function(cGSEAcoreOutput, alpha = 0.05, directoryPath, pvalAdjMeth
     # resumPlot2Table[[condi]] = as.data.frame(resumPlot2[[condi]])
     # resumPlot2Table[[condi]][,c("x.data","y.data","dir","rank","sig", "gsSize")] = apply(resumPlot2[[condi]][,c("x.data","y.data","dir","rank","sig", "gsSize")],2, as.numeric)
   }
+  if (length(colnames(cGSEAcoreOutput$contrast)) > 1){
+    output$comparison = comparisonSummaryData(resumPlot2)
+  }
 
-  output$comparison = comparisonSummaryData(resumPlot2)
   # output$comparisonTable = comparisonSummaryData(resumPlot2Table)
 
   output$time = time
@@ -602,4 +608,3 @@ prepareData = function(cGSEAcoreOutput, alpha = 0.05, directoryPath, pvalAdjMeth
   return(output)
 
 }
-