8
|
1 #!/usr/bin/env Rscript
|
|
2 # setup R error handling to go to stderr
|
|
3 options( show.error.messages=F, error = function () { cat( geterrmessage(), file=stderr() ); q( "no", 1, F ) } )
|
|
4
|
|
5 # we need that to not crash galaxy with an UTF8 error on German LC settings.
|
|
6 loc <- Sys.setlocale("LC_MESSAGES", "en_US.UTF-8")
|
|
7
|
|
8 library("optparse")
|
|
9
|
|
10 ##### Read options
|
|
11 option_list=list(
|
|
12 make_option("--input",type="character",default="NULL",help="rdata object containing eset object"),
|
|
13 make_option("--conditions",type="character",default=NULL,help="Text file summarizing conditions of the experiment (required)"),
|
|
14 make_option("--normalization",type="character",default=NULL,help="log2 transformation"),
|
|
15 make_option("--annotations",type="character",default="NULL",help="rdata object containing eset object"),
|
|
16 make_option("--rdataoutput",type="character",default="NULL",help="output rdata object containing eset object"),
|
|
17 make_option("--htmloutput",type="character",default=NULL,help="Output html report"),
|
|
18 make_option("--htmloutputpath",type="character",default="NULL",help="Path of output html report"),
|
|
19 make_option("--htmltemplate",type="character",default="NULL",help="html template)")
|
|
20
|
|
21
|
|
22 );
|
|
23
|
|
24
|
|
25 opt_parser = OptionParser(option_list=option_list);
|
|
26 opt = parse_args(opt_parser);
|
2
|
27
|
8
|
28 if(is.null(opt$input)){
|
|
29 print_help(opt_parser)
|
|
30 stop("input required.", call.=FALSE)
|
|
31 }
|
|
32
|
|
33 if(is.null(opt$conditions)){
|
|
34 print_help(opt_parser)
|
|
35 stop("conditions input required.", call.=FALSE)
|
|
36 }
|
|
37
|
|
38
|
|
39 #loading libraries
|
|
40 suppressPackageStartupMessages(require(GEOquery))
|
2
|
41
|
8
|
42 suppressPackageStartupMessages(require(Biobase))
|
|
43 suppressPackageStartupMessages(require(GEOquery))
|
|
44 suppressPackageStartupMessages(require(GEOmetadb))
|
|
45 suppressPackageStartupMessages(require(limma))
|
|
46 suppressPackageStartupMessages(require(jsonlite))
|
|
47 suppressPackageStartupMessages(require(affy))
|
|
48 suppressPackageStartupMessages(require(dplyr))
|
|
49 suppressPackageStartupMessages(require(affyPLM))
|
|
50
|
|
51 dataFile=opt$input
|
|
52 normalization=opt$normalization
|
|
53 conditionsFile=opt$conditions
|
|
54 annotation=opt$annotations
|
|
55 result_export_eset=opt$rdataoutput
|
|
56 result=opt$htmloutput
|
|
57 result.path=opt$htmloutputpath
|
|
58 result.template=opt$htmltemplate
|
2
|
59
|
|
60 dir.create(result.path, showWarnings = TRUE, recursive = FALSE)
|
|
61
|
8
|
62 data=as.matrix(read.table(file = dataFile,row.names=1,header=TRUE))
|
2
|
63 conditions=read.table(file=conditionsFile,sep = "\t",row.names=1)
|
|
64 htmlfile=readChar(result.template, file.info(result.template)$size)
|
|
65
|
|
66 colnames(conditions)=c("source_name_ch1","description")
|
|
67 phenodata<-new("AnnotatedDataFrame",data=conditions)
|
|
68
|
8
|
69 head(data)
|
|
70 conditions
|
|
71
|
2
|
72 eset=ExpressionSet(assayData=data,phenoData=phenodata,annotation=annotation)
|
|
73
|
|
74 if (normalization == "quantile") {
|
8
|
75 eset <- normalize.ExpressionSet.quantiles(eset, transfn="log2")
|
2
|
76 } else if (normalization == "log2") {
|
|
77 exprs(eset) = log2(exprs(eset))
|
|
78 }
|
|
79
|
|
80 boxplotnorm="boxplotnorm.png"
|
|
81 png(boxplotnorm,width=800,height = 400)
|
|
82 par(mar=c(7,5,1,1))
|
|
83 boxplot(data.frame(exprs(eset)),las=2,outline=FALSE)
|
|
84 dev.off()
|
|
85 htmlfile=gsub(x=htmlfile,pattern = "###BOXPLOTNORM###",replacement = boxplotnorm, fixed = TRUE)
|
|
86 file.copy(boxplotnorm,result.path)
|
|
87
|
|
88 plotMAnorm="plotMAnorm.png"
|
|
89 nblines=length(colnames(data))%/%3 + as.numeric((length(colnames(data))%%3)!=0)
|
|
90 png(plotMAnorm,width=800,height =300*nblines )
|
|
91 par(mfrow=c(nblines,3))
|
8
|
92 ##for (i in 1:length(colnames(data))){
|
2
|
93 MAplot(eset)
|
|
94 #}
|
|
95
|
|
96 dev.off()
|
|
97 htmlfile=gsub(x=htmlfile,pattern = "###PLOTMANORM###",replacement = plotMAnorm, fixed = TRUE)
|
|
98 file.copy(plotMAnorm,result.path)
|
|
99 #write.table(tolower(c(condition1Name,condition2Name)),quote = FALSE,col.names = FALSE, row.names=FALSE,file=result_export_conditions)
|
|
100 #saveConditions=c(condition1Name,condition2Name)
|
|
101 save(eset,file=result_export_eset)
|
|
102 write(htmlfile,result)
|
|
103
|
|
104 #l=list()
|
|
105 #for(i in 1:length(esets))
|
|
106 #{
|
|
107 # l[[paste("study",i,sep="")]]<-res[[i]]
|
|
108 #}
|
|
109 #l[["Meta"]]=res[[length(res)-1]]
|
|
110 #showVenn(res,file.path(temp.files.path,"venn.png"))
|
|
111 #writeLines(c("<h2>Venn diagram</h2>"),file.conn)
|
|
112 #writeLines(c("<img src='venn.png'><br/><br/>"),file.conn)
|
|
113 #writeLines(c("</body></html>"),file.conn)
|
|
114 #close(file.conn) |