diff GO-enrich.R @ 10:d951677a50d4 draft

planemo upload commit 4ba1ebe7b3f5e3fabf78b5fed7ed0b92e2cbf9e5-dirty
author proteore
date Fri, 28 Jun 2019 05:08:48 -0400
parents 2f67202ffdb3
children f6107b8ae8f8
line wrap: on
line diff
--- a/GO-enrich.R	Wed Feb 27 03:39:16 2019 -0500
+++ b/GO-enrich.R	Fri Jun 28 05:08:48 2019 -0400
@@ -44,7 +44,7 @@
   return (width)
 }
 
-repartition.GO <- function(geneid, orgdb, ontology, level=3, readable=TRUE) {
+repartition_GO <- function(geneid, orgdb, ontology, level=3, readable=TRUE) {
   ggo<-groupGO(gene=geneid, 
                OrgDb = orgdb, 
                ont=ontology, 
@@ -66,7 +66,7 @@
 }
 
 # GO over-representation test
-enrich.GO <- function(geneid, universe, orgdb, ontology, pval_cutoff, qval_cutoff,plot) {
+enrich_GO <- function(geneid, universe, orgdb, ontology, pval_cutoff, qval_cutoff,plot) {
   ego<-enrichGO(gene=geneid,
                 universe=universe,
                 OrgDb=orgdb,
@@ -107,6 +107,15 @@
   }
 }
 
+clean_ids <- function(ids){
+  ids = gsub(" ","",ids)
+  ids = ids[which(ids!="")]
+  ids = ids[which(ids!="NA")]
+  ids = ids[!is.na(ids)]
+ 
+  return(ids) 
+}
+
 check_ids <- function(vector,type) {
   uniprot_pattern = "^([OPQ][0-9][A-Z0-9]{3}[0-9]|[A-NR-Z][0-9]([A-Z][A-Z0-9]{2}[0-9]){1,2})$"
   entrez_id = "^([0-9]+|[A-Z]{1,2}_[0-9]+|[A-Z]{1,2}_[A-Z]{1,4}[0-9]+)$"
@@ -117,7 +126,7 @@
   }
 }
 
-clusterProfiler = function() {
+get_args <- function(){
   args <- commandArgs(TRUE)
   if(length(args)<1) {
     args <- c("--help")
@@ -153,10 +162,18 @@
   args <- as.list(as.character(argsDF$V2))
   names(args) <- argsDF$V1
   
+  return(args)
+}
+
+
+main <- function() {
+  
+  #get args from command
+  args <- get_args()
+
   #save(args,file="/home/dchristiany/proteore_project/ProteoRE/tools/cluster_profiler/args.Rda")
   #load("/home/dchristiany/proteore_project/ProteoRE/tools/cluster_profiler/args.Rda")
   
-  
   go_represent=str2bool(args$go_represent)
   go_enrich=str2bool(args$go_enrich)
   if (go_enrich){
@@ -179,7 +196,7 @@
   id_type = args$id_type
   
   if (input_type == "text") {
-    input = strsplit(args$input, "[ \t\n]+")[[1]]
+    input = unlist(strsplit(strsplit(args$input, "[ \t\n]+")[[1]],";"))
   } else if (input_type == "file") {
     filename = args$input
     ncol = args$ncol
@@ -193,7 +210,7 @@
     file = read_file(filename, header)              # Extract Protein IDs list
     input =  unlist(sapply(as.character(file[,ncol]),function(x) rapply(strsplit(x,";"),c),USE.NAMES = FALSE))
   }
-  
+  input = clean_ids(input)
   
   ## Get input gene list from input IDs
   #ID format Conversion 
@@ -222,7 +239,7 @@
     if (!is.null(args$universe_type)) {
       universe_type = args$universe_type
       if (universe_type == "text") {
-        universe = strsplit(args$universe, "[ \t\n]+")[[1]]
+        universe = unlist(strsplit(strsplit(args$input, "[ \t\n]+")[[1]],";"))
       } else if (universe_type == "file") {
         universe_filename = args$universe
         universe_ncol = args$uncol
@@ -238,6 +255,7 @@
         # Extract Protein IDs list
         universe <- unlist(sapply(universe_file[,universe_ncol], function(x) rapply(strsplit(x,";"),c),USE.NAMES = FALSE))
       }
+      universe = clean_ids(input)
       universe_id_type = args$universe_id_type
       ##to initialize
       if (universe_id_type=="Uniprot" & any(check_ids(universe,"uniprot"))) {
@@ -265,14 +283,14 @@
   ##enrichGO : GO over-representation test
   for (onto in ontology) {
     if (go_represent) {
-      ggo<-repartition.GO(gene, orgdb, onto, level, readable=TRUE)
+      ggo<-repartition_GO(gene, orgdb, onto, level, readable=TRUE)
       if (is.list(ggo)){ggo <- as.data.frame(apply(ggo, c(1,2), function(x) gsub("^$|^ $", NA, x)))}  #convert "" and " " to NA
       output_path = paste("cluster_profiler_GGO_",onto,".tsv",sep="")
       write.table(ggo, output_path, sep="\t", row.names = FALSE, quote = FALSE )
     }
 
     if (go_enrich) {
-      ego<-enrich.GO(gene, universe_gene, orgdb, onto, pval_cutoff, qval_cutoff,plot)
+      ego<-enrich_GO(gene, universe_gene, orgdb, onto, pval_cutoff, qval_cutoff,plot)
       if (is.list(ego)){ego <- as.data.frame(apply(ego, c(1,2), function(x) gsub("^$|^ $", NA, x)))}  #convert "" and " " to NA
       output_path = paste("cluster_profiler_EGO_",onto,".tsv",sep="")
       write.table(ego, output_path, sep="\t", row.names = FALSE, quote = FALSE )
@@ -280,4 +298,6 @@
   }
 }
 
-clusterProfiler()
+if(!interactive()) {
+  main()
+}