Mercurial > repos > sblanck > smagexp

comparison Analyse.R @ 6:a2b8c2aabeb0 draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Uploaded
author sblanck
date Wed, 12 Apr 2017 03:45:25 -0400
parents 93451f832736
children f3c021bdc000
comparison
equal deleted inserted replaced
5:5d25fe2fcab2 6:a2b8c2aabeb0
1 library(Biobase) 1 #!/usr/bin/env Rscript
2 library(GEOquery) 2 # setup R error handling to go to stderr
3 library(GEOmetadb) 3 options( show.error.messages=F, error = function () { cat( geterrmessage(), file=stderr() ); q( "no", 1, F ) } )
4 library(limma) 4
5 library(jsonlite) 5 # we need that to not crash galaxy with an UTF8 error on German LC settings.
6 library(affy) 6 loc <- Sys.setlocale("LC_MESSAGES", "en_US.UTF-8")
7 library(dplyr) 7
8 library("optparse")
9
10 ##### Read options
11 option_list=list(
12 make_option("--rdatainput",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("--selectcondition1",type="character",default=NULL,help="log2 transformation"),
15 # make_option("--condition1",type="character",default=NULL,help="A table containing the expression data"),
16 make_option("--selectcondition2",type="character",default="NULL",help="rdata object containing eset object"),
17 # make_option("--condition2",type="character",default=NULL,help="Text file summarizing conditions of the experiment (required)"),
18 make_option("--nbresult",type="character",default=NULL,help="number of result displayed results"),
19 make_option("--rdataoutput",type="character",default="NULL",help="output rdata object containing eset object"),
20 make_option("--htmloutput",type="character",default=NULL,help="Output html report"),
21 make_option("--htmloutputpath",type="character",default="NULL",help="Path of output html report"),
22 make_option("--tabularoutput",type="character",default=NULL,help="Output text file"),
23 make_option("--htmltemplate",type="character",default=NULL,help="html template)")
24
25
26 );
27
28 opt_parser = OptionParser(option_list=option_list);
29 opt = parse_args(opt_parser);
30
31 if(is.null(opt$rdatainput)){
32 print_help(opt_parser)
33 stop("rData input required.", call.=FALSE)
34 }
35
36 if(is.null(opt$conditions)){
37 print_help(opt_parser)
38 stop("conditions input required.", call.=FALSE)
39 }
40
41
42 #loading libraries
43 suppressPackageStartupMessages(require(GEOquery))
44
45 suppressPackageStartupMessages(require(Biobase))
46 suppressPackageStartupMessages(require(GEOquery))
47 suppressPackageStartupMessages(require(GEOmetadb))
48 suppressPackageStartupMessages(require(limma))
49 suppressPackageStartupMessages(require(jsonlite))
50 suppressPackageStartupMessages(require(affy))
51 suppressPackageStartupMessages(require(dplyr))
8 52
9 cargs<-commandArgs() 53 cargs<-commandArgs()
10 cargs<-cargs[(which(cargs=="--args")+1):length(cargs)] 54 cargs<-cargs[(which(cargs=="--args")+1):length(cargs)]
11 nbargs=length(cargs) 55 nbargs=length(cargs)
12 56
13 load(cargs[[nbargs-13]]) 57 load(opt$rdatainput)
14 targetFile=cargs[[nbargs-12]] 58 targetFile=opt$conditions
15 condition1Name=cargs[[nbargs-11]] 59 condition1Name=opt$selectcondition1
16 condition1=cargs[[nbargs-10]] 60 #condition1=opt$condition1
17 condition2Name=cargs[[nbargs-9]] 61 condition2Name=opt$selectcondition2
18 condition2=cargs[[nbargs-8]] 62 #condition2=opt$condition2
19 nbresult=cargs[[nbargs-7]] 63 nbresult=opt$nbresult
20 result_export_eset=cargs[[nbargs-6]] 64 result_export_eset=opt$rdataoutput
21 result=cargs[[nbargs-5]] 65 result=opt$htmloutput
22 result.path=cargs[[nbargs-4]] 66 result.path=opt$htmloutputpath
23 result.tabular=cargs[[nbargs-3]] 67 result.tabular=opt$tabularoutput
24 result.template=cargs[[nbargs-2]] 68 result.template=opt$htmltemplate
25 69
26 #file.copy(targetFile,"./targetFile.txt") 70 #file.copy(targetFile,"./targetFile.txt")
27 71
28 condition1_tmp <- strsplit(condition1,",") 72 targets <- read.table(targetFile,sep="\t",stringsAsFactors=FALSE)
29 condition1 <-unlist(condition1_tmp)
30 73
31 condition2_tmp <- strsplit(condition2,",") 74 #condition1_tmp <- strsplit(condition1,",")
32 condition2<-unlist(condition2_tmp) 75 condition1 <-targets[which(targets$V2==condition1Name),1]
76
77 #condition2_tmp <- strsplit(condition2,",")
78 #condition2<-unlist(condition2_tmp)
79 condition2 <-targets[which(targets$V2==condition2Name),1]
33 80
34 conditions=c(condition1,condition2) 81 conditions=c(condition1,condition2)
35 82
36 #nbresult=1000 83 #nbresult=1000
37 dir.create(result.path, showWarnings = TRUE, recursive = FALSE) 84 dir.create(result.path, showWarnings = TRUE, recursive = FALSE)
38 85
39 targets <- read.table(targetFile,sep="\t") 86 eset=eset[,which(rownames(eset@phenoData@data) %in% conditions)]
40 87
41 eset=eset[,which(rownames(eset@phenoData@data) %in% conditions)] 88 rownames(eset@phenoData@data)
89 which(rownames(eset@phenoData@data) %in% conditions)
42 #condition1Name=make.names(condition1Name) 90 #condition1Name=make.names(condition1Name)
43 #condition2Name=make.names(condition2Name) 91 #condition2Name=make.names(condition2Name)
44 #condition1Name=gsub("_","",condition1Name) 92 #condition1Name=gsub("_","",condition1Name)
45 #condition2Name=gsub("_","",condition2Name) 93 #condition2Name=gsub("_","",condition2Name)
46 #condition1Name 94 #condition1Name
47 #condition2Name 95 #condition2Name
48
49 96
50 eset@phenoData@data$source_name_ch1="" 97 eset@phenoData@data$source_name_ch1=""
51 eset@phenoData@data$source_name_ch1[which(rownames(eset@phenoData@data) %in% condition1)]=condition1Name 98 eset@phenoData@data$source_name_ch1[which(rownames(eset@phenoData@data) %in% condition1)]=condition1Name
52 eset@phenoData@data$source_name_ch1[which(rownames(eset@phenoData@data) %in% condition2)]=condition2Name 99 eset@phenoData@data$source_name_ch1[which(rownames(eset@phenoData@data) %in% condition2)]=condition2Name
53 #condition1Name 100 #condition1Name
54 #condition2Name 101 #condition2Name
55 102
56 condNames=paste0("G",as.numeric(as.character(pData(eset)["source_name_ch1"][,1])!=condition1Name)) 103 condNames=paste0("G",as.numeric(as.character(pData(eset)["source_name_ch1"][,1])!=condition1Name))
57 #condNames=make.names(targets[,2]) 104 #condNames=make.names(targets[,2])
58 #condNames=gsub("_","",condNames) 105 #condNames=gsub("_","",condNames)
59
60 f <- as.factor(condNames) 106 f <- as.factor(condNames)
61 #eset$description <- factors 107 #eset$description <- factors
62 design <- model.matrix(~ 0+f) 108 design <- model.matrix(~ 0+f)
63 design 109 design
64 110