...
|
...
|
@@ -23,66 +23,56 @@ threshold="all", lp=NA, full=FALSE)
|
23
|
23
|
|
24
|
24
|
# find the A with the highest magnitude
|
25
|
25
|
Arowmax <- t(apply(Amatrix, 1, function(x) x/max(x)))
|
26
|
|
- #Arowmax <- t(apply(Amatrix, 1, function(x){
|
27
|
|
- # if (mean(x) == 0){ return(rep(0, times=length(x))) }
|
28
|
|
- # else { return(x/max(x)) }
|
29
|
|
- #})) ## this prevents NA values from creeping due to rows with zero means
|
30
|
|
-
|
31
|
|
- pmax<-apply(Amatrix, 1, max)
|
|
26
|
+
|
32
|
27
|
# determine which genes are most associated with each pattern
|
|
28
|
+ sstat<-matrix(NA, nrow=nrow(Amatrix), ncol=ncol(Amatrix),dimnames=dimnames(Amatrix))
|
33
|
29
|
ssranks<-matrix(NA, nrow=nrow(Amatrix), ncol=ncol(Amatrix),dimnames=dimnames(Amatrix))#list()
|
34
|
30
|
ssgenes<-matrix(NA, nrow=nrow(Amatrix), ncol=ncol(Amatrix),dimnames=NULL)
|
35
|
31
|
nP=dim(Amatrix)[2]
|
36
|
32
|
if(!is.na(lp))
|
37
|
33
|
{
|
38
|
|
- if(length(lp)!=dim(Amatrix)[2])
|
39
|
|
- {
|
|
34
|
+ if(length(lp)!=dim(Amatrix)[2]){
|
40
|
35
|
warning("lp length must equal the number of columns of the Amatrix")
|
41
|
36
|
}
|
42
|
|
- for (i in 1:nP)
|
43
|
|
- {
|
44
|
|
- sstat <- apply(Arowmax, 1, function(x) sqrt(t(x-lp)%*%(x-lp)))
|
45
|
|
- ssranks[order(sstat),i] <- 1:length(sstat)
|
46
|
|
- ssgenes[,i]<-names(sort(sstat,decreasing=FALSE,na.last=TRUE))
|
|
37
|
+ for (i in 1:nP){
|
|
38
|
+ sstat[,i] <- apply(Arowmax, 1, function(x) sqrt(t(x-lp)%*%(x-lp)))
|
|
39
|
+ ssranks[order(sstat[,i]),i] <- 1:dim(sstat)[1]
|
|
40
|
+ ssgenes[,i]<-names(sort(sstat[,i],decreasing=FALSE,na.last=TRUE))
|
47
|
41
|
}
|
48
|
42
|
}
|
49
|
43
|
else
|
50
|
44
|
{
|
51
|
|
- for(i in 1:nP)
|
52
|
|
- {
|
|
45
|
+ for(i in 1:nP){
|
53
|
46
|
lp <- rep(0,dim(Amatrix)[2])
|
54
|
47
|
lp[i] <- 1
|
55
|
|
- sstat <- apply(Arowmax, 1, function(x) sqrt(t(x-lp)%*%(x-lp)))
|
56
|
|
- ssranks[order(sstat),i] <- 1:length(sstat)
|
57
|
|
- ssgenes[,i]<-names(sort(sstat,decreasing=FALSE,na.last=TRUE))
|
|
48
|
+ sstat[,i] <- unlist(apply(Arowmax, 1, function(x) sqrt(t(x-lp)%*%(x-lp))))
|
|
49
|
+ ssranks[order(sstat[,i]),i] <- 1:dim(sstat)[1]
|
|
50
|
+ ssgenes[,i]<-names(sort(sstat[,i],decreasing=FALSE,na.last=TRUE))
|
58
|
51
|
}
|
59
|
52
|
}
|
60
|
53
|
|
61
|
|
- if(threshold=="cut")
|
62
|
|
- {
|
|
54
|
+ if(threshold=="cut"){
|
63
|
55
|
geneThresh <- sapply(1:nP,function(x) min(which(ssranks[ssgenes[,x],x] > apply(ssranks[ssgenes[,x],],1,min))))
|
64
|
56
|
ssgenes.th <- sapply(1:nP,function(x) ssgenes[1:geneThresh[x],x])
|
65
|
|
- #geneThresh <- apply(sweep(ssranks,1,t(apply(ssranks, 1, min)),"-"),2,function(x) which(x==0))
|
66
|
|
- #ssgenes.th <- lapply(geneThresh,names)
|
67
|
57
|
}
|
68
|
|
- else if (threshold=="all")
|
69
|
|
- {
|
70
|
|
- pIndx<-apply(ssranks,1,which.min)
|
|
58
|
+ else if (threshold=="all"){
|
|
59
|
+ pIndx<-unlist(apply(sstat,1,which.min))
|
71
|
60
|
gBYp <- list()
|
72
|
61
|
for(i in sort(unique(pIndx))){
|
73
|
|
- gBYp[[i]]<-names(pIndx[pIndx==i])
|
|
62
|
+ #gBYp[[i]]<-names(pIndx[pIndx==i])
|
|
63
|
+ gBYp[[i]]<-sapply(strsplit(names(pIndx[pIndx==i]),"[.]"),function(x) x[[1]][1])
|
|
64
|
+
|
74
|
65
|
}
|
75
|
66
|
ssgenes.th <- lapply(1:max(sort(unique(pIndx))), function(x) {
|
76
|
67
|
ssgenes[which(ssgenes[,x] %in% gBYp[[x]]),x]
|
77
|
68
|
})
|
78
|
69
|
}
|
79
|
|
- else
|
80
|
|
- {
|
|
70
|
+ else{
|
81
|
71
|
stop("Threshold arguement not viable option")
|
82
|
72
|
}
|
83
|
73
|
|
84
|
74
|
if (full)
|
85
|
|
- return(list("PatternMarkers"=ssgenes.th,"PatternRanks"=ssranks))
|
|
75
|
+ return(list("PatternMarkers"=ssgenes.th,"PatternRanks"=ssranks,"PatternMarkerScores"=sstat))
|
86
|
76
|
else
|
87
|
77
|
return("PatternMarkers"=ssgenes.th)
|
88
|
78
|
}
|