HW 3
NastaranTaabodi -dew050
2026-07-14
library(ISLR2)
data(Weekly)
head(Weekly)## Year Lag1 Lag2 Lag3 Lag4 Lag5 Volume Today Direction
## 1 1990 0.816 1.572 -3.936 -0.229 -3.484 0.1549760 -0.270 Down
## 2 1990 -0.270 0.816 1.572 -3.936 -0.229 0.1485740 -2.576 Down
## 3 1990 -2.576 -0.270 0.816 1.572 -3.936 0.1598375 3.514 Up
## 4 1990 3.514 -2.576 -0.270 0.816 1.572 0.1616300 0.712 Up
## 5 1990 0.712 3.514 -2.576 -0.270 0.816 0.1537280 1.178 Up
## 6 1990 1.178 0.712 3.514 -2.576 -0.270 0.1544440 -1.372 Downsummary(Weekly)## Year Lag1 Lag2 Lag3
## Min. :1990 Min. :-18.1950 Min. :-18.1950 Min. :-18.1950
## 1st Qu.:1995 1st Qu.: -1.1540 1st Qu.: -1.1540 1st Qu.: -1.1580
## Median :2000 Median : 0.2410 Median : 0.2410 Median : 0.2410
## Mean :2000 Mean : 0.1506 Mean : 0.1511 Mean : 0.1472
## 3rd Qu.:2005 3rd Qu.: 1.4050 3rd Qu.: 1.4090 3rd Qu.: 1.4090
## Max. :2010 Max. : 12.0260 Max. : 12.0260 Max. : 12.0260
## Lag4 Lag5 Volume Today
## Min. :-18.1950 Min. :-18.1950 Min. :0.08747 Min. :-18.1950
## 1st Qu.: -1.1580 1st Qu.: -1.1660 1st Qu.:0.33202 1st Qu.: -1.1540
## Median : 0.2380 Median : 0.2340 Median :1.00268 Median : 0.2410
## Mean : 0.1458 Mean : 0.1399 Mean :1.57462 Mean : 0.1499
## 3rd Qu.: 1.4090 3rd Qu.: 1.4050 3rd Qu.:2.05373 3rd Qu.: 1.4050
## Max. : 12.0260 Max. : 12.0260 Max. :9.32821 Max. : 12.0260
## Direction
## Down:484
## Up :605
##
##
##
## cor(Weekly[, -9])## Year Lag1 Lag2 Lag3 Lag4
## Year 1.00000000 -0.032289274 -0.03339001 -0.03000649 -0.031127923
## Lag1 -0.03228927 1.000000000 -0.07485305 0.05863568 -0.071273876
## Lag2 -0.03339001 -0.074853051 1.00000000 -0.07572091 0.058381535
## Lag3 -0.03000649 0.058635682 -0.07572091 1.00000000 -0.075395865
## Lag4 -0.03112792 -0.071273876 0.05838153 -0.07539587 1.000000000
## Lag5 -0.03051910 -0.008183096 -0.07249948 0.06065717 -0.075675027
## Volume 0.84194162 -0.064951313 -0.08551314 -0.06928771 -0.061074617
## Today -0.03245989 -0.075031842 0.05916672 -0.07124364 -0.007825873
## Lag5 Volume Today
## Year -0.030519101 0.84194162 -0.032459894
## Lag1 -0.008183096 -0.06495131 -0.075031842
## Lag2 -0.072499482 -0.08551314 0.059166717
## Lag3 0.060657175 -0.06928771 -0.071243639
## Lag4 -0.075675027 -0.06107462 -0.007825873
## Lag5 1.000000000 -0.05851741 0.011012698
## Volume -0.058517414 1.00000000 -0.033077783
## Today 0.011012698 -0.03307778 1.000000000plot(Weekly$Volume, type = "l", xlab = "Index (time order)", ylab = "Volume",
main = "Weekly Trading Volume Over Time")
glm.fit <- glm(Direction ~ Lag1 + Lag2 + Lag3 + Lag4 + Lag5 + Volume,
data = Weekly, family = binomial)
summary(glm.fit)##
## Call:
## glm(formula = Direction ~ Lag1 + Lag2 + Lag3 + Lag4 + Lag5 +
## Volume, family = binomial, data = Weekly)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.26686 0.08593 3.106 0.0019 **
## Lag1 -0.04127 0.02641 -1.563 0.1181
## Lag2 0.05844 0.02686 2.175 0.0296 *
## Lag3 -0.01606 0.02666 -0.602 0.5469
## Lag4 -0.02779 0.02646 -1.050 0.2937
## Lag5 -0.01447 0.02638 -0.549 0.5833
## Volume -0.02274 0.03690 -0.616 0.5377
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1496.2 on 1088 degrees of freedom
## Residual deviance: 1486.4 on 1082 degrees of freedom
## AIC: 1500.4
##
## Number of Fisher Scoring iterations: 4#Among the six predictors, only Lag2 is statistically significant (p = 0.0296), with a positive coefficient (0.0584) suggesting that a positive return two weeks prior is associated with a higher probability of the market going up this week. All other predictors (Lag1, Lag3, Lag4, Lag5, Volume) are not statistically significant at the 0.05 level
glm.probs <- predict(glm.fit, type = "response")
glm.pred <- rep("Down", length(glm.probs))
glm.pred[glm.probs > 0.5] <- "Up"
table(Predicted = glm.pred, Actual = Weekly$Direction)## Actual
## Predicted Down Up
## Down 54 48
## Up 430 557mean(glm.pred == Weekly$Direction)## [1] 0.5610652train <- (Weekly$Year <= 2008)
Weekly.test <- Weekly[!train, ]
Direction.test <- Weekly$Direction[!train]
glm.fit.d <- glm(Direction ~ Lag2, data = Weekly, family = binomial, subset = train)
summary(glm.fit.d)##
## Call:
## glm(formula = Direction ~ Lag2, family = binomial, data = Weekly,
## subset = train)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.20326 0.06428 3.162 0.00157 **
## Lag2 0.05810 0.02870 2.024 0.04298 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1354.7 on 984 degrees of freedom
## Residual deviance: 1350.5 on 983 degrees of freedom
## AIC: 1354.5
##
## Number of Fisher Scoring iterations: 4glm.probs.d <- predict(glm.fit.d, Weekly.test, type = "response")
glm.pred.d <- rep("Down", length(glm.probs.d))
glm.pred.d[glm.probs.d > 0.5] <- "Up"
table(Predicted = glm.pred.d, Actual = Direction.test)## Actual
## Predicted Down Up
## Down 9 5
## Up 34 56mean(glm.pred.d == Direction.test)## [1] 0.625#62.5% test accuracy, with the model still favoring “Up” predictions but performing much better than the full model.
library(MASS)##
## Attaching package: 'MASS'## The following object is masked from 'package:ISLR2':
##
## Bostonlda.fit <- lda(Direction ~ Lag2, data = Weekly, subset = train)
lda.pred <- predict(lda.fit, Weekly.test)
table(Predicted = lda.pred$class, Actual = Direction.test)## Actual
## Predicted Down Up
## Down 9 5
## Up 34 56mean(lda.pred$class == Direction.test)## [1] 0.625qda.fit <- qda(Direction ~ Lag2, data = Weekly, subset = train)
qda.pred <- predict(qda.fit, Weekly.test)
table(Predicted = qda.pred$class, Actual = Direction.test)## Actual
## Predicted Down Up
## Down 0 0
## Up 43 61mean(qda.pred$class == Direction.test)## [1] 0.5865385#QDA predicts “Up” for every test observation, giving 58.65% accuracy (equal to the proportion of actual Up weeks
library(class)
train.X <- as.matrix(Weekly$Lag2[train])
test.X <- as.matrix(Weekly$Lag2[!train])
train.Direction <- Weekly$Direction[train]
set.seed(1)
knn.pred <- knn(train.X, test.X, train.Direction, k = 1)
table(Predicted = knn.pred, Actual = Direction.test)## Actual
## Predicted Down Up
## Down 21 30
## Up 22 31mean(knn.pred == Direction.test)## [1] 0.5#Accuracy = 0.50 (50%), very close to but not identical to the Python result (49.04%) — this small difference comes from how R and Python each break ties randomly when multiple training points are equally close to a test point. Like before, KNN with K=1 performs poorly here, close to a coin flip, since it’s highly sensitive to noise in the weak Lag2-Direction relationship
library(e1071)
nb.fit <- naiveBayes(Direction ~ Lag2, data = Weekly, subset = train)
nb.pred <- predict(nb.fit, Weekly.test)
table(Predicted = nb.pred, Actual = Direction.test)## Actual
## Predicted Down Up
## Down 0 0
## Up 43 61mean(nb.pred == Direction.test)## [1] 0.5865385#Accuracy = 0.5865 (58.65%), matching the Python result exactly — Naive Bayes predicts “Up” for every test observation, just like QDA, effectively defaulting to the majority class
#part (i): Which method performs best?Logistic regression and LDA tie for best performance (62.5%), followed by QDA and Naive Bayes (58.65%, both defaulting to majority-class predictions), with KNN (K=1) performing worst (50%) due to its sensitivity to noise. This matches the pattern found in the Python analysis
for (k in c(1, 3, 5, 10, 15, 20, 30, 50)) {
set.seed(1)
knn.pred.k <- knn(train.X, test.X, train.Direction, k = k)
acc <- mean(knn.pred.k == Direction.test)
cat("K =", k, ": Accuracy =", round(acc, 4), "\n")
}## K = 1 : Accuracy = 0.5
## K = 3 : Accuracy = 0.5481
## K = 5 : Accuracy = 0.5385
## K = 10 : Accuracy = 0.5481
## K = 15 : Accuracy = 0.5865
## K = 20 : Accuracy = 0.5865
## K = 30 : Accuracy = 0.5481
## K = 50 : Accuracy = 0.5577#Accuracy generally improves as K increases from 1 (50%) up to a peak of 58.65% at K=15 and K=20, similar to the pattern seen in Python (where K=20 peaked at 59.6%). Even at its best, KNN doesn’t beat logistic regression/LDA (62.5%), confirming that a simple linear boundary using Lag2 works better than KNN’s flexible, local approach on this weak-signal data
question 14:
data(Auto)
median.mpg <- median(Auto$mpg)
Auto$mpg01 <- as.numeric(Auto$mpg > median.mpg)
head(Auto)## mpg cylinders displacement horsepower weight acceleration year origin
## 1 18 8 307 130 3504 12.0 70 1
## 2 15 8 350 165 3693 11.5 70 1
## 3 18 8 318 150 3436 11.0 70 1
## 4 16 8 304 150 3433 12.0 70 1
## 5 17 8 302 140 3449 10.5 70 1
## 6 15 8 429 198 4341 10.0 70 1
## name mpg01
## 1 chevrolet chevelle malibu 0
## 2 buick skylark 320 0
## 3 plymouth satellite 0
## 4 amc rebel sst 0
## 5 ford torino 0
## 6 ford galaxie 500 0table(Auto$mpg01)##
## 0 1
## 196 196par(mfrow = c(2, 3))
boxplot(cylinders ~ mpg01, data = Auto, main = "Cylinders")
boxplot(displacement ~ mpg01, data = Auto, main = "Displacement")
boxplot(horsepower ~ mpg01, data = Auto, main = "Horsepower")
boxplot(weight ~ mpg01, data = Auto, main = "Weight")
boxplot(acceleration ~ mpg01, data = Auto, main = "Acceleration")
boxplot(year ~ mpg01, data = Auto, main = "Year")
#These match the Python results exactly — cylinders, displacement,
horsepower, and weight all show strong, clear separation between mpg01
groups (mpg01=1 cars have fewer cylinders, lower displacement, less
horsepower, lower weight), while acceleration and year show weaker
separation. We’ll use cylinders, displacement, horsepower, and weight as
our predictors going forward
set.seed(1)
train.idx <- sample(1:nrow(Auto), size = 0.7 * nrow(Auto))
Auto.train <- Auto[train.idx, ]
Auto.test <- Auto[-train.idx, ]
mpg01.test <- Auto$mpg01[-train.idx]
dim(Auto.train)## [1] 274 10dim(Auto.test)## [1] 118 10#274 training rows, 118 test rows, matching the 70/30 split from Python
lda.fit.auto <- lda(mpg01 ~ cylinders + displacement + horsepower + weight,
data = Auto.train)
lda.pred.auto <- predict(lda.fit.auto, Auto.test)
table(Predicted = lda.pred.auto$class, Actual = mpg01.test)## Actual
## Predicted 0 1
## 0 50 3
## 1 11 54mean(lda.pred.auto$class == mpg01.test)## [1] 0.8813559Accuracy = 0.8814 , Test Error = 0.1186 (11.9%). This is somewhat worse than the Python result (8.5% error) — likely due to a different random train/test split (R’s sample() and Python’s train_test_split() produce different splits even with similar seeds/methods). Still, LDA performs reasonably well, with a bit more false positives (11) than false negatives (3)
qda.fit.auto <- qda(mpg01 ~ cylinders + displacement + horsepower + weight,
data = Auto.train)
qda.pred.auto <- predict(qda.fit.auto, Auto.test)
table(Predicted = qda.pred.auto$class, Actual = mpg01.test)## Actual
## Predicted 0 1
## 0 52 5
## 1 9 52mean(qda.pred.auto$class == mpg01.test)## [1] 0.8813559#Accuracy = 0.8814 → Test Error = 0.1186 (11.9%), identical accuracy to LDA but with a slightly different error pattern (9 false positives, 5 false negatives vs. LDA’s 11/3). Unlike the Python run (where QDA outperformed LDA), here they tie exactly — again likely due to the different train/test split
glm.fit.auto <- glm(mpg01 ~ cylinders + displacement + horsepower + weight,
data = Auto.train, family = binomial)
glm.probs.auto <- predict(glm.fit.auto, Auto.test, type = "response")
glm.pred.auto <- rep(0, length(glm.probs.auto))
glm.pred.auto[glm.probs.auto > 0.5] <- 1
table(Predicted = glm.pred.auto, Actual = mpg01.test)## Actual
## Predicted 0 1
## 0 53 3
## 1 8 54mean(glm.pred.auto == mpg01.test)## [1] 0.9067797#Accuracy = 0.9068 → Test Error = 0.0932 (9.3%). This is the best result so far among LDA/QDA/logistic regression on this R split, with a fairly balanced error pattern (8 false positives, 3 false negatives)
nb.fit.auto <- naiveBayes(mpg01 ~ cylinders + displacement + horsepower + weight,
data = Auto.train)
nb.pred.auto <- predict(nb.fit.auto, Auto.test)
table(Predicted = nb.pred.auto, Actual = mpg01.test)## Actual
## Predicted 0 1
## 0 52 4
## 1 9 53mean(nb.pred.auto == mpg01.test)## [1] 0.8898305#Accuracy = 0.8898 → Test Error = 0.1102 (11.0%), landing between LDA/QDA (11.9%) and logistic regression (9.3%). A reasonably balanced error pattern (9 false positives, 4 false negatives
train.X.auto <- scale(Auto.train[, c("cylinders", "displacement", "horsepower", "weight")])
test.X.auto <- scale(Auto.test[, c("cylinders", "displacement", "horsepower", "weight")],
center = attr(train.X.auto, "scaled:center"),
scale = attr(train.X.auto, "scaled:scale"))
train.mpg01 <- Auto.train$mpg01
for (k in c(1, 3, 5, 10, 15, 20, 30)) {
set.seed(1)
knn.pred.auto <- knn(train.X.auto, test.X.auto, train.mpg01, k = k)
acc <- mean(knn.pred.auto == mpg01.test)
cat("K =", k, ": Accuracy =", round(acc, 4), "\n")
}## K = 1 : Accuracy = 0.8644
## K = 3 : Accuracy = 0.8729
## K = 5 : Accuracy = 0.8644
## K = 10 : Accuracy = 0.8644
## K = 15 : Accuracy = 0.8898
## K = 20 : Accuracy = 0.8898
## K = 30 : Accuracy = 0.8898#Best KNN accuracy is 0.8898 (11.0% test error), achieved at K=15, 20, and 30 — tying with Naive Bayes but still slightly worse than logistic regression’s 9.3% error. Unlike K=1 sensitivity issues seen in Question 13, here K=1 already performs reasonably well (86.4%), and accuracy stabilizes/improves slightly with larger K
##Logistic regression performs best on this train/test split, though all methods perform reasonably well (all above 88% accuracy), confirming that cylinders, displacement, horsepower, and weight are strong predictors of mpg01 regardless of method. This differs slightly from the Python run (where QDA won), which is expected given the different random train/test split — a good reminder that these small differences in performance between methods aren’t hugely meaningful with a single train/test split; results can shift somewhat with different splits
##question 16 :
Boston <- ISLR2::Boston
median.crim <- median(Boston$crim)
Boston$crim01 <- as.numeric(Boston$crim > median.crim)
table(Boston$crim01)##
## 0 1
## 253 253par(mfrow = c(3, 4), mar = c(4, 4, 2, 1))
boxplot(zn ~ crim01, data = Boston, main = "zn")
boxplot(indus ~ crim01, data = Boston, main = "indus")
boxplot(chas ~ crim01, data = Boston, main = "chas")
boxplot(nox ~ crim01, data = Boston, main = "nox")
boxplot(rm ~ crim01, data = Boston, main = "rm")
boxplot(age ~ crim01, data = Boston, main = "age")
boxplot(dis ~ crim01, data = Boston, main = "dis")
boxplot(rad ~ crim01, data = Boston, main = "rad")
boxplot(tax ~ crim01, data = Boston, main = "tax")
boxplot(ptratio ~ crim01, data = Boston, main = "ptratio")
boxplot(lstat ~ crim01, data = Boston, main = "lstat")
boxplot(medv ~ crim01, data = Boston, main = "medv")
#Great, that worked. These match the Python results — nox, age, dis, rad, tax, and indus show the strongest separation between low-crime and high-crime suburbs (higher crime associated with more industrial land, more pollution, older housing, closer proximity to employment centers, greater highway access, and higher taxes), while chas, rm, ptratio, lstat, and medv show weaker separation. We’ll use nox, age, dis, rad, tax, indus as our predictors going forward.
set.seed(1)
train.idx.b <- sample(1:nrow(Boston), size = 0.7 * nrow(Boston))
Boston.train <- Boston[train.idx.b, ]
Boston.test <- Boston[-train.idx.b, ]
crim01.test <- Boston$crim01[-train.idx.b]
dim(Boston.train)## [1] 354 14dim(Boston.test)## [1] 152 14#354 training rows, 152 test rows, matching the split size from Python
glm.fit.b <- glm(crim01 ~ nox + age + dis + rad + tax + indus,
data = Boston.train, family = binomial)
glm.probs.b <- predict(glm.fit.b, Boston.test, type = "response")
glm.pred.b <- rep(0, length(glm.probs.b))
glm.pred.b[glm.probs.b > 0.5] <- 1
table(Predicted = glm.pred.b, Actual = crim01.test)## Actual
## Predicted 0 1
## 0 58 6
## 1 15 73mean(glm.pred.b == crim01.test)## [1] 0.8618421Accuracy = 0.8618 → Test Error = 0.1382 (13.8%). Reasonably good performance, with more false negatives (15 — high-crime suburbs misclassified as low-crime) than false positives (6), similar to the mild bias seen in the Python run
lda.fit.b <- lda(crim01 ~ nox + age + dis + rad + tax + indus,
data = Boston.train)
lda.pred.b <- predict(lda.fit.b, Boston.test)
table(Predicted = lda.pred.b$class, Actual = crim01.test)## Actual
## Predicted 0 1
## 0 71 20
## 1 2 59mean(lda.pred.b$class == crim01.test)## [1] 0.8552632#Accuracy = 0.8553 → Test Error = 0.1447 (14.5%), close to logistic regression but slightly worse, and much better than what we saw in the Python run (where LDA had 17.8% error). Same pattern though — more false negatives (20) than false positives (2), suggesting LDA still tends to miss actual high-crime suburbs
nb.fit.b <- naiveBayes(crim01 ~ nox + age + dis + rad + tax + indus,
data = Boston.train)
nb.pred.b <- predict(nb.fit.b, Boston.test)
table(Predicted = nb.pred.b, Actual = crim01.test)## Actual
## Predicted 0 1
## 0 65 19
## 1 8 60mean(nb.pred.b == crim01.test)## [1] 0.8223684#Accuracy = 0.8224 → Test Error = 0.1776 (17.8%), worse than both logistic regression (13.8%) and LDA (14.5%), consistent with the pattern from Python where Naive Bayes underperformed due to violated independence/normality assumptions among correlated predictors like tax, rad, and indus
train.X.b <- scale(Boston.train[, c("nox", "age", "dis", "rad", "tax", "indus")])
test.X.b <- scale(Boston.test[, c("nox", "age", "dis", "rad", "tax", "indus")],
center = attr(train.X.b, "scaled:center"),
scale = attr(train.X.b, "scaled:scale"))
train.crim01 <- Boston.train$crim01
for (k in c(1, 3, 5, 10, 15, 20)) {
set.seed(1)
knn.pred.b <- knn(train.X.b, test.X.b, train.crim01, k = k)
acc <- mean(knn.pred.b == crim01.test)
cat("K =", k, ": Accuracy =", round(acc, 4), "\n")
}## K = 1 : Accuracy = 0.9079
## K = 3 : Accuracy = 0.9276
## K = 5 : Accuracy = 0.9276
## K = 10 : Accuracy = 0.9079
## K = 15 : Accuracy = 0.9145
## K = 20 : Accuracy = 0.9013#Best KNN accuracy is 0.9276 (7.24% test error), achieved at K=3 and K=5 — the best result of any method in this R analysis of Question 16, consistent with the Python findings where KNN also performed best. This confirms that the relationship between these predictors and crime classification is highly local/non-linear, which KNN captures better than the linear or normality-based methods. ##. Overall conclusion: KNN clearly outperforms the parametric methods (logistic regression, LDA, Naive Bayes) on this dataset, suggesting the boundary between high-crime and low-crime suburbs is non-linear and doesn’t fit well with methods assuming linear boundaries or normal/independent predictor distributions. Naive Bayes performs worst, likely due to violated independence assumptions among correlated predictors like tax, rad, and indus.