STA-6543 Predictive Modeling
Emrah Akbas
2026-07-24
Assignment 8
ISLR Chapter 09 (page 398): 5, 7, 8
Q5. Non-Linear Decision Boundaries with Logistic Regression and SVMs
We have seen that we can fit an SVM with a non-linear kernel in order to perform classification using a non-linear decision boundary. We will now see that we can also obtain a non-linear decision boundary by performing logistic regression using non-linear transformations of the features.
(a) Generate a data set with n = 500 and p = 2, such that the observations belong to two classes with a quadratic decision boundary between them. For instance, you can do this as follows:
> x1 <- runif(500) - 0.5
> x2 <- runif(500) - 0.5
> y <- 1 * (x1^2 - x2^2 > 0)(b) Plot the observations, colored according to their class labels. Your plot should display X1 on the x-axis, and X2 on the y-axis.
library(ggplot2)
ggplot(df, aes(x = X1, y = X2, color = Class)) +
geom_point(alpha = 0.8, size = 2) +
labs(x = "X1", y = "X2") +
theme_minimal()
(c) Fit a logistic regression model to the data, using X1 and X2 as predictors.
##
## Call:
## glm(formula = Class ~ x1 + x2, family = binomial, data = df)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.02449 0.08988 0.272 0.7853
## x1 0.23573 0.31281 0.754 0.4511
## x2 -0.52372 0.30367 -1.725 0.0846 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 693.12 on 499 degrees of freedom
## Residual deviance: 689.57 on 497 degrees of freedom
## AIC: 695.57
##
## Number of Fisher Scoring iterations: 3(d) Apply this model to the training data in order to obtain a predicted class label for each training observation. Plot the observations, colored according to the predicted class labels. The decision boundary should be linear.
df$pred_prob <- predict(log_model, newdata = df, type = "response")
df$pred_class <- as.factor(ifelse(df$pred_prob > 0.5, 1, 0))
ggplot(df, aes(x = X1, y = X2, color = pred_class)) +
geom_point(alpha = 0.8, size = 2) +
labs(title = "Logistic Regression", x = "X1", y = "X2", color = "Predicted Class") +
theme_minimal()
(e) Now fit a logistic regression model to the data using non-linear functions of X1 and X2 as predictors (e.g. X21 , X1×X2, log(X2),and so forth).
##
## Call:
## glm(formula = Class ~ I(X1^2) + I(X2^2), family = binomial, data = df)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.752e-02 3.592e+02 0.000 1.000
## I(X1^2) 5.101e+04 7.197e+05 0.071 0.943
## I(X2^2) -5.102e+04 7.199e+05 -0.071 0.943
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 6.9312e+02 on 499 degrees of freedom
## Residual deviance: 2.4094e-05 on 497 degrees of freedom
## AIC: 6
##
## Number of Fisher Scoring iterations: 25(f) Apply this model to the training data in order to obtain a predicted class label for each training observation. Plot the observations, colored according to the predicted class labels. The decision boundary should be obviously non-linear. If it is not, then repeat (a)-(e) until you come up with an example in which the predicted class labels are obviously non-linear.
df$pred_prob_nl <- predict(nl_log_model, newdata = df, type = "response")
df$pred_class_nl <- as.factor(ifelse(df$pred_prob_nl > 0.5, 1, 0))
ggplot(df, aes(x = X1, y = X2, color = pred_class_nl)) + geom_point(alpha = 0.8, size = 2) +
labs(title = "Logistic Regression (Quadratic)", x = "X1", y = "X2", color = "Predicted Class") +
theme_minimal()
(g) Fit a support vector classifier to the data with X1 and X2 as predictors. Obtain a class prediction for each training observation. Plot the observations, colored according to the predicted class labels.
##
## Attaching package: 'e1071'## The following object is masked from 'package:ggplot2':
##
## elementsvc_model <- svm(Class ~ X1 + X2, data = df, kernel = "linear", cost = 1)
df$pred_svc <- predict(svc_model, newdata = df)
ggplot(df, aes(x = X1, y = X2, color = pred_svc)) +
geom_point(alpha = 0.8, size = 2) + labs(title = "Support Vector Classifier", x = "X1", y = "X2", color = "Predicted Class") +
theme_minimal()
(h) Fit a SVM using a non-linear kernel to the data. Obtain a class prediction for each training observation. Plot the observations, colored according to the predicted class labels.
svm_model <- svm(Class ~ X1 + X2, data = df, kernel = "radial", gamma = 1, cost = 1)
df$pred_svm <- predict(svm_model, newdata = df)
ggplot(df, aes(x = X1, y = X2, color = pred_svm)) +
geom_point(alpha = 0.8, size = 2) + labs(title = "SVM", x = "X1", y = "X2", color = "Predicted Class") +
theme_minimal()
(i) Comment on your results.
Linear models like basic logistic regression and the linear SVC completely fail here because they try to draw a straight line through a boundary that is naturally curved. Since they can’t bend, they end up predicting almost everything as one class. On the other hand, both the nonlinear logistic regression and the radial SVM do an amazing job capturing that curved, quadratic shape. The big difference is that with logistic regression, I had to manually figure out and add those squared features, whereas the SVM handles all that complex mapping automatically behind the scenes.
Q7. Support Vector Approaches on the Auto Dataset
In this problem, you will use support vector approaches in order to predict whether a given car gets high or low gas mileage based on the Auto data set.
(a) Create a binary variable that takes on a 1 for cars with gas mileage above the median, and a 0 for cars with gas mileage below the median.
library(dplyr)
median_mpg <- median(Auto$mpg)
Auto_clean <- Auto %>%
mutate(mpg_binary = as.factor(ifelse(mpg > median(mpg), 1, 0))) %>%
select(-mpg)(b) Fit a support vector classifier to the data with various values of cost, in order to predict whether a car gets high or low gas mileage. Report the cross-validation errors associated with different values of this parameter. Comment on your results. Note you will need to fit the classifier without the gas mileage variable to produce sensible results.
set.seed(210)
svc_auto <- tune(svm, mpg_binary ~ ., data = Auto_clean,
kernel = "linear",
ranges = list(cost = c(0.001, 0.01, 0.1, 1, 5, 10, 100)))
summary(svc_auto)##
## Parameter tuning of 'svm':
##
## - sampling method: 10-fold cross validation
##
## - best parameters:
## cost
## 0.01
##
## - best performance: 0.08948718
##
## - Detailed performance results:
## cost error dispersion
## 1 1e-03 0.12762821 0.05283998
## 2 1e-02 0.08948718 0.06068316
## 3 1e-01 0.09711538 0.07032950
## 4 1e+00 0.09185897 0.05683925
## 5 5e+00 0.10211538 0.06835740
## 6 1e+01 0.09955128 0.06443224
## 7 1e+02 0.11474359 0.05550705The linear model works best when the cost is set right at 0.01, giving us our lowest error rate of about 8.9%. If you drop the cost lower than that, the model underfits and the error jumps up to nearly 13%. If you push the cost higher, it starts to overfit, causing the error rate to bounce around between 9% and 11.5% across the rest of the values.
(c) Now repeat (b), this time using SVMs with radial and polynomial basis kernels, with different values of gamma and degree and cost. Comment on your results.
set.seed(210)
svm_radial <- tune(svm, mpg_binary ~ ., data = Auto_clean, kernel = "radial",
ranges = list(cost = c(0.1, 1, 10, 100),
gamma = c(0.01, 0.1, 1, 5)))
summary(svm_radial)##
## Parameter tuning of 'svm':
##
## - sampling method: 10-fold cross validation
##
## - best parameters:
## cost gamma
## 10 0.1
##
## - best performance: 0.07416667
##
## - Detailed performance results:
## cost gamma error dispersion
## 1 0.1 0.01 0.11237179 0.04735888
## 2 1.0 0.01 0.08948718 0.06068316
## 3 10.0 0.01 0.08685897 0.06532122
## 4 100.0 0.01 0.09705128 0.06808378
## 5 0.1 0.10 0.08948718 0.05565999
## 6 1.0 0.10 0.08948718 0.06419310
## 7 10.0 0.10 0.07416667 0.06215859
## 8 100.0 0.10 0.09698718 0.05636819
## 9 0.1 1.00 0.55865385 0.04971543
## 10 1.0 1.00 0.07923077 0.05966338
## 11 10.0 1.00 0.07910256 0.05704384
## 12 100.0 1.00 0.07910256 0.05704384
## 13 0.1 5.00 0.55865385 0.04971543
## 14 1.0 5.00 0.49237179 0.06830555
## 15 10.0 5.00 0.48974359 0.07644916
## 16 100.0 5.00 0.48974359 0.07644916set.seed(210)
svm_poly <- tune(svm, mpg_binary ~ ., data = Auto_clean, kernel = "polynomial",
ranges = list(cost = c(0.1, 1, 10, 100),
degree = c(2, 3, 4)))
summary(svm_poly)##
## Parameter tuning of 'svm':
##
## - sampling method: 10-fold cross validation
##
## - best parameters:
## cost degree
## 100 2
##
## - best performance: 0.3112179
##
## - Detailed performance results:
## cost degree error dispersion
## 1 0.1 2 0.5586538 0.04971543
## 2 1.0 2 0.5586538 0.04971543
## 3 10.0 2 0.5100000 0.11956315
## 4 100.0 2 0.3112179 0.08596398
## 5 0.1 3 0.5586538 0.04971543
## 6 1.0 3 0.5586538 0.04971543
## 7 10.0 3 0.5586538 0.04971543
## 8 100.0 3 0.4316026 0.08983006
## 9 0.1 4 0.5586538 0.04971543
## 10 1.0 4 0.5586538 0.04971543
## 11 10.0 4 0.5586538 0.04971543
## 12 100.0 4 0.5586538 0.04971543Our svm_radial model wins by hitting the lowest overall error rate at 0.07923077 when the cost is 1.0 and gamma is 1.00. It performs really well across most settings, though if you set a low cost of 0.1 with a high gamma of 1.0, it completely breaks down to a 55.9% error rate. On the other side, our svm_poly model absolutely struggles with this dataset. In almost every single configuration across degrees 2, 3, and 4, it completely stalls at that exact same 55.9% error.
(d) Make some plots to back up your assertions in (b) and (c). 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)where svmfit contains your fitted model and dat is a data frame containing your data, you can type
> plot(svmfit , dat , x1 ∼ x4)in order to plot just the first and fourth variables. However, you must replace x1 and x4 with the correct variable names. To find out more, type ?plot.svm.
best_linear <- svc_auto$best.model
best_radial <- svm_radial$best.model
best_poly <- svm_poly$best.model
slice_values <- list(
cylinders = median(Auto_clean$cylinders),
displacement = median(Auto_clean$displacement),
acceleration = median(Auto_clean$acceleration),
year = median(Auto_clean$year),
origin = median(as.numeric(Auto_clean$origin))
)
plot(best_linear, data = Auto_clean, formula = horsepower ~ weight,
slice = slice_values, main = "Linear SVC (HP vs Weight)")
plot(best_radial, data = Auto_clean, formula = horsepower ~ weight,
slice = slice_values, main = "Radial SVM (HP vs Weight)")
plot(best_poly, data = Auto_clean, formula = horsepower ~ weight,
slice = slice_values, main = "Polynomial SVM (HP vs Weight)")
Q8. Kernel Performance Comparison on the OJ Dataset
This problem involves the OJ data set which is part of the ISLR2 package.
## 'data.frame': 1070 obs. of 18 variables:
## $ Purchase : Factor w/ 2 levels "CH","MM": 1 1 1 2 1 1 1 1 1 1 ...
## $ WeekofPurchase: num 237 239 245 227 228 230 232 234 235 238 ...
## $ StoreID : num 1 1 1 1 7 7 7 7 7 7 ...
## $ PriceCH : num 1.75 1.75 1.86 1.69 1.69 1.69 1.69 1.75 1.75 1.75 ...
## $ PriceMM : num 1.99 1.99 2.09 1.69 1.69 1.99 1.99 1.99 1.99 1.99 ...
## $ DiscCH : num 0 0 0.17 0 0 0 0 0 0 0 ...
## $ DiscMM : num 0 0.3 0 0 0 0 0.4 0.4 0.4 0.4 ...
## $ SpecialCH : num 0 0 0 0 0 0 1 1 0 0 ...
## $ SpecialMM : num 0 1 0 0 0 1 1 0 0 0 ...
## $ LoyalCH : num 0.5 0.6 0.68 0.4 0.957 ...
## $ SalePriceMM : num 1.99 1.69 2.09 1.69 1.69 1.99 1.59 1.59 1.59 1.59 ...
## $ SalePriceCH : num 1.75 1.75 1.69 1.69 1.69 1.69 1.69 1.75 1.75 1.75 ...
## $ PriceDiff : num 0.24 -0.06 0.4 0 0 0.3 -0.1 -0.16 -0.16 -0.16 ...
## $ Store7 : Factor w/ 2 levels "No","Yes": 1 1 1 1 2 2 2 2 2 2 ...
## $ PctDiscMM : num 0 0.151 0 0 0 ...
## $ PctDiscCH : num 0 0 0.0914 0 0 ...
## $ ListPriceDiff : num 0.24 0.24 0.23 0 0 0.3 0.3 0.24 0.24 0.24 ...
## $ STORE : num 1 1 1 1 0 0 0 0 0 0 ...(a) Create a training set containing a random sample of 800 observations, and a test set containing the remaining observations.
set.seed(210)
train_index_oj <- sample(1:nrow(OJ), 800)
train_set_oj<- OJ[train_index_oj,]
test_set_oj <- OJ[-train_index_oj, ]
dim(train_set_oj)## [1] 800 18## [1] 270 18(b) Fit a support vector classifier to the training data using cost = 0.01, with Purchase as the response and the other variables as predictors. Use the summary() function to produce summary statistics, and describe the results obtained.
svm_model_oj <- svm(Purchase ~ ., kernel="linear", data = train_set_oj, cost = 0.01)
summary(svm_model_oj)##
## Call:
## svm(formula = Purchase ~ ., data = train_set_oj, kernel = "linear",
## cost = 0.01)
##
##
## Parameters:
## SVM-Type: C-classification
## SVM-Kernel: linear
## cost: 0.01
##
## Number of Support Vectors: 441
##
## ( 221 220 )
##
##
## Number of Classes: 2
##
## Levels:
## CH MMThis output shows a Support Vector Classifier model trained to categorize customer purchases as CH or MM. It uses a linear kernel with a low cost of 0.01, which creates a rigid, wider margin that tolerates more training errors to prevent overfitting. To establish this decision boundary, the algorithm relies on 441 support vectors, split almost perfectly down the middle with 221 points from one class and 220 from the other.
(c) What are the training and test error rates?
train_pred <-predict(svm_model_oj, train_set_oj)
train_error <- mean(train_pred !=train_set_oj$Purchase)
cat("Training Error:", train_error)## Training Error: 0.17125test_pred <-predict(svm_model_oj, test_set_oj)
test_error <- mean(test_pred !=test_set_oj$Purchase)
cat("Test Error:", test_error)## Test Error: 0.1666667(d) Use the tune() function to select an optimal cost. Consider values in the range 0.01 to 10.
set.seed(10)
tune_linear <- tune(svm, Purchase ~ ., data = train_set_oj, kernel = "linear",
ranges = list(cost = c(0.01, 0.1, 0.5, 1, 5, 10)))(e) Compute the training and test error rates using this new value for cost.
best_linear <- tune_linear$best.model
pred_train_lin <- predict(best_linear, train_set_oj)
mean(pred_train_lin != train_set_oj$Purchase)## [1] 0.16375## [1] 0.1592593(f) Repeat parts (b) through (e) using a support vector machine with a radial kernel. Use the default value for gamma.
svm_radial_oj <- svm(Purchase ~ ., data = train_set_oj, kernel = "radial", cost = 0.01)
mean(predict(svm_radial_oj, train_set_oj) != train_set_oj$Purchase)## [1] 0.39375## [1] 0.3777778tune_radial <- tune(svm, Purchase ~ ., data = train_set_oj, kernel = "radial",
ranges = list(cost = c(0.01, 0.1, 0.5, 1, 5, 10)))
best_radial <- tune_radial$best.model
mean(predict(best_radial, train_set_oj) != train_set_oj$Purchase)## [1] 0.15625## [1] 0.1518519(g) Repeat parts (b) through (e) using a support vector machine with a polynomial kernel. Set degree = 2.
svm_poly_oj <- svm(Purchase ~ ., data = train_set_oj, kernel = "polynomial", degree = 2, cost = 0.01)
mean(predict(svm_poly_oj, train_set_oj) != train_set_oj$Purchase)## [1] 0.39375## [1] 0.3777778tune_poly <- tune(svm, Purchase ~ ., data = train_set_oj, kernel = "polynomial", degree = 2,
ranges = list(cost = c(0.01, 0.1, 0.5, 1, 5, 10)))
best_poly <- tune_poly$best.model
mean(predict(best_poly, train_set_oj) != train_set_oj$Purchase)## [1] 0.15375## [1] 0.162963(h) Overall, which approach seems to give the best results on this data?
## Warning: package 'DT' was built under R version 4.5.3err_svm_model_oj <- mean(predict(svm_model_oj, test_set_oj) != test_set_oj$Purchase)
err_tune_linear <- mean(predict(tune_linear$best.model, test_set_oj) != test_set_oj$Purchase)
err_svm_radial_oj <- mean(predict(svm_radial_oj, test_set_oj) != test_set_oj$Purchase)
err_tune_radial <- mean(predict(tune_radial$best.model, test_set_oj) != test_set_oj$Purchase)
err_svm_poly_oj <- mean(predict(svm_poly_oj, test_set_oj) != test_set_oj$Purchase)
err_tune_poly <- mean(predict(tune_poly$best.model, test_set_oj) != test_set_oj$Purchase)
summary_df <- data.frame(
Model = c("Linear (Untuned)", "Linear (Tuned)",
"Radial (Untuned)", "Radial (Tuned)",
"Poly d=2 (Untuned)", "Poly d=2 (Tuned)"),
Test_Error = c(err_svm_model_oj, err_tune_linear,
err_svm_radial_oj, err_tune_radial,
err_svm_poly_oj, err_tune_poly)
)
datatable(summary_df, rownames = FALSE) %>%
formatRound(columns = 'Test_Error', digits = 4)Based on the test error rates, the Tuned Radial SVM (`) performed best, achieving the lowest misclassification rate (15.19%), closely followed by the Tuned Linear model (15.93%). Tuning the models was critical for performance; when left untuned with a low cost of 0.01, both the Radial and Polynomial models severely underperformed with a high error rate of 37.78%. Overall, the Tuned Radial approach provides the most accurate and reliable predictions for this dataset.