#chapter 9 page 368: 5, 7, 8
# 5
set.seed(1)
x1=runif(500)-0.5
x2=runif(500)-0.5
y=1*(x1^2-x2^2>0)
set.seed(1)
x1=runif(500)-0.5
x2=runif(500)-0.5
y=1*(x1^2-x2^2>0)
plot(x1,x2,xlab="X1",ylab="X2",col=(3-y),pch=(2-y))

logit.fit=glm(y~x1+x2,family="binomial")
summary(logit.fit)
##
## Call:
## glm(formula = y ~ x1 + x2, family = "binomial")
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.179 -1.139 -1.112 1.206 1.257
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.087260 0.089579 -0.974 0.330
## x1 0.196199 0.316864 0.619 0.536
## x2 -0.002854 0.305712 -0.009 0.993
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 692.18 on 499 degrees of freedom
## Residual deviance: 691.79 on 497 degrees of freedom
## AIC: 697.79
##
## Number of Fisher Scoring iterations: 3
data=data.frame(x1=x1,x2=x2,y=y)
probs=predict(logit.fit,data,type="response")
preds=rep(0,500)
preds[probs>0.47]=1
plot(data[preds==1,]$x1,data[preds==1,]$x2,col=(3-1),pch=(2-1),xlab="X1",ylab="X2")
points(data[preds==0,]$x1,data[preds==0,]$x2,col=(3-0),pch=(2-0))

logitnl.fit=glm(y~poly(x1,2)+poly(x2,2)+I(x1*x2),family="binomial")
## Warning: glm.fit: algorithm did not converge
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
summary(logitnl.fit)
##
## Call:
## glm(formula = y ~ poly(x1, 2) + poly(x2, 2) + I(x1 * x2), family = "binomial")
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -8.240e-04 -2.000e-08 -2.000e-08 2.000e-08 1.163e-03
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -102.2 4302.0 -0.024 0.981
## poly(x1, 2)1 2715.3 141109.5 0.019 0.985
## poly(x1, 2)2 27218.5 842987.2 0.032 0.974
## poly(x2, 2)1 -279.7 97160.4 -0.003 0.998
## poly(x2, 2)2 -28693.0 875451.3 -0.033 0.974
## I(x1 * x2) -206.4 41802.8 -0.005 0.996
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 6.9218e+02 on 499 degrees of freedom
## Residual deviance: 3.5810e-06 on 494 degrees of freedom
## AIC: 12
##
## Number of Fisher Scoring iterations: 25
probs=predict(logitnl.fit,data,type="response")
preds=rep(0,500)
preds[probs>0.47]=1
plot(data[preds==1,]$x1,data[preds==1,]$x2,col=(3-1),pch=(2-1),xlab="X1",ylab="X2")
points(data[preds==0,]$x1,data[preds==0,]$x2,col=(3-0),pch=(2-0))

library(e1071)
## Warning: package 'e1071' was built under R version 4.1.3
data$y=as.factor(data$y)
svm.fit=svm(y~x1+x2,data,kernel="linear",cost=0.01)
preds=predict(svm.fit,data)
plot(data[preds==0,]$x1,data[preds==0,]$x2,col=(3-0),pch=(2-0),xlab="X1",ylab="X2")
points(data[preds==1,]$x1,data[preds==1,]$x2,col=(3-1),pch=(2-1))

data$y=as.factor(data$y)
svmnl.fit=svm(y~x1+x2,data,kernel="radial",gamma=1)
preds=predict(svmnl.fit,data)
plot(data[preds==0,]$x1,data[preds==0,]$x2,col=(3-0),pch=(2-0),xlab="X1",ylab="X2")
points(data[preds==1,]$x1,data[preds==1,]$x2,col=(3-1),pch=(2-1))

# 7
library(ISLR)
## Warning: package 'ISLR' was built under R version 4.1.3
var=ifelse(Auto$mpg>median(Auto$mpg),1,0)
Auto$mpglevel=as.factor(var)
set.seed(1)
tune.out=tune(svm,mpglevel~.,data=Auto,kernel="linear",ranges=list(cost=c(0.01,0.1,1,5,10,100,1000)))
summary(tune.out)
##
## Parameter tuning of 'svm':
##
## - sampling method: 10-fold cross validation
##
## - best parameters:
## cost
## 1
##
## - best performance: 0.01025641
##
## - Detailed performance results:
## cost error dispersion
## 1 1e-02 0.07653846 0.03617137
## 2 1e-01 0.04596154 0.03378238
## 3 1e+00 0.01025641 0.01792836
## 4 5e+00 0.02051282 0.02648194
## 5 1e+01 0.02051282 0.02648194
## 6 1e+02 0.03076923 0.03151981
## 7 1e+03 0.03076923 0.03151981
set.seed(1)
tune.out=tune(svm,mpglevel~.,data=Auto,kernel="polynomial",ranges=list(cost=c(0.01,0.1,1,5,10,100),degree=c(2,3,4)))
summary(tune.out)
##
## Parameter tuning of 'svm':
##
## - sampling method: 10-fold cross validation
##
## - best parameters:
## cost degree
## 100 2
##
## - best performance: 0.3013462
##
## - Detailed performance results:
## cost degree error dispersion
## 1 1e-02 2 0.5511538 0.04366593
## 2 1e-01 2 0.5511538 0.04366593
## 3 1e+00 2 0.5511538 0.04366593
## 4 5e+00 2 0.5511538 0.04366593
## 5 1e+01 2 0.5130128 0.08963366
## 6 1e+02 2 0.3013462 0.09961961
## 7 1e-02 3 0.5511538 0.04366593
## 8 1e-01 3 0.5511538 0.04366593
## 9 1e+00 3 0.5511538 0.04366593
## 10 5e+00 3 0.5511538 0.04366593
## 11 1e+01 3 0.5511538 0.04366593
## 12 1e+02 3 0.3446154 0.09821588
## 13 1e-02 4 0.5511538 0.04366593
## 14 1e-01 4 0.5511538 0.04366593
## 15 1e+00 4 0.5511538 0.04366593
## 16 5e+00 4 0.5511538 0.04366593
## 17 1e+01 4 0.5511538 0.04366593
## 18 1e+02 4 0.5511538 0.04366593
set.seed(1)
tune.out=tune(svm,mpglevel~.,data=Auto,kernel="radial",ranges=list(cost=c(0.01,0.1,1,5,10,100),gamma=c(0.01,0.1,1,5,10,100)))
summary(tune.out)
##
## Parameter tuning of 'svm':
##
## - sampling method: 10-fold cross validation
##
## - best parameters:
## cost gamma
## 100 0.01
##
## - best performance: 0.01282051
##
## - Detailed performance results:
## cost gamma error dispersion
## 1 1e-02 1e-02 0.55115385 0.04366593
## 2 1e-01 1e-02 0.08929487 0.04382379
## 3 1e+00 1e-02 0.07403846 0.03522110
## 4 5e+00 1e-02 0.04852564 0.03303346
## 5 1e+01 1e-02 0.02557692 0.02093679
## 6 1e+02 1e-02 0.01282051 0.01813094
## 7 1e-02 1e-01 0.21711538 0.09865227
## 8 1e-01 1e-01 0.07903846 0.03874545
## 9 1e+00 1e-01 0.05371795 0.03525162
## 10 5e+00 1e-01 0.02820513 0.03299190
## 11 1e+01 1e-01 0.03076923 0.03375798
## 12 1e+02 1e-01 0.03583333 0.02759051
## 13 1e-02 1e+00 0.55115385 0.04366593
## 14 1e-01 1e+00 0.55115385 0.04366593
## 15 1e+00 1e+00 0.06384615 0.04375618
## 16 5e+00 1e+00 0.05884615 0.04020934
## 17 1e+01 1e+00 0.05884615 0.04020934
## 18 1e+02 1e+00 0.05884615 0.04020934
## 19 1e-02 5e+00 0.55115385 0.04366593
## 20 1e-01 5e+00 0.55115385 0.04366593
## 21 1e+00 5e+00 0.49493590 0.04724924
## 22 5e+00 5e+00 0.48217949 0.05470903
## 23 1e+01 5e+00 0.48217949 0.05470903
## 24 1e+02 5e+00 0.48217949 0.05470903
## 25 1e-02 1e+01 0.55115385 0.04366593
## 26 1e-01 1e+01 0.55115385 0.04366593
## 27 1e+00 1e+01 0.51794872 0.05063697
## 28 5e+00 1e+01 0.51794872 0.04917316
## 29 1e+01 1e+01 0.51794872 0.04917316
## 30 1e+02 1e+01 0.51794872 0.04917316
## 31 1e-02 1e+02 0.55115385 0.04366593
## 32 1e-01 1e+02 0.55115385 0.04366593
## 33 1e+00 1e+02 0.55115385 0.04366593
## 34 5e+00 1e+02 0.55115385 0.04366593
## 35 1e+01 1e+02 0.55115385 0.04366593
## 36 1e+02 1e+02 0.55115385 0.04366593
svm.linear=svm(mpglevel~.,data=Auto,kernel="linear",cost=1)
svm.poly=svm(mpglevel~.,data=Auto,kernel="polynomial",cost=100,degree=2)
svm.radial=svm(mpglevel~.,data=Auto,kernel="radial",cost=100,gamma=0.01)
# 8
set.seed(1)
train=sample(nrow(OJ),800)
OJ.train=OJ[train,]
OJ.test=OJ[-train,]
svm.linear=svm(Purchase~.,data=OJ.train,kernel="linear",cost=0.01)
summary(svm.linear)
##
## Call:
## svm(formula = Purchase ~ ., data = OJ.train, kernel = "linear", cost = 0.01)
##
##
## Parameters:
## SVM-Type: C-classification
## SVM-Kernel: linear
## cost: 0.01
##
## Number of Support Vectors: 435
##
## ( 219 216 )
##
##
## Number of Classes: 2
##
## Levels:
## CH MM
train.pred=predict(svm.linear,OJ.train)
table(OJ.train$Purchase,train.pred)
## train.pred
## CH MM
## CH 420 65
## MM 75 240
train.pred=predict(svm.linear,OJ.test)
svm.linear=svm(Purchase~.,kernel="linear",data=OJ.train,cost=tune.out$best.parameter$cost)
train.pred=predict(svm.linear,OJ.train)
table(OJ.train$Purchase,train.pred)
## train.pred
## CH MM
## CH 423 62
## MM 70 245
train.pred=predict(svm.linear,OJ.test)
svm.radial=svm(Purchase~.,kernel="radial",data=OJ.train)
summary(svm.radial)
##
## Call:
## svm(formula = Purchase ~ ., data = OJ.train, kernel = "radial")
##
##
## Parameters:
## SVM-Type: C-classification
## SVM-Kernel: radial
## cost: 1
##
## Number of Support Vectors: 373
##
## ( 188 185 )
##
##
## Number of Classes: 2
##
## Levels:
## CH MM
train.pred=predict(svm.radial,OJ.train)
table(OJ.train$Purchase,train.pred)
## train.pred
## CH MM
## CH 441 44
## MM 77 238
test.pred=predict(svm.radial,OJ.test)
table(OJ.test$Purchase,test.pred)
## test.pred
## CH MM
## CH 151 17
## MM 33 69
svm.poly=svm(Purchase~.,kernel="polynomial",data=OJ.train,degree=2)
summary(svm.poly)
##
## Call:
## svm(formula = Purchase ~ ., data = OJ.train, kernel = "polynomial",
## degree = 2)
##
##
## Parameters:
## SVM-Type: C-classification
## SVM-Kernel: polynomial
## cost: 1
## degree: 2
## coef.0: 0
##
## Number of Support Vectors: 447
##
## ( 225 222 )
##
##
## Number of Classes: 2
##
## Levels:
## CH MM
train.pred=predict(svm.poly,OJ.train)
table(OJ.train$Purchase,train.pred)
## train.pred
## CH MM
## CH 449 36
## MM 110 205
test.pred=predict(svm.poly,OJ.test)
table(OJ.test$Purchase,test.pred)
## test.pred
## CH MM
## CH 153 15
## MM 45 57