x1=runif (500) - 0.5
x2=runif (500) - 0.5
y=1*(x1^2-x2^2 > 0)
plot(x1[y == 0], x2[y == 0], col = "red", xlab = "X1", ylab = "X2")
points(x1[y == 1], x2[y == 1], )
lm.fit = glm(y ~ x1 + x2, family = binomial)
summary(lm.fit)
##
## Call:
## glm(formula = y ~ x1 + x2, family = binomial)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.234 -1.134 -1.074 1.198 1.278
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.06401 0.08985 -0.712 0.476
## x1 -0.03169 0.32131 -0.099 0.921
## x2 -0.37902 0.30517 -1.242 0.214
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 692.64 on 499 degrees of freedom
## Residual deviance: 691.09 on 497 degrees of freedom
## AIC: 697.09
##
## Number of Fisher Scoring iterations: 3
data = data.frame(x1 = x1, x2 = x2, y = y)
lm.prob = predict(lm.fit, data, type = "response")
lm.pred = ifelse(lm.prob > 0.50, 1, 0)
data.pos = data[lm.pred == 1, ]
data.neg = data[lm.pred == 0, ]
plot(data.pos$x1, data.pos$x2, xlab = "X1", ylab = "X2")
points(data.neg$x1, data.neg$x2, col = "red")
lm.fit = glm(y ~ poly(x1, 2) + poly(x2, 2) + I(x1 * x2), data = data, family = binomial)
## Warning: glm.fit: algorithm did not converge
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
lm.prob = predict(lm.fit, data, type = "response")
lm.pred = ifelse(lm.prob > 0.5, 1, 0)
data.pos = data[lm.pred == 1, ]
data.neg = data[lm.pred == 0, ]
plot(data.pos$x1, data.pos$x2, xlab = "X1", ylab = "X2")
points(data.neg$x1, data.neg$x2, col = "red")
####Hint: In the lab, we used the plot() function for svm objects only in cases with p = 2. When p > 2, you can use the plot() function to create plots displaying pairs of variables at a time. Essentially, instead of typing plot(svmfit , dat)
train = sample(dim(OJ)[1], 800)
OJ.train = OJ[train, ]
OJ.test = OJ[-train, ]
library(e1071)
svm.linear = svm(Purchase ~ ., kernel = "linear", data = OJ.train, 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: 438
##
## ( 219 219 )
##
##
## 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 441 58
## MM 87 214
(87 + 58)/(441 + 58 + 87 + 214)
## [1] 0.18125
test.pred = predict(svm.linear, OJ.test)
table(OJ.test$Purchase, test.pred)
## test.pred
## CH MM
## CH 140 14
## MM 26 90
(14 + 26)/(140 + 14 + 26 + 90)
## [1] 0.1481481
The training error rate is 18.1% and test error rate is about 14.8%.
tune.out = tune(svm, Purchase ~ ., data = OJ.train, kernel = "linear", ranges = list(cost = 10^seq(-2, 1, by = 0.25)))
summary(tune.out)
##
## Parameter tuning of 'svm':
##
## - sampling method: 10-fold cross validation
##
## - best parameters:
## cost
## 0.03162278
##
## - best performance: 0.1725
##
## - Detailed performance results:
## cost error dispersion
## 1 0.01000000 0.18375 0.04084609
## 2 0.01778279 0.18000 0.04133199
## 3 0.03162278 0.17250 0.04322101
## 4 0.05623413 0.17750 0.04281744
## 5 0.10000000 0.17375 0.04387878
## 6 0.17782794 0.17500 0.04249183
## 7 0.31622777 0.17375 0.04387878
## 8 0.56234133 0.17750 0.04031129
## 9 1.00000000 0.17375 0.04619178
## 10 1.77827941 0.17500 0.04166667
## 11 3.16227766 0.17250 0.04241004
## 12 5.62341325 0.18250 0.04417453
## 13 10.00000000 0.18250 0.04794383
svm.linear = svm(Purchase ~ ., kernel = "linear", data = OJ.train, cost = tune.out$best.parameters$cost)
train.pred = predict(svm.linear, OJ.train)
table(OJ.train$Purchase, train.pred)
## train.pred
## CH MM
## CH 443 56
## MM 80 221
(56 + 80)/(443 + 56 + 80 + 221)
## [1] 0.17
test.pred = predict(svm.linear, OJ.test)
table(OJ.test$Purchase, test.pred)
## test.pred
## CH MM
## CH 139 15
## MM 25 91
(15 + 25)/(139 + 15 + 25 + 91)
## [1] 0.1481481
The training error decreases to 17% but test error slightly increases to 14.8% by using best cost.