--- a/enrichment_v3.R	Fri Sep 21 05:35:46 2018 -0400
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,380 +0,0 @@
-# enrichment_v3.R
-# Usage : Rscript --vanilla enrichment_v3.R --inputtype tabfile (or
-# copypaste) --input file.txt --ontology "BP/CC/MF" --option option (e.g
-# : classic/elim...) --threshold threshold --correction correction --textoutput
-# text --barplotoutput barplot
-# --dotplotoutput dotplot --column column --geneuniver human 
-# e.g : Rscript --vanilla enrichment_v3.R --inputtype tabfile --input file.txt
-# --ontology BP --option classic --threshold 1e-15 --correction holm
-# --textoutput TRUE
-# --barplotoutput TRUE --dotplotoutput TRUE --column c1 --geneuniverse
-# org.Hs.eg.db
-# INPUT :
-# - type of input. Can be ids separated by a blank space (copypast), or a text
-# file (tabfile)
-#	- file with at least one column of ensembl ids 
-#	- gene ontology category : Biological Process (BP), Cellular Component (CC), Molecular Function (MF)
-#	- test option (relative to topGO algorithms) : elim, weight01, parentchild, or no option (classic)
-#	- threshold for enriched GO term pvalues (e.g : 1e-15)     
-#	- correction for multiple testing (see p.adjust options : holm, hochberg, hommel, bonferroni, BH, BY,fdr,none
-#	- outputs wanted in this order text, barplot, dotplot with boolean value (e.g
-#	: TRUE TRUE TRUE ).
-#	Declare the output not wanted as none
-#	- column containing the ensembl ids if the input file is a tabfile
-# - gene universe reference for the user chosen specie
-# - header : if the input is a text file, does this text file have a header
-# (TRUE/FALSE)
-#
-# OUTPUT :
-#	- outputs commanded by the user named respectively result.tsv for the text
-#	results file, barplot.png for the barplot image file and dotplot.png for the
-#	dotplot image file 
-
-options(warn=-1)  #TURN OFF WARNINGS !!!!!!
-
-# loading topGO library
-suppressMessages(library(topGO))
-
-# Read file and return file content as data.frame
-readfile = function(filename, header) {
-  if (header == "true") {
-    # Read only first line of the file as header:
-    headers <- read.table(filename, nrows = 1, header = FALSE, sep = "\t", stringsAsFactors = FALSE, fill = TRUE, na.strings=c("", "NA"), blank.lines.skip = TRUE, quote = "")
-    #Read the data of the files (skipping the first row)
-    file <- read.table(filename, skip = 1, header = FALSE, sep = "\t", stringsAsFactors = FALSE, fill = TRUE, na.strings=c("", "NA"), blank.lines.skip = TRUE, quote = "")
-    # Remove empty rows
-    file <- file[!apply(is.na(file) | file == "", 1, all), , drop=FALSE]
-    #And assign the header to the data
-    names(file) <- headers
-  }
-  else {
-    file <- read.table(filename, header = FALSE, sep = "\t", stringsAsFactors = FALSE, fill = TRUE, na.strings=c("", "NA"), blank.lines.skip = TRUE, quote = "")
-    # Remove empty rows
-    file <- file[!apply(is.na(file) | file == "", 1, all), , drop=FALSE]
-  }
-  return(file)
-} 
-
-check_ens_ids <- function(vector) {
-  ens_pattern = "^(ENS[A-Z]+[0-9]{11}|[A-Z]{3}[0-9]{3}[A-Za-z](-[A-Za-z])?|CG[0-9]+|[A-Z0-9]+\\.[0-9]+|YM[A-Z][0-9]{3}[a-z][0-9])$"
-  return(grepl(ens_pattern,vector))
-}
-
-'%!in%' <- function(x,y)!('%in%'(x,y))
-
-
-# Parse command line arguments
-
-args = commandArgs(trailingOnly = TRUE)
-
-# create a list of the arguments from the command line, separated by a blank space
-hh <- paste(unlist(args),collapse=' ')
-
-# delete the first element of the list which is always a blank space
-listoptions <- unlist(strsplit(hh,'--'))[-1]
-
-# for each input, split the arguments with blank space as separator, unlist,
-# and delete the first element which is the input name (e.g --inputtype) 
-options.args <- sapply(listoptions,function(x){
-         unlist(strsplit(x, '[ \t\n]+'))[-1]
-        })
-# same as the step above, except that only the names are kept
-options.names <- sapply(listoptions,function(x){
-  option <-  unlist(strsplit(x, '[ \t\n]+'))[1]
-})
-names(options.args) <- unlist(options.names)
-
-
-if (length(options.args) != 12) {
-    stop("Not enough/Too many arguments", call. = FALSE)
-}
-
-#save(options.args,file="/home/dchristiany/proteore_project/ProteoRE/tools/topGO/args.Rda")
-#load("/home/dchristiany/proteore_project/ProteoRE/tools/topGO/args.Rda")
-
-
-typeinput = options.args[1]
-listfile = options.args[2]
-onto = as.character(options.args[3])
-option = as.character(options.args[4])
-correction = as.character(options.args[6])
-threshold = as.numeric(options.args[5])
-text = as.character(options.args[7])
-barplot = as.character(options.args[8])
-dotplot = as.character(options.args[9])
-column = as.numeric(gsub("c","",options.args[10]))
-geneuniverse = as.character(options.args[11])
-header = as.character(options.args[12])
-
-if (typeinput=="copypaste"){
-  sample = as.data.frame(unlist(listfile))
-  sample = sample[,column]
-}
-if (typeinput=="tabfile"){
-  
-  if (header=="TRUE"){
-    sample = readfile(listfile, "true")  
-  }else{
-    sample = readfile(listfile, "false")
-  }
-  sample = sample[,column]
-}
-
-#check of ENS ids
-if (! any(check_ens_ids(sample))){
-  print("no ensembl gene ids found in your ids list, please check your IDs in input or the selected column of your input file")
-  stop()
-}
-
-# Launch enrichment analysis and return result data from the analysis or the null
-# object if the enrichment could not be done.
-goEnrichment = function(geneuniverse,sample,onto){
-
-	# get all the GO terms of the corresponding ontology (BP/CC/MF) and all their
-  # associated ensembl ids according to the org package
-	xx = annFUN.org(onto,mapping=geneuniverse,ID="ensembl")
-	allGenes = unique(unlist(xx))
-	# check if the genes given by the user can be found in the org package (gene
-  # universe), that is in
-	# allGenes 
-	if (length(intersect(sample,allGenes))==0){
-	
-	    print("None of the input ids can be found in the org package data, enrichment analysis cannot be realized. \n The inputs ids probably have no associated GO terms.")
-      return(c(NULL,NULL))
-	
-	}
-	
-	geneList = factor(as.integer(allGenes %in% sample)) 
-	names(geneList) <- allGenes
-	
-	
-	#topGO enrichment 
-	
-	
-	# Creation of a topGOdata object
-	# It will contain : the list of genes of interest, the GO annotations and the GO hierarchy
-	# Parameters : 
-	# ontology : character string specifying the ontology of interest (BP, CC, MF)
-	# allGenes : named vector of type numeric or factor 
-	# annot : tells topGO how to map genes to GO annotations.
-	# argument not used here : nodeSize : at which minimal number of GO annotations
-	# do we consider a gene  
-	 
-	myGOdata = new("topGOdata", description="SEA with TopGO", ontology=onto, allGenes=geneList,  annot = annFUN.org, mapping=geneuniverse,ID="ensembl")
-	
-	
-	# Performing enrichment tests
-	result <- runTest(myGOdata, algorithm=option, statistic="fisher")
-  return(c(result,myGOdata))	
-}
-
-# Some libraries such as GOsummaries won't be able to treat the values such as
-# "< 1e-30" produced by topGO. As such it is important to delete the < char
-# with the deleteInfChar function. Nevertheless the user will have access to the original results in the text output.
-deleteInfChar = function(values){
-
-	lines = grep("<",values)
-	if (length(lines)!=0){
-		for (line in lines){
-		  values[line]=gsub("<","",values[line])
-		}
-	}
-	return(values)
-}
-
-corrMultipleTesting = function(result, myGOdata,correction,threshold){
-	
-	# adjust for multiple testing
-	if (correction!="none"){	
-	  # GenTable : transforms the result object into a list. Filters can be applied
-	  # (e.g : with the topNodes argument, to get for instance only the n first
-	  # GO terms with the lowest pvalues), but as we want to  apply a correction we
-	  # take all the GO terms, no matter their pvalues 
-	  allRes <- GenTable(myGOdata, test = result, orderBy = "result", ranksOf = "result",topNodes=length(attributes(result)$score))
-    # Some pvalues given by topGO are not numeric (e.g : "<1e-30). As such, these
-	  # values are converted to 1e-30 to be able to correct the pvalues 
-    pvaluestmp = deleteInfChar(allRes$test)
-	
-	  # the correction is done from the modified pvalues  
-	  allRes$qvalues = p.adjust(pvaluestmp, method = as.character(correction), n = length(pvaluestmp))
-	  allRes = as.data.frame(allRes)
-
-	  # Rename the test column by pvalues, so that is more explicit
-	  nb = which(names(allRes) %in% c("test"))
-	 
-    names(allRes)[nb] = "pvalues"
-	
-	  allRes = allRes[which(as.numeric(allRes$pvalues) <= threshold),]
-	  if (length(allRes$pvalues)==0){
-		  print("Threshold was too stringent, no GO term found with pvalue equal or lesser than the threshold value")
-      return(NULL)
-	  }
-	  allRes = allRes[order(allRes$qvalues),]
-	}
-	
-	if (correction=="none"){
-	  # get all the go terms under user threshold 
-	  mysummary <- summary(attributes(result)$score <= threshold)
-	  numsignif <- as.integer(mysummary[[3]])
-	  # get all significant nodes 
-	  allRes <- GenTable(myGOdata, test = result, orderBy = "result", ranksOf = "result",topNodes=numsignif)
-	
-
-	  allRes = as.data.frame(allRes)
-	  # Rename the test column by pvalues, so that is more explicit
-	  nb = which(names(allRes) %in% c("test"))
-	  names(allRes)[nb] = "pvalues"
-	  if (numsignif==0){
-	
-		  print("Threshold was too stringent, no GO term found with pvalue equal or lesser than the threshold value")
-      return(NULL)
-	  }
-	
-	 allRes = allRes[order(allRes$pvalues),] 
-	} 
-
-  return(allRes)  
-}
-
-# roundValues will simplify the results by rounding down the values. For instance 1.1e-17 becomes 1e-17
-roundValues = function(values){
-	for (line in 1:length(values)){
-    		values[line]=as.numeric(gsub(".*e","1e",as.character(values[line])))
-  }
-  return(values)
-}
-
-createDotPlot = function(data, onto){
-  
-	values  = deleteInfChar(data$pvalues)
-  values = roundValues(values)
-  values = as.numeric(values)
-	
-  geneRatio = data$Significant/data$Annotated
-	goTerms = data$Term
-	count = data$Significant
-  
-	labely = paste("GO terms",onto,sep=" ")
-	ggplot(data,aes(x=geneRatio,y=goTerms, color=values,size=count)) +geom_point( ) + scale_colour_gradientn(colours=c("red","violet","blue")) + xlab("Gene Ratio") + ylab(labely) + labs(color="p-values\n" ) 
-  ggsave("dotplot.png", device = "png", dpi = 320, limitsize = TRUE, width = 15, height = 15, units="cm")
-}
-
-createBarPlot = function(data, onto){
-
-  
-	values  = deleteInfChar(data$pvalues)
-  values = roundValues(values)
-
-  values = as.numeric(values)
-  goTerms = data$Term
-	count = data$Significant
-	
-	labely = paste("GO terms",onto,sep=" ")
-  ggplot(data, aes(x=goTerms, y=count,fill=values,scale(scale = 0.5))) + ylab("Gene count") + xlab(labely) +geom_bar(stat="identity") + scale_fill_gradientn(colours=c("red","violet","blue")) + coord_flip() + labs(fill="p-values\n") 
-  ggsave("barplot.png", device = "png", dpi = 320, limitsize = TRUE, width = 15, height = 15, units="cm")
-}
-
-
-# Produce the different outputs
-createOutputs = function(result, cut_result,text, barplot, dotplot, onto){
-
-
-  if (is.null(result)){
-    
-	  if (text=="TRUE"){
-
-      err_msg = "None of the input ids can be found in the org package data, enrichment analysis cannot be realized. \n The inputs ids probably either have no associated GO terms or are not ENSG identifiers (e.g : ENSG00000012048)."
-      write.table(err_msg, file='result.csv', quote=FALSE, sep='\t', col.names = T, row.names = F)
-	  
-    }
-
-	  if (barplot=="TRUE"){
-
-	    png(filename="barplot.png")
-      plot.new()
-      #text(0,0,err_msg)
-	    dev.off()
-    }
-	
-	  if (dotplot=="TRUE"){
-	
-	    png(filename="dotplot.png")
-      plot.new()
-      #text(0,0,err_msg)
-	    dev.off()
-	
-	  }
-    return(TRUE)
-  }
-
-	
-  if (is.null(cut_result)){
-
-
-	  if (text=="TRUE"){
-
-		  err_msg = "Threshold was too stringent, no GO term found with pvalue equal or lesser than the threshold value."
-      write.table(err_msg, file='result.csv', quote=FALSE, sep='\t', col.names = T, row.names = F)
-	  
-    }
-
-	  if (barplot=="TRUE"){
-
-	    png(filename="barplot.png")
-      plot.new()
-      text(0,0,err_msg)
-	    dev.off()
-    }
-	
-	  if (dotplot=="TRUE"){
-	
-	    png(filename="dotplot.png")
-      plot.new()
-      text(0,0,err_msg)
-	    dev.off()
-	
-	  }
-    return(TRUE)
-
-
-
-  }
-
-	if (text=="TRUE"){
-		write.table(cut_result, file='result.csv', quote=FALSE, sep='\t', col.names = T, row.names = F)
-	}
-	
-	if (barplot=="TRUE"){
-	
-		createBarPlot(cut_result, onto)
-	}
-	
-	if (dotplot=="TRUE"){
-	
-		createDotPlot(cut_result, onto)
-	}
-}
-
-
-
-# Load R library ggplot2 to plot graphs
-suppressMessages(library(ggplot2))
-
-# Launch enrichment analysis
-allresult = suppressMessages(goEnrichment(geneuniverse,sample,onto))
-result = allresult[1][[1]]
-myGOdata = allresult[2][[1]]
-if (!is.null(result)){
-
-	# Adjust the result with a multiple testing correction or not and with the user
-	# p-value cutoff
-  cut_result = corrMultipleTesting(result,myGOdata, correction,threshold)
-}else{
-
-  cut_result=NULL
-
-}
-
-
-createOutputs(result, cut_result,text, barplot, dotplot, onto)
-

--- a/topGO.xml	Fri Sep 21 05:35:46 2018 -0400
+++ b/topGO.xml	Mon Oct 15 11:30:04 2018 -0400
@@ -1,7 +1,5 @@
-<tool id="topGO" name="topGO" version="2018.09.21">
-    <description>
-    Enrichment analysis for Gene Ontology
-    </description>
+<tool id="topGO" name="Enrichment analysis for Gene Ontology" version="2018.10.12">
+    <description>(Human, Mouse, Rat) (topGO)</description>
     <requirements>
         <requirement type="package" version="3.4.1">R</requirement>
         <requirement type="package" version="3.0.0">r-ggplot2</requirement>
@@ -20,48 +18,31 @@
         <exit_code range="1:" />
     </stdio>
     <command><![CDATA[
-    
-  #if $inputtype.filetype == "file_all": 
-  Rscript --vanilla $__tool_directory__/enrichment_v3.R 
-  --inputtype tabfile 
-  --input '$inputtype.genelist' 
-  --ontology '$ontocat' 
-  --option '$option' 
-  --threshold '$threshold' 
-  --correction '$correction' 
-  --textoutput '$condtext.textoutput' 
-  --barplotoutput '$condbar.barplotoutput' 
-  --dotplotoutput '$conddot.dotplotoutput' 
-  --column '$inputtype.column' 
-  --geneuniverse '$geneuniverse' 
-  --header '$inputtype.header'
-  #end if
-    
-    
-  #if $inputtype.filetype == "copy_paste": 
-  Rscript --vanilla $__tool_directory__/enrichment_v3.R 
-  --inputtype copypaste 
-  --input '$inputtype.genelist' 
-  --ontology '$ontocat' 
-  --option '$option' 
-  --threshold '$threshold' 
-  --correction '$correction' 
-  --textoutput '$condtext.textoutput' 
-  --barplotoutput '$condbar.barplotoutput' 
-  --dotplotoutput '$conddot.dotplotoutput' 
-  --column c1 
-  --geneuniverse '$geneuniverse' 
-  --header None
+       
+  Rscript --vanilla $__tool_directory__/topGO_enrichment.R 
+  --inputtype="$inputtype.filetype"
+  --input='$inputtype.genelist'
+
+  #if $inputtype.filetype == "file"
+    --column='$inputtype.column'
+    --header='$inputtype.header'
   #end if
 
-
+  --ontology='$ontocat'
+  --option='$option'
+  --threshold='$threshold'
+  --correction='$correction'
+  --textoutput='$textoutput'
+  --barplotoutput='$barplot'
+  --dotplotoutput='$dotplot'
+  --geneuniverse='$geneuniverse' 
 
     ]]></command>
 
     <inputs>
   <conditional name="inputtype">
     <param name="filetype" type="select" label="Select your type of input file" help="The identifiers must be Ensembl gene IDs (e.g : ENSG00000139618). If it is not the case, please use the ID Mapping tool.">
-      <option value="file_all" selected="true">Input file containing your identifiers</option>
+      <option value="file" selected="true">Input file containing your identifiers</option>
       <option value="copy_paste">Copy/paste your list of IDs</option> 
     </param>
     <when value="copy_paste">
@@ -76,85 +57,58 @@
         </sanitizer>
       </param>
     </when>
-    <when value="file_all">
+    <when value="file">
       <param name="genelist" type="data" format="txt,tabular" label="Choose an input file" help="This file must imperatively have 1 column filled with IDs consistent with the database that will be used. Please use the MappingIDs component if this is not the case."/>
       <param name="column" type="text" label="Please specify the column where your Ensembl IDs are (e.g : Enter 'c1' for column n°1..)" value="c1"/> 
-    
-      <param name="header" type="select" label="Does your file have a header?" multiple="false" optional="false"> 
- 	<option value="TRUE" selected="true">Yes</option>
-        <option value="FALSE" selected="false">No</option>
-      </param>
+      <param name="header" type="boolean" checked="true" truevalue="true" falsevalue="false" label="Does your file have a header?" /> 
     </when>
   </conditional>
-      <param name="geneuniverse" type="select" label="Select a specie">
-        <!--option value="org.At.tair.db" >Arabidopsis</option-->
-        <option value="org.Ce.eg.db" >Worm (C. elegans)</option>
-        <option value="org.Dm.eg.db" >Fly (D. melanogaster)</option>
-        <option value="org.Hs.eg.db" selected="true">Human (H. sapiens)</option>
-        <option value="org.Mm.eg.db" >Mouse (M. musculus)</option>
-        <option value="org.Sc.sgd.db" >Yeast (S. cerevisiae)</option>
-      </param>
-
-      <param name="ontocat" type="select" label="Ontology category">
-              <option value="BP" >Biological Process</option>
-              <option value="CC" >Cellular Component</option>
-              <option value="MF" >Molecular Function</option>
-        </param>
-
-        <param name="option" type="select" label="Choose the topGO option for your analysis">
-                <option value="classic" >Classic fisher test</option>
-                <option value="elim" selected="true">Elim</option>
-                <option value="weight01" >Weight01</option>
-                <option value="parentchild" >ParentChild</option>
-          </param>
-        <param name="threshold" type="text" label="Enter the p-value threshold level under the form 1e-level wanted (e.g : 1e-3)" value="1e-3"/>
-	<param name="correction" label="Choose a correction for multiple testing" type="select">
-             <option value="none" >None</option>
-             <option value="holm">Holm correction</option>
-             <option value="hochberg" >Hochberg correction</option>
-             <option value="hommel" >Hommel correction</option>
-             <option value="bonferroni" >Bonferroni correction</option>
-             <option value="BH" selected="true">Benjamini and Hochberg</option>
-             <option value="BY" >Benjamini and Yekutieli</option>
-             <option value="fdr" >FDR</option>
-           </param>
-            <conditional name="condtext">
-            <param name="textoutput" type="select" label="Generate a text file for results">
-              <option value="TRUE">Yes</option>
-              <option value="FALSE">No</option>
-            </param>
-            <when value="TRUE"/>
-            <when value="FALSE"/>
-          </conditional>
-          <conditional name="condbar">
-            <param name="barplotoutput" type="select" label="Generate a barplot of over-represented GO terms">
-              <option value="TRUE">Yes</option>
-              <option value="FALSE">No</option>
-            </param>
-            <when value="TRUE"/>
-            <when value="FALSE"/>
-          </conditional>
-          <conditional name="conddot">
-              <param name="dotplotoutput" type="select" label="Generate a dotplot of over-represented GO terms">
-              <option value="TRUE">Yes</option>
-              <option value="FALSE">No</option>
-            </param>
-            <when value="TRUE"/>
-            <when value="FALSE"/>
-          </conditional>
+    <param name="geneuniverse" type="select" label="Select a specie">
+      <!--option value="org.At.tair.db" >Arabidopsis</option-->
+      <option value="org.Ce.eg.db" >Worm (C. elegans)</option>
+      <option value="org.Dm.eg.db" >Fly (D. melanogaster)</option>
+      <option value="org.Hs.eg.db" selected="true">Human (H. sapiens)</option>
+      <option value="org.Mm.eg.db" >Mouse (M. musculus)</option>
+      <option value="org.Sc.sgd.db" >Yeast (S. cerevisiae)</option>
+    </param>
+    <param name="ontocat" type="select" label="Ontology category">
+      <option value="BP" >Biological Process</option>
+      <option value="CC" >Cellular Component</option>
+      <option value="MF" >Molecular Function</option>
+    </param>
+    <param name="option" type="select" label="Choose the topGO option for your analysis">
+      <option value="classic" >Classic fisher test</option>
+      <option value="elim" selected="true">Elim</option>
+      <option value="weight01" >Weight01</option>
+      <option value="parentchild" >ParentChild</option>
+    </param>
+    <param name="threshold" type="text" label="Enter the p-value threshold level under the form 1e-level wanted (e.g : 1e-3)" value="1e-3"/>
+	  <param name="correction" label="Choose a correction for multiple testing" type="select">
+      <option value="none" >None</option>
+      <option value="holm">Holm correction</option>
+      <option value="hochberg" >Hochberg correction</option>
+      <option value="hommel" >Hommel correction</option>
+      <option value="bonferroni" >Bonferroni correction</option>
+      <option value="BH" selected="true">Benjamini and Hochberg</option>
+      <option value="BY" >Benjamini and Yekutieli</option>
+      <option value="fdr" >FDR</option>
+    </param>
+    <param name="textoutput" type="boolean" checked="true" truevalue="true" falsevalue="false" label="Generate a text file for results" />
+    <param name="barplot" type="boolean" checked="true" truevalue="true" falsevalue="false" label="Generate a barplot of over-represented GO terms" />
+    <param name="dotplot" type="boolean" checked="true" truevalue="true" falsevalue="false" label="Generate a dotplot of over-represented GO terms" />
     </inputs>
     <outputs>
 
-      <data name="outputtext" format="tabular" label="Text output for topGO analysis $ontocat category" from_work_dir="result.csv">
-        <filter>condtext['textoutput']=="TRUE"</filter>
+      <data name="outputtext" format="tsv" label="Text output for topGO analysis $ontocat category" from_work_dir="result">
+        <filter>textoutput</filter>
       </data>
 
       <data name="outputbarplot" format="png" label="Barplot output for topGO analysis $ontocat category" from_work_dir="barplot.png">
-        <filter>condbar['barplotoutput']=="TRUE"</filter>
+        <filter>barplot</filter>
       </data>
       
       <data name="outputdotplot" format="png" label="Dotplot output for topGO analysis $ontocat category" from_work_dir="dotplot.png">
-        <filter>conddot['dotplotoutput']=="TRUE"</filter>
+        <filter>dotplot</filter>
       </data>
       
    </outputs>

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/topGO_enrichment.R	Mon Oct 15 11:30:04 2018 -0400
@@ -0,0 +1,320 @@
+options(warn=-1)  #TURN OFF WARNINGS !!!!!!
+
+suppressMessages(library(ggplot2))
+suppressMessages(library(topGO))
+
+get_args <- function(){
+  
+  ## Collect arguments
+  args <- commandArgs(TRUE)
+  
+  ## Default setting when no arguments passed
+  if(length(args) < 1) {
+    args <- c("--help")
+  }
+  
+  ## Help section
+  if("--help" %in% args) {
+    cat("Pathview R script
+    Arguments:
+      --help                  Print this test
+      --input_type             
+      --onto
+      --option
+      --correction
+      --threshold
+      --text
+      --barplot
+      --dotplot
+      --column
+      --geneuniverse
+      --header
+
+      Example:
+      Rscript --vanilla enrichment_v3.R --inputtype=tabfile (or copypaste) --input=file.txt --ontology='BP/CC/MF' --option=option
+      (e.g : classic/elim...) --threshold=threshold --correction=correction --textoutput=text --barplotoutput=barplot --dotplotoutput=dotplot 
+      --column=column --geneuniver=human \n\n")
+    
+    q(save="no")
+  }
+  
+  parseArgs <- function(x) strsplit(sub("^--", "", x), "=")
+  argsDF <- as.data.frame(do.call("rbind", parseArgs(args)))
+  args <- as.list(as.character(argsDF$V2))
+  names(args) <- argsDF$V1
+  
+  return(args)
+}
+
+read_file <- function(path,header){
+  file <- try(read.table(path,header=header, sep="\t",stringsAsFactors = FALSE, quote=""),silent=TRUE)
+  if (inherits(file,"try-error")){
+    stop("File not found !")
+  }else{
+    return(file)
+  }
+}
+
+check_ens_ids <- function(vector) {
+  ens_pattern = "^(ENS[A-Z]+[0-9]{11}|[A-Z]{3}[0-9]{3}[A-Za-z](-[A-Za-z])?|CG[0-9]+|[A-Z0-9]+\\.[0-9]+|YM[A-Z][0-9]{3}[a-z][0-9])$"
+  return(grepl(ens_pattern,vector))
+}
+
+str2bool <- function(x){
+  if (any(is.element(c("t","true"),tolower(x)))){
+    return (TRUE)
+  }else if (any(is.element(c("f","false"),tolower(x)))){
+    return (FALSE)
+  }else{
+    return(NULL)
+  }
+}
+
+# Some libraries such as GOsummaries won't be able to treat the values such as
+# "< 1e-30" produced by topGO. As such it is important to delete the < char
+# with the deleteInfChar function. Nevertheless the user will have access to the original results in the text output.
+deleteInfChar = function(values){
+  
+  lines = grep("<",values)
+  if (length(lines)!=0){
+    for (line in lines){
+      values[line]=gsub("<","",values[line])
+    }
+  }
+  return(values)
+}
+
+corrMultipleTesting = function(result, myGOdata,correction,threshold){
+  
+  # adjust for multiple testing
+  if (correction!="none"){	
+    # GenTable : transforms the result object into a list. Filters can be applied
+    # (e.g : with the topNodes argument, to get for instance only the n first
+    # GO terms with the lowest pvalues), but as we want to  apply a correction we
+    # take all the GO terms, no matter their pvalues 
+    allRes <- GenTable(myGOdata, test = result, orderBy = "result", ranksOf = "result",topNodes=length(attributes(result)$score))
+    # Some pvalues given by topGO are not numeric (e.g : "<1e-30). As such, these
+    # values are converted to 1e-30 to be able to correct the pvalues 
+    pvaluestmp = deleteInfChar(allRes$test)
+    
+    # the correction is done from the modified pvalues  
+    allRes$qvalues = p.adjust(pvaluestmp, method = as.character(correction), n = length(pvaluestmp))
+    allRes = as.data.frame(allRes)
+    
+    # Rename the test column by pvalues, so that is more explicit
+    nb = which(names(allRes) %in% c("test"))
+    
+    names(allRes)[nb] = "pvalues"
+    
+    allRes = allRes[which(as.numeric(allRes$pvalues) <= threshold),]
+    if (length(allRes$pvalues)==0){
+      print("Threshold was too stringent, no GO term found with pvalue equal or lesser than the threshold value")
+      return(NULL)
+    }
+    allRes = allRes[order(allRes$qvalues),]
+  }
+  
+  if (correction=="none"){
+    # get all the go terms under user threshold 
+    mysummary <- summary(attributes(result)$score <= threshold)
+    numsignif <- as.integer(mysummary[[3]])
+    # get all significant nodes 
+    allRes <- GenTable(myGOdata, test = result, orderBy = "result", ranksOf = "result",topNodes=numsignif)
+    
+    
+    allRes = as.data.frame(allRes)
+    # Rename the test column by pvalues, so that is more explicit
+    nb = which(names(allRes) %in% c("test"))
+    names(allRes)[nb] = "pvalues"
+    if (numsignif==0){
+      
+      print("Threshold was too stringent, no GO term found with pvalue equal or lesser than the threshold value")
+      return(NULL)
+    }
+    
+    allRes = allRes[order(allRes$pvalues),] 
+  } 
+  
+  return(allRes)  
+}
+
+# roundValues will simplify the results by rounding down the values. For instance 1.1e-17 becomes 1e-17
+roundValues = function(values){
+  for (line in 1:length(values)){
+    values[line]=as.numeric(gsub(".*e","1e",as.character(values[line])))
+  }
+  return(values)
+}
+
+createDotPlot = function(data, onto){
+  
+  values  = deleteInfChar(data$pvalues)
+  values = roundValues(values)
+  values = as.numeric(values)
+  
+  geneRatio = data$Significant/data$Annotated
+  goTerms = data$Term
+  count = data$Significant
+  
+  labely = paste("GO terms",onto,sep=" ")
+  ggplot(data,aes(x=geneRatio,y=goTerms, color=values,size=count)) +geom_point( ) + scale_colour_gradientn(colours=c("red","violet","blue")) + xlab("Gene Ratio") + ylab(labely) + labs(color="p-values\n" ) 
+  ggsave("dotplot.png", device = "png", dpi = 320, limitsize = TRUE, width = 15, height = 15, units="cm")
+}
+
+createBarPlot = function(data, onto){
+  
+  
+  values  = deleteInfChar(data$pvalues)
+  values = roundValues(values)
+  values = as.numeric(values)
+  
+  goTerms = data$Term
+  count = data$Significant
+  
+  labely = paste("GO terms",onto,sep=" ")
+  ggplot(data, aes(x=goTerms, y=count,fill=values,scale(scale = 0.5))) + ylab("Gene count") + xlab(labely) +geom_bar(stat="identity") + scale_fill_gradientn(colours=c("red","violet","blue")) + coord_flip() + labs(fill="p-values\n") 
+  ggsave("barplot.png", device = "png", dpi = 320, limitsize = TRUE, width = 15, height = 15, units="cm")
+}
+
+# Produce the different outputs
+createOutputs = function(result, cut_result,text, barplot, dotplot, onto){
+  
+  if (is.null(result)){
+    if (text){
+      err_msg = "None of the input ids can be found in the org package data, enrichment analysis cannot be realized. \n The inputs ids probably either have no associated GO terms or are not ENSG identifiers (e.g : ENSG00000012048)."
+      write.table(err_msg, file='result', quote=FALSE, sep='\t', col.names = T, row.names = F)
+    }
+    if (barplot){
+      png(filename="barplot.png")
+      plot.new()
+      #text(0,0,err_msg)
+      dev.off()
+    }
+    if (dotplot){
+      png(filename="dotplot.png")
+      plot.new()
+      #text(0,0,err_msg)
+      dev.off()
+    }
+    opt <- options(show.error.messages=FALSE)
+    on.exit(options(opt))
+    stop("null result")
+  }
+  
+  if (is.null(cut_result)){
+    if (text){
+      err_msg = "Threshold was too stringent, no GO term found with pvalue equal or lesser than the threshold value."
+      write.table(err_msg, file='result', quote=FALSE, sep='\t', col.names = T, row.names = F)
+    }
+    if (barplot){
+      png(filename="barplot.png")
+      plot.new()
+      text(0,0,err_msg)
+      dev.off()
+    }
+    if (dotplot){
+      png(filename="dotplot.png")
+      plot.new()
+      text(0,0,err_msg)
+      dev.off()
+    }
+    opt <- options(show.error.messages=FALSE)
+    on.exit(options(opt))
+    stop("null cut_result")
+  }
+  
+  if (text){
+    write.table(cut_result, file='result', quote=FALSE, sep='\t', col.names = T, row.names = F)
+  }
+  
+  if (barplot){
+    createBarPlot(cut_result, onto)
+  }
+  
+  if (dotplot){
+    createDotPlot(cut_result, onto)
+  }
+}
+
+# Launch enrichment analysis and return result data from the analysis or the null
+# object if the enrichment could not be done.
+goEnrichment = function(geneuniverse,sample,onto){
+  
+  # get all the GO terms of the corresponding ontology (BP/CC/MF) and all their
+  # associated ensembl ids according to the org package
+  xx = annFUN.org(onto,mapping=geneuniverse,ID="ensembl")
+  allGenes = unique(unlist(xx))
+  # check if the genes given by the user can be found in the org package (gene
+  # universe), that is in
+  # allGenes 
+  if (length(intersect(sample,allGenes))==0){
+    
+    print("None of the input ids can be found in the org package data, enrichment analysis cannot be realized. \n The inputs ids probably have no associated GO terms.")
+    return(c(NULL,NULL))
+    
+  }
+  
+  geneList = factor(as.integer(allGenes %in% sample)) 
+  names(geneList) <- allGenes
+  
+  #topGO enrichment 
+  
+  # Creation of a topGOdata object
+  # It will contain : the list of genes of interest, the GO annotations and the GO hierarchy
+  # Parameters : 
+  # ontology : character string specifying the ontology of interest (BP, CC, MF)
+  # allGenes : named vector of type numeric or factor 
+  # annot : tells topGO how to map genes to GO annotations.
+  # argument not used here : nodeSize : at which minimal number of GO annotations
+  # do we consider a gene  
+  
+  myGOdata = new("topGOdata", description="SEA with TopGO", ontology=onto, allGenes=geneList,  annot = annFUN.org, mapping=geneuniverse,ID="ensembl")
+  
+  # Performing enrichment tests
+  result <- runTest(myGOdata, algorithm=option, statistic="fisher")
+  return(c(result,myGOdata))	
+}
+
+args <- get_args()
+
+#save(args,file="/home/dchristiany/proteore_project/ProteoRE/tools/topGO/args.rda")
+#load("/home/dchristiany/proteore_project/ProteoRE/tools/topGO/args.rda")
+
+input_type = args$inputtype
+input = args$input
+onto = args$ontology
+option = args$option
+correction = args$correction
+threshold = as.numeric(args$threshold)
+text = str2bool(args$textoutput)
+barplot = str2bool(args$barplotoutput)
+dotplot = str2bool(args$dotplotoutput)
+column = as.numeric(gsub("c","",args$column))
+geneuniverse = args$geneuniverse
+header = str2bool(args$header)
+
+#get input
+if (input_type=="copy_paste"){
+  sample <- unlist(strsplit(input,","))
+} else if (input_type=="file"){
+  tab=read_file(input,header)
+  sample = tab[,column]
+}
+
+#check of ENS ids
+if (! any(check_ens_ids(sample))){
+  print("no ensembl gene ids found in your ids list, please check your IDs in input or the selected column of your input file")
+  stop()
+}
+
+# Launch enrichment analysis
+allresult = suppressMessages(goEnrichment(geneuniverse,sample,onto))
+result = allresult[1][[1]]
+myGOdata = allresult[2][[1]]
+if (!is.null(result)){
+  cut_result = corrMultipleTesting(result,myGOdata, correction,threshold) # Adjust the result with a multiple testing correction or not and with the user, p-value cutoff
+}else{
+  cut_result=NULL
+}
+
+createOutputs(result, cut_result,text, barplot, dotplot, onto)