Bài tập 1_ tuần 10

Tên: Mai Huy

MSSV: 43.01.104.065

Số thứ tự: 08

Câu a) Remove the observations for whom the salary information is unknown, and then log-transform the salaries.

# Load thư viện Mass
library(MASS)
# Load thư viện ISLR
library(ISLR)

# Bỏ những quan sát mà thông tin bị trống và cập nhập lại tập dữ liệu Hitlers mới
Hitters <- na.omit(Hitters)
# Trả về cột salary sau khi được log
Hitters$Salary <- log(Hitters$Salary)

Câu b) Create a training set consisting of the first 200 observations, and a test set consisting of the remaining observations.

# Tập train gồm 200 quan sát đầu
train <- 1:200
# Lấy 200 quan sát đầu trong tập dữ liệu Hitters làm tập train
Hitters.train <- Hitters[train, ]
# Lấy các quan sát còn lại ngoài tập train làm tập test
Hitters.test <- Hitters[-train, ]

câu c) Perform boosting on the training set with 1,000 trees for a range of values of the shrinkage parameter A. Produce a plot with different shrinkage values on the x-axis and the corresponding training set MSE on the y-axis.

# Load thư viện gbm
library(gbm)
# set.seed dùng để tái tạo những vector random giống nhau theo tương ứng với giá trị được đưa vào hàm seed
set.seed(1)
# Tạo một dãy số với số đầu là -10, số cuối là -0.2, giá trị mỗi bước nhảy là 0.1 
pows <- seq(-10, -0.2, by = 0.1)
# lambdas là một vector chứa giá trị 10^x với x là các giá trị trong vector pows
lambdas <- 10^pows
# Tạo ra một mảng trống với độ dài bằng vector lambdas
train.err <- rep(NA, length(lambdas))

for (i in 1:length(lambdas)) {
    boost.hitters <- gbm(Salary ~ ., data = Hitters.train, distribution = "gaussian", n.trees = 1000, shrinkage = lambdas[i])
    pred.train <- predict(boost.hitters, Hitters.train, n.trees = 1000)
    train.err[i] <- mean((pred.train - Hitters.train$Salary)^2)
}
plot(lambdas, train.err, type = "b", xlab = "Shrinkage values", ylab = "Training MSE")

Với mỗi giá trị trong vector lambdas, ta gán cho mô hình một giá trị learning rate khác nhau, sau đó tiến hành dự đoán và tính MSE trên tập huấn luyện, lưu lại từng giá trị MSE ứng với mỗi learning rate đó trong mảng trống train.err

Sau đó ta biểu diễn trên biểu đồ với x là giá trị learning của mô hình, y là giá trị MSE trên tập huấn luyện, ta thấy với learning rate khoảng gần 0.65 thì mô hình cho kết quả MSE trên tập train là nhỏ nhất

Câu d) Produce a plot with different shrinkage values on the x-axis and the corresponding test set MSE on the y-axis.

# set.seed dùng để tái tạo những vector random giống nhau theo tương ứng với giá trị được đưa vào hàm seed
set.seed(1)
# Tạo ra một mảng trống với độ dài bằng vector lambdas để lưu các giá trị MSE tương ứng mỗi learning rate khác nhau
test.err <- rep(NA, length(lambdas))
for (i in 1:length(lambdas)) {
    boost.hitters <- gbm(Salary ~ ., data = Hitters.train, distribution = "gaussian", n.trees = 1000, shrinkage = lambdas[i])
    yhat <- predict(boost.hitters, Hitters.test, n.trees = 1000)
    test.err[i] <- mean((yhat - Hitters.test$Salary)^2)
}
plot(lambdas, test.err, type = "b", xlab = "Shrinkage values", ylab = "Test MSE")

Với mỗi giá trị trong vector lambdas, ta gán cho mô hình một giá trị learning rate khác nhau, sau đó tiến hành dự đoán và tính MSE trên tập test, lưu lại từng giá trị MSE ứng với mỗi learning rate đó trong mảng trống test.err

Sau đó ta biểu diễn trên biểu đồ với x là giá trị learning của mô hình, y là giá trị MSE trên tập test

# Giá trị MSE nhỏ nhất trên tập test
min(test.err)

[1] 0.2540265

# Giá trị learning khiến cho MSE trên tập test là nhỏ nhất 
lambdas[which.min(test.err)]

[1] 0.07943282

Câu e) Compare the test MSE of boosting to the test MSE that results from applying two of the regression approaches seen in Chapters 3 and 6.

# Load thư viện glmnet
library(glmnet)

# Tạo mô hình hồi quy bội với Salary là biến đầu ra và các biến lại là biến đầu vào, tập huấn luyện là Hitters.train
fit1 <- lm(Salary ~ ., data = Hitters.train)
# Tiến hành dự đoán trên tập test
pred1 <- predict(fit1, Hitters.test)
# Tính Mean Squared Error (MSE) giữa giá trị dự đoán và giá trị thật
mean((pred1 - Hitters.test$Salary)^2)

[1] 0.4917959

MSE khi dùng hồi quy bội là 0.4917959

# Tạo 1 một ma trận mô hình với biến đầu ra là Salary, biến đầu vào là các biến loại, tập dữ liệu là tập Hitters.train
x <- model.matrix(Salary ~ ., data = Hitters.train)
# Tạo 1 một ma trận mô hình với biến đầu ra là Salary, biến đầu vào là các biến loại, tập dữ liệu là tập Hitters.test
x.test <- model.matrix(Salary ~ ., data = Hitters.test)
# Lấy biến Salary trên tập Hitters.train
y <- Hitters.train$Salary
# Tạo mô hình ridge regression với x là biến đầu vào, đầu ra là y, alpha =0 nhằm sử dụng ridge penalty
fit2 <- glmnet(x, y, alpha = 0)
# Tiến hành dự đoán trên tập x.test với s=0.01 là giá trị áp dụng cho ridge penalty
pred2 <- predict(fit2, s = 0.01, newx = x.test)
# Tính Mean Squared Error (MSE) giữa giá trị dự đoán và giá trị thật\
mean((pred2 - Hitters.test$Salary)^2)

[1] 0.4570283

MSE khi dùng ridge regression là 0.4570283

Ta thấy là MSE khi dùng mô hình boosting thấp hơn cả linear regression và ridge regression

Câu f)Which variables appear to be the most important predictors in the boosted model?

# Xây dụng mô hình boosting với biến đầu ra là Salary và các biến còn lại là biến đầu vào, mô hình dùng phân phối gaussian, số lượng cây = 1000, learning rate bằng giá trị mà cho ra MSE nhỏ nhất trên tập test
boost.hitters <- gbm(Salary ~ ., data = Hitters.train, distribution = "gaussian", n.trees = 1000, shrinkage = lambdas[which.min(test.err)])
# Hiển thị biểu dô mức quan trọng các biến
summary(boost.hitters)

Ta thấy biến CAtBat là biến quan trọng nhất là mức độ ảnh hưởng là 22.9336528

Câu g) Now apply bagging to the training set. What is the test set MSE for this approach?

# set.seed dùng để tái tạo những vector random giống nhau theo tương ứng với giá trị được đưa vào hàm seed
set.seed(1)
# Tạo mô hình randomforest để dự đoan giá trị Salary là biến đầu ra và lấy tất cả các biến còn lại làm biến đầu vào trong tập Hitters.train
bag.hitters <- randomForest(Salary ~ ., data = Hitters.train, mtry = 19, ntree = 500)

• mtry =13: có 19 yếu tố đầu vào được sử dụng cho mỗi lúc phân chia cây
• ntree = 500: số lượng cây là 500

# Đưa ra dự đoán trên tập test ,biến yhat.bag trả về các giá trị dự đoán trên tập test
yhat.bag <- predict(bag.hitters, newdata = Hitters.test)
# Tính Mean Squared Error (MSE) giữa giá trị dự đoán và giá trị thật
mean((yhat.bag - Hitters.test$Salary)^2)

[1] 0.2299324

Ta thấy MSE của mô hình bagging là 0.2299324, nhỏ hơn so với mô hình boosting(0.2540265)