utilizaremos la funcion svm para la clasificacion vectorial
set.seed(1)
x=matrix(rnorm (20*2), ncol=2)
y=c(rep(-1,10), rep(1,10))
x[y==1,]=x[y==1,] + 1
plot(x, col=(3-y))
crearemos un data frame con factor codificando a y, la siguiente grafica mostrará donde se parte el espacio
dat=data.frame(x=x, y=as.factor(y))
library(e1071)
svmfit=svm(y∼., data=dat , kernel ="linear", cost=10,scale=FALSE)
plot(svmfit , dat)
Notemos que los valores estan clasificados entre -1 y 1 para cada una de las casificaciones, el vector de indices nos indican las identidades del proceso
svmfit$index
[1] 1 2 5 7 14 16 17
en la descripcion del svm vemos que le costo es de 10 y clasifica 4 para la primera clase y 3 para la segunda (el grupo en parentisis)
summary(svmfit)
Call:
svm(formula = y ~ ., data = dat, kernel = "linear", cost = 10, scale = FALSE)
Parameters:
SVM-Type: C-classification
SVM-Kernel: linear
cost: 10
Number of Support Vectors: 7
( 4 3 )
Number of Classes: 2
Levels:
-1 1
Ahora bien como reduciremos el costo, veremos que se incrementa la cantidad de entidades a operar
svmfit=svm(y∼., data=dat , kernel ="linear", cost =0.1,scale=FALSE)
plot(svmfit , dat)
svmfit$index
[1] 1 2 3 4 5 7 9 10 12 13 14 15 16 17 18 20
la libreria e1071 contiene la funcion tune que correra un cross validation tomemos encuenta que le estamos indicando un kernel lineal al procedimiento
library(e1071)
set.seed(1)
tune.out=tune(svm ,y∼.,data=dat ,kernel ="linear",
ranges=list(cost=c(0.001, 0.01, 0.1, 1,5,10,100) ))
summary (tune.out)
Parameter tuning of ‘svm’:
- sampling method: 10-fold cross validation
- best parameters:
- best performance: 0.1
- Detailed performance results:
NA
la respuesta esta considerando el modelo del mejor costo con valor de .1 , para extraer el modelo
bestmod=tune.out$best.model
summary(bestmod)
Call:
best.tune(method = svm, train.x = y ~ ., data = dat, ranges = list(cost = c(0.001, 0.01, 0.1, 1, 5,
10, 100)), kernel = "linear")
Parameters:
SVM-Type: C-classification
SVM-Kernel: linear
cost: 0.1
Number of Support Vectors: 16
( 8 8 )
Number of Classes: 2
Levels:
-1 1
haremos el tes data basado en los factores esperados, tomemos encuentra que las categorias tendra su matriz de confusion par avalidar el metodo
xtest=matrix(rnorm (20*2) , ncol=2)
ytest=sample (c(-1,1), 20, rep=TRUE)
xtest[ytest==1,]= xtest[ytest==1,] + 1
testdat=data.frame(x=xtest , y=as.factor(ytest))
ypred=predict (bestmod ,testdat)
table(predict =ypred , truth=testdat$y )
truth
predict -1 1
-1 10 5
1 3 2
utilizando svm sobre el modelo de costo ya previamente estimado
svmfit=svm(y∼., data=dat , kernel ="linear", cost =.01,scale=FALSE)
ypred=predict (svmfit ,testdat )
table(predict =ypred , truth=testdat$y )
truth
predict -1 1
-1 11 7
1 2 0
x[y==1,]=x[y==1,]+0.5
plot(x, col=(y+5)/2, pch =19)
ahora bien la linea es a duras penas separable por loq ue operaremos un costo pobre
dat=data.frame(x=x,y=as.factor(y))
svmfit=svm(y∼., data=dat , kernel ="linear", cost=1e5)
summary(svmfit)
Call:
svm(formula = y ~ ., data = dat, kernel = "linear", cost = 1e+05)
Parameters:
SVM-Type: C-classification
SVM-Kernel: linear
cost: 1e+05
Number of Support Vectors: 3
( 1 2 )
Number of Classes: 2
Levels:
-1 1
svmfit=svm(y∼., data=dat , kernel ="linear", cost=1)
summary(svmfit)
Call:
svm(formula = y ~ ., data = dat, kernel = "linear", cost = 1)
Parameters:
SVM-Type: C-classification
SVM-Kernel: linear
cost: 1
Number of Support Vectors: 7
( 4 3 )
Number of Classes: 2
Levels:
-1 1
plot(svmfit ,dat)
como podemos observar usando cost=1 no clasificamos correctamente las observaciones del training, por consiguiente una reduccion al costo da una pqueña mejora
con el objetivo de encajar mejor el modelo de svm utilizaremos un kernel no lineal , para nuestro caso las opciones podran ser polynomial o radial
set.seed(1)
x=matrix(rnorm (200*2) , ncol=2)
x[1:100,]=x[1:100,]+2
x[101:150 ,]=x[101:150,]-2
y=c(rep(1,150) ,rep(2,50))
dat=data.frame(x=x,y=as.factor(y))
plot(x, col=y)
la data sera segmentada al azar dentro del training y el testing , para el caso de clasificaciones radiales, nuestra vairable a realizar tweaking sera gamma
train=sample (200,100)
svmfit=svm(y∼., data=dat[train ,], kernel ="radial", gamma=1,cost=1)
plot(svmfit , dat[train ,])
summary(svmfit)
Call:
svm(formula = y ~ ., data = dat[train, ], kernel = "radial", gamma = 1, cost = 1)
Parameters:
SVM-Type: C-classification
SVM-Kernel: radial
cost: 1
Number of Support Vectors: 39
( 22 17 )
Number of Classes: 2
Levels:
1 2
entonces si reducimos la cantidad de costo podemos reducir el error asociado al training
svmfit=svm(y∼., data=dat[train ,], kernel ="radial",gamma=1,cost=1e5)
plot(svmfit ,dat[train ,])
como vimos en el ejercicio anterior tune nos permitira mejorar nuestro clasificador vectorial
set.seed(1)
tune.out=tune(svm , y∼., data=dat[train ,], kernel ="radial",
ranges=list(cost=c(0.1,1,10,100,1000),
gamma=c(0.5,1,2,3,4) ))
summary (tune.out)
Parameter tuning of ‘svm’:
- sampling method: 10-fold cross validation
- best parameters:
- best performance: 0.11
- Detailed performance results:
NA
table(true=dat[-train ,"y"], pred=predict (tune.out$best.model ,newdata =dat[-train ,]))
pred
true 1 2
1 72 5
2 4 19
bajo este modelo solo 11% esta mal clasificado
Dado un vector que contiene un clasificador numerico para cada bservacion y un vector que contiene la clase ROC estara aplicado por SVM
Empezamos definiendo nuestra funcion
library(ROCR)
rocplot =function (pred , truth , ...){
predob = prediction (pred , truth)
perf = performance (predob , "tpr", "fpr")
plot(perf ,...)}
tomemos encuenta que
X = (X1, X2,…,Xp)T tomara la forma de la combinacion lineal m βˆ0 + βˆ1X1 + βˆ2X2 + … + βˆpXp.
svmfit.opt=svm(y∼., data=dat[train ,], kernel ="radial",
gamma=2, cost=1, decision.values =T)
fitted =attributes (predict (svmfit.opt ,dat[train ,], decision.values=TRUE))$decision.values
par(mfrow=c(1,2))
rocplot(fitted ,dat[train ,"y"], main="Training Data")
tomemos encuenta que ROC curves estan en la data del traingin y nosotros estamos interesados mas en la prediccion por lo que es necesario hacer ajustes al proceso de la definicion del modelo
svmfit.flex=svm(y∼., data=dat[train ,], kernel ="radial",
gamma=50, cost=1, decision.values =T)
fitted=attributes (predict (svmfit.flex ,dat[train ,], decision.values=T))$decision.values
#rocplot(fitted ,dat[train ,"y"],add=T,col="red ")
fitted =attributes(predict (svmfit.opt ,dat[-train ,], decision.values=T))$decision.values
rocplot (fitted ,dat[-train ,"y"], main="Test Data")
fitted=attributes(predict (svmfit.flex ,dat[- train ,], decision.values=T))$decision.values
rocplot (fitted ,dat[-train ,"y"],add=T,col="red")