proteore_topgo: enrichment

annotate enrichment_v3.R @ 9:70c0c8757f5f draft

planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty

author	proteore
date	Fri, 21 Sep 2018 05:32:38 -0400
parents	ddaa0c318d65
children

rev	line source
0 472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	1 # enrichment_v3.R
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	2 # Usage : Rscript --vanilla enrichment_v3.R --inputtype tabfile (or
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	3 # copypaste) --input file.txt --ontology "BP/CC/MF" --option option (e.g
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	4 # : classic/elim...) --threshold threshold --correction correction --textoutput
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	5 # text --barplotoutput barplot
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	6 # --dotplotoutput dotplot --column column --geneuniver human
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	7 # e.g : Rscript --vanilla enrichment_v3.R --inputtype tabfile --input file.txt
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	8 # --ontology BP --option classic --threshold 1e-15 --correction holm
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	9 # --textoutput TRUE
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	10 # --barplotoutput TRUE --dotplotoutput TRUE --column c1 --geneuniverse
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	11 # org.Hs.eg.db
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	12 # INPUT :
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	13 # - type of input. Can be ids separated by a blank space (copypast), or a text
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	14 # file (tabfile)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	15 # - file with at least one column of ensembl ids
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	16 # - gene ontology category : Biological Process (BP), Cellular Component (CC), Molecular Function (MF)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	17 # - test option (relative to topGO algorithms) : elim, weight01, parentchild, or no option (classic)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	18 # - threshold for enriched GO term pvalues (e.g : 1e-15)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	19 # - correction for multiple testing (see p.adjust options : holm, hochberg, hommel, bonferroni, BH, BY,fdr,none
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	20 # - outputs wanted in this order text, barplot, dotplot with boolean value (e.g
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	21 # : TRUE TRUE TRUE ).
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	22 # Declare the output not wanted as none
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	23 # - column containing the ensembl ids if the input file is a tabfile
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	24 # - gene universe reference for the user chosen specie
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	25 # - header : if the input is a text file, does this text file have a header
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	26 # (TRUE/FALSE)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	27 #
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	28 # OUTPUT :
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	29 # - outputs commanded by the user named respectively result.tsv for the text
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	30 # results file, barplot.png for the barplot image file and dotplot.png for the
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	31 # dotplot image file
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	32
9 70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	33 options(warn=-1) #TURN OFF WARNINGS !!!!!!
0 472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	34
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	35 # loading topGO library
9 70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	36 suppressMessages(library(topGO))
8 ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	37
ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	38 # Read file and return file content as data.frame
ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	39 readfile = function(filename, header) {
ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	40 if (header == "true") {
ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	41 # Read only first line of the file as header:
ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	42 headers <- read.table(filename, nrows = 1, header = FALSE, sep = "\t", stringsAsFactors = FALSE, fill = TRUE, na.strings=c("", "NA"), blank.lines.skip = TRUE, quote = "")
ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	43 #Read the data of the files (skipping the first row)
ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	44 file <- read.table(filename, skip = 1, header = FALSE, sep = "\t", stringsAsFactors = FALSE, fill = TRUE, na.strings=c("", "NA"), blank.lines.skip = TRUE, quote = "")
ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	45 # Remove empty rows
ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	46 file <- file[!apply(is.na(file) \| file == "", 1, all), , drop=FALSE]
ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	47 #And assign the header to the data
ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	48 names(file) <- headers
ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	49 }
ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	50 else {
ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	51 file <- read.table(filename, header = FALSE, sep = "\t", stringsAsFactors = FALSE, fill = TRUE, na.strings=c("", "NA"), blank.lines.skip = TRUE, quote = "")
ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	52 # Remove empty rows
ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	53 file <- file[!apply(is.na(file) \| file == "", 1, all), , drop=FALSE]
ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	54 }
ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	55 return(file)
ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	56 }
0 472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	57
9 70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	58 check_ens_ids <- function(vector) {
70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	59 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])$"
70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	60 return(grepl(ens_pattern,vector))
70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	61 }
70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	62
0 472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	63 '%!in%' <- function(x,y)!('%in%'(x,y))
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	64
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	65
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	66 # Parse command line arguments
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	67
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	68 args = commandArgs(trailingOnly = TRUE)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	69
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	70 # create a list of the arguments from the command line, separated by a blank space
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	71 hh <- paste(unlist(args),collapse=' ')
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	72
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	73 # delete the first element of the list which is always a blank space
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	74 listoptions <- unlist(strsplit(hh,'--'))[-1]
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	75
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	76 # for each input, split the arguments with blank space as separator, unlist,
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	77 # and delete the first element which is the input name (e.g --inputtype)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	78 options.args <- sapply(listoptions,function(x){
8 ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	79 unlist(strsplit(x, '[ \t\n]+'))[-1]
0 472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	80 })
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	81 # same as the step above, except that only the names are kept
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	82 options.names <- sapply(listoptions,function(x){
8 ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	83 option <- unlist(strsplit(x, '[ \t\n]+'))[1]
0 472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	84 })
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	85 names(options.args) <- unlist(options.names)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	86
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	87
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	88 if (length(options.args) != 12) {
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	89 stop("Not enough/Too many arguments", call. = FALSE)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	90 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	91
9 70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	92 #save(options.args,file="/home/dchristiany/proteore_project/ProteoRE/tools/topGO/args.Rda")
70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	93 #load("/home/dchristiany/proteore_project/ProteoRE/tools/topGO/args.Rda")
70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	94
70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	95
0 472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	96 typeinput = options.args[1]
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	97 listfile = options.args[2]
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	98 onto = as.character(options.args[3])
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	99 option = as.character(options.args[4])
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	100 correction = as.character(options.args[6])
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	101 threshold = as.numeric(options.args[5])
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	102 text = as.character(options.args[7])
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	103 barplot = as.character(options.args[8])
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	104 dotplot = as.character(options.args[9])
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	105 column = as.numeric(gsub("c","",options.args[10]))
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	106 geneuniverse = as.character(options.args[11])
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	107 header = as.character(options.args[12])
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	108
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	109 if (typeinput=="copypaste"){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	110 sample = as.data.frame(unlist(listfile))
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	111 sample = sample[,column]
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	112 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	113 if (typeinput=="tabfile"){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	114
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	115 if (header=="TRUE"){
8 ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	116 sample = readfile(listfile, "true")
0 472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	117 }else{
8 ddaa0c318d65 planemo upload commit d23909db127a655212b79d6826033889812f59b5-dirty proteore parents: 5 diff changeset	118 sample = readfile(listfile, "false")
0 472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	119 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	120 sample = sample[,column]
9 70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	121 }
0 472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	122
9 70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	123 #check of ENS ids
70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	124 if (! any(check_ens_ids(sample))){
70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	125 print("no ensembl gene ids found in your ids list, please check your IDs in input or the selected column of your input file")
70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	126 stop()
0 472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	127 }
9 70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	128
0 472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	129 # Launch enrichment analysis and return result data from the analysis or the null
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	130 # object if the enrichment could not be done.
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	131 goEnrichment = function(geneuniverse,sample,onto){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	132
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	133 # get all the GO terms of the corresponding ontology (BP/CC/MF) and all their
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	134 # associated ensembl ids according to the org package
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	135 xx = annFUN.org(onto,mapping=geneuniverse,ID="ensembl")
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	136 allGenes = unique(unlist(xx))
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	137 # check if the genes given by the user can be found in the org package (gene
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	138 # universe), that is in
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	139 # allGenes
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	140 if (length(intersect(sample,allGenes))==0){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	141
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	142 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.")
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	143 return(c(NULL,NULL))
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	144
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	145 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	146
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	147 geneList = factor(as.integer(allGenes %in% sample))
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	148 names(geneList) <- allGenes
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	149
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	150
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	151 #topGO enrichment
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	152
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	153
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	154 # Creation of a topGOdata object
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	155 # It will contain : the list of genes of interest, the GO annotations and the GO hierarchy
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	156 # Parameters :
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	157 # ontology : character string specifying the ontology of interest (BP, CC, MF)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	158 # allGenes : named vector of type numeric or factor
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	159 # annot : tells topGO how to map genes to GO annotations.
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	160 # argument not used here : nodeSize : at which minimal number of GO annotations
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	161 # do we consider a gene
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	162
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	163 myGOdata = new("topGOdata", description="SEA with TopGO", ontology=onto, allGenes=geneList, annot = annFUN.org, mapping=geneuniverse,ID="ensembl")
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	164
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	165
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	166 # Performing enrichment tests
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	167 result <- runTest(myGOdata, algorithm=option, statistic="fisher")
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	168 return(c(result,myGOdata))
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	169 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	170
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	171 # Some libraries such as GOsummaries won't be able to treat the values such as
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	172 # "< 1e-30" produced by topGO. As such it is important to delete the < char
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	173 # with the deleteInfChar function. Nevertheless the user will have access to the original results in the text output.
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	174 deleteInfChar = function(values){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	175
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	176 lines = grep("<",values)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	177 if (length(lines)!=0){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	178 for (line in lines){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	179 values[line]=gsub("<","",values[line])
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	180 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	181 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	182 return(values)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	183 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	184
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	185 corrMultipleTesting = function(result, myGOdata,correction,threshold){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	186
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	187 # adjust for multiple testing
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	188 if (correction!="none"){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	189 # GenTable : transforms the result object into a list. Filters can be applied
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	190 # (e.g : with the topNodes argument, to get for instance only the n first
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	191 # GO terms with the lowest pvalues), but as we want to apply a correction we
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	192 # take all the GO terms, no matter their pvalues
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	193 allRes <- GenTable(myGOdata, test = result, orderBy = "result", ranksOf = "result",topNodes=length(attributes(result)$score))
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	194 # Some pvalues given by topGO are not numeric (e.g : "<1e-30). As such, these
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	195 # values are converted to 1e-30 to be able to correct the pvalues
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	196 pvaluestmp = deleteInfChar(allRes$test)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	197
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	198 # the correction is done from the modified pvalues
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	199 allRes$qvalues = p.adjust(pvaluestmp, method = as.character(correction), n = length(pvaluestmp))
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	200 allRes = as.data.frame(allRes)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	201
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	202 # Rename the test column by pvalues, so that is more explicit
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	203 nb = which(names(allRes) %in% c("test"))
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	204
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	205 names(allRes)[nb] = "pvalues"
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	206
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	207 allRes = allRes[which(as.numeric(allRes$pvalues) <= threshold),]
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	208 if (length(allRes$pvalues)==0){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	209 print("Threshold was too stringent, no GO term found with pvalue equal or lesser than the threshold value")
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	210 return(NULL)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	211 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	212 allRes = allRes[order(allRes$qvalues),]
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	213 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	214
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	215 if (correction=="none"){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	216 # get all the go terms under user threshold
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	217 mysummary <- summary(attributes(result)$score <= threshold)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	218 numsignif <- as.integer(mysummary[[3]])
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	219 # get all significant nodes
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	220 allRes <- GenTable(myGOdata, test = result, orderBy = "result", ranksOf = "result",topNodes=numsignif)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	221
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	222
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	223 allRes = as.data.frame(allRes)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	224 # Rename the test column by pvalues, so that is more explicit
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	225 nb = which(names(allRes) %in% c("test"))
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	226 names(allRes)[nb] = "pvalues"
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	227 if (numsignif==0){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	228
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	229 print("Threshold was too stringent, no GO term found with pvalue equal or lesser than the threshold value")
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	230 return(NULL)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	231 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	232
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	233 allRes = allRes[order(allRes$pvalues),]
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	234 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	235
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	236 return(allRes)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	237 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	238
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	239 # roundValues will simplify the results by rounding down the values. For instance 1.1e-17 becomes 1e-17
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	240 roundValues = function(values){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	241 for (line in 1:length(values)){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	242 values[line]=as.numeric(gsub(".*e","1e",as.character(values[line])))
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	243 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	244 return(values)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	245 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	246
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	247 createDotPlot = function(data, onto){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	248
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	249 values = deleteInfChar(data$pvalues)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	250 values = roundValues(values)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	251 values = as.numeric(values)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	252
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	253 geneRatio = data$Significant/data$Annotated
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	254 goTerms = data$Term
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	255 count = data$Significant
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	256
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	257 labely = paste("GO terms",onto,sep=" ")
9 70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	258 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" )
70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	259 ggsave("dotplot.png", device = "png", dpi = 320, limitsize = TRUE, width = 15, height = 15, units="cm")
0 472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	260 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	261
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	262 createBarPlot = function(data, onto){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	263
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	264
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	265 values = deleteInfChar(data$pvalues)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	266 values = roundValues(values)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	267
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	268 values = as.numeric(values)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	269 goTerms = data$Term
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	270 count = data$Significant
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	271
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	272 labely = paste("GO terms",onto,sep=" ")
9 70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	273 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")
70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	274 ggsave("barplot.png", device = "png", dpi = 320, limitsize = TRUE, width = 15, height = 15, units="cm")
0 472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	275 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	276
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	277
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	278 # Produce the different outputs
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	279 createOutputs = function(result, cut_result,text, barplot, dotplot, onto){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	280
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	281
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	282 if (is.null(result)){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	283
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	284 if (text=="TRUE"){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	285
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	286 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)."
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	287 write.table(err_msg, file='result.csv', quote=FALSE, sep='\t', col.names = T, row.names = F)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	288
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	289 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	290
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	291 if (barplot=="TRUE"){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	292
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	293 png(filename="barplot.png")
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	294 plot.new()
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	295 #text(0,0,err_msg)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	296 dev.off()
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	297 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	298
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	299 if (dotplot=="TRUE"){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	300
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	301 png(filename="dotplot.png")
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	302 plot.new()
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	303 #text(0,0,err_msg)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	304 dev.off()
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	305
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	306 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	307 return(TRUE)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	308 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	309
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	310
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	311 if (is.null(cut_result)){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	312
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	313
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	314 if (text=="TRUE"){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	315
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	316 err_msg = "Threshold was too stringent, no GO term found with pvalue equal or lesser than the threshold value."
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	317 write.table(err_msg, file='result.csv', quote=FALSE, sep='\t', col.names = T, row.names = F)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	318
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	319 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	320
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	321 if (barplot=="TRUE"){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	322
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	323 png(filename="barplot.png")
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	324 plot.new()
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	325 text(0,0,err_msg)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	326 dev.off()
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	327 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	328
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	329 if (dotplot=="TRUE"){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	330
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	331 png(filename="dotplot.png")
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	332 plot.new()
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	333 text(0,0,err_msg)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	334 dev.off()
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	335
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	336 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	337 return(TRUE)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	338
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	339
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	340
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	341 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	342
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	343 if (text=="TRUE"){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	344 write.table(cut_result, file='result.csv', quote=FALSE, sep='\t', col.names = T, row.names = F)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	345 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	346
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	347 if (barplot=="TRUE"){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	348
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	349 createBarPlot(cut_result, onto)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	350 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	351
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	352 if (dotplot=="TRUE"){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	353
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	354 createDotPlot(cut_result, onto)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	355 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	356 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	357
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	358
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	359
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	360 # Load R library ggplot2 to plot graphs
9 70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	361 suppressMessages(library(ggplot2))
0 472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	362
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	363 # Launch enrichment analysis
9 70c0c8757f5f planemo upload commit 9d3e0b226140b566fc529fd0ffe7aa9e8388c6e5-dirty proteore parents: 8 diff changeset	364 allresult = suppressMessages(goEnrichment(geneuniverse,sample,onto))
0 472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	365 result = allresult[1][[1]]
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	366 myGOdata = allresult[2][[1]]
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	367 if (!is.null(result)){
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	368
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	369 # Adjust the result with a multiple testing correction or not and with the user
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	370 # p-value cutoff
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	371 cut_result = corrMultipleTesting(result,myGOdata, correction,threshold)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	372 }else{
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	373
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	374 cut_result=NULL
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	375
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	376 }
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	377
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	378
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	379 createOutputs(result, cut_result,text, barplot, dotplot, onto)
472ad7da3d92 planemo upload commit 9f9e8b72a239e58b5087b2b3737262c25cc2671e-dirty proteore parents: diff changeset	380

Mercurial > repos > proteore > proteore_topgo

annotate enrichment_v3.R @ 9:70c0c8757f5f draft