In [75]:
sessionInfo()
In [8]:
# we load a few packages
library(ggplot2)
library(corrplot)
library(reshape2)
library(caret)
library(MASS) # for LDA
library(klaR) # for pls
library(pls)
library(e1071)# for pls
library(pROC) # for pls
library(vegan)
require(scales)
require(gridExtra);
In [77]:
options(repr.plot.width=4, repr.plot.height=4)
In [78]:
#rm(list=ls(all=TRUE))
#dataset <- read.csv('/Users/jul/git/pipe-generate-dataset/dataset.csv', header=FALSE)
#dataClass <- read.csv('/Users/jul/git/pipe-generate-dataset/dataClass.csv', header=FALSE)
dataset <- read.csv2("/Users/jul/git/pipe-generate-dataset/dataset.csv",header=TRUE,sep=",",dec=".")
dataClass <- read.csv2("/Users/jul/git/pipe-generate-dataset/classMatrix.csv",header=TRUE,sep=",",dec=".")
In [79]:
dim(dataset)
In [80]:
colnames(dataset) <- seq(1:768)
In [81]:
data <- data.frame("x"=I(dataset), "y"=as.factor(dataClass[,1]))
In [82]:
filter <- data$x[1,] > 1e9
In [89]:
options(repr.plot.width=8, repr.plot.height=4)
matplot(t(stdev(data$x,1)), type='l')
In [102]:
dim(data$x)
In [88]:
plot(as.numeric(data$x[1, filter]), type='l', ylim=range(data$x[1,]))
In [68]:
pca <- prcomp(data$x[, data$x[1,] > 1e9], scale=FALSE)
In [69]:
options(repr.plot.width=8, repr.plot.height=4)
pc <- c(1,2)
COLOR <- c(2:4)
ev <- c(round(pca$sdev[pc[1]]/sum(pca$sdev)*100,0),
round(pca$sdev[pc[2]]/sum(pca$sdev)*100,0),
round(pca$sdev[pc[3]]/sum(pca$sdev)*100,0))
plot(pca$x[,pc[1]], pca$x[,pc[2]], col=COLOR[data$y], cex=0.7,
xlab=paste0("PC ", pc[1], " (", ev[pc[1]], "%)"),
ylab=paste0("PC ", pc[2], " (", ev[pc[2]], "%)"))
legend("topright", legend=levels(data$y), col=COLOR, pch=1)
In [62]:
op <- par(mar = c(4,4,1,1))
plot(pca$rotation[,1], type='l', ylim=range(pca$rotation[,1:3]))
lines(pca$rotation[,2], type='l', col=2)
lines(pca$rotation[,3], type='l', col=3)
par(op)
In [39]:
plot(pca$rotation[,1], pca$rotation[,2])
In [23]:
options(repr.plot.width=5, repr.plot.height=4)
qplot(x=Var1, y=Var2, data=melt(cor(data$x)), fill=value, geom="tile")
In [35]:
options(repr.plot.width=8, repr.plot.height=4)
layout(matrix(c(1,2,3), 1, 3, byrow = TRUE),
widths=c(3,3,3),
heights=c(1))
op <- par(mar = c(4,4,0,0))
plot(pca$rotation[,1],pca$rotation[,2]);
plot(pca$rotation[,1],pca$rotation[,3]);
plot(pca$rotation[,2],pca$rotation[,3])
par(op)
In [45]:
X<-as.matrix(dataset)
Y<-as.numeric(dataClass[,1])
In [46]:
nf = dim(X)[1]
T = c()
P = c()
C = c()
W = c()
In [47]:
for (j in 1:nf) {
w = (t(X) %*% Y) %*% solve(t(Y) %*% Y)
w1 = t(w) %*% w
w2 = abs(sqrt(w1))
w = w %*% solve(w2)
t = (X %*% w)
t1 = t(t) %*% t
c = t(Y) %*% t %*% solve(t1)
c1 = t(c) %*% c
u = Y %*% c %*% solve(c1)
u1 = t(u) %*% u
u2 = abs(sqrt(u1))
p = (t(X) %*% t) %*% solve(t1)
X = X - (t)%*%t(p)
T = matrix(c(T, t))
P = matrix(c(P, p))
C = matrix(c(C, c))
W = matrix(c(W, w))
}
In [48]:
T = matrix(T, ncol = nf)
T = scale(T, scale = TRUE, center = TRUE)
P = matrix(P, ncol = nf)
C = matrix(C, ncol = nf)
W = matrix(W, ncol = nf)
In [49]:
plot(T[,1], T[,2], col = factor(dataClass[,1]))
#text(T[,1], T[,2], dataClass[,1], pos = 2)
In [ ]:
In [50]:
nf=15
T = c()
P = c()
C = c()
W = c()
Tortho = c()
Portho = c()
Wortho = c()
Cortho = c()
In [51]:
for (j in 1:nf) {
w = (t(X) %*% Y) %*% solve(t(Y) %*% Y)
w=w%*%solve(norm(w))
t = (X %*% w) %*% solve(t(w) %*% w)
c = t(Y) %*% t %*% solve(t(t) %*% t)
u = Y %*% c %*% solve(t(c) %*% c)
p = (t(X) %*% t) %*% solve(t(t) %*% t)
wortho = p - w
wortho= p - t(((t(w)%*%p)%*% solve(t(w)%*%w))%*% t(w))
wortho= wortho%*%solve(norm(wortho))
tortho = X %*% wortho %*% solve(t(wortho) %*% wortho)
portho = t(X) %*% tortho %*% solve(t(tortho) %*% tortho)
cortho = t(Y) %*% tortho %*% solve(t(tortho) %*% tortho)
X = X - tortho %*% t(portho)
T = matrix(c(T, t))
P = matrix(c(P, p))
C = matrix(c(C, c))
W = matrix(c(W, w))
Tortho = matrix(c(Tortho, tortho))
Portho = matrix(c(Portho, portho))
Wortho = matrix(c(Wortho, wortho))
Cortho = matrix(c(Cortho, cortho))
}
In [52]:
T = matrix(T, ncol = nf)
T = scale(T, scale = TRUE, center = TRUE)
P = matrix(P, ncol = nf)
C = matrix(C, ncol = nf)
W = matrix(W, ncol = nf)
Tortho = matrix(Tortho, ncol = nf)
Tortho = scale(Tortho, scale = TRUE, center = TRUE)
Portho = matrix(Portho, ncol = nf)
Wortho = matrix(Wortho, ncol = nf)
Xortho = Tortho %*% t(Portho)
In [54]:
plot(T[,1], Tortho[, 1], pch = 19, col=factor(dataClass[,1]))
#text(T[, 1], Tortho[,1],dataClass[,1], cex=1, pos=2)
In [ ]:
In [18]:
opls = function(data = data, nf) {
X <- data$X
Y <- as.numeric(data$Y)
nf <- 15
T = c()
P = c()
C = c()
W = c()
Tortho = c()
Portho = c()
Wortho = c()
Cortho = c()
for (j in 1:nf) {
w = (t(X) %*% Y) %*% solve(t(Y) %*% Y)
w=w%*%solve(norm(w))
t = (X %*% w) %*% solve(t(w) %*% w)
c = t(Y) %*% t %*% solve(t(t) %*% t)
u = Y %*% c %*% solve(t(c) %*% c)
p = (t(X) %*% t) %*% solve(t(t) %*% t)
wortho = p - w
wortho= p - t(((t(w)%*%p)%*% solve(t(w)%*%w))%*% t(w))
wortho= wortho%*%solve(norm(wortho))
tortho = X %*% wortho %*% solve(t(wortho) %*% wortho)
portho = t(X) %*% tortho %*% solve(t(tortho) %*% tortho)
cortho = t(Y) %*% tortho %*% solve(t(tortho) %*% tortho)
X = X - tortho %*% t(portho)
T = matrix(c(T, t))
P = matrix(c(P, p))
C = matrix(c(C, c))
W = matrix(c(W, w))
Tortho = matrix(c(Tortho, tortho))
Portho = matrix(c(Portho, portho))
Wortho = matrix(c(Wortho, wortho))
Cortho = matrix(c(Cortho, cortho))
}
T = matrix(T, ncol = nf)
T = scale(T, scale = TRUE, center = TRUE)
P = matrix(P, ncol = nf)
C = matrix(C, ncol = nf)
W = matrix(W, ncol = nf)
Tortho = matrix(Tortho, ncol = nf)
Tortho = scale(Tortho, scale = TRUE, center = TRUE)
Portho = matrix(Portho, ncol = nf)
Wortho = matrix(Wortho, ncol = nf)
Xortho = Tortho %*% t(Portho)
oplsModel <- data.frame(I(T), I(Tortho))
return(oplsModel)
}
In [19]:
m<-matrix(sample(6e6:7e6,140), nrow= 20,ncol =7, byrow="TRUE")
species<-c("1", "1", "1","1", "1", "1","1", "1", "1","1",
"-1","-1","-1","-1","-1","-1","-1","-1","-1","-1")
data <- data.frame(X = I(m), Y = species)
In [20]:
oplsModel <- opls(data, 15)
In [21]:
plot(oplsModel$T[,1], oplsModel$Tortho[, 1], pch = 19, col=factor(data$Y))
text(oplsModel$T[, 1], oplsModel$Tortho[,1],data$Y, cex=1, pos=2)
In [22]:
names(oplsModel)
In [ ]: