Exercise 5: Non-linear decision boundary for logistic regression vs. SVM

(a) Generate data with a quadratic decision boundary

set.seed(1)
x1 <- runif(500) - 0.5
x2 <- runif(500) - 0.5
y  <- 1 * (x1^2 - x2^2 > 0)   # 1/0 class label
dat <- data.frame(x1 = x1, x2 = x2, y = as.factor(y))

(b) Plot the observations

plot(x1, x2, col = (y + 2), pch = 19,
     main = "5(b): True classes (quadratic boundary)",
     xlab = "X1", ylab = "X2")

(c) Logistic regression

glm.fit <- glm(y ~ x1 + x2, data = dat, family = "binomial")
summary(glm.fit)
## 
## Call:
## glm(formula = y ~ x1 + x2, family = "binomial", data = dat)
## 
## 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

(d) Predicted labels from linear logistic model, plotted

glm.probs <- predict(glm.fit, type = "response")
glm.pred  <- ifelse(glm.probs > 0.5, 1, 0)

plot(x1, x2, col = (glm.pred + 2), pch = 19,
     main = "5(d): Logistic regression (linear terms) predictions",
     xlab = "X1", ylab = "X2")

(e) Logistic regression with non-linear transformations

glm.fit2 <- glm(y ~ poly(x1, 2) + poly(x2, 2) + I(x1 * x2),
                data = dat, family = "binomial")
## Warning: glm.fit: algorithm did not converge
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
summary(glm.fit2)
## 
## Call:
## glm(formula = y ~ poly(x1, 2) + poly(x2, 2) + I(x1 * x2), family = "binomial", 
##     data = dat)
## 
## 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

(f) Predicted labels from non-linear logistic model, plotted

glm.probs2 <- predict(glm.fit2, type = "response")
glm.pred2  <- ifelse(glm.probs2 > 0.5, 1, 0)

plot(x1, x2, col = (glm.pred2 + 2), pch = 19,
     main = "5(f): Logistic regression (non-linear terms) predictions",
     xlab = "X1", ylab = "X2")

(g) Support vector classifier

svm.linear <- svm(y ~ x1 + x2, data = dat, kernel = "linear", cost = 0.1)
svm.linear.pred <- predict(svm.linear, dat)

plot(x1, x2, col = (as.numeric(as.character(svm.linear.pred)) + 2), pch = 19,
     main = "5(g): SVM linear kernel predictions",
     xlab = "X1", ylab = "X2")

(h) SVM with a non-linear kernel

svm.radial <- svm(y ~ x1 + x2, data = dat, kernel = "radial", gamma = 1, cost = 1)
svm.radial.pred <- predict(svm.radial, dat)

plot(x1, x2, col = (as.numeric(as.character(svm.radial.pred)) + 2), pch = 19,
     main = "5(h): SVM radial kernel predictions",
     xlab = "X1", ylab = "X2")

(i) Comment on results

Linear logistic regression (c):None of the coefficients for x1 or x2 are significant (p-values 0.536 and 0.993), and the null vs. residual deviance barely moved (692.18 to 691.79). This is expected — the true boundary (x1^2 - x2^2 = 0) is a shape a linear-in-X1,X2 model cannot represent, so the fit is essentially useless at separating the classes.

Non-linear logistic regression (e): Residual deviance collapsed to essentially zero (3.58e-06) with AIC dropping to 12 — a massive improvement. However, the huge standard errors (e.g. a coefficient of 27218.5 with SE 842987.2) and the “algorithm did not converge” / “fitted probabilities numerically 0 or 1” warnings indicate perfect separation: with the quadratic and interaction terms, the model can draw an almost perfect boundary between the classes, so perfectly that the MLE diverges. The fitted class predictions are still excellent even though the coefficients themselves are unstable and uninterpretable.

Linear SVM (g): Every point gets predicted as the same class. With cost = 0.1, the optimizer finds that no straight line separates the classes meaningfully better than guessing — because the true boundary has zero linear separability. Any straight line cuts through both classes roughly equally, so the trivial single-class solution wins.

Radial SVM (h): A genuine non-linear boundary emerges — green regions on the left/right, red running through the middle in a curved band, roughly tracking the true |x1| = |x2| split. It isn’t perfect (visible misclassification near the boundary), but it clearly recovers the real structure that the linear SVM completely missed.

Overall: This example shows that when the true decision boundary is non-linear, both feature engineering and kernel choice matter. The linear logistic regression and linear-kernel SVM — both restricted to straight-line boundaries — fail completely. Both non-linear approaches (logistic regression with quadratic/interaction terms, and the radial-kernel SVM) successfully recover the shape of the true boundary, whether the non-linearity is introduced explicitly via engineered polynomial features or implicitly via a non-linear kernel. Note also the asymmetry in how well each fits: the logistic regression achieves perfect separation (and the overfitting risk that implies), while the radial SVM produces a softer, more moderate boundary that may generalize better despite being visually “less clean.”

Exercise 7: Predicting high/low gas mileage in the Auto data set

(a) Binary variable: 1 = above median mpg, 0 = below median

Auto2 <- Auto
mpg.median <- median(Auto2$mpg)
Auto2$mpglevel <- as.factor(ifelse(Auto2$mpg > mpg.median, 1, 0))
Auto2$mpg <- NULL   # remove mpg itself so it isn't used as a predictor

(b) Support vector classifier, tune over cost

set.seed(1)
tune.linear <- tune(svm, mpglevel ~ ., data = Auto2, kernel = "linear",
                    ranges = list(cost = c(0.01, 0.1, 1, 5, 10, 100)))
summary(tune.linear)
## 
## Parameter tuning of 'svm':
## 
## - sampling method: 10-fold cross validation 
## 
## - best parameters:
##  cost
##   0.1
## 
## - best performance: 0.08673077 
## 
## - Detailed performance results:
##    cost      error dispersion
## 1 1e-02 0.08923077 0.04698309
## 2 1e-01 0.08673077 0.04040897
## 3 1e+00 0.09961538 0.04923181
## 4 5e+00 0.11230769 0.05826857
## 5 1e+01 0.11237179 0.05701890
## 6 1e+02 0.11750000 0.06208951

(c) SVM with radial and polynomial kernels

set.seed(1)
tune.radial <- tune(svm, mpglevel ~ ., data = Auto2, kernel = "radial",
                     ranges = list(cost = c(0.1, 1, 5, 10, 100),
                                   gamma = c(0.01, 0.1, 1, 5)))
summary(tune.radial)
## 
## Parameter tuning of 'svm':
## 
## - sampling method: 10-fold cross validation 
## 
## - best parameters:
##  cost gamma
##    10     1
## 
## - best performance: 0.07897436 
## 
## - Detailed performance results:
##     cost gamma      error dispersion
## 1    0.1  0.01 0.11224359 0.03836937
## 2    1.0  0.01 0.08673077 0.04551036
## 3    5.0  0.01 0.08673077 0.04040897
## 4   10.0  0.01 0.08673077 0.03855882
## 5  100.0  0.01 0.09692308 0.05742483
## 6    0.1  0.10 0.08666667 0.04193895
## 7    1.0  0.10 0.08923077 0.04376306
## 8    5.0  0.10 0.08423077 0.04689205
## 9   10.0  0.10 0.08416667 0.05256241
## 10 100.0  0.10 0.10211538 0.04535762
## 11   0.1  1.00 0.55115385 0.04366593
## 12   1.0  1.00 0.07903846 0.04891067
## 13   5.0  1.00 0.08147436 0.04910668
## 14  10.0  1.00 0.07897436 0.04869339
## 15 100.0  1.00 0.07897436 0.04869339
## 16   0.1  5.00 0.55115385 0.04366593
## 17   1.0  5.00 0.48967949 0.05080301
## 18   5.0  5.00 0.48211538 0.05914633
## 19  10.0  5.00 0.48211538 0.05914633
## 20 100.0  5.00 0.48211538 0.05914633
set.seed(1)
tune.poly <- tune(svm, mpglevel ~ ., data = Auto2, kernel = "polynomial",
                   ranges = list(cost = c(0.1, 1, 5, 10, 100),
                                 degree = c(2, 3, 4)))
summary(tune.poly)
## 
## Parameter tuning of 'svm':
## 
## - sampling method: 10-fold cross validation 
## 
## - best parameters:
##  cost degree
##   100      2
## 
## - best performance: 0.3167308 
## 
## - Detailed performance results:
##     cost degree     error dispersion
## 1    0.1      2 0.5511538 0.04366593
## 2    1.0      2 0.5511538 0.04366593
## 3    5.0      2 0.5511538 0.04366593
## 4   10.0      2 0.5206410 0.08505283
## 5  100.0      2 0.3167308 0.09410274
## 6    0.1      3 0.5511538 0.04366593
## 7    1.0      3 0.5511538 0.04366593
## 8    5.0      3 0.5511538 0.04366593
## 9   10.0      3 0.5511538 0.04366593
## 10 100.0      3 0.4032692 0.10793388
## 11   0.1      4 0.5511538 0.04366593
## 12   1.0      4 0.5511538 0.04366593
## 13   5.0      4 0.5511538 0.04366593
## 14  10.0      4 0.5511538 0.04366593
## 15 100.0      4 0.5511538 0.04366593

Interpretation of (b) and (c):

The best linear-kernel result is cost = 0.1 with CV error 8.67%.

For the radial kernel, the best combination is cost = 10, gamma = 1, giving CV error 7.90% — slightly better than linear. For low gamma (0.01–0.1), radial performance is very close to the linear kernel’s regardless of cost, since a small gamma makes the radial kernel behave nearly linearly. As gamma increases, performance depends much more sharply on cost: at gamma = 1, error is good for cost >= 1 (~7.9–8.1%), but at gamma = 5, error is very poor across the board (48–55%) — the boundary becomes so localized around individual training points that it stops generalizing (severe overfitting), regardless of cost. Gamma is the parameter driving the qualitative behavior here, with cost fine-tuning within that.

The polynomial kernel performs much worse overall. Even at its best (cost = 100, degree = 2), CV error is 31.67% — nearly 4x worse than radial or linear. Most of the grid (low-to-moderate cost, any degree) sits at ~55% error, which is essentially the trivial “predict one class for everyone” baseline. Only very high cost values pull performance down, and higher degree (3, 4) doesn’t help — degree 2 is best, and degree 4 never escapes the baseline error at any cost tried.

Ranking: radial (7.90%) < linear (8.67%) < polynomial (31.67%). The modest gap between linear and radial suggests the true relationship between the predictors and mpg-level is close to linear (or needs little curvature), while the polynomial kernel’s poor showing suggests it needs a very different tuning grid to be competitive, or its basis functions just don’t match the shape of this boundary well.

(d) Plots to visualize decision boundaries

svm.linear.best <- tune.linear$best.model
svm.radial.best <- tune.radial$best.model
svm.poly.best   <- tune.poly$best.model

# Example usage (uncomment to view individual plots):
# plot(svm.linear.best, Auto2, horsepower ~ weight)
# plot(svm.radial.best, Auto2, horsepower ~ weight)
# plot(svm.poly.best,   Auto2, horsepower ~ weight)

Exercise 8: OJ data set — Purchase classification

(a) Training/test split: 800 training observations

set.seed(1)
train.idx <- sample(1:nrow(OJ), 800)
OJ.train <- OJ[train.idx, ]
OJ.test  <- OJ[-train.idx, ]

(b) Support vector classifier, cost = 0.01

svm.oj.linear <- svm(Purchase ~ ., data = OJ.train, kernel = "linear", cost = 0.01)
summary(svm.oj.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

(c) Training and test error rates

train.pred.linear <- predict(svm.oj.linear, OJ.train)
test.pred.linear  <- predict(svm.oj.linear, OJ.test)

table(OJ.train$Purchase, train.pred.linear)
##     train.pred.linear
##       CH  MM
##   CH 420  65
##   MM  75 240
table(OJ.test$Purchase, test.pred.linear)
##     test.pred.linear
##       CH  MM
##   CH 153  15
##   MM  33  69
train.err.linear <- mean(train.pred.linear != OJ.train$Purchase)
test.err.linear  <- mean(test.pred.linear != OJ.test$Purchase)
train.err.linear
## [1] 0.175
test.err.linear
## [1] 0.1777778

(d) Tune cost from 0.01 to 10

set.seed(1)
tune.oj.linear <- tune(svm, Purchase ~ ., data = OJ.train, kernel = "linear",
                        ranges = list(cost = 10^seq(-2, 1, length.out = 10)))
summary(tune.oj.linear)
## 
## Parameter tuning of 'svm':
## 
## - sampling method: 10-fold cross validation 
## 
## - best parameters:
##       cost
##  0.4641589
## 
## - best performance: 0.16875 
## 
## - Detailed performance results:
##           cost   error dispersion
## 1   0.01000000 0.17625 0.02853482
## 2   0.02154435 0.17625 0.02972676
## 3   0.04641589 0.17500 0.02568506
## 4   0.10000000 0.17250 0.03162278
## 5   0.21544347 0.17250 0.02751262
## 6   0.46415888 0.16875 0.02651650
## 7   1.00000000 0.17500 0.02946278
## 8   2.15443469 0.17125 0.03064696
## 9   4.64158883 0.17125 0.03175973
## 10 10.00000000 0.17375 0.03197764
best.cost.linear <- tune.oj.linear$best.parameters$cost

(e) Training/test error with optimal cost

svm.oj.linear.best <- svm(Purchase ~ ., data = OJ.train, kernel = "linear",
                          cost = best.cost.linear)
train.pred.linear.best <- predict(svm.oj.linear.best, OJ.train)
test.pred.linear.best  <- predict(svm.oj.linear.best, OJ.test)

train.err.linear.best <- mean(train.pred.linear.best != OJ.train$Purchase)
test.err.linear.best  <- mean(test.pred.linear.best != OJ.test$Purchase)
train.err.linear.best
## [1] 0.165
test.err.linear.best
## [1] 0.1555556

(f) Radial kernel

svm.oj.radial <- svm(Purchase ~ ., data = OJ.train, kernel = "radial", cost = 0.01)
summary(svm.oj.radial)
## 
## Call:
## svm(formula = Purchase ~ ., data = OJ.train, kernel = "radial", cost = 0.01)
## 
## 
## Parameters:
##    SVM-Type:  C-classification 
##  SVM-Kernel:  radial 
##        cost:  0.01 
## 
## Number of Support Vectors:  634
## 
##  ( 319 315 )
## 
## 
## Number of Classes:  2 
## 
## Levels: 
##  CH MM
train.err.radial0 <- mean(predict(svm.oj.radial, OJ.train) != OJ.train$Purchase)
test.err.radial0  <- mean(predict(svm.oj.radial, OJ.test)  != OJ.test$Purchase)
train.err.radial0
## [1] 0.39375
test.err.radial0
## [1] 0.3777778
set.seed(1)
tune.oj.radial <- tune(svm, Purchase ~ ., data = OJ.train, kernel = "radial",
                        ranges = list(cost = 10^seq(-2, 1, length.out = 10)))
summary(tune.oj.radial)
## 
## Parameter tuning of 'svm':
## 
## - sampling method: 10-fold cross validation 
## 
## - best parameters:
##       cost
##  0.4641589
## 
## - best performance: 0.17125 
## 
## - Detailed performance results:
##           cost   error dispersion
## 1   0.01000000 0.39375 0.04007372
## 2   0.02154435 0.39375 0.04007372
## 3   0.04641589 0.20875 0.04041881
## 4   0.10000000 0.18625 0.02853482
## 5   0.21544347 0.18250 0.03016160
## 6   0.46415888 0.17125 0.02045490
## 7   1.00000000 0.17125 0.02128673
## 8   2.15443469 0.17875 0.02128673
## 9   4.64158883 0.18125 0.02144923
## 10 10.00000000 0.18625 0.02853482
best.cost.radial <- tune.oj.radial$best.parameters$cost

svm.oj.radial.best <- svm(Purchase ~ ., data = OJ.train, kernel = "radial",
                          cost = best.cost.radial)
train.err.radial.best <- mean(predict(svm.oj.radial.best, OJ.train) != OJ.train$Purchase)
test.err.radial.best  <- mean(predict(svm.oj.radial.best, OJ.test)  != OJ.test$Purchase)
train.err.radial.best
## [1] 0.1475
test.err.radial.best
## [1] 0.1777778

(g) Polynomial kernel, degree = 2

svm.oj.poly <- svm(Purchase ~ ., data = OJ.train, kernel = "polynomial",
                   degree = 2, cost = 0.01)
summary(svm.oj.poly)
## 
## Call:
## svm(formula = Purchase ~ ., data = OJ.train, kernel = "polynomial", 
##     degree = 2, cost = 0.01)
## 
## 
## Parameters:
##    SVM-Type:  C-classification 
##  SVM-Kernel:  polynomial 
##        cost:  0.01 
##      degree:  2 
##      coef.0:  0 
## 
## Number of Support Vectors:  636
## 
##  ( 321 315 )
## 
## 
## Number of Classes:  2 
## 
## Levels: 
##  CH MM
train.err.poly0 <- mean(predict(svm.oj.poly, OJ.train) != OJ.train$Purchase)
test.err.poly0  <- mean(predict(svm.oj.poly, OJ.test)  != OJ.test$Purchase)
train.err.poly0
## [1] 0.3725
test.err.poly0
## [1] 0.3666667
set.seed(1)
tune.oj.poly <- tune(svm, Purchase ~ ., data = OJ.train, kernel = "polynomial",
                     degree = 2,
                     ranges = list(cost = 10^seq(-2, 1, length.out = 10)))
summary(tune.oj.poly)
## 
## Parameter tuning of 'svm':
## 
## - sampling method: 10-fold cross validation 
## 
## - best parameters:
##      cost
##  2.154435
## 
## - best performance: 0.1775 
## 
## - Detailed performance results:
##           cost   error dispersion
## 1   0.01000000 0.39125 0.04210189
## 2   0.02154435 0.36625 0.03537988
## 3   0.04641589 0.35125 0.03839216
## 4   0.10000000 0.32125 0.05001736
## 5   0.21544347 0.21875 0.03830162
## 6   0.46415888 0.20625 0.04497299
## 7   1.00000000 0.20250 0.04116363
## 8   2.15443469 0.17750 0.04158325
## 9   4.64158883 0.18375 0.03387579
## 10 10.00000000 0.18125 0.02779513
best.cost.poly <- tune.oj.poly$best.parameters$cost

svm.oj.poly.best <- svm(Purchase ~ ., data = OJ.train, kernel = "polynomial",
                        degree = 2, cost = best.cost.poly)
train.err.poly.best <- mean(predict(svm.oj.poly.best, OJ.train) != OJ.train$Purchase)
test.err.poly.best  <- mean(predict(svm.oj.poly.best, OJ.test)  != OJ.test$Purchase)
train.err.poly.best
## [1] 0.15875
test.err.poly.best
## [1] 0.2111111

(h) Compare all six error rates

results <- data.frame(
  Model = c("Linear (cost=0.01)", "Linear (tuned)",
            "Radial (cost=0.01)", "Radial (tuned)",
            "Poly deg2 (cost=0.01)", "Poly deg2 (tuned)"),
  TrainError = c(train.err.linear, train.err.linear.best,
                 train.err.radial0, train.err.radial.best,
                 train.err.poly0, train.err.poly.best),
  TestError = c(test.err.linear, test.err.linear.best,
                test.err.radial0, test.err.radial.best,
                test.err.poly0, test.err.poly.best)
)
print(results)
##                   Model TrainError TestError
## 1    Linear (cost=0.01)    0.17500 0.1777778
## 2        Linear (tuned)    0.16500 0.1555556
## 3    Radial (cost=0.01)    0.39375 0.3777778
## 4        Radial (tuned)    0.14750 0.1777778
## 5 Poly deg2 (cost=0.01)    0.37250 0.3666667
## 6     Poly deg2 (tuned)    0.15875 0.2111111

Interpretation of (h):

At the default cost = 0.01, the linear kernel already does reasonably (17.5% / 17.8% train/test), while radial and polynomial are both badly underfit (~37–39% error) — with cost this low, the margin is so wide/soft that both non-linear kernels can’t yet build a useful boundary. This shows cost matters enormously for these kernels specifically.

After tuning, all three models land in a similar training-error range (14.75%–16.5%), but they diverge on test error: linear generalizes best (15.6%), radial is respectable but worse (17.8%), and polynomial is the weakest (21.1%).

Overfitting signal: radial has the lowest training error of all six models (14.75%) but not the lowest test error — its test error (17.8%) is actually identical to the untuned linear model’s. That gap between train and test error is a classic sign that the radial kernel fits some training-specific noise that doesn’t carry over to new data. The polynomial kernel shows this pattern even more strongly (15.9% train vs. 21.1% test — the largest train/test gap of the three tuned models).

Answer: Overall, the tuned support vector classifier with a linear kernel gives the best results on this data, achieving the lowest test error (15.6%) among all six models compared. While the tuned radial kernel achieves a slightly lower training error (14.75% vs. 16.5%), it does not generalize as well, and its test error (17.8%) equals the untuned linear model’s test error — meaning the radial kernel’s extra flexibility isn’t translating into better real-world performance here. The polynomial kernel (degree 2) performs worst overall on the test set (21.1%), despite having a similar training error to the other two tuned models, suggesting it overfits the most relative to what it gains in flexibility.

This suggests the true relationship between the predictors and Purchase in the OJ data is close to linear, or at least does not require the additional curvature that radial or polynomial kernels provide. More flexible models are not automatically betterm cross-validated test performance should guide model selection, and here the simplest kernel wins. Since tune() used 10-fold CV within the training set only, the held-out test set gives an independent check, and it confirms the same ranking — adding confidence that linear really is the better choice here rather than an artifact of the particular CV folds.