Mercurial > repos > proteore > proteore_topgo

annotate topGO_enrichment.R @ 15:159596ec807d draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
planemo upload commit b8671ffe2e12dc6612b971a3e6e1dc71496aefd0-dirty
author proteore
date Fri, 24 Jan 2020 10:45:15 -0500
parents 8eaa43ba1bfc
children 7f1ce70f0f09
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
10
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
1 options(warn=-1) #TURN OFF WARNINGS !!!!!!
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
2
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
3 suppressMessages(library(ggplot2))
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
4 suppressMessages(library(topGO))
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
5
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
6 get_args <- function(){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
7
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
8 ## Collect arguments
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
9 args <- commandArgs(TRUE)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
10
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
11 ## Default setting when no arguments passed
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
12 if(length(args) < 1) {
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
13 args <- c("--help")
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
14 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
15
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
16 ## Help section
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
17 if("--help" %in% args) {
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
18 cat("Pathview R script
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
19 Arguments:
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
20 --help Print this test
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
21 --input_type
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
22 --onto
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
23 --option
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
24 --correction
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
25 --threshold
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
26 --text
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
27 --plot
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
28 --column
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
29 --geneuniverse
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
30 --header
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
31
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
32 Example:
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
33 Rscript --vanilla enrichment_v3.R --inputtype=tabfile (or copypaste) --input=file.txt --ontology='BP/CC/MF' --option=option
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
34 (e.g : classic/elim...) --threshold=threshold --correction=correction --textoutput=text --barplotoutput=barplot --dotplotoutput=dotplot
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
35 --column=column --geneuniver=human \n\n")
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
36
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
37 q(save="no")
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
38 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
39
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
40 parseArgs <- function(x) strsplit(sub("^--", "", x), "=")
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
41 argsDF <- as.data.frame(do.call("rbind", parseArgs(args)))
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
42 args <- as.list(as.character(argsDF$V2))
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
43 names(args) <- argsDF$V1
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
44
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
45 return(args)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
46 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
47
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
48 read_file <- function(path,header){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
49 file <- try(read.csv(path,header=header, sep="\t",stringsAsFactors = FALSE, quote="\"", check.names = F),silent=TRUE)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
50 if (inherits(file,"try-error")){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
51 stop("File not found !")
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
52 }else{
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
53 return(file)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
54 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
55 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
56
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
57 get_list_from_cp <-function(list){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
58 list = gsub(";"," ",list)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
59 list = strsplit(list, "[ \t\n]+")[[1]]
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
60 list = list[list != ""] #remove empty entry
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
61 list = gsub("-.+", "", list) #Remove isoform accession number (e.g. "-2")
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
62 return(list)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
63 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
64
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
65 check_ens_ids <- function(vector) {
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
66 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])$"
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
67 return(grepl(ens_pattern,vector))
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
68 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
69
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
70 str2bool <- function(x){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
71 if (any(is.element(c("t","true"),tolower(x)))){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
72 return (TRUE)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
73 }else if (any(is.element(c("f","false"),tolower(x)))){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
74 return (FALSE)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
75 }else{
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
76 return(NULL)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
77 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
78 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
79
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
80 # Some libraries such as GOsummaries won't be able to treat the values such as
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
81 # "< 1e-30" produced by topGO. As such it is important to delete the < char
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
82 # with the deleteInfChar function. Nevertheless the user will have access to the original results in the text output.
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
83 deleteInfChar = function(values){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
84
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
85 lines = grep("<",values)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
86 if (length(lines)!=0){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
87 for (line in lines){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
88 values[line]=gsub("<","",values[line])
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
89 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
90 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
91 return(values)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
92 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
93
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
94 corrMultipleTesting = function(result, myGOdata,correction,threshold){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
95
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
96 # adjust for multiple testing
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
97 if (correction!="none"){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
98 # GenTable : transforms the result object into a list. Filters can be applied
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
99 # (e.g : with the topNodes argument, to get for instance only the n first
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
100 # GO terms with the lowest pvalues), but as we want to apply a correction we
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
101 # take all the GO terms, no matter their pvalues
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
102 allRes <- GenTable(myGOdata, test = result, orderBy = "result", ranksOf = "result",topNodes=length(attributes(result)$score))
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
103 # Some pvalues given by topGO are not numeric (e.g : "<1e-30). As such, these
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
104 # values are converted to 1e-30 to be able to correct the pvalues
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
105 pvaluestmp = deleteInfChar(allRes$test)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
106
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
107 # the correction is done from the modified pvalues
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
108 allRes$qvalues = p.adjust(pvaluestmp, method = as.character(correction), n = length(pvaluestmp))
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
109 allRes = as.data.frame(allRes)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
110
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
111 # Rename the test column by pvalues, so that is more explicit
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
112 nb = which(names(allRes) %in% c("test"))
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
113
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
114 names(allRes)[nb] = "pvalues"
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
115
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
116 allRes = allRes[which(as.numeric(allRes$pvalues) <= threshold),]
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
117 if (length(allRes$pvalues)==0){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
118 print("Threshold was too stringent, no GO term found with pvalue equal or lesser than the threshold value")
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
119 return(NULL)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
120 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
121 allRes = allRes[order(allRes$qvalues),]
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
122 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
123
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
124 if (correction=="none"){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
125 # get all the go terms under user threshold
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
126 mysummary <- summary(attributes(result)$score <= threshold)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
127 numsignif <- as.integer(mysummary[[3]])
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
128 # get all significant nodes
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
129 allRes <- GenTable(myGOdata, test = result, orderBy = "result", ranksOf = "result",topNodes=numsignif)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
130
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
131
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
132 allRes = as.data.frame(allRes)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
133 # Rename the test column by pvalues, so that is more explicit
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
134 nb = which(names(allRes) %in% c("test"))
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
135 names(allRes)[nb] = "pvalues"
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
136 if (numsignif==0){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
137
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
138 print("Threshold was too stringent, no GO term found with pvalue equal or lesser than the threshold value")
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
139 return(NULL)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
140 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
141
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
142 allRes = allRes[order(allRes$pvalues),]
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
143 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
144
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
145 return(allRes)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
146 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
147
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
148 # roundValues will simplify the results by rounding down the values. For instance 1.1e-17 becomes 1e-17
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
149 roundValues = function(values){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
150 for (line in 1:length(values)){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
151 values[line]=as.numeric(gsub(".*e","1e",as.character(values[line])))
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
152 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
153 return(values)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
154 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
155
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
156 createDotPlot = function(data, onto){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
157
11
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
158 values = deleteInfChar(data$pvalues)
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
159 values = roundValues(values)
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
160 values = as.numeric(values)
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
161
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
162 geneRatio = data$Significant/data$Annotated
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
163 goTerms = data$Term
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
164 count = data$Significant
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
165
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
166 labely = paste("GO terms",onto,sep=" ")
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
167 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" )
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
168 ggsave("dotplot.png", device = "png", dpi = 320, limitsize = TRUE, width = 15, height = 15, units="cm")
10
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
169 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
170
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
171 createBarPlot = function(data, onto){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
172
11
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
173 values = deleteInfChar(data$pvalues)
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
174 values = roundValues(values)
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
175 values = as.numeric(values)
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
176
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
177 goTerms = data$Term
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
178 count = data$Significant
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
179
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
180 labely = paste("GO terms",onto,sep=" ")
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
181 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")
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
182 ggsave("barplot.png", device = "png", dpi = 320, limitsize = TRUE, width = 15, height = 15, units="cm")
10
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
183 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
184
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
185 # Produce the different outputs
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
186 createOutputs = function(result, cut_result,text, barplot, dotplot, onto){
11
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
187
10
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
188 if (is.null(result)){
11
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
189 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)."
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
190 write.table(err_msg, file='result', quote=FALSE, sep='\t', col.names = F, row.names = F)
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
191 }else if (is.null(cut_result)){
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
192 err_msg = "Threshold was too stringent, no GO term found with pvalue equal or lesser than the threshold value."
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
193 write.table(err_msg, file='result.tsv', quote=FALSE, sep='\t', col.names = F, row.names = F)
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
194 }else {
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
195 write.table(cut_result, file='result.tsv', quote=FALSE, sep='\t', col.names = T, row.names = F)
10
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
196
11
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
197 if (barplot){createBarPlot(cut_result, onto)}
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
198 if (dotplot){createDotPlot(cut_result, onto)}
10
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
199 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
200 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
201
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
202 # Launch enrichment analysis and return result data from the analysis or the null
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
203 # object if the enrichment could not be done.
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
204 goEnrichment = function(geneuniverse,sample,background_sample,onto){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
205
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
206 if (is.null(background_sample)){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
207 xx = annFUN.org(onto,mapping=geneuniverse,ID="ensembl") # get all the GO terms of the corresponding ontology (BP/CC/MF) and all their associated ensembl ids according to the org package
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
208 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
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
209 } else {
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
210 allGenes = background_sample
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
211 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
212
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
213 if (length(intersect(sample,allGenes))==0){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
214 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.")
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
215 return(c(NULL,NULL))
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
216 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
217
12
8eaa43ba1bfc planemo upload commit 4ba1ebe7b3f5e3fabf78b5fed7ed0b92e2cbf9e5-dirty
proteore
parents: 11
diff changeset
218 geneList = factor(as.integer(allGenes %in% sample)) #duplicated ids in sample count only for one
10
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
219 if (length(levels(geneList)) == 1 ){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
220 stop("All or none of the background genes are found in tested genes dataset, enrichment analysis can't be done")
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
221 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
222 names(geneList) <- allGenes
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
223
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
224 #topGO enrichment
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
225
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
226 # Creation of a topGOdata object
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
227 # It will contain : the list of genes of interest, the GO annotations and the GO hierarchy
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
228 # Parameters :
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
229 # ontology : character string specifying the ontology of interest (BP, CC, MF)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
230 # allGenes : named vector of type numeric or factor
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
231 # annot : tells topGO how to map genes to GO annotations.
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
232 # argument not used here : nodeSize : at which minimal number of GO annotations
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
233 # do we consider a gene
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
234
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
235 myGOdata = new("topGOdata", description="SEA with TopGO", ontology=onto, allGenes=geneList, annot = annFUN.org, mapping=geneuniverse,ID="ensembl")
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
236
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
237 # Performing enrichment tests
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
238 result <- runTest(myGOdata, algorithm=option, statistic="fisher")
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
239 return(c(result,myGOdata))
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
240 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
241
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
242 args <- get_args()
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
243
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
244 #save(args,file="/home/dchristiany/proteore_project/ProteoRE/tools/topGO/args.rda")
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
245 #load("/home/dchristiany/proteore_project/ProteoRE/tools/topGO/args.rda")
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
246
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
247 input_type = args$inputtype
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
248 input = args$input
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
249 onto = args$ontology
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
250 option = args$option
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
251 correction = args$correction
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
252 threshold = as.numeric(args$threshold)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
253 text = str2bool(args$textoutput)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
254 barplot = "barplot" %in% unlist(strsplit(args$plot,","))
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
255 dotplot = "dotplot" %in% unlist(strsplit(args$plot,","))
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
256 column = as.numeric(gsub("c","",args$column))
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
257 geneuniverse = args$geneuniverse
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
258 header = str2bool(args$header)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
259 background = str2bool(args$background)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
260 if (background){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
261 background_genes = args$background_genes
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
262 background_input_type = args$background_input_type
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
263 background_header = str2bool(args$background_header)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
264 background_column = as.numeric(gsub("c","",args$background_column))
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
265 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
266
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
267 #get input
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
268 if (input_type=="copy_paste"){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
269 sample <- get_list_from_cp(input)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
270 } else if (input_type=="file"){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
271 tab=read_file(input,header)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
272 sample = trimws(unlist(strsplit(tab[,column],";")))
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
273 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
274
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
275 #check of ENS ids
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
276 if (! any(check_ens_ids(sample))){
11
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
277 stop("no ensembl gene ids found in your ids list, please check your IDs in input or the selected column of your input file")
10
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
278 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
279
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
280 #get input if background genes
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
281 if (background){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
282 if (background_input_type=="copy_paste"){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
283 background_sample <- get_list_from_cp(background_genes)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
284 } else if (background_input_type=="file"){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
285 background_tab=read_file(background_genes,background_header)
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
286 background_sample = unique(trimws(unlist(strsplit(background_tab[,background_column],";"))))
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
287 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
288 #check of ENS ids
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
289 if (! any(check_ens_ids(background_sample))){
11
fa2e27165d5d planemo upload commit 4efc56eb769fbceb66c64181441ff8781d523454-dirty
proteore
parents: 10
diff changeset
290 stop("no ensembl gene ids found in your background ids list, please check your IDs in input or the selected column of your input file")
10
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
291 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
292 } else {
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
293 background_sample=NULL
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
294 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
295
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
296 # Launch enrichment analysis
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
297 allresult = suppressMessages(goEnrichment(geneuniverse,sample,background_sample,onto))
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
298 result = allresult[1][[1]]
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
299 myGOdata = allresult[2][[1]]
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
300 if (!is.null(result)){
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
301 cut_result = corrMultipleTesting(result,myGOdata, correction,threshold) # Adjust the result with a multiple testing correction or not and with the user, p-value cutoff
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
302 }else{
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
303 cut_result=NULL
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
304 }
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
305
e3430084c996 planemo upload commit ad5f1c5a1a71d7fa2bc8bac408856aa80b0fc2a3
proteore
parents:
diff changeset
306 createOutputs(result, cut_result,text, barplot, dotplot, onto)