Mercurial > repos > mmonsoor > probmetab

annotate export.class.table-color-graph.R @ 4:52b222a626b0 draft default tip

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
planemo upload commit 00684d80f032fee5bd1cb86e05a477fcdcb1c3fc
author lecorguille
date Fri, 07 Apr 2017 09:11:22 -0400
parents e13ec2c3fabe
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
0
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
1 #' export.class.table
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
2 #'
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
3 #' Builds a matrix with the probability for all mass to candidate compounds
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
4 #' assignments, by averaging the number of assignments obtained by the gibbs sampler algorithm
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
5 #' or ordering the compound candidates with the likelihood matrix.
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
6 #'
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
7 #' @param gibbsL a list of attributions and probabilities from gibbs.samp function.
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
8 #' @param reactionM data.frame with compound annotation information.
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
9 #' @param molIon non redundant ion annotation.
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
10 #' @param probM optionally to gibbsL, a matrix of likelihoods.
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
11 #' @param html logical, indicating whether a html file should be generated. This parameter uses the raw data to plot EICs and may be time consuming.
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
12 #' @param filename html file name, the default is "test".
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
13 #' @param burnIn how many samples of the gibbs sampler should be discarded.
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
14 #' @param linkPattern which pattern should be linked to compound id, for now we have
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
15 #' implemented "kegg", "pubchem" and "chebi" patterns.
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
16 #' @param m.test statistical test to compare mean differences. This option
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
17 #' is only available to single acquisition mode analysis, with options
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
18 #' "t.test" and "anova".
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
19 #' @param class1 if the m.test is "t.test" first class to compare in the test,
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
20 #' according with xcmsSet phenoData.
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
21 #' @param class2 if the m.test is "t.test" second class to compare in the test,
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
22 #' according with xcmsSet phenoData.
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
23 #' @param norm logical, if TRUE performs median normalization from (Anal. Chem. 2011, 83, 5864-5872).
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
24 #' @param DB data.frame table used to search compounds, with the field name to be incorporated in the html table.
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
25 #' @param prob how to calculate the probability to attribute a mass to a compound.
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
26 #' Default is "count", which divide the number of times each identity was
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
27 #' was attributed by the number of samples. Optionally the user could
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
28 #' choose to use the mean of the probabilities of the identity, "mean".
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
29 #' @return A list with a matrix "classTable" with attributions and probabilities and
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
30 #' indexes of selected masses from xcms peak table.
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
31 #'
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
32 #' @export
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
33
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
34 export.class.table <- function(gibbsL=NULL, reactionM, molIon=NULL, probM=NULL, html=FALSE, filename="test", burnIn=3000, linkPattern="kegg", m.test="none", class1=NULL, class2=NULL, norm=FALSE, DB, prob="count", colorplot=FALSE, addLink=NULL) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
35
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
36 plotEIC <- function (xcmsObject, figidx, pngidx, colorplot, mode=NULL) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
37 dir.create(paste(filename,"_fig",sep=""))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
38 gt<-groups(xcmsObject)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
39 if(colorplot==TRUE){
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
40 gt2 <- gt[figidx,]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
41 rgt <- gt2[,c("rtmin","rtmax")]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
42 rgt[,1] <- rgt[,1]-100
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
43 rgt[,2] <- rgt[,2]+100
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
44 #require(doMC)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
45 #registerDoMC()
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
46 #system.time(
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
47 #foreach(i=1:nrow(gt2)) %dopar% {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
48 for(i in 1:nrow(gt2)){
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
49 groupidx1 <- which(gt[,"rtmed"] > rgt[i,1] & gt[,"rtmed"] < rgt[i,2] & gt[,"mzmed"]> gt2[i,"mzmed"]-0.1 & gt[,"mzmed"]< gt2[i,"mzmed"]+0.1)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
50 eiccor <- getEIC(xcmsObject, groupidx = groupidx1)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
51 png(paste(filename, "_fig/", sprintf("%003d", i), ".png", sep=""))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
52 plot(eiccor, xcmsObject, groupidx = 1)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
53 dev.off()
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
54 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
55 } else {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
56 gt <- gt[figidx,]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
57 rgt <- gt[,c("rtmin","rtmax")]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
58 rgt[,1] <- rgt[,1]-100
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
59 rgt[,2] <- rgt[,2]+100
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
60
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
61 eics <- getEIC(xcmsObject, mzrange=gt, rtrange =rgt, groupidx = 1:nrow(gt))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
62 png(file.path(paste(filename, "_fig/%003d.png", sep="")), height=768, width=1024)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
63 #png(file.path(paste(filename, "_fig/", pngidx, sep="")), h=768, w=1024)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
64 plot(eics, xcmsObject)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
65 dev.off()
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
66 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
67 if(!is.null(mode)) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
68 pngs <- dir(paste(filename, "_fig/", sep=""))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
69 if(length(grep("pos|neg" , pngs))) pngs <- pngs[-grep("pos|neg" , pngs)]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
70 opng <- as.numeric(sub(".png","", pngs))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
71 pngs <- pngs[order(opng)]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
72 name1 <- paste(filename, "_fig/",pngs, sep="")
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
73 name2 <- paste(filename, "_fig/",pngidx, mode, ".png", sep="")
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
74 for(i in 1:length(name1)) file.rename(name1[i], name2[i])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
75 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
76
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
77 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
78 allion <- molIon$molIon[molIon$molIon[,"isotope"]==0,]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
79 ReactMatrix <- reactionM[reactionM[,5]!="unknown",]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
80 x <- apply(unique(ReactMatrix[,c(2, 3)]), 2, as.numeric) # Have to look for all pairs
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
81 y <- as.numeric(ReactMatrix[,4])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
82 prob_mean_ma <- matrix(0, nrow = length(y), ncol = nrow(x))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
83 # z_average <- matrix(0, nrow = length(y), ncol = length(x))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
84
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
85 if (!is.null(gibbsL)){
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
86 prob_table <- gibbsL$prob_table[,-c(1:burnIn)]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
87 class_table <- gibbsL$class_table[,-c(1:burnIn)]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
88 #indList <- tapply(1:nrow(ReactMatrix), as.numeric(ReactMatrix[,1]), function(x) x)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
89 coords <- tapply(1:nrow(ReactMatrix), ReactMatrix[,"molIonID"], function(x) x)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
90 coords2 <- unlist(lapply(coords, function(x) rep(x[1], length(x))))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
91 indList <- coords[order(unique(coords2))]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
92 fillMatrix <- function(j,i) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
93 idp <- which(class_table[i,] == j)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
94 if(prob=="count") prob_mean_ma[j,i] <<- length(idp)/ncol(class_table)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
95 if(prob=="mean") prob_mean_ma[j,i] <<- mean(prob_table[i,idp])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
96 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
97
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
98
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
99 for ( i in 1:nrow(x) ) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
100
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
101 sapply(indList[[i]], fillMatrix, i)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
102 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
103 if(sum(prob_mean_ma=="NaN")) prob_mean_ma[prob_mean_ma=="NaN"] <- 0
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
104 # for ( i in 1:nrow(x) ) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
105 # for ( j in 1:length(y) ) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
106 # idp <- which(class_table[i,] == j)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
107 # prob_mean_ma[j,i] <- mean(prob_table[i,idp])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
108 # # this is an alternative way to calculate the probabilities, should try latter, and compare results
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
109 # #prob_mean_ma[j,i] <- length(idp)/ncol(class_table)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
110 # if ( prob_mean_ma[j,i] == "NaN" ) prob_mean_ma[j,i] <- 0
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
111 # }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
112 # # do I still need this matrix?
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
113 # k=which(prob_mean_ma[,i]==max(prob_mean_ma[,i]))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
114 # z_average[k[1],i]=1
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
115 # }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
116 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
117 else {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
118 prob_mean_ma <- probM
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
119 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
120 # think about natural probabilities
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
121 # prob_mean_ma[prob_mean_ma[,1]!=0,1]/sum(prob_mean_ma[prob_mean_ma[,1]!=0,1])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
122 prob_mean_ma <- apply(prob_mean_ma, 2, function(x){ x[x!=0] <- x[x!=0]/sum(x[x!=0]); return(x)} )
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
123
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
124 # create a dir to figures
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
125 lpattern <- function(type){
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
126 switch(type,
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
127 kegg = "http://www.genome.jp/dbget-bin/www_bget?",
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
128 chebi = "http://www.ebi.ac.uk/chebi/searchId.do;EFB7DFF9E88306BBCD6AB78B32664A85?chebiId=",
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
129 pubchem = "http://www.ncbi.nlm.nih.gov/pccompound/?term="
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
130 )
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
131 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
132 linkURL <- lpattern(linkPattern)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
133 fig <- paste("file://", getwd(), paste("/",filename,"_fig/",sep=""), sep="")
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
134 if(!is.null(molIon$cameraobj)) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
135 figidx <- c("")
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
136 coords <- gsub("(^\\d)","X\\1",rownames(molIon$cameraobj@xcmsSet@phenoData))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
137 # experimental! Which set of characters????
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
138 coords <- gsub("-|\\,|~","\\.",coords)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
139 coords <- gsub("\\s+","\\.",coords)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
140 peaklist <- getPeaklist(molIon$cameraobj)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
141 rpeaklist <- peaklist[,c("mz","rt","isotopes","adduct","pcgroup")]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
142 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
143 else {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
144 figidx <- c("","")
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
145 coordsP <- gsub("(^\\d)","X\\1",rownames(molIon$pos@xcmsSet@phenoData))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
146 # experimental! Which set of characters????
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
147 coordsP <- gsub("-|\\,|~","\\.",coordsP)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
148 coordsP <- gsub("\\s+","\\.",coordsP)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
149 coordsN <- gsub("(^\\d)","X\\1",rownames(molIon$neg@xcmsSet@phenoData))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
150 # experimental! Which set of characters????
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
151 coordsN <- gsub("-|\\,|~","\\.",coordsN)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
152 coordsN <- gsub("\\s+","\\.",coordsN)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
153 coords <- coordsP
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
154 if(length(coordsP)!=length(coordsN)) cat("\n Warning: The number of samples are different\n")
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
155
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
156 peaklistP <- getPeaklist(molIon$pos)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
157 rpeaklistP <- peaklistP[,c("mz","rt","isotopes","adduct","pcgroup")]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
158 peaklistN <- getPeaklist(molIon$neg)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
159 rpeaklistN <- peaklistN[,c("mz","rt","isotopes","adduct","pcgroup")]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
160 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
161
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
162 # if(sum(is.na(peaklist))) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
163 # cat("\nWarning: NAs Found in peaklist\n\nSubstituting for \"ones\"\n")
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
164 # na.ids <- which(is.na(peaklist),arr.ind=TRUE)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
165 # for(l in 1:nrow(na.ids)){
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
166 # peaklist[na.ids[l,][1], na.ids[l,][2]] <- 1
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
167 # }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
168 # }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
169 #
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
170
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
171 ans <- matrix("", nrow=1, ncol=7+length(coords))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
172 unq <- unique(ReactMatrix[,2:3])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
173 for (i in 1:nrow(unq)) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
174 coord <- which(ReactMatrix[,2]==unq[i,1] & ReactMatrix[,3]==unq[i,2])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
175 coord2 <- which(allion[,2]==unq[i,1] & allion[,1]==unq[i,2])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
176 # idx2 <- unique(which(allion[,1] %in% reactionM[reactionM[,5]=="unknown",2]))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
177 # work with the higher intensities for a given ion annotation, not necessarily the right one
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
178
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
179 if(!is.null(molIon$cameraobj)) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
180 idx <- as.vector(unlist(sapply(allion[coord2,"trace"],
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
181 function(x) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
182 x <- as.matrix(x)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
183 raw <- strsplit(x,";")[[1]]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
184 mraw <- apply(peaklist[raw, coords], 1, mean)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
185 raw[which.max(mraw)]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
186 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
187
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
188 )
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
189 )
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
190 )
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
191
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
192 idx <- unique(idx)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
193 figidx <- append(figidx,idx)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
194 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
195 else {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
196 idx <- c()
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
197
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
198 for(l in 1:nrow( allion[coord2,c("trace","comb")])) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
199 x <- as.matrix(allion[coord2,c("trace","comb")][l,])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
200 raw <- strsplit(x[1],";")[[1]]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
201 if(x[2]!="neg"){
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
202 mraw <- apply(peaklistP[raw, coordsP], 1, mean, na.rm=TRUE)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
203 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
204 else {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
205
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
206 mraw <- apply(peaklistN[raw, coordsN], 1, mean, na.rm=TRUE)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
207 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
208 idx <- c(idx, raw[which.max(mraw)])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
209 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
210
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
211
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
212 idx <- unique(idx)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
213 figidx <- rbind(figidx,c(idx,allion[coord2,"comb"][1]))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
214 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
215 #figidx <- append(figidx,strsplit(allion[coord2,5], ";")[[1]][1])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
216 ans1 <- matrix("", nrow=length(coord), ncol=7+length(coords))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
217 ans1[,2]<-as.matrix(ReactMatrix[coord,5])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
218 prob <- as.matrix(prob_mean_ma[coord, i]) # need to change and compare a pair of mass/rt
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
219 # number figs
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
220 if ( i >= 100 ) { ans1[1,6]=i }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
221 else { if ( i >= 10 ) { ans1[1,6]=paste(0,i, sep="") } else { ans1[1,6]=paste("00",i, sep="") } }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
222
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
223 if (sum(prob)>0) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
224 #prob <- prob/sum(prob)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
225 o <- order(prob, decreasing=TRUE)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
226 ans1[,-6] <- ans1[o,-6]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
227 ans1 <- matrix(ans1, nrow=length(o))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
228 ans1[1,1] <- ReactMatrix[coord[1],3]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
229 #ans1[,3] <- round(prob/min(prob[prob!=0]), 3)[o]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
230 ans1[,3] <- round(prob, 3)[o]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
231 if (length(prob[prob!=0])>1) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
232 entropy <- -sum(prob[prob!=0]*log(prob[prob!=0], length(prob[prob!=0])))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
233 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
234 else { entropy <- 0
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
235 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
236 ans1[1,4] <- round(entropy, 3)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
237 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
238 else {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
239 ans1[1,1] <- ReactMatrix[coord[1],3]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
240 ans1[1,3] <- "undef"
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
241 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
242
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
243 if(!is.null(molIon$cameraobj)) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
244 ans1[1,7] <- apply(rpeaklist[idx,], 1, function(x) paste(x[c(1,2,3,4)], collapse="#"))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
245 ans1[1,8:ncol(ans1)] <- as.matrix(peaklist[idx, coords])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
246 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
247 else {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
248 if(allion[coord2,"comb"]=="pos"|allion[coord2,"comb"]=="both") {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
249 ans1[1,7] <- apply(rpeaklistP[idx,], 1, function(x) paste(x[c(1,2,3,4)], collapse="#"))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
250 ans1[1,8:ncol(ans1)] <- as.matrix(peaklistP[idx, coordsP])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
251 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
252 else {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
253 ans1[1,7] <- apply(rpeaklistN[idx,], 1, function(x) paste(x[c(1,2,3,4)], collapse="#"))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
254 ans1[1,8:ncol(ans1)] <- as.matrix(peaklistN[idx, coordsN])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
255 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
256 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
257 ans <- rbind(ans, as.matrix(ans1))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
258 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
259 ans <- ans[-1,]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
260 # this option should change according with the bank
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
261 if(html) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
262 nid <- unlist(sapply(ans[,2], function(x) which(DB$id==x)))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
263 #ans[,2] <- as.character(DB$name[nid])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
264 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
265 unk <- reactionM[reactionM[,5]=="unknown",]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
266 ans1 <- matrix("", nrow=nrow(unk), ncol=7+length(coords))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
267 ans1[,1] <- unk[,3]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
268 ans1[,2] <- unk[,5]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
269 for(j in 1:nrow(ans1)) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
270 i <- j + max(as.numeric(ans[,6]),na.rm=TRUE)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
271 if ( i >= 100 ) { ans1[j,6]=i }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
272 else { if ( i >= 10 ) { ans1[j,6]=paste(0,i, sep="") } else { ans1[j,6]=paste("00",i, sep="") } }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
273 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
274 # this step try to recover ids of ion annotation for masses without database annotation
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
275 idx2 <- c(); #for(m in 1:nrow(allion)) if(sum(allion[m,2]==as.numeric(unk[,2])) & sum(allion[m,1]==as.numeric(unk[,3]))) idx2 <- append(idx2, m)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
276 # temp changes made 03/03/2014 have to check carefuly
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
277 lidx <- lapply(1:nrow(allion), function(m) which(allion[m,2]==unk[,2] & allion[m,1]==unk[,3]))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
278 idx2 <- which(lapply(lidx, length)>0)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
279
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
280 if(!is.null(molIon$cameraobj)) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
281 idx <- as.vector(unlist(sapply(allion[idx2,"trace"],
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
282 function(x) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
283 x <- as.matrix(x)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
284 raw <- strsplit(x,";")[[1]]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
285 mraw <- apply(peaklist[raw, coords], 1, mean)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
286 raw[which.max(mraw)]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
287 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
288
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
289 )
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
290 )
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
291 )
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
292 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
293
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
294 else {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
295 # don't know what happened here with apply
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
296 idx <- c()
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
297
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
298 for(i in 1:nrow( allion[idx2,c("trace","comb")])) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
299 x <- as.matrix(allion[idx2,c("trace","comb")][i,])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
300 raw <- strsplit(x[1],";")[[1]]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
301 if(x[2]!="neg"){
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
302 mraw <- apply(peaklistP[raw, coordsP], 1, mean, na.rm=TRUE)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
303 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
304 else {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
305
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
306 mraw <- apply(peaklistN[raw, coordsN], 1, mean, na.rm=TRUE)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
307 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
308 idx <- c(idx, raw[which.max(mraw)])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
309 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
310
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
311
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
312
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
313 tmpidx <- cbind(idx,allion[idx2,"comb"])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
314 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
315 if(!is.null(molIon$cameraobj)) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
316 ans1[,7] <- apply(rpeaklist[idx,], 1, function(x) paste(x[c(1,2,3,4)], collapse="#"))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
317 ans1[,8:ncol(ans1)] <- as.matrix(peaklist[idx, coords])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
318 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
319 else {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
320 idxP <- tmpidx[tmpidx[,2]!="neg",1]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
321 ans1[1:length(idxP),7] <- apply(rpeaklistP[idxP,], 1, function(x) paste(x[c(1,2,3,4)], collapse="#"))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
322 ans1[1:length(idxP),8:ncol(ans1)] <- as.matrix(peaklistP[idxP, coordsP])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
323 idxN <- tmpidx[tmpidx[,2]=="neg",1]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
324 ans1[(length(idxP)+1):nrow(ans1),7] <- apply(rpeaklistN[idxN,], 1, function(x) paste(x[c(1,2,3,4)], collapse="#"))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
325 ans1[(length(idxP)+1):nrow(ans1),8:ncol(ans1)] <- as.matrix(peaklistN[idxN, coordsN])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
326 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
327 ans <- rbind(ans, as.matrix(ans1))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
328
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
329 if(!is.null(molIon$cameraobj)) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
330 figidx <- c(figidx,idx)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
331 figidx <- as.numeric(figidx[-1])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
332 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
333 else {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
334 figidx <- rbind(figidx,tmpidx)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
335 allidx <- figidx[-1,]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
336 allidx <- cbind(allidx, ans[ans[,6]!="",6])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
337 figidx <- as.numeric(figidx[-1,1])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
338 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
339
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
340
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
341 if(m.test=="none") {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
342 testname <- "none"
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
343 #testname <- "Formula"
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
344 #ans[ans[,2]!="unknown",][,5] <- as.character(DB$formula[nid])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
345 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
346 if(m.test=="t.test") {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
347 normalize.medFC <- function(mat) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
348 # Perform median fold change normalisation
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
349 # X - data set [Variables & Samples]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
350 medSam <- apply(mat, 1, median)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
351 medSam[which(medSam==0)] <- 0.0001
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
352 mat <- apply(mat, 2, function(mat, medSam){
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
353 medFDiSmpl <- mat/medSam
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
354 vec<-mat/median(medFDiSmpl)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
355 return(vec)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
356 }, medSam)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
357 return (mat)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
358 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
359 # this piece of code was copied from xcms
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
360 pval <- function(X, classlabel, teststat) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
361
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
362 n1 <- rowSums(!is.na(X[,classlabel == 0]))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
363 n2 <- rowSums(!is.na(X[,classlabel == 1]))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
364 A <- apply(X[,classlabel == 0], 1, sd, na.rm=TRUE)^2/n1 ## sd(t(X[,classlabel == 0]), na.rm = TRUE)^2/n1
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
365 B <- apply(X[,classlabel == 1], 1, sd, na.rm=TRUE)^2/n2 ## sd(t(X[,classlabel == 1]), na.rm = TRUE)^2/n2
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
366 df <- (A+B)^2/(A^2/(n1-1)+B^2/(n2-1))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
367
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
368 pvalue <- 2 * (1 - pt(abs(teststat), df))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
369 invisible(pvalue)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
370 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
371
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
372 c1 <- grep(class1, molIon$cameraobj@xcmsSet@phenoData[,1])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
373 c2 <- grep(class2, molIon$cameraobj@xcmsSet@phenoData[,1])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
374 testclab <- c(rep(0, length(c1)), rep(1, length(c2)))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
375 testval <- groupval(molIon$cameraobj@xcmsSet, "medret", "into")
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
376 if(norm) testval <- normalize.medFC(testval)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
377 tstat <- mt.teststat(testval, testclab)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
378 pvalue <- pval(testval, testclab, tstat)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
379
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
380 #
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
381 # rport <- diffreport(molIon$cameraobj@xcmsSet, class1=class1, class2= class2, sortpval=FALSE)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
382 # ans[ans[,6]!="",5] <- rport[figidx, "pvalue"]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
383 ans[ans[,6]!="",5] <- pvalue[figidx]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
384 testname <- "t.test p-value"
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
385 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
386 if(m.test=="anova"){
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
387 class <- molIon$cameraobj@xcmsSet@phenoData
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
388 getPvalue <- function(dataidx) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
389 aov.data <- data.frame(resp=as.numeric(peaklist[dataidx,coords]), class=class)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
390 anova(aov(resp~class, aov.data))$Pr[1]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
391 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
392 testname <- "anova p-value"
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
393 ans[ans[,6]!="",5] <- sapply(figidx, getPvalue)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
394 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
395
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
396 header <- matrix(c("Proposed Mass","Most probable Compound","Probability","Entropy", testname,"EIC-plot", "Ion annotation",coords), nrow=1 , ncol=7+length(coords) )
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
397 ans <- rbind(header, ans)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
398
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
399
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
400 # additional field
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
401 # ans <- cbind(ans[,1:2], ans[,2], ans[,3:ncol(ans)])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
402 #ans[ans[,3]!="unknown",][-1,3] <- as.character(DB$sbml.id[nid])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
403
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
404 if(html) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
405 #require(hwriter)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
406 ansb <- ans
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
407 ans[ans[,2]!="unknown",][-1,2] <- as.character(DB$name[nid])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
408 if(linkPattern=="pubchem") ansb <- ans
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
409
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
410 hyper=matrix(paste(linkURL, ansb[-1,2], sep=""),ncol=1 )
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
411 if(!is.null(molIon$cameraobj)) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
412 hyper1=matrix(paste(fig, ans[-1,6],".png", sep=""),ncol=1 )
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
413 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
414 else {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
415 hyper1 <- ans[-1,6]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
416 hyper1[ans[-1,6]!=""][allidx[,2]!="neg"] <- paste(hyper1[ans[-1,6]!=""][allidx[,2]!="neg"], "pos", sep="")
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
417 hyper1[ans[-1,6]!=""][allidx[,2]=="neg"] <- paste(hyper1[ans[-1,6]!=""][allidx[,2]=="neg"], "neg", sep="")
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
418 hyper1=matrix(paste(fig, hyper1,".png", sep=""),ncol=1 )
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
419 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
420 p=openPage(paste(filename,".html",sep=""))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
421 ans2 <- ans[,1:7]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
422 link <- cbind(matrix(NA,nrow(ans2),1),rbind(NA,hyper),matrix(NA,nrow(ans2),3),rbind(NA,hyper1),matrix(NA,nrow(ans2),1))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
423 # This block is responsible to add as many columns as the user
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
424 # wants
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
425 if(!is.null(addLink)){
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
426 for(l in 1:length(addLink$link)) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
427 link <- cbind(link, rbind(NA, addLink[[l]]))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
428 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
429 for(a in 1:length(addLink$ans)) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
430 ans2 <- cbind(ans2,addLink$ans[[a]])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
431 #colnames(ans2)[7+a] <- addLink$ans[[a]][1]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
432 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
433 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
434 hwrite(ans2, p,row.bgcolor='#ffdc98', link=link )
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
435 closePage(p)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
436 if(!is.null(molIon$cameraobj)) {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
437 plotEIC(molIon$cameraobj@xcmsSet, figidx, ans[ans[,6]!="",6][-1], colorplot=colorplot)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
438 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
439 else {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
440 dataidxP <- as.numeric(allidx[allidx[,2]!="neg",1])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
441 pngidxP <- allidx[allidx[,2]!="neg",3]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
442 plotEIC(molIon$pos@xcmsSet, dataidxP, pngidxP, "pos", colorplot=colorplot)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
443 dataidxN <- as.numeric(allidx[allidx[,2]=="neg",1])
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
444 pngidxN <- allidx[allidx[,2]=="neg",3]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
445 plotEIC(molIon$neg@xcmsSet, dataidxN, pngidxN, "neg", colorplot=colorplot)
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
446 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
447
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
448 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
449 else {
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
450 ansb <- ans
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
451 }
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
452 colnames(ansb) <- ansb[1,]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
453 ansb <- ansb[-1,]
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
454 return(list(classTable=ansb, figidx=figidx))
e13ec2c3fabe planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
mmonsoor
parents:
diff changeset
455 }