Mercurial > repos > sblanck > smagexp

annotate MetaRNASeq.R @ 30:4048dae6ace8 draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit 904ac44909a062b118757205861451a3007158b6
author sblanck
date Thu, 21 Jun 2018 08:03:02 -0400
parents 102e4ef8fed2
children 8191cc63589c
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
14
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
1 #!/usr/bin/env Rscript
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
2 # setup R error handling to go to stderr
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
3 options( show.error.messages=F, error = function () { cat( geterrmessage(), file=stderr() ); q( "no", 1, F ) } )
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
4
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
5 # we need that to not crash galaxy with an UTF8 error on German LC settings.
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
6 loc <- Sys.setlocale("LC_MESSAGES", "en_US.UTF-8")
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
7
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
8 library("optparse")
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
9
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
10 ##### Read options
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
11 option_list=list(
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
12 make_option("--input",type="character",default="NULL",help="list of rdata objects containing eset objects"),
29
102e4ef8fed2 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit 808ffb34e11a74cd2396e32868ea56f7ce8d83ff
sblanck
parents: 28
diff changeset
13 make_option("--fdr",type="character",default=NULL,help="Adjusted p-value threshold to be declared differentially expressed"),
14
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
14 make_option("--result",type="character",default=NULL,help="text file containing result of the meta-analysis"),
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
15 make_option("--htmloutput",type="character",default=NULL,help="Output html report"),
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
16 make_option("--htmloutputpath",type="character",default="NULL",help="Path of output html report"),
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
17 make_option("--htmltemplate",type="character",default=NULL,help="html template)")
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
18 );
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
19
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
20 opt_parser = OptionParser(option_list=option_list);
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
21 opt = parse_args(opt_parser);
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
22
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
23 if(is.null(opt$input)){
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
24 print_help(opt_parser)
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
25 stop("input required.", call.=FALSE)
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
26 }
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
27
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
28 #loading libraries
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
29
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
30 suppressPackageStartupMessages(require(affy))
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
31 suppressPackageStartupMessages(require(annaffy))
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
32 suppressPackageStartupMessages(require(VennDiagram))
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
33 suppressPackageStartupMessages(require(GEOquery))
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
34
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
35 listInput <- trimws( unlist( strsplit(trimws(opt$input), ",") ) )
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
36
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
37 listfiles=vector()
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
38 listfilenames=vector()
29
102e4ef8fed2 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit 808ffb34e11a74cd2396e32868ea56f7ce8d83ff
sblanck
parents: 28
diff changeset
39 nbreplicates=vector()
14
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
40
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
41 for (i in 1:length(listInput))
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
42 {
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
43 inputFileInfo <- unlist( strsplit( listInput[i], ';' ) )
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
44 listfiles=c(listfiles,inputFileInfo[1])
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
45 listfilenames=c(listfilenames,inputFileInfo[2])
29
102e4ef8fed2 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit 808ffb34e11a74cd2396e32868ea56f7ce8d83ff
sblanck
parents: 28
diff changeset
46 nbreplicates[i]=inputFileInfo[3]
14
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
47 }
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
48
29
102e4ef8fed2 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit 808ffb34e11a74cd2396e32868ea56f7ce8d83ff
sblanck
parents: 28
diff changeset
49
14
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
50 outputfile <- opt$result
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
51 result.html = opt$htmloutput
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
52 html.files.path=opt$htmloutputpath
e5a94bc69bd6 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit da1436252b5db4e90c39c95140558f0f93544160
sblanck
parents: 2
diff changeset
53 result.template=opt$htmltemplate
2
93451f832736 Uploaded
sblanck
parents:
diff changeset
54
30
4048dae6ace8 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit 904ac44909a062b118757205861451a3007158b6
sblanck
parents: 29
diff changeset
55 alpha=as.numeric(opt$fdr)
2
93451f832736 Uploaded
sblanck
parents:
diff changeset
56
93451f832736 Uploaded
sblanck
parents:
diff changeset
57 #print(comparison)
93451f832736 Uploaded
sblanck
parents:
diff changeset
58
93451f832736 Uploaded
sblanck
parents:
diff changeset
59 listData=lapply(listfiles,read.table)
93451f832736 Uploaded
sblanck
parents:
diff changeset
60 orderData=lapply(listData, function(x) x[order(x[1]), ])
93451f832736 Uploaded
sblanck
parents:
diff changeset
61 rawpval=lapply(orderData,function(x) x[6])
93451f832736 Uploaded
sblanck
parents:
diff changeset
62 rawpval=lapply(rawpval, function(x) as.numeric(unlist(x)))
93451f832736 Uploaded
sblanck
parents:
diff changeset
63
93451f832736 Uploaded
sblanck
parents:
diff changeset
64 DE=list()
93451f832736 Uploaded
sblanck
parents:
diff changeset
65 DE=lapply(orderData, function(x) ifelse(x[7]<=0.05,1,0))
93451f832736 Uploaded
sblanck
parents:
diff changeset
66
93451f832736 Uploaded
sblanck
parents:
diff changeset
67 FC=list()
93451f832736 Uploaded
sblanck
parents:
diff changeset
68 FC=lapply(orderData, function(x) x[3])
93451f832736 Uploaded
sblanck
parents:
diff changeset
69
93451f832736 Uploaded
sblanck
parents:
diff changeset
70 DE=as.data.frame(DE)
93451f832736 Uploaded
sblanck
parents:
diff changeset
71 colnames(DE)=listfilenames
93451f832736 Uploaded
sblanck
parents:
diff changeset
72 FC=as.data.frame(FC)
93451f832736 Uploaded
sblanck
parents:
diff changeset
73 colnames(FC)=listfilenames
93451f832736 Uploaded
sblanck
parents:
diff changeset
74 # the comparison must only have two values and the conds must
93451f832736 Uploaded
sblanck
parents:
diff changeset
75 # be a vector from those values, at least one of each.
93451f832736 Uploaded
sblanck
parents:
diff changeset
76
93451f832736 Uploaded
sblanck
parents:
diff changeset
77 #if (length(comparison) != 2) {
93451f832736 Uploaded
sblanck
parents:
diff changeset
78 # stop("Comparison type must be a tuple: ", cargs[length(cargs) - 8])
93451f832736 Uploaded
sblanck
parents:
diff changeset
79 #}
93451f832736 Uploaded
sblanck
parents:
diff changeset
80
93451f832736 Uploaded
sblanck
parents:
diff changeset
81 sink("/dev/null")
93451f832736 Uploaded
sblanck
parents:
diff changeset
82 dir.create(html.files.path, recursive=TRUE)
93451f832736 Uploaded
sblanck
parents:
diff changeset
83 #library(DESeq)
93451f832736 Uploaded
sblanck
parents:
diff changeset
84 #library(HTSFilter)
93451f832736 Uploaded
sblanck
parents:
diff changeset
85
93451f832736 Uploaded
sblanck
parents:
diff changeset
86 #DE=list()
93451f832736 Uploaded
sblanck
parents:
diff changeset
87 #FC=list()
93451f832736 Uploaded
sblanck
parents:
diff changeset
88 #i=1
93451f832736 Uploaded
sblanck
parents:
diff changeset
89
93451f832736 Uploaded
sblanck
parents:
diff changeset
90 # Open the html output file
93451f832736 Uploaded
sblanck
parents:
diff changeset
91 #file.conn = file(diag.html, open="w")
93451f832736 Uploaded
sblanck
parents:
diff changeset
92
93451f832736 Uploaded
sblanck
parents:
diff changeset
93 #writeLines( c("<html><body>"), file.conn)
93451f832736 Uploaded
sblanck
parents:
diff changeset
94
93451f832736 Uploaded
sblanck
parents:
diff changeset
95 # Perform deseq analysis on each study
93451f832736 Uploaded
sblanck
parents:
diff changeset
96 #for(i in 1:length(listfiles))
93451f832736 Uploaded
sblanck
parents:
diff changeset
97 #{
93451f832736 Uploaded
sblanck
parents:
diff changeset
98 # f=listfiles[i]
93451f832736 Uploaded
sblanck
parents:
diff changeset
99 # fname=listfilenames[i]
93451f832736 Uploaded
sblanck
parents:
diff changeset
100 # study_name=unlist(strsplit(fname,"[.]"))[1]
93451f832736 Uploaded
sblanck
parents:
diff changeset
101 # print(paste0("study.name ",study_name))
93451f832736 Uploaded
sblanck
parents:
diff changeset
102 # d <- read.table(f, sep=" ", header=TRUE, row.names=1)
93451f832736 Uploaded
sblanck
parents:
diff changeset
103 # conds<-sapply(strsplit(colnames(d),"[.]"),FUN=function(x) x[1])
93451f832736 Uploaded
sblanck
parents:
diff changeset
104 # if (length(unique(conds)) != 2) {
93451f832736 Uploaded
sblanck
parents:
diff changeset
105 # warning(as.data.frame(strsplit(colnames(d),"[.]")))
93451f832736 Uploaded
sblanck
parents:
diff changeset
106 # stop("You can only have two columns types: ", paste(conds,collapse=" "))
93451f832736 Uploaded
sblanck
parents:
diff changeset
107 # }
93451f832736 Uploaded
sblanck
parents:
diff changeset
108 # if (!identical(sort(comparison), sort(unique(conds)))) {
93451f832736 Uploaded
sblanck
parents:
diff changeset
109 # stop("Column types must use the two names from Comparison type, and vice versa. Must have at least one of each in the Column types.\nColumn types: ", cargs[2], "\n", "Comparison type: ", cargs[3])
93451f832736 Uploaded
sblanck
parents:
diff changeset
110 # }
93451f832736 Uploaded
sblanck
parents:
diff changeset
111 # if (length(d) != length(conds)) {
93451f832736 Uploaded
sblanck
parents:
diff changeset
112 # stop("Number of total sample columns in counts file must correspond to the columns types field. E.g. if column types is 'kidney,kidney,liver,liver' then number of sample columns in counts file must be 4 as well.")
93451f832736 Uploaded
sblanck
parents:
diff changeset
113 # }
93451f832736 Uploaded
sblanck
parents:
diff changeset
114 #
93451f832736 Uploaded
sblanck
parents:
diff changeset
115 # cds <- newCountDataSet(d, conds)
93451f832736 Uploaded
sblanck
parents:
diff changeset
116 # cds <- estimateSizeFactors(cds)
93451f832736 Uploaded
sblanck
parents:
diff changeset
117 #
93451f832736 Uploaded
sblanck
parents:
diff changeset
118 # cdsBlind <- estimateDispersions( cds, method="blind" )
93451f832736 Uploaded
sblanck
parents:
diff changeset
119 #
93451f832736 Uploaded
sblanck
parents:
diff changeset
120 # if (length(conds) != 2) {
93451f832736 Uploaded
sblanck
parents:
diff changeset
121 # cds <- estimateDispersions( cds )
93451f832736 Uploaded
sblanck
parents:
diff changeset
122 # norep = FALSE
93451f832736 Uploaded
sblanck
parents:
diff changeset
123 # }
93451f832736 Uploaded
sblanck
parents:
diff changeset
124 #
93451f832736 Uploaded
sblanck
parents:
diff changeset
125 # if (length(conds) == 2) {
93451f832736 Uploaded
sblanck
parents:
diff changeset
126 # cds <- estimateDispersions( cds, method=method, sharingMode=mod, fitType="parametric" )
93451f832736 Uploaded
sblanck
parents:
diff changeset
127 # norep = TRUE
93451f832736 Uploaded
sblanck
parents:
diff changeset
128 # }
93451f832736 Uploaded
sblanck
parents:
diff changeset
129 #
93451f832736 Uploaded
sblanck
parents:
diff changeset
130 # filter<-HTSFilter(cds, plot=FALSE)
93451f832736 Uploaded
sblanck
parents:
diff changeset
131 # cds.filter<-filter$filteredData
93451f832736 Uploaded
sblanck
parents:
diff changeset
132 # on.index<-which(filter$on==1)
93451f832736 Uploaded
sblanck
parents:
diff changeset
133 #
93451f832736 Uploaded
sblanck
parents:
diff changeset
134 # res<-as.data.frame(matrix(NA,nrow=nrow(cds),ncol=ncol(cds)))
93451f832736 Uploaded
sblanck
parents:
diff changeset
135 # nbT <- nbinomTest(cds.filter, comparison[1], comparison[2])
93451f832736 Uploaded
sblanck
parents:
diff changeset
136 # colnames(res)<-colnames(nbT)
93451f832736 Uploaded
sblanck
parents:
diff changeset
137 # res[on.index,]<-nbT
93451f832736 Uploaded
sblanck
parents:
diff changeset
138 # #write.table(res[order(res$padj), ], file=outputfile, quote=FALSE, row.names=FALSE, sep="\t")
93451f832736 Uploaded
sblanck
parents:
diff changeset
139 #
93451f832736 Uploaded
sblanck
parents:
diff changeset
140 #
93451f832736 Uploaded
sblanck
parents:
diff changeset
141 # temp.pval.plot = file.path( html.files.path, paste("PvalHist",i,".png",sep=""))
93451f832736 Uploaded
sblanck
parents:
diff changeset
142 # png( temp.pval.plot, width=500, height=500 )
93451f832736 Uploaded
sblanck
parents:
diff changeset
143 # hist(res$pval, breaks=100, col="skyblue", border="slateblue", main="")
93451f832736 Uploaded
sblanck
parents:
diff changeset
144 # dev.off()
93451f832736 Uploaded
sblanck
parents:
diff changeset
145 #
93451f832736 Uploaded
sblanck
parents:
diff changeset
146 # writeLines( c("<h2>P-value histogram for ",study_name,"</h2>"), file.conn)
93451f832736 Uploaded
sblanck
parents:
diff changeset
147 # writeLines( c("<img src='PvalHist",i,".png'><br/><br/>"), file.conn)
93451f832736 Uploaded
sblanck
parents:
diff changeset
148 #
93451f832736 Uploaded
sblanck
parents:
diff changeset
149 # #on enregistre la p-value
93451f832736 Uploaded
sblanck
parents:
diff changeset
150 # rawpval[[study_name]]<-res$pval
93451f832736 Uploaded
sblanck
parents:
diff changeset
151 # DE[[study_name]]<-ifelse(res$padj<=alpha,1,0)
93451f832736 Uploaded
sblanck
parents:
diff changeset
152 # FC[[study_name]]<-res$log2FoldChange
93451f832736 Uploaded
sblanck
parents:
diff changeset
153 #
93451f832736 Uploaded
sblanck
parents:
diff changeset
154 # i=i+1
93451f832736 Uploaded
sblanck
parents:
diff changeset
155 #}
93451f832736 Uploaded
sblanck
parents:
diff changeset
156
93451f832736 Uploaded
sblanck
parents:
diff changeset
157
93451f832736 Uploaded
sblanck
parents:
diff changeset
158 # combinations
93451f832736 Uploaded
sblanck
parents:
diff changeset
159 library(metaRNASeq)
28
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
160 library(UpSetR)
2
93451f832736 Uploaded
sblanck
parents:
diff changeset
161 fishcomb<-fishercomb(rawpval, BHth=alpha)
93451f832736 Uploaded
sblanck
parents:
diff changeset
162 warning(length(rawpval))
28
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
163
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
164 #nrep=c(length(listFiles))
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
165
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
166
29
102e4ef8fed2 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit 808ffb34e11a74cd2396e32868ea56f7ce8d83ff
sblanck
parents: 28
diff changeset
167 invnormcomb<-invnorm(rawpval, nrep=nbreplicates, BHth=alpha)
2
93451f832736 Uploaded
sblanck
parents:
diff changeset
168 #DE[["fishercomb"]]<-ifelse(fishcomb$adjpval<=alpha,1,0)
93451f832736 Uploaded
sblanck
parents:
diff changeset
169 #DE[["invnormcomb"]]<-ifelse(invnormcomb$adjpval<=alpha,1,0)
93451f832736 Uploaded
sblanck
parents:
diff changeset
170
93451f832736 Uploaded
sblanck
parents:
diff changeset
171 signsFC<-mapply(FC,FUN=function(x) sign(x))
93451f832736 Uploaded
sblanck
parents:
diff changeset
172 sumsigns<-apply(signsFC,1,sum)
93451f832736 Uploaded
sblanck
parents:
diff changeset
173 commonsgnFC<-ifelse(abs(sumsigns)==dim(signsFC)[2],sign(sumsigns),0)
93451f832736 Uploaded
sblanck
parents:
diff changeset
174
93451f832736 Uploaded
sblanck
parents:
diff changeset
175 DEresults <- data.frame(DE=DE,"DE.fishercomb"=ifelse(fishcomb$adjpval<=alpha,1,0),"DE.invnorm"=ifelse(invnormcomb$adjpval<=alpha,1,0))
93451f832736 Uploaded
sblanck
parents:
diff changeset
176
93451f832736 Uploaded
sblanck
parents:
diff changeset
177 unionDE <- unique(c(fishcomb$DEindices,invnormcomb$DEindices))
93451f832736 Uploaded
sblanck
parents:
diff changeset
178 FC.selecDE <- data.frame(DEresults[unionDE,],FC[unionDE,],signFC=commonsgnFC[unionDE])
93451f832736 Uploaded
sblanck
parents:
diff changeset
179 keepDE <- FC.selecDE[which(abs(FC.selecDE$signFC)==1),]
93451f832736 Uploaded
sblanck
parents:
diff changeset
180
93451f832736 Uploaded
sblanck
parents:
diff changeset
181 fishcomb_de <- rownames(keepDE)[which(keepDE[,"DE.fishercomb"]==1)]
93451f832736 Uploaded
sblanck
parents:
diff changeset
182 invnorm_de <- rownames(keepDE)[which(keepDE[,"DE.invnorm"]==1)]
93451f832736 Uploaded
sblanck
parents:
diff changeset
183 indstudy_de = list()
93451f832736 Uploaded
sblanck
parents:
diff changeset
184 for (i in 1:length(listfiles)) {
93451f832736 Uploaded
sblanck
parents:
diff changeset
185 currentIndstudy_de = rownames(keepDE)[which(keepDE[,i]==1)]
93451f832736 Uploaded
sblanck
parents:
diff changeset
186 indstudy_de[[listfilenames[i]]]=currentIndstudy_de
93451f832736 Uploaded
sblanck
parents:
diff changeset
187 }
93451f832736 Uploaded
sblanck
parents:
diff changeset
188
93451f832736 Uploaded
sblanck
parents:
diff changeset
189 IDDIRRfishcomb=IDD.IRR(fishcomb_de,indstudy_de)
93451f832736 Uploaded
sblanck
parents:
diff changeset
190 IDDIRRinvnorm=IDD.IRR(invnorm_de,indstudy_de)
93451f832736 Uploaded
sblanck
parents:
diff changeset
191
93451f832736 Uploaded
sblanck
parents:
diff changeset
192 #conflits<-data.frame(ID=listData[[1]][rownames(DEresults),1],Fishercomb=DEresults[["DE.fishercomb"]],Invnormcomb=DEresults[["DE.invnorm"]],sign=commonsgnFC)
93451f832736 Uploaded
sblanck
parents:
diff changeset
193 conflits<-data.frame(ID=listData[[1]][rownames(DEresults),1],DE=DEresults,FC=FC,signFC=commonsgnFC)
93451f832736 Uploaded
sblanck
parents:
diff changeset
194 #write DE outputfile
93451f832736 Uploaded
sblanck
parents:
diff changeset
195 write.table(conflits, outputfile,sep="\t",,row.names=FALSE)
93451f832736 Uploaded
sblanck
parents:
diff changeset
196 library(VennDiagram)
28
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
197 DE_num=apply(keepDE[,1:(length(listfiles)+2)], 2, FUN=function(x) which(x==1))
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
198 #DE_num=apply(DEresults, 2, FUN=function(x) which(x==1))
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
199 if (length(listfiles<=2)) {
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
200 title="VENN DIAGRAM"
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
201 width=500
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
202 venn.plot<-venn.diagram(x=as.list(DE_num),filename=NULL, col="black", fill=1:length(DE_num)+1,alpha=0.6)
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
203 temp.venn.plot = file.path( html.files.path, paste("venn.png"))
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
204 png(temp.venn.plot,width=width,height=500)
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
205 grid.draw(venn.plot)
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
206 dev.off()
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
207 } else {
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
208 title="UPSETR DIAGRAM"
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
209 width=1000
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
210 png(temp.venn.plot,width=width,height=500)
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
211 upset(fromList(as.list(DE_num)), order.by = "freq")
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
212 dev.off()
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
213
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
214 }
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
215
2
93451f832736 Uploaded
sblanck
parents:
diff changeset
216
93451f832736 Uploaded
sblanck
parents:
diff changeset
217 library(jsonlite)
93451f832736 Uploaded
sblanck
parents:
diff changeset
218 matrixConflits=as.matrix(conflits)
93451f832736 Uploaded
sblanck
parents:
diff changeset
219 datajson=toJSON(matrixConflits,pretty = TRUE)
93451f832736 Uploaded
sblanck
parents:
diff changeset
220 summaryFishcombjson=toJSON(as.matrix(t(IDDIRRfishcomb)),pretty = TRUE)
93451f832736 Uploaded
sblanck
parents:
diff changeset
221 summaryinvnormjson=toJSON(as.matrix(t(IDDIRRinvnorm)),pretty = TRUE)
93451f832736 Uploaded
sblanck
parents:
diff changeset
222
93451f832736 Uploaded
sblanck
parents:
diff changeset
223
93451f832736 Uploaded
sblanck
parents:
diff changeset
224 #vennsplit=strsplit(result.venn,split="/")[[1]]
93451f832736 Uploaded
sblanck
parents:
diff changeset
225 #venn=paste0("./",vennsplit[length(vennsplit)])
93451f832736 Uploaded
sblanck
parents:
diff changeset
226
93451f832736 Uploaded
sblanck
parents:
diff changeset
227
93451f832736 Uploaded
sblanck
parents:
diff changeset
228 vennFilename="venn.png"
93451f832736 Uploaded
sblanck
parents:
diff changeset
229 vennFile=file.path(html.files.path,vennFilename)
93451f832736 Uploaded
sblanck
parents:
diff changeset
230 htmlfile=readChar(result.template, file.info(result.template)$size)
93451f832736 Uploaded
sblanck
parents:
diff changeset
231 htmlfile=gsub(x=htmlfile,pattern = "###DATAJSON###",replacement = datajson, fixed = TRUE)
93451f832736 Uploaded
sblanck
parents:
diff changeset
232 htmlfile=gsub(x=htmlfile,pattern = "###FISHSUMMARYJSON###",replacement = summaryFishcombjson, fixed = TRUE)
93451f832736 Uploaded
sblanck
parents:
diff changeset
233 htmlfile=gsub(x=htmlfile,pattern = "###INVSUMMARYJSON###",replacement = summaryinvnormjson, fixed = TRUE)
93451f832736 Uploaded
sblanck
parents:
diff changeset
234 htmlfile=gsub(x=htmlfile,pattern = "###VENN###",replacement = vennFilename, fixed = TRUE)
28
49ce6fbe11de planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit cb90fdea2a64908e301355aef89de403e107764d
sblanck
parents: 24
diff changeset
235 htmlfile=gsub(x=htmlfile,pattern = "###WIDTH##",replacement = as.character(width), fixed = TRUE)
29
102e4ef8fed2 planemo upload for repository https://github.com/sblanck/smagexp/tree/master/smagexp_tools commit 808ffb34e11a74cd2396e32868ea56f7ce8d83ff
sblanck
parents: 28
diff changeset
236 htmlfile=gsub(x=htmlfile,pattern = "###TITLE##",replacement = title, fixed = TRUE)
2
93451f832736 Uploaded
sblanck
parents:
diff changeset
237 write(htmlfile,result.html)
93451f832736 Uploaded
sblanck
parents:
diff changeset
238
93451f832736 Uploaded
sblanck
parents:
diff changeset
239 #library(VennDiagram)
93451f832736 Uploaded
sblanck
parents:
diff changeset
240 #flog.threshold(ERROR)
93451f832736 Uploaded
sblanck
parents:
diff changeset
241 #
93451f832736 Uploaded
sblanck
parents:
diff changeset
242 ##venn.plot<-venn.diagram(x = c(res[c(1:(length(res)-3))],meta=list(res$Meta)),filename = v, col = "black", fill = c(1:(length(res)-2)), margin=0.05, alpha = 0.6,imagetype = "png")
93451f832736 Uploaded
sblanck
parents:
diff changeset
243 #dir.create(result.path, showWarnings = TRUE, recursive = FALSE)
93451f832736 Uploaded
sblanck
parents:
diff changeset
244 #
93451f832736 Uploaded
sblanck
parents:
diff changeset
245 #showVenn<-function(liste,file)
93451f832736 Uploaded
sblanck
parents:
diff changeset
246 #{
93451f832736 Uploaded
sblanck
parents:
diff changeset
247 # venn.plot<-venn.diagram(x = liste,
93451f832736 Uploaded
sblanck
parents:
diff changeset
248 # filename = vennFilename, col = "black",
93451f832736 Uploaded
sblanck
parents:
diff changeset
249 # fill = 1:length(liste)+1,
93451f832736 Uploaded
sblanck
parents:
diff changeset
250 # margin=0.05, alpha = 0.6,imagetype = "png")
93451f832736 Uploaded
sblanck
parents:
diff changeset
251 ## png(file);
93451f832736 Uploaded
sblanck
parents:
diff changeset
252 ## grid.draw(venn.plot);
93451f832736 Uploaded
sblanck
parents:
diff changeset
253 ## dev.off();
93451f832736 Uploaded
sblanck
parents:
diff changeset
254 #
93451f832736 Uploaded
sblanck
parents:
diff changeset
255 #}
93451f832736 Uploaded
sblanck
parents:
diff changeset
256 #
93451f832736 Uploaded
sblanck
parents:
diff changeset
257 #l=list()
93451f832736 Uploaded
sblanck
parents:
diff changeset
258 #for(i in 1:length(esets))
93451f832736 Uploaded
sblanck
parents:
diff changeset
259 #{
93451f832736 Uploaded
sblanck
parents:
diff changeset
260 # l[[paste("study",i,sep="")]]<-res[[i]]
93451f832736 Uploaded
sblanck
parents:
diff changeset
261 #}
93451f832736 Uploaded
sblanck
parents:
diff changeset
262 #l[["Meta"]]=res[[length(res)-1]]
93451f832736 Uploaded
sblanck
parents:
diff changeset
263 #showVenn(l,vennFile)
93451f832736 Uploaded
sblanck
parents:
diff changeset
264 #file.copy(vennFilename,result.path)
93451f832736 Uploaded
sblanck
parents:
diff changeset
265
93451f832736 Uploaded
sblanck
parents:
diff changeset
266
93451f832736 Uploaded
sblanck
parents:
diff changeset
267 #writeLines( c("<h2>Venn Plot</h2>"), file.conn)
93451f832736 Uploaded
sblanck
parents:
diff changeset
268 #writeLines( c("<img src='venn.png'><br/><br/>"), file.conn)
93451f832736 Uploaded
sblanck
parents:
diff changeset
269 #writeLines( c("</body></html>"), file.conn)
93451f832736 Uploaded
sblanck
parents:
diff changeset
270 #close(file.conn)
93451f832736 Uploaded
sblanck
parents:
diff changeset
271 #print("passe6")
93451f832736 Uploaded
sblanck
parents:
diff changeset
272 #sink(NULL)