#install.packages("GGally")
#install.packages("e1071")
#install.packages("modelr")
library(MASS)
library(ISLR)
library(GGally)
## Loading required package: ggplot2
library(class)
library(e1071)
##
## Attaching package: 'e1071'
## The following object is masked from 'package:ggplot2':
##
## element
library(modelr)
This question should be answered using the Weekly data set, which is part of the ISLR2 package. This data is similar in nature to the Smarket data from this chapter’s lab, except that it contains 1,089 weekly returns for 21 years, from the beginning of 1990 to the end of 2010.
# loading in dataset
market_data <- Weekly
summary(market_data)
## 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
##
##
##
##
Based on the summary of the data set, Lags 1 through 5 and today have the same minimum and maximum return. We can also see that according to the data, there has been more up (positive) return then down (negative) returns.
ggpairs(
market_data,
ggplot2::aes(colour = Direction),
progress = FALSE
)
## `stat_bin()` using `bins = 30`. Pick better value `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value `binwidth`.
Based on the correlation scatterplot, we can see that Year and Volume are statistically significant and have a positive correlation of 0.842. Another thing to note is that the graph shows an upward trend in the points, meaning that as the years increase so does the volume of shares traded.
market_glm <- glm(Direction ~ Lag1 + Lag2 + Lag3 + Lag4 + Lag5 + Volume,
data = market_data,
family = "binomial")
summary(market_glm)
##
## Call:
## glm(formula = Direction ~ Lag1 + Lag2 + Lag3 + Lag4 + Lag5 +
## Volume, family = "binomial", data = market_data)
##
## 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
From the summary of the logistic regression, there is only one variable that is statistically significant: Lag2.
market_preds <- predict(market_glm , type = "response")
market_glm.pred <- rep("Down", 1089)
market_glm.pred[market_preds > .5] = "Up"
table(market_glm.pred, market_data$Direction)
##
## market_glm.pred Down Up
## Down 54 48
## Up 430 557
mean(market_glm.pred == market_data$Direction)
## [1] 0.5610652
From the confusion matrix we can see that the logistic regression model was able to correctly predict the direction of the market 56% of the time. In the confusion matrix, it shows that the model predicted the market to be up 987 weeks out of the 1,089 weeks in the data set. This could be any issue as it is showing the model has a high false positive rate. On the other hand, the model correctly predicted the market being up when it actually was up 92% of the time (557 / (557 + 48)).
market_train <- (market_data$Year < 2009)
market_test <- market_data[!market_train,]
direction_test <- market_data$Direction[!market_train]
train_glm <- glm(Direction ~ Lag2,
data = market_data,
family = binomial,
subset = market_train)
train_glm_pred <- predict(train_glm, market_test, type = "response")
train_glm.pred <- rep("Down", 104)
train_glm.pred[train_glm_pred > .5] = "Up"
table(train_glm.pred, direction_test)
## direction_test
## train_glm.pred Down Up
## Down 9 5
## Up 34 56
mean(train_glm.pred == direction_test)
## [1] 0.625
market_lda <- lda(Direction ~ Lag2,
data = market_data,
subset = market_train)
market_lda
## Call:
## lda(Direction ~ Lag2, data = market_data, subset = market_train)
##
## Prior probabilities of groups:
## Down Up
## 0.4477157 0.5522843
##
## Group means:
## Lag2
## Down -0.03568254
## Up 0.26036581
##
## Coefficients of linear discriminants:
## LD1
## Lag2 0.4414162
lda_pred <- predict(market_lda, market_test)
lda_class <- lda_pred$class
table(lda_class, direction_test)
## direction_test
## lda_class Down Up
## Down 9 5
## Up 34 56
mean(lda_class == direction_test)
## [1] 0.625
market_qda <- qda(Direction ~ Lag2,
data = market_data,
subset = market_train)
market_qda
## Call:
## qda(Direction ~ Lag2, data = market_data, subset = market_train)
##
## Prior probabilities of groups:
## Down Up
## 0.4477157 0.5522843
##
## Group means:
## Lag2
## Down -0.03568254
## Up 0.26036581
qda_preds <- predict(market_qda, market_test)
qda_class <- qda_preds$class
table(qda_class, direction_test)
## direction_test
## qda_class Down Up
## Down 0 0
## Up 43 61
mean(qda_class == direction_test)
## [1] 0.5865385
train_knn <- as.matrix(market_data$Lag2[market_train])
test_knn <- as.matrix(market_test$Lag2)
train_Direction <- market_data$Direction[market_train]
set.seed(10)
knn_pred <- knn(train = train_knn,
test = test_knn,
cl = train_Direction,
k = 1)
table(knn_pred, direction_test)
## direction_test
## knn_pred Down Up
## Down 21 29
## Up 22 32
mean(knn_pred == direction_test)
## [1] 0.5096154
nb_fit <- naiveBayes(Direction ~ Lag2,
data = market_data,
subset = market_train)
nb_fit
##
## Naive Bayes Classifier for Discrete Predictors
##
## Call:
## naiveBayes.default(x = X, y = Y, laplace = laplace)
##
## A-priori probabilities:
## Y
## Down Up
## 0.4477157 0.5522843
##
## Conditional probabilities:
## Lag2
## Y [,1] [,2]
## Down -0.03568254 2.199504
## Up 0.26036581 2.317485
nb_preds <- predict(nb_fit, market_test)
table(nb_preds, direction_test)
## direction_test
## nb_preds Down Up
## Down 0 0
## Up 43 61
mean(nb_preds == direction_test)
## [1] 0.5865385
| Test Name | Accuracy |
|---|---|
| Logistic Regression | 0.625 |
| LDA | 0.625 |
| QDA | 0.587 |
| KNN with K = 1 | 0.510 |
| Naive Bayes | 0.587 |
Of these models, logistic regression and LDA had the bust accuracy in prediction the direction of the market.
full_lda <- lda(Direction ~ Lag1 + Lag2 + Lag3 + Lag4 + Lag5 + Volume,
data = market_data,
subset = market_train)
full_lda_pred <- predict(full_lda, market_test)
full_lda_class <- full_lda_pred$class
table(full_lda_class, direction_test)
## direction_test
## full_lda_class Down Up
## Down 31 44
## Up 12 17
mean(full_lda_class == direction_test)
## [1] 0.4615385
interaction_qda <- qda(Direction ~ Lag2:Volume,
data = market_data,
subset = market_train)
interaction_qda_preds <- predict(interaction_qda, market_test)
interaction_qda_class <- interaction_qda_preds$class
table(interaction_qda_class, direction_test)
## direction_test
## interaction_qda_class Down Up
## Down 24 28
## Up 19 33
mean(interaction_qda_class == direction_test)
## [1] 0.5480769
train_knn_5 <- as.matrix(market_data$Lag2[market_train])
test_knn_5 <- as.matrix(market_test$Lag2)
train_Direction_5 <- market_data$Direction[market_train]
set.seed(10)
knn_pred_5 <- knn(train = train_knn_5,
test = test_knn_5,
cl = train_Direction_5,
k = 5)
table(knn_pred_5, direction_test)
## direction_test
## knn_pred_5 Down Up
## Down 16 20
## Up 27 41
mean(knn_pred_5 == direction_test)
## [1] 0.5480769
nb_fit_more_pred <- naiveBayes(Direction ~ Lag2 + Volume + Year,
data = market_data,
subset = market_train)
nb_more_preds <- predict(nb_fit_more_pred, market_test)
table(nb_more_preds, direction_test)
## direction_test
## nb_more_preds Down Up
## Down 43 59
## Up 0 2
mean(nb_more_preds == direction_test)
## [1] 0.4326923
Updated Accuracies:
| Test Name | Accuracy |
|---|---|
| Logistic Regression | 0.625 |
| LDA | 0.461 |
| QDA | 0.548 |
| KNN with K = 5 | 0.548 |
| Naive Bayes | 0.433 |
After trying different variations of interactions, adding more predictors, and changing the K value on the LDA, QDA, KNN, and Naive Bayes model; the model performing the best is a tie between QDA and KNN.
In this problem, you will develop a model to predict whether a given car gets high or low gas mileage based on the Auto data set.
# called in data set
auto <- read.csv("Auto.csv", header = TRUE)
# removed any rows where horsepower was ?
clean_auto <- auto[auto$horsepower != "?",]
# converting horsepower to a numeric value not that all ? have been removed
clean_auto$horsepower <- as.numeric(clean_auto$horsepower)
# found median of mpg
mpg_median <- median(auto$mpg)
# added new column to auto data set. Anything over the median is listed as 1
auto$mpg01 <- ifelse(auto$mpg > mpg_median, 1, 0)
auto$mpg01 <- as.factor(auto$mpg01)
plot(auto)
auto |>
ggplot(
aes(
x = mpg01,
y = displacement,
fill = mpg01
)
) +
geom_boxplot()
auto |>
ggplot(
aes(
x = mpg01,
y = weight,
fill = mpg01
)
) +
geom_boxplot()
auto |>
ggplot(
aes(
x = mpg01,
y = acceleration,
fill = mpg01
)
) +
geom_boxplot()
Based on the scatter plot, it seems that displacement, weight, and
acceleration would be good predictors for mpg01. According to the box
plots, of the three displacement and weight show the most distance
between 0 and 1 making them the most useful predictors.
set.seed(42)
partitions <- resample_partition(auto, c(train = 0.75, test = 0.25))
auto_train <- as.data.frame(partitions$train)
auto_test <- as.data.frame(partitions$test)
auto_lda <- lda(mpg01 ~ displacement + weight + acceleration,
data = auto_train)
auto_lda_pred <- predict(auto_lda, auto_test)
auto_lda_class <- auto_lda_pred$class
# Test Error of the model
mean(auto_lda_class != auto_test$mpg01)
## [1] 0.11
auto_qda <- qda(mpg01 ~ displacement + weight + acceleration,
data = auto_train)
auto_qda_pred <- predict(auto_qda, auto_test)
auto_qda_class <- auto_qda_pred$class
# Test Error of the model
mean(auto_qda_class != auto_test$mpg01)
## [1] 0.1
auto_glm <- glm(mpg01 ~ displacement + weight + acceleration,
data = auto_train,
family = "binomial")
auto_glm_preds <- predict(auto_glm, newdata = auto_test, type = "response")
auto_glm.preds <- rep(0, nrow(auto_test))
auto_glm.preds[auto_glm_preds > .5] = 1
mean(auto_glm.preds != auto_test$mpg01)
## [1] 0.14
auto_nb <- naiveBayes(mpg01 ~ displacement + weight + acceleration,
data = auto_train)
auto_nb_preds <- predict(auto_nb, auto_test)
mean(auto_nb_preds != auto_test$mpg01)
## [1] 0.1
preds <- c("displacement", "weight", "acceleration")
auto_train_x <- auto_train[, preds]
auto_test_x <- auto_test[, preds]
auto_knn_1 <- knn(train = auto_train_x,
test = auto_test_x,
cl = auto_train$mpg01,
k = 1)
auto_knn_3 <- knn(train = auto_train_x,
test = auto_test_x,
cl = auto_train$mpg01,
k = 3)
auto_knn_5 <- knn(train = auto_train_x,
test = auto_test_x,
cl = auto_train$mpg01,
k = 5)
auto_knn_10 <- knn(train = auto_train_x,
test = auto_test_x,
cl = auto_train$mpg01,
k = 10)
mean(auto_knn_1 != auto_test$mpg01)
## [1] 0.14
mean(auto_knn_3 != auto_test$mpg01)
## [1] 0.12
mean(auto_knn_5 != auto_test$mpg01)
## [1] 0.13
mean(auto_knn_10 != auto_test$mpg01)
## [1] 0.16
| KNN Value | Test Erro |
|---|---|
| K = 1 | 0.14 |
| K = 3 | 0.12 |
| K = 5 | 0.13 |
| K = 10 | 0.15 |
Of the KNN tests ran, when K = 3 it had the lowest test error.
Using the Boston data set, fit classification models in order to predict whether a given census tract has a crime rate above or below the median. Explore logistic regression, LDA, naive Bayes, and KNN models using various subsets of the predictors. Describe your findings.
Hint: You will have to create the response variable yourself, using the variables that are contained in the Boston data set.
boston_data <- Boston
median_crime <- median(boston_data$crim)
boston_data$high_crime <- ifelse(boston_data$crim > median_crime, 1, 0)
full_boston_glm <- glm(high_crime ~ ., data = boston_data)
summary(full_boston_glm)
##
## Call:
## glm(formula = high_crime ~ ., data = boston_data)
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.4811899 0.3596766 -4.118 4.48e-05 ***
## crim 0.0011350 0.0022288 0.509 0.610817
## zn -0.0013810 0.0009315 -1.483 0.138837
## indus 0.0030875 0.0041258 0.748 0.454618
## chas -0.0145084 0.0583661 -0.249 0.803793
## nox 2.0025143 0.2618145 7.649 1.08e-13 ***
## rm 0.0213894 0.0303113 0.706 0.480736
## age 0.0027612 0.0008862 3.116 0.001941 **
## dis 0.0111019 0.0141047 0.787 0.431598
## rad 0.0171097 0.0045460 3.764 0.000188 ***
## tax -0.0002067 0.0002550 -0.810 0.418098
## ptratio 0.0125743 0.0092372 1.361 0.174055
## black -0.0002421 0.0001824 -1.328 0.184948
## lstat 0.0035227 0.0037541 0.938 0.348518
## medv 0.0094707 0.0030244 3.131 0.001843 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for gaussian family taken to be 0.1013345)
##
## Null deviance: 126.500 on 505 degrees of freedom
## Residual deviance: 49.755 on 491 degrees of freedom
## AIC: 294.34
##
## Number of Fisher Scoring iterations: 2
# creating test and train
set.seed(42)
boston_partitions <- resample_partition(boston_data, c(train = 0.75, test = 0.25))
boston_train <- as.data.frame(boston_partitions$train)
boston_test <- as.data.frame(boston_partitions$test)
# Logistic Regression
boston_glm <- glm(high_crime ~ nox + age + rad + medv,
data = boston_train,
family = "binomial")
boston_glm_preds <- predict(boston_glm, newdata = boston_test, type = "response")
boston_glm.preds <- rep(0, nrow(boston_test))
boston_glm.preds[boston_glm_preds > .5] = 1
mean(boston_glm.preds == boston_test$high_crime)
## [1] 0.8818898
# LDA
boston_lda <- lda(high_crime ~ nox + age + rad + medv,
data = boston_train)
boston_lda_pred <- predict(boston_lda, boston_test)
boston_lda_class <- boston_lda_pred$class
# Test Error of the model
mean(boston_lda_class == boston_test$high_crime)
## [1] 0.8897638
# Naive Bayes
boston_nb <- naiveBayes(high_crime ~ nox + age + rad + medv,
data = boston_train)
boston_nb_preds <- predict(boston_nb, boston_test)
mean(boston_nb_preds == boston_test$high_crime)
## [1] 0.8897638
boston_preds <- c("nox", "age", "rad", "medv")
boston_train_x <- boston_train[, boston_preds]
boston_test_x <- boston_test[, boston_preds]
boston_knn_1 <- knn(train = boston_train_x,
test = boston_test_x,
cl = boston_train$high_crime,
k = 1)
boston_knn_3 <- knn(train = boston_train_x,
test = boston_test_x,
cl = boston_train$high_crime,
k = 3)
boston_knn_5 <- knn(train = boston_train_x,
test = boston_test_x,
cl = boston_train$high_crime,
k = 5)
boston_knn_10 <- knn(train = boston_train_x,
test = boston_test_x,
cl = boston_train$high_crime,
k = 10)
mean(boston_knn_1 == boston_test$high_crime)
## [1] 0.7795276
mean(boston_knn_3 == boston_test$high_crime)
## [1] 0.8503937
mean(boston_knn_5 == boston_test$high_crime)
## [1] 0.8740157
mean(boston_knn_10 == boston_test$high_crime)
## [1] 0.8897638