Nothing
######################################################################################################################
#
# Software that implements
#
# GAUSSIAN PROCESS REGRESSION AND CLASSIFICATION
#
# Copyright (c) 2005 - 2007 by Carl Edward Rasmussen and Chris Williams
#
# Permission is granted for anyone to copy, use, or modify these programs for
# purposes of research or education, provided this copyright notice is retained,
# and note is made of any changes that have been made.
#
# These programs are distributed without any warranty, express or
# implied. As these programs were written for research purposes only, they
# have not been tested to the degree that would be advisable in any
# important application. All use of these programs is entirely at the
# user's own risk.
#
# The code and associated documentation are available from
#
# http://www.GaussianProcess.org/gpml/code
#
######################################################################################################################
#
# Function: covNoise
# This function is an R version of MATLAB implementation (gpml) by Carl Edward Rasmussen and Chris Williams.
# This R version is derived from "gpr" package (http://cran.r-project.org/web/packages/gpr/)
#
######################################################################################################################
covNoise = function(loghyper= NULL , x = NULL , z = NULL , testset.covariances= FALSE){
toprint=FALSE
if (is.null(loghyper) )
{
return(1)
} # report number of parameters
s2 = exp(2*loghyper) # noise variance
if (is.null(z))
{ # compute covariance matrix
A = s2* diag(dim(x)[1]) #I change the code as here it was giving error by %*% (in source code too)
B=0
}else if (testset.covariances== TRUE ) # compute test set covariances
{
A = s2
B = 0 # zeros cross covariance by independence
}else if (testset.covariances== FALSE )
{ # compute derivative matrix
A = 2*s2 * diag(dim(x)[1])
B= 0
}
result= list(A, B)
return (result)
}
######################################################################################################################
#
# Software that implements
#
# GAUSSIAN PROCESS REGRESSION AND CLASSIFICATION
#
# Copyright (c) 2005 - 2007 by Carl Edward Rasmussen and Chris Williams
#
# Permission is granted for anyone to copy, use, or modify these programs for
# purposes of research or education, provided this copyright notice is retained,
# and note is made of any changes that have been made.
#
# These programs are distributed without any warranty, express or
# implied. As these programs were written for research purposes only, they
# have not been tested to the degree that would be advisable in any
# important application. All use of these programs is entirely at the
# user's own risk.
#
# The code and associated documentation are available from
#
# http://www.GaussianProcess.org/gpml/code
#
######################################################################################################################
#
# Function: sq_dist
# This function is an R version of MATLAB implementation (gpml) by Carl Edward Rasmussen and Chris Williams.
# This R version is derived from "gpr" package (http://cran.r-project.org/web/packages/gpr/)
#
######################################################################################################################
sq_dist = function(A , B = NULL){
mym= 0
resut = 0
if( is.null(B)) {
B=A
}
toprint=FALSE
if(dim(A)[1] != dim( B)[1] ){
print("Error: column lengths must agree in sq_dist")
return
}
if(length(A)==1 && length(B)==1){
A = as.vector(A)
B = as.vector(B)
}
# result= ((as.matrix(dist(cbind (t(A[]),t(B[])))))^2)/2
n=dim(A)[2]
m=dim(B)[2]
C= array(0, dim=c(n,m))
if(m==1) { #special case, R automatically turns a column matrix into a one row matrix
for( d in 1:dim(A)[1]){
C= C+t((B[rep(d,n),] - ( t(A[rep(d,m),]) ))^2 )
}
}else{
for( d in 1:dim(A)[1]){
#C = C + (repmat(b(d,:), n, 1) - repmat(a(d,:)', 1, m)).^2
C = C + (B[rep(d,n),] - ( t(A[rep(d,m),]) ))^2
}
}
return (C)
}
######################################################################################################################
#
# Software that implements
#
# GAUSSIAN PROCESS REGRESSION AND CLASSIFICATION
#
# Copyright (c) 2005 - 2007 by Carl Edward Rasmussen and Chris Williams
#
# Permission is granted for anyone to copy, use, or modify these programs for
# purposes of research or education, provided this copyright notice is retained,
# and note is made of any changes that have been made.
#
# These programs are distributed without any warranty, express or
# implied. As these programs were written for research purposes only, they
# have not been tested to the degree that would be advisable in any
# important application. All use of these programs is entirely at the
# user's own risk.
#
# The code and associated documentation are available from
#
# http://www.GaussianProcess.org/gpml/code
#
######################################################################################################################
#
# Function: covSEiso
# This function is an R version of MATLAB implementation (gpml) by Carl Edward Rasmussen and Chris Williams.
# This R version is derived from "gpr" package (http://cran.r-project.org/web/packages/gpr/)
#
######################################################################################################################
covSEiso =function(loghyper = NULL , x = NULL , z = NULL , testset.covariances= FALSE){
A=B=0
toprint=FALSE
if (is.null(loghyper) )
{
return(2)
} # report number of parameters
n = dim(x)[1]
D = dim(x)[2]
ell = exp(loghyper[1]) # characteristic length scale
sf2 = exp(2*loghyper[2]) # signal variance
A=B= array()
if (is.null(z)){
# as.matrix(dist(cbind(x,y)))
A = sf2*exp(-sq_dist(t(x)/ell)/2)
B=0
}else if (testset.covariances== TRUE ) # compute test set covariances
{
if (is.null(dim(z) ))
{
matlab.dim = length(z)
}else{
matlab.dim = dim(z)[1]
}
A=sf2*rep(1, matlab.dim)
B = sf2*exp(-sq_dist(t(x)/ell,t(z)/ell)/2)
}else if (testset.covariances== FALSE )
{ #compute derivative matrix
if (length(z)== 1 && z == 1)#dim returns null if there is one number, different from matlab
{
#-> # first parameter
A = sf2*exp(-sq_dist(t(x)/ell)/2)*sq_dist(t(x)/ell) # A = sf2*exp(-sq_dist(x'/ell)/2).*sq_dist(x'/ell);
}else{ #second parameter
A = 2*sf2*exp(-sq_dist(t(x)/ell)/2)
}
B= 0
}
result = list(A,B)
return (result)
}
######################################################################################################################
#
# Software that implements
#
# GAUSSIAN PROCESS REGRESSION AND CLASSIFICATION
#
# Copyright (c) 2005 - 2007 by Carl Edward Rasmussen and Chris Williams
#
# Permission is granted for anyone to copy, use, or modify these programs for
# purposes of research or education, provided this copyright notice is retained,
# and note is made of any changes that have been made.
#
# These programs are distributed without any warranty, express or
# implied. As these programs were written for research purposes only, they
# have not been tested to the degree that would be advisable in any
# important application. All use of these programs is entirely at the
# user's own risk.
#
# The code and associated documentation are available from
#
# http://www.GaussianProcess.org/gpml/code
#
######################################################################################################################
#
# Function: covSum
# This function is an R version of MATLAB implementation (gpml) by Carl Edward Rasmussen and Chris Williams.
# This R version is derived from "gpr" package (http://cran.r-project.org/web/packages/gpr/)
#
######################################################################################################################
covSum= function(covfuncsum , logtheta = NULL, x = NULL, z = NULL, testset.covariances= FALSE){
A=B=0
toprint= FALSE
temp1=strsplit(covfuncsum , ",")
covfuncs1 =temp1[[1]][1]
covfuncs2 =temp1[[1]][2]
covarray = as.vector(c(covfuncs1, covfuncs2 ))
j=as.array( c(eval(call(covfuncs1) ) , eval(call(covfuncs2 ))) ) #j stores number of parameters: j=[2 1]
v = NULL
for (i in 1:length(j)){
v= cbind( v, array(rep(i,j[i] ), dim= c(1,j[i])) )
}
#v is parameter mapper array.
#v= [1,1,2] ,means for the fist 2 parameter of loghyper must be taken for the fist function
# and third parameter of loghyper must be used for the second function
if (is.null(logtheta))
{ # report number of parameters
A = eval(call(covfuncs1) ) + eval(call(covfuncs2 ))
B=0
}else{
n = dim(x)[1]
D = dim(x)[2]
# v = array() # v vector indicates to which covariance parameters belong
if(is.null(z)){ #nargin case==3 # compute covariance matrix
A1 = eval(call(covfuncs1, logtheta[v==1], x)) [[1]]
A2= eval(call(covfuncs2, logtheta[v==2], x))[[1]]
A=A1+A2
B=0
}else{ #nargin case==4
if (testset.covariances == TRUE){ #nargout==2
ans1 = eval(call(covfuncs1, logtheta[v==1], x, z, testset.covariances))
ans2= eval(call(covfuncs2, logtheta[v==2], x, z, testset.covariances))
A = ans1[[1]] [1]+ ans2[[1]]
B = ans1[[2] ]+ans2[[2]]
}else if (testset.covariances == FALSE){
i=v[z]
j=sum(v[1:z]==i)
f = covarray[i]
#print(paste("logtheta[v==i] ,f,z,i,j:",logtheta[v==i] ,f,z,i,j,sep="; "))
A= eval(call(f, logtheta[v==i], x, j, testset.covariances)) [[1]]
B=0
}
}
}
result = list(A,B)
return( result)
}
######################################################################################################################
#
# Software that implements
#
# GAUSSIAN PROCESS REGRESSION AND CLASSIFICATION
#
# Copyright (c) 2005 - 2007 by Carl Edward Rasmussen and Chris Williams
#
# Permission is granted for anyone to copy, use, or modify these programs for
# purposes of research or education, provided this copyright notice is retained,
# and note is made of any changes that have been made.
#
# These programs are distributed without any warranty, express or
# implied. As these programs were written for research purposes only, they
# have not been tested to the degree that would be advisable in any
# important application. All use of these programs is entirely at the
# user's own risk.
#
# The code and associated documentation are available from
#
# http://www.GaussianProcess.org/gpml/code
#
######################################################################################################################
#
# Function: gpr
# This function is an R version of MATLAB implementation (gpml) by Carl Edward Rasmussen and Chris Williams.
# This R version is derived from "gpr" package (http://cran.r-project.org/web/packages/gpr/)
#
######################################################################################################################
gpr = function(logtheta, covfunc.gpr, x, y, xstar = NULL, partial.derivatives = FALSE){
#toprint = TRUE
#if(toprint)print("in gprt to compare with matlab \n")
#if(toprint)print(paste(logtheta[1],covfunc.gpr,x[1,1],y[1],xstar[1,1],sep="; "))
n= dim(x)[1]
D = dim(x)[2]
tempF1 = strsplit(covfunc.gpr, ",")
covfunc1 = tempF1[[1]][1]
covfunc23 = paste(tempF1[[1]][2] , tempF1[[1]][3], sep=",")
if (length(tempF1[[1]])==3)
{
mo= (eval(call( covfunc1, covfunc23)))
if (mo[1] != dim(logtheta)[1])
{
print("Error: Number of parameters do not agree with covariance function")
return (-1)
}
K = eval(call(covfunc1, covfunc23, logtheta , x ))[[1]]
}
else
{
K = eval(call(covfunc.gpr, logtheta, x)) [[1]]
}
L = t(chol(K))
alpha = solve(t(L),solve(L,y))
if (is.null(xstar)) #nargin ==4
{
out1 = 0.5*t(y)%*%alpha + sum(log(diag(L))) + 0.5*n*log(2*pi)
out2=0
if( partial.derivatives == TRUE)
{
out2 = rep(0,length(logtheta)) # set the size of the derivative vector
W = solve( t(L) , (solve(L, diag(n)) ) ) - ( alpha %*% t(alpha) ) # precompute for convenience
for (i in 1:length(out2))
{
if (length(tempF1[[1]])==1)
{
out2[i] = sum(W*eval(call(covfunc.gpr, logtheta, x, i , FALSE)) [[1]] )/2
}
if (length(tempF1[[1]])==3)
{
out2[i] = sum(W*eval(call(covfunc1, covfunc23, logtheta , x, i , FALSE))[[1]])/2
}
}
}# end partial derivative
}
else if(!is.null(xstar))
{
K = eval(call(covfunc1, covfunc23, logtheta , x, xstar ,TRUE ))
Kss= K[[1]]
Kstar = K[[2]]
out1 = t(Kstar) %*% alpha # predicted means
out2 = 0
if(partial.derivatives == TRUE)
{
if( !is.null(dim(x)))
{ #can not say if it is null, because if it is an array it gives an error
v = solve(L,Kstar)#dim v: 27,51
out2 = Kss - t(colSums(v * v))
}
}# end partial derivative
}
mout1 = out1
mout2 = t(out2)
result = list(mout1, mout2)
return(result)
}
######################################################################################################################
#
# Software that implements
#
# GAUSSIAN PROCESS REGRESSION AND CLASSIFICATION
#
# Copyright (c) 2005 - 2007 by Carl Edward Rasmussen and Chris Williams
#
# Permission is granted for anyone to copy, use, or modify these programs for
# purposes of research or education, provided this copyright notice is retained,
# and note is made of any changes that have been made.
#
# These programs are distributed without any warranty, express or
# implied. As these programs were written for research purposes only, they
# have not been tested to the degree that would be advisable in any
# important application. All use of these programs is entirely at the
# user's own risk.
#
# The code and associated documentation are available from
#
# http://www.GaussianProcess.org/gpml/code
#
######################################################################################################################
#
# Function: minimize
# This function is an R version of MATLAB implementation (gpml) by Carl Edward Rasmussen and Chris Williams.
# This R version is derived from "gpr" package (http://cran.r-project.org/web/packages/gpr/)
#
######################################################################################################################
minimize = function (X, f, .length, covfunc, x, y){
INT = 0.1
EXT = 3.0
MAX = 20
RATIO = 10
SIG = 0.1
RHO = SIG/2
if(is.array( .length)){
if (max(dim( .length)) == 2)
{
red= .length[2]
.length= .length[1]
}
} else{
red=1
}
if ( .length>0)
{
S='Linesearch'
} else{
S='Function evaluation'
}
i.m= 0 # zero the run length counter
ls_failed = 0 # no previous line search has failed
f_out =eval(call (f, X, covfunc , x, y, NULL , TRUE))
#print(paste("fout", f_out,sep="|"));
f0= f_out[1][[1]]
df0 = t(f_out[2][[1]]) #out put is a colum, R makes it a row when put it in list to output, I conver it back here
fX = f0
i.m = i.m + ( .length<0) # count epochs?!
s = -df0
d0 = -t(s)%*%s # initial search direction (steepest) and slope
x3 = red/(1-d0) # initial step is red/(|s|+1)
mainloop = TRUE
while( i.m < abs( .length) && mainloop)
{ # while not finished
i.m = i.m + (.length > 0) # count iterations?!
X0 = X
F0 = f0
dF0 = df0 # make a copy of current values
if ( .length>0)
M = MAX
else
M = min(MAX, - .length- i.m)
whilerun = TRUE
f3_=c(0)
df3_=c(0)
while (whilerun ==TRUE) # keep extrapolating as long as necessary
{
x2 = 0
f2 = f0
d2 = d0
f3 = f0
df3 = df0
success = FALSE
while (success == FALSE && M > 0)
{
M = M - 1
i.m = i.m + (.length<0) #count epochs?!
options(show.error.messages = FALSE)
f_out2 =eval(call (f , X+ x3[1]*s, covfunc , x, y, NULL, TRUE))
f3= f_out2[1][[ 1 ]][[ 1 ]]
df3 = t(f_out2[2][[1]])
f3_=rbind(f3_,f3,0);
df3_=rbind(df3_,df3,0);
if (is.na(f3) || is.infinite(f3) || is.nan(f3) || any(is.nan(df3) || is.na(df3) || is.infinite(df3)) )
{
cat( " ")
x3 = (x2+x3)/2
}
else
{
success = TRUE
}
options(show.error.messages = TRUE)
}
if (f3 < F0)
{
X0 = X+x3[1]*s
F0 = f3
dF0 = df3
} # keep best values
d3 = t(df3)%*%s # new slope
if( d3 > SIG*d0 || f3 > f0+x3*RHO*d0 || M == 0 )# are we done extrapolating?
{
whilerun = FALSE
break
}
x1 = x2
f1 = f2
d1 = d2 # move point 2 to point 1
x2 = x3
f2 = f3
d2 = d3 # move point 3 to point 2
A = 6*(f1-f2)+3*(d2+d1)*(x2-x1) # make cubic extrapolation
B = 3*(f2-f1)-(2*d1+d2)*(x2-x1)
x3 = x1-d1*(x2-x1)^2/(B+sqrt(abs(B*B-A*d1*(x2-x1)))) # num. error possible, ok!
if ( (B*B-A*d1*(x2-x1) < 0)[1] || is.nan(x3) || is.infinite(x3) || x3 < 0) # num prob | wrong sign?
{
x3 = x2*EXT # extrapolate maximum amount
}else if( x3 > x2*EXT) # new point beyond extrapolation limit?
{
x3 = x2*EXT # extrapolate maximum amount
}else if (x3 < x2+INT*(x2-x1) ) # new point too close to previous point?
{
x3 = x2+INT*(x2-x1)
}
#x3= round(x3*10000)/10000 #this is just to make same answers with matlab code
}
#print(f3_)
#print(df3_)
#ANSWER <- readline("continue? ")
#if (substr(ANSWER, 1, 1) == "n")
# stop();
while ((abs(d3) > -SIG*d0 || f3 > f0+x3*RHO*d0) && M > 0 ) # keep interpolating
{
if( d3 > 0 || f3 > f0+x3*RHO*d0) # choose subinterval
{
x4 = x3
f4 = f3
d4 = d3 #move point 3 to point 4
}else{
x2 = x3
f2 = f3
d2 = d3 # move point 3 to point 2
}
if (f4 > f0 ){
x3 = x2-(0.5*d2*(x4-x2)^2)/(f4-f2-d2*(x4-x2)) # quadratic interpolation
}else{
A = 6*(f2-f4)/(x4-x2)+3*(d4+d2) # cubic interpolation
B = 3*(f4-f2)-(2*d2+d4)*(x4-x2)
x3 = x2+(sqrt(B*B-A*d2*(x4-x2)^2)-B)/A # num. error possible, ok!
}
if (is.nan(x3) || is.infinite(x3) )
{
x3 = (x2+x4)/2 # if we had a numerical problem then bisect
}
x3 = max(min(x3, x4-INT*(x4-x2)),x2+INT*(x4-x2)) # don't accept too close
f_out3 =eval(call (f , X+ x3*s, covfunc , x, y, NULL, TRUE))
f3 = f_out3[1] [[1]][[1]]
df3 = t(f_out3[2][[1]])
if (f3 < F0)
{
x3=x3[[1]] # arrays of one in matlab converst to number, in R it is not
X0 = X+x3*s
F0 = f3
dF0 = df3
}# keep best values
M = M - 1
i.m = i.m + (.length<0) # count epochs?!
d3 = t(df3)%*%s # new slope
}#end while #end interpolation
if (abs(d3) < -SIG*d0 && f3 < f0+x3*RHO*d0 ) # if line search succeeded
{
x3=x3[[1]]
X = X+x3*s
f0 = f3
fX= t(cbind (t(fX), f0))
s=( ( (t(df3)%*%df3- t(df0)%*%(df3))[1] )/ ( (t(df0)%*%(df0))[[1]] ) *s ) - df3
df0 = df3 # swap derivatives
d3 = d0
d0 = t(df0)%*%s
if (d0 > 0 )
{ # new slope must be negative
s = -df0
d0 = -t(s)%*%s # otherwise use steepest direction
}
x3 = x3 * min(RATIO, d3/(d0- ( 2^(-1022) ) )) # slope ratio but max RATIO
ls_failed = 0 # this line search did not fail
} else{
X = X0
f0 = F0
df0 = dF0 # restore best point so far
if (ls_failed || i.m > abs(.length) ) # line search failed twice in a row
{
mainloop = 0 #break; #or we ran out of time, so we give up
break
}
s = -df0
d0 = -t(s)%*%s #try steepest
x3 = 1/(1-d0)
ls_failed = 1 # this line search failed
}#end if
}#end first while
return (list(X)) #i.m is number of epochs
}
########################################################
# #
# Function: normalize #
# Author: Mohsen Ahmadi #
# Email: mohsen_ahmadi989@yahoo.com #
# Date: Feb 2014 #
# #
########################################################
normalize <- function(x){
mu = colMeans(x);
x_norm = sweep(x, length(dim(x)) , mu, FUN="-")
if (class(x) == "data.frame")
{
sigma = sapply(x_norm, sd)
}
else
{
sigma = sd(x_norm)
}
sigma[sigma==0]=1
norm = sweep(x_norm, length(dim(x)), sigma, FUN="/" )
return (list(norm, mu, sigma))
}
########################################################
# #
# Function: one_nn_tc #
# Author: Mohsen Ahmadi #
# Email: mohsen_ahmadi989@yahoo.com #
# Date: Feb 2014 #
# #
########################################################
one_nn_tc = function(X1,X2){
index=0
tc=-1
for (k in 1:dim(X2)[1])
{
K=tanimoto_coef(X1,X2[k,])
#CheckTc(K, X1, X2[k,])
if ( K>tc)
{
tc=K
index=k
}
}
return( index )
}
########################################################
# #
# Function: tanimoto_coef #
# Author: Mohsen Ahmadi #
# Email: mohsen_ahmadi989@yahoo.com #
# Date: Feb 2014 #
# #
########################################################
tanimoto_coef = function(X, Y){
X = as.matrix(X)#R automatically converts matrixes to data frame, we have to convert them back
Y = as.matrix(Y)
res = ( t(X) %*% Y ) / (( t(X) %*% X ) + ( t(Y) %*% Y ) - ( t(X) %*% Y ))
return(res)
}
########################################################
# #
# Function: trainChemPC #
# Author: Mohsen Ahmadi #
# Email: mohsen_ahmadi989@yahoo.com #
# Date: Feb 2014 #
# #
########################################################
trainChemPC = function(xTrain = NULL, yTrain = NULL){
options(warn=-1);
if (is.null(xTrain) || is.null(yTrain))
{
print("Error : Input parameters can not be null.")
return;
}
yTrain=log10(yTrain);
xNorm = normalize(xTrain);
yNorm = normalize(yTrain);
x = xNorm[[1]];
y = yNorm[[1]];
#rho = cor(x, y, method='spearman');
p_values = c();
for(j in 1:dim(x)[2])
{
p = cor.test(as.matrix(x[, j]), y, method='spearman')$p.value
p_values = rbind(p_values, p);
}
adj_pVals = p.adjust(p_values, method = "BH");
hVals = as.integer(adj_pVals <= 0.05);
if (sum(hVals) > 0)
{
attr = grep(1, hVals);
}
else
{
min_p = min(adj_pVals);
attr = grep(TRUE, adj_pVals==min_p);
}
x = as.matrix(x[, c(attr)]);
covfunc="covSum,covSEiso,covNoise"
loghyper= array(c(-1,-1,-1), dim=c(3,1))
loghyper = minimize(loghyper, 'gpr', -100, covfunc, x, y)
return (loghyper[[1]])
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.