
#Example 1


library(archetypes)
source('LstepArchetypesMod.R')

library(DDoutlier)
library(MASS)
set.seed(4321)

X=c()
Sigma1=matrix(c(0.1,0,0,0.1),2,2)
X=rbind(X,mvrnorm(n = 60, c(0,0), Sigma1))


Sigma1=matrix(c(0.1,0,0,0.1),2,2)
X=rbind(X,mvrnorm(n = 35, c(5,3), Sigma1))

x=runif(5,min=0,max=3)
y=runif(5,min=1,max=5)

X=rbind(X,cbind(x,y))

label=as.factor(c(rep(0,95),rep(1,5)))
plot(X)

#k-nn
pe=KNN_SUM(X, k=5) 

library(AUC)
rr=roc(pe, label)
auc(rr)

b=boxplot(pe)

oo=which(is.element(pe,b$out))


sum(label== as.integer(is.element(pe,b$out)))


table(Predictions = as.integer(is.element(pe,b$out)), TrueLabels = label)



library(archetypes)
nada=6
norep=20
set.seed(4321)
ai=LstepArchetypesMod(X, 1:nada, nrep=norep)

 screeplot(ai) #elbow 2;



e=2
nada=6
i=5
pei=rep(0,dim(X)[1])
for (j in e:nada){

a=bestModel(ai[[j]])

pei=pei+KNN_SUM(a$alphas, k=i) 
}

rr=roc(pei, label)
auc(rr)

pei

b=boxplot(pei)

oo=which(is.element(pei,b$out))


sum(label== as.integer(is.element(pei,b$out)))


table(Predictions = as.integer(is.element(pei,b$out)), TrueLabels = label)


library(vcd)
a=bestModel(ai[[3]])

ternaryplot(a$alphas,dimnames=1:3,col=c(rep('black',95),rep('red',5)),cex=0.65,pch=c(rep(19,95),rep(17,5)),main='')


library(ggplot2)
Xd=as.data.frame(X)
ggplot(Xd, aes(x, y)) + geom_point(colour = "red", size = pei)
#ggplot(Xd, aes(x, y)) + geom_point(colour = "red", size = pei,shape=label)
ccc=rep('green',length(label))
ccc[label==1]='red'
 ggplot(Xd, aes(x, y)) + geom_point(colour = ccc, size = pei,shape=as.numeric(label)+19) +theme_bw()

Xd2=rbind(Xd,a$archetypes)
 ggplot(Xd2, aes(x, y)) + geom_point(colour = c(ccc,rep('black',3)), size = c(pei,rep(1,3)),shape=c(as.numeric(label)+19,rep(4,3))) +theme_bw()


######################

#Example 2



x=c(seq(2,20,2),5)
y=c(seq(1,10,1),5)


label=as.factor(c(rep(0,10),1))

X=cbind(x,y)

 hpts <- chull(X)
#> hpts
#[1]  1 11 10


library(archetypes)
nada=10
norep=20
set.seed(1)
ai=LstepArchetypesMod(X, 1:nada, nrep=norep)

 screeplot(ai) #elbow 2;


#ensemble
e=2
nada=3

auci=c()

for (i in 1:6)
{
pei=rep(0,dim(X)[1])
for (j in e:nada){

a=bestModel(ai[[j]])

pei=pei+KNN_SUM(a$alphas, k=i) 

}
auci=c(auci,auc(roc(pei, label)))

}
summary(auci) 
sd(auci)

#auci= 1


i=3

pei=rep(0,dim(X)[1])
for (j in e:nada){

a=bestModel(ai[[j]])

pei=pei+KNN_SUM(a$alphas, k=i) 
}

auc(roc(pei, label)) #1
rank(pei)


library(ggplot2)
Xd=as.data.frame(X)
ccc=rep('green',length(label))
ccc[label==1]='red'
ggplot(Xd, aes(x, y)) + geom_point(colour = ccc, size = pei*2,shape=as.numeric(label)+19) +theme_bw()

#Knn
pei=c()
i=3
pei=KNN_SUM(X, k=i) 
Xd=as.data.frame(X)
ggplot(Xd, aes(x, y)) + geom_point(colour = ccc, size = pei*.5,shape=as.numeric(label)+19) +theme_bw()

auc(roc(pei, label))
#0.3
rank(pei)


#LOF

pei=c()
i=3
pei=LOF(X, k=i) 
Xd=as.data.frame(X)
ggplot(Xd, aes(x, y)) + geom_point(colour = ccc, size = pei*5,shape=as.numeric(label)+19) +theme_bw()
auc(roc(pei, label)) #0
rank(pei)



# COF
pei=c()
i=3
pei=COF(X, k=i) 
Xd=as.data.frame(X)
ggplot(Xd, aes(x, y)) + geom_point(colour = ccc, size = pei*5,shape=as.numeric(label)+19) +theme_bw()
auc(roc(pei, label)) #0.5
 rank(pei)



#LOOP
pei=c()
i=3
pei=LOOP(X, k=i) 
Xd=as.data.frame(X)
ggplot(Xd, aes(x, y)) + geom_point(colour = ccc, size = pei*20,shape=as.numeric(label)+19) +theme_bw()

auc(roc(pei, label)) #0.7
rank(pei)


#PCA robust

pc.rob <- princomp(X, covmat = MASS::cov.rob(X))
i=3
pei=KNN_SUM(pc.rob$scores, k=i) 
Xd=as.data.frame(X)
ggplot(Xd, aes(x, y)) + geom_point(colour = ccc, size = pei,shape=as.numeric(label)+19) +theme_bw()
auc(roc(pei, label)) #0.3
rank(pei)


#RDOS

i=3
pei=RDOS(X, k=i) 
Xd=as.data.frame(X)
ggplot(Xd, aes(x, y)) + geom_point(colour = ccc, size = pei*5,shape=as.numeric(label)+19) +theme_bw()
auc(roc(pei, label)) #0.6
 rank(pei)



#hbos
library('reticulate')
os <- import("os")
source_python('hbos.py')
hbos = HBOS()
hr = hbos$fit_predict(X)
 hr = hbos$fit_predict(data.frame(X))
Xd=as.data.frame(X)
ggplot(Xd, aes(x, y)) + geom_point(colour = ccc, size = hr*5,shape=as.numeric(label)+19) +theme_bw()
auc(roc(hr, label)) #0.2 
 rank(hr)

