Modify QuSAGE function to not create cor plots and to handle multiple sets with the same name

r.scripts/sc-TissueMapper.R

@@ -633,6 +633,115 @@ scQuSAGE <- function(sc10x,gs,res.use=0.1,ds=25000,nm="Lin",folder="lin"){
 }
 
 
+scQuSAGE.nocorplot <- function(sc10x,gs,res.use=0.1,ds=25000,nm="Lin",folder="lin"){
+    #Runs QuSAGE
+    
+    #Inputs:
+    #sc10x = Seruat object
+    #gs = geneset to use for correlation
+    #res.use = calculated PC to use
+    #nm = name of test
+    #folder = folder for output
+    
+    #Outputs:
+    #results[1] = Seurat object
+    #results[2] = correlation table
+    #results[3] = correlation results
+    
+    #Downsample and store expression matrix for QuSAGE
+    sc10x <- SetAllIdent(object=sc10x,id=paste0("res",res.use))
+    number.clusters <- length(unique(sc10x@ident))
+    
+    labels <- paste0("Cluster_",as.vector(factor(sc10x@ident)))
+    
+    cell.sample <- NULL
+    for (i in 1:number.clusters){
+      cell <- names(sc10x@ident[sc10x@ident==i])
+      if (length(cell)>ds & ds!=0){
+        rnd <- sample(1:length(cell),ds)
+        cell <- cell[rnd]
+      }
+      cell.sample <- c(cell.sample,cell)
+    }
+    data <- sc10x@data[,colnames(sc10x@data) %in% cell.sample]
+    data <- as.data.frame(as.matrix(data))
+    labels <- labels[colnames(sc10x@data) %in% cell.sample]
+    
+    #Make labels for QuSAGE
+    clust <- list()
+    clust.comp <- list()
+    for (i in 1:number.clusters){
+      t <- labels
+      t[!(t %in% paste0("Cluster_",i))] <- "REST"
+      clust[i] <- list(i=t)
+      rm(t)
+      clust.comp[i] <- paste0("Cluster_",i,"-REST")
+    }
+    
+    #Run QuSAGE
+    for (i in 1:number.clusters){
+      assign(paste0("results.",i),qusage(data,unlist(clust[i]),unlist(clust.comp[i]),gs))
+    }
+    
+  #Generate ID table
+  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
+    results.cor <- rbind(results.cor,qs)
+  }
+  results.cor <- results.cor[,-3]
+  rownames(results.cor) <- NULL
+  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){
+    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,data.frame(pathway.name="Unknown",log.fold.change=0,FDR=0,Cluster=i))
+    }
+  }
+  rownames(results.clust.id) <- NULL
+  
+  #Save ident
+  merge.cluster <- NULL
+  for (i in 1:number.clusters){
+    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=c(unique(names(gs)),"Unknown"))
+  sc10x <- StashIdent(object=sc10x,save.name=nm)
+  
+  postscript(paste0("./analysis/tSNE/",folder,"/tSNE_",nm,".eps"))
+  plot <- TSNEPlot(object=sc10x,pt.size=2.5,do.return=TRUE,vector.friendly=FALSE)
+  plot <- plot+theme(axis.text.x=element_text(size=20),axis.text.y=element_text(size=20),axis.title.x=element_text(size=20),axis.title.y=element_text(size=20),legend.text=element_text(size=20))
+  plot <- plot+guides(colour=guide_legend(override.aes=list(size=10)))
+  plot(plot)
+  dev.off()
+  
+  results <- list(
+    sc10x=sc10x,
+    results.cor=results.cor,
+    results.clust.id=results.clust.id
+  )
+  names(results)=c("sc10x",paste0("results.cor.",nm),paste0("results.clust.",nm,".id"))
+  return(results)
+}
+
+
 scQuSAGEsm <- function(sc10x,gs,ds=25000,nm="Lin",folder="lin"){
   #Runs QuSAGE