#question chapter 7
p <- seq(0, 1, by = 0.0001)
Gini <- 2 * p * (1 - p)
ClassError <- 1 - pmax(p, 1 - p)
Entropy <- -(ifelse(p == 0, 0, p * log(p)) +
ifelse(p == 1, 0, (1 - p) * log(1 - p)))
plot(p, Entropy,
type = "l",
lwd = 3,
col = "red",
xlab = expression(p[m1]),
ylab = "Impurity Measure",
main = "Comparison of Classification Impurity Measures",
ylim = c(0, max(Entropy)))
lines(p, ClassError, col = "darkgreen", lwd = 3)
lines(p, Gini, col = "blue", lwd = 3)
legend("top",
legend = c("Entropy", "Classification Error", "Gini Index"),
col = c("red", "darkgreen", "blue"),
lwd = 3,
bty = "n")
library(ISLR2)
## Warning: package 'ISLR2' was built under R version 4.5.3
library(tree)
## Warning: package 'tree' was built under R version 4.5.3
library(randomForest)
## Warning: package 'randomForest' was built under R version 4.5.3
## randomForest 4.7-1.2
## Type rfNews() to see new features/changes/bug fixes.
library(BART)
## Warning: package 'BART' was built under R version 4.5.3
## Loading required package: nlme
## Loading required package: survival
data("Carseats")
set.seed(1)
train_id <- sample(
1:nrow(Carseats),
size = 200
)
carseats_train <- Carseats[train_id, ]
carseats_test <- Carseats[-train_id, ]
nrow(carseats_train)
## [1] 200
nrow(carseats_test)
## [1] 200
set.seed(1)
tree_fit <- tree(
Sales ~ .,
data = carseats_train
)
summary(tree_fit)
##
## Regression tree:
## tree(formula = Sales ~ ., data = carseats_train)
## Variables actually used in tree construction:
## [1] "ShelveLoc" "Price" "Age" "Advertising" "CompPrice"
## [6] "US"
## Number of terminal nodes: 18
## Residual mean deviance: 2.167 = 394.3 / 182
## Distribution of residuals:
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -3.88200 -0.88200 -0.08712 0.00000 0.89590 4.09900
plot(tree_fit)
text(
tree_fit,
pretty = 0,
cex = 0.8
)
tree_pred <- predict(
tree_fit,
newdata = carseats_test
)
tree_mse <- mean(
(carseats_test$Sales - tree_pred)^2
)
tree_mse
## [1] 4.922039
set.seed(1)
cv_tree <- cv.tree(tree_fit)
cv_tree
## $size
## [1] 18 17 16 15 14 13 12 11 10 8 7 6 5 4 3 2 1
##
## $dev
## [1] 984.3936 1031.3372 1036.0021 1027.2166 1027.2166 1055.8168 1044.6955
## [8] 1061.0899 1061.0899 1225.5973 1221.3487 1219.0219 1231.6886 1337.3952
## [15] 1300.0524 1338.3702 1605.0221
##
## $k
## [1] -Inf 16.99544 20.56322 25.01730 25.57104 28.01938 30.36962
## [8] 31.56747 31.80816 40.75445 44.44673 52.57126 76.21881 99.59459
## [15] 116.69889 159.79501 337.60153
##
## $method
## [1] "deviance"
##
## attr(,"class")
## [1] "prune" "tree.sequence"
plot(
cv_tree$size,
cv_tree$dev,
type = "b",
pch = 19,
xlab = "Number of Terminal Nodes",
ylab = "Cross-Validated Deviance",
main = "Cross-Validation for Regression Tree Size"
)
best_size <- cv_tree$size[
which.min(cv_tree$dev)
]
best_size
## [1] 18
pruned_tree <- prune.tree(
tree_fit,
best = best_size
)
summary(pruned_tree)
##
## Regression tree:
## tree(formula = Sales ~ ., data = carseats_train)
## Variables actually used in tree construction:
## [1] "ShelveLoc" "Price" "Age" "Advertising" "CompPrice"
## [6] "US"
## Number of terminal nodes: 18
## Residual mean deviance: 2.167 = 394.3 / 182
## Distribution of residuals:
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -3.88200 -0.88200 -0.08712 0.00000 0.89590 4.09900
plot(pruned_tree)
text(
pruned_tree,
pretty = 0,
cex = 0.8
)
pruned_pred <- predict(
pruned_tree,
newdata = carseats_test
)
pruned_mse <- mean(
(carseats_test$Sales - pruned_pred)^2
)
pruned_mse
## [1] 4.922039
tree_comparison <- data.frame(
Model = c(
"Original regression tree",
"Pruned regression tree"
),
Test_MSE = c(
tree_mse,
pruned_mse
)
)
tree_comparison
## Model Test_MSE
## 1 Original regression tree 4.922039
## 2 Pruned regression tree 4.922039
set.seed(1)
rf_fit <- randomForest(
Sales ~ .,
data = carseats_train,
mtry = 3,
ntree = 500,
importance = TRUE
)
rf_fit
##
## Call:
## randomForest(formula = Sales ~ ., data = carseats_train, mtry = 3, ntree = 500, importance = TRUE)
## Type of random forest: regression
## Number of trees: 500
## No. of variables tried at each split: 3
##
## Mean of squared residuals: 3.363781
## % Var explained: 57.22
rf_pred <- predict(
rf_fit,
newdata = carseats_test
)
rf_mse <- mean(
(carseats_test$Sales - rf_pred)^2
)
rf_mse
## [1] 2.960559
#question 9 chapter 8
library(ISLR2)
library(tree)
data(OJ)
set.seed(1)
train_id <- sample(
1:nrow(OJ),
size = 800
)
OJ_train <- OJ[train_id, ]
OJ_test <- OJ[-train_id, ]
dim(OJ_train)
## [1] 800 18
dim(OJ_test)
## [1] 270 18
oj_tree <- tree(
Purchase ~ .,
data = OJ_train
)
summary(oj_tree)
##
## Classification tree:
## tree(formula = Purchase ~ ., data = OJ_train)
## Variables actually used in tree construction:
## [1] "LoyalCH" "PriceDiff" "SpecialCH" "ListPriceDiff"
## [5] "PctDiscMM"
## Number of terminal nodes: 9
## Residual mean deviance: 0.7432 = 587.8 / 791
## Misclassification error rate: 0.1588 = 127 / 800
oj_tree
## node), split, n, deviance, yval, (yprob)
## * denotes terminal node
##
## 1) root 800 1073.00 CH ( 0.60625 0.39375 )
## 2) LoyalCH < 0.5036 365 441.60 MM ( 0.29315 0.70685 )
## 4) LoyalCH < 0.280875 177 140.50 MM ( 0.13559 0.86441 )
## 8) LoyalCH < 0.0356415 59 10.14 MM ( 0.01695 0.98305 ) *
## 9) LoyalCH > 0.0356415 118 116.40 MM ( 0.19492 0.80508 ) *
## 5) LoyalCH > 0.280875 188 258.00 MM ( 0.44149 0.55851 )
## 10) PriceDiff < 0.05 79 84.79 MM ( 0.22785 0.77215 )
## 20) SpecialCH < 0.5 64 51.98 MM ( 0.14062 0.85938 ) *
## 21) SpecialCH > 0.5 15 20.19 CH ( 0.60000 0.40000 ) *
## 11) PriceDiff > 0.05 109 147.00 CH ( 0.59633 0.40367 ) *
## 3) LoyalCH > 0.5036 435 337.90 CH ( 0.86897 0.13103 )
## 6) LoyalCH < 0.764572 174 201.00 CH ( 0.73563 0.26437 )
## 12) ListPriceDiff < 0.235 72 99.81 MM ( 0.50000 0.50000 )
## 24) PctDiscMM < 0.196196 55 73.14 CH ( 0.61818 0.38182 ) *
## 25) PctDiscMM > 0.196196 17 12.32 MM ( 0.11765 0.88235 ) *
## 13) ListPriceDiff > 0.235 102 65.43 CH ( 0.90196 0.09804 ) *
## 7) LoyalCH > 0.764572 261 91.20 CH ( 0.95785 0.04215 ) *
plot(oj_tree)
text(
oj_tree,
pretty = 0,
cex = 0.7
)
oj_test_pred <- predict(
oj_tree,
newdata = OJ_test,
type = "class"
)
test_confusion <- table(
Predicted = oj_test_pred,
Actual = OJ_test$Purchase
)
test_confusion
## Actual
## Predicted CH MM
## CH 160 38
## MM 8 64
test_error <- mean(
oj_test_pred != OJ_test$Purchase
)
test_error
## [1] 0.1703704