Change QuSAGE ID labeling to highest positive enrichment score, reguardless of FDR

r.scripts/sc-TissueMapper.R

@@ -489,26 +489,32 @@ scQuSAGE <- function(sc10x,gs,res.use=0.1,ds=25000,nm="Lin",folder="lin"){
   }
   
   #Generate ID table
-  for (i in 1:number.clusters){
-    qs <- qsTable(get(paste0("results.",i)),number=length(gs))
-    qs <- qs[(qs$FDR<0.05) & (qs$log.fold.change>0),-3]
-    if (nrow(qs)==0){
-      qs <- NULL
-      qs$pathway.name <- "Unknown"
-      qs$log.fold.change <- 0
-      qs$FDR <- 0
-      qs <- as.data.frame(qs)
-    }
+  results.cor <- NULL
+  results.cor <- qsTable(results.1)
+  results.cor$Cluster <- 1
+  for (i in 2:number.clusters){
+    qs <- qsTable(get(paste0("results.",i)))
     qs$Cluster <- i
-    if (i==1){
-      results.cor <- qs
-    } else {
-      results.cor <- rbind(results.cor,qs)
-  }}
+    results.cor <- rbind(results.cor,qs)
+  }
+  results.cor <- results.cor[,-3]
   rownames(results.cor) <- NULL
-  results.clust.id <- results.cor[results.cor[,4]==1,][which.max(results.cor[results.cor[,4]==1,][,2]),]
+  results.clust.id <- NULL
+  if (max(results.cor[results.cor[,4]==1,][,2])>=0){
+    results.clust.id <- results.cor[results.cor[,4]==1,][which.max(results.cor[results.cor[,4]==1,][,2]),]
+  } else {
+    results.clust.id$pathway.name <- "Unknown"
+    results.clust.id$log.fold.change <- 0
+    results.clust.id$FDR <- 0
+    results.clust.id$Cluster <- 1
+    results.clust.id <- as.data.frame(results.clust.id)
+  }
   for (i in 2:number.clusters){
-    results.clust.id <- rbind(results.clust.id,results.cor[results.cor[,4]==i,][which.max(results.cor[results.cor[,4]==i,][,2]),])
+    if (max(results.cor[results.cor[,4]==i,][,2])>=0){
+      results.clust.id <- rbind(results.clust.id,results.cor[results.cor[,4]==i,][which.max(results.cor[results.cor[,4]==i,][,2]),])
+    } else {
+      results.clust.id <- rbind(results.clust.id,c("Unknown",0,0,i))
+    }
   }
   rownames(results.clust.id) <- NULL
   
@@ -559,8 +565,11 @@ scQuSAGE <- function(sc10x,gs,res.use=0.1,ds=25000,nm="Lin",folder="lin"){
   #Save ident
   merge.cluster <- NULL
   for (i in 1:number.clusters){
-    merge.cluster <- c(merge.cluster,as.character(results.clust.id[results.clust.id$Cluster==i,1]))
-  }
+    if (max(qsTable(get(paste0("results.",i)))[,2])>=0){
+      merge.cluster <- c(merge.cluster,as.character(qsTable((get(paste0("results.",i))))[qsTable((get(paste0("results.",i))))[,2]==max(qsTable((get(paste0("results.",i))))[,2]) & qsTable((get(paste0("results.",i))))[,2]>0,1]))
+    } else {
+      merge.cluster <- c(merge.cluster,"Unknown")
+    }}
   sc10x@ident <- plyr::mapvalues(x=sc10x@ident,from=1:number.clusters,to=merge.cluster)
   sc10x@ident <- factor(sc10x@ident,levels=names(gs))
   sc10x <- StashIdent(object=sc10x,save.name=nm)
@@ -828,10 +837,11 @@ scCCA <-  function(sc10x.1,sc10x.2,nm.1="D17",nm.2="D27",cc=FALSE){
   plot2 <- VlnPlot(sc10x,features.plot="ACC2",group.by="patient",do.return=TRUE,point.size.use=0.05)
   plot_grid(plot1,plot2)
   dev.off()
-
+  
   results <- list(
     sc10x=sc10x,
     genes.hvg.cca=genes.hvg.Comb
   )
   return(results)
 }
+