Mercurial > repos > ning > ebseq

comparison EBGeneMultiCondTest.R @ 13:4d6d52fafb86 draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Uploaded
author ning
date Fri, 15 Mar 2013 18:23:26 -0400
parents
children
comparison
equal deleted inserted replaced
12:1108020e2d84 13:4d6d52fafb86
1 sink(file="/tmp/none")
2 sink("/dev/null")
3 options(warn=-1)
4 options(echo=F)
5
6 invisible("EBSeq")
7 suppressMessages(library("EBSeq"))
8
9 args <- commandArgs(trailingOnly = T)
10 inputfile <- args[1]
11 WhetherSampleName <- args[2]
12 CondIn <- args[3]
13 PatternFile <- args[4]
14 outputfile <- args[5]
15 MAP.out<-args[6]
16 Sizesout <-args[7]
17
18 #write.table(PatternFile,outputfile,quote=F,col.names=T,row.names=T,sep = "\t")
19
20 #write.table(args,outputfile,quote=F,col.names=T,row.names=T,sep = "\t")
21
22 Conditions=strsplit(CondIn,split=",")[[1]]
23
24 if(WhetherSampleName=="y"){
25 ReadIn=read.table(inputfile,stringsAsFactors=F,header=T, sep="\t")
26 Names=names(ReadIn)[-1]
27 }
28 if(WhetherSampleName=="n"){
29 ReadIn=read.table(inputfile,stringsAsFactors=F,header=F, sep="\t")
30 Names=paste0("S",1:length(Conditions))
31 }
32
33 PatternIn=read.table(PatternFile,stringsAsFactors=F,header=T,sep="\t")
34
35 if(class(ReadIn[[1]])=="character"){
36 GeneMat=do.call(cbind,ReadIn[-1])
37 rownames(GeneMat)=ReadIn[[1]]
38 colnames(GeneMat)=Names
39 }
40 if(class(ReadIn[[1]])=="numeric"){
41 GeneMat=data.matrix(ReadIn)
42 colnames(GeneMat)=Names
43 }
44
45 Patterns=data.matrix(PatternIn)
46
47
48 Sizes=MedianNorm(GeneMat)
49 EBOut=EBMultiTest(Data=GeneMat,Conditions=as.factor(Conditions),
50 AllParti=Patterns,sizeFactors=Sizes, maxround=5)
51 PPout=GetMultiPP(EBOut)
52 MultiPP=PPout$PP
53 MultiMAP=PPout$MAP
54 Data.norm=round(GetNormalizedMat(GeneMat, Sizes),2)
55
56 Mat=cbind(MultiMAP,Data.norm[names(MultiMAP),])
57
58 colnames(Mat)=
59 c("MAP",Names)
60 options(warn=-1)
61
62 write.table(round(MultiPP,2),file=outputfile,quote=F,col.names=T,row.names=T,sep = "\t")
63 write.table(Mat,file=MAP.out ,quote=F,col.names=T,row.names=T,sep = "\t")
64 write.table(Sizes,file=Sizesout,quote=F,col.names=F,row.names=F,sep = "\t")
65