Khai báo thư viện
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(ggplot2)
library(corrplot)
## corrplot 0.95 loaded
library(MASS)
##
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
##
## select
library(lmtest)
## Loading required package: zoo
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
library(car)
## Loading required package: carData
##
## Attaching package: 'car'
## The following object is masked from 'package:dplyr':
##
## recode
data <- read.csv("California_housing.csv")
head(data)
str(data)
## 'data.frame': 20640 obs. of 10 variables:
## $ longitude : num -122 -122 -122 -122 -122 ...
## $ latitude : num 37.9 37.9 37.9 37.9 37.9 ...
## $ housing_median_age: num 41 21 52 52 52 52 52 52 42 52 ...
## $ total_rooms : num 880 7099 1467 1274 1627 ...
## $ total_bedrooms : num 129 1106 190 235 280 ...
## $ population : num 322 2401 496 558 565 ...
## $ households : num 126 1138 177 219 259 ...
## $ median_income : num 8.33 8.3 7.26 5.64 3.85 ...
## $ median_house_value: num 452600 358500 352100 341300 342200 ...
## $ ocean_proximity : chr "NEAR BAY" "NEAR BAY" "NEAR BAY" "NEAR BAY" ...
Bộ dữ liệu có 9 biến định lượng và 1 biến định tính là ocean_proximit. Kiểm tra xem biến định tính này có bao nhiêu biến độc lập.
length(unique(data$ocean_proximity))
## [1] 5
unique(data$ocean_proximity)
## [1] "NEAR BAY" "<1H OCEAN" "INLAND" "NEAR OCEAN" "ISLAND"
Biến ocean proximity có 5 giá trị riêng biệt, biến đổi biến này thành kiểu dữ liệu factor để hồi quy.
data$ocean_proximity <- as.factor(data$ocean_proximity)
str(data)
## 'data.frame': 20640 obs. of 10 variables:
## $ longitude : num -122 -122 -122 -122 -122 ...
## $ latitude : num 37.9 37.9 37.9 37.9 37.9 ...
## $ housing_median_age: num 41 21 52 52 52 52 52 52 42 52 ...
## $ total_rooms : num 880 7099 1467 1274 1627 ...
## $ total_bedrooms : num 129 1106 190 235 280 ...
## $ population : num 322 2401 496 558 565 ...
## $ households : num 126 1138 177 219 259 ...
## $ median_income : num 8.33 8.3 7.26 5.64 3.85 ...
## $ median_house_value: num 452600 358500 352100 341300 342200 ...
## $ ocean_proximity : Factor w/ 5 levels "<1H OCEAN","INLAND",..: 4 4 4 4 4 4 4 4 4 4 ...
Kiểm tra giá trị khuyết trong 20.640 quang trắc:
colSums(is.na(data))
## longitude latitude housing_median_age total_rooms
## 0 0 0 0
## total_bedrooms population households median_income
## 207 0 0 0
## median_house_value ocean_proximity
## 0 0
Biến total_bedrooms có 207 diểm dữ liệu bị khuyết, chiếm 1% bộ dữ liệu, nhóm chọn bỏ các quang trắc này.
data_dropna <- na.omit(data)
colSums(!is.na(data_dropna))
## longitude latitude housing_median_age total_rooms
## 20433 20433 20433 20433
## total_bedrooms population households median_income
## 20433 20433 20433 20433
## median_house_value ocean_proximity
## 20433 20433
###Thống kê mô tả của các biến
summary(data_dropna)
## longitude latitude housing_median_age total_rooms
## Min. :-124.3 Min. :32.54 Min. : 1.00 Min. : 2
## 1st Qu.:-121.8 1st Qu.:33.93 1st Qu.:18.00 1st Qu.: 1450
## Median :-118.5 Median :34.26 Median :29.00 Median : 2127
## Mean :-119.6 Mean :35.63 Mean :28.63 Mean : 2637
## 3rd Qu.:-118.0 3rd Qu.:37.72 3rd Qu.:37.00 3rd Qu.: 3143
## Max. :-114.3 Max. :41.95 Max. :52.00 Max. :39320
## total_bedrooms population households median_income
## Min. : 1.0 Min. : 3 Min. : 1.0 Min. : 0.4999
## 1st Qu.: 296.0 1st Qu.: 787 1st Qu.: 280.0 1st Qu.: 2.5637
## Median : 435.0 Median : 1166 Median : 409.0 Median : 3.5365
## Mean : 537.9 Mean : 1425 Mean : 499.4 Mean : 3.8712
## 3rd Qu.: 647.0 3rd Qu.: 1722 3rd Qu.: 604.0 3rd Qu.: 4.7440
## Max. :6445.0 Max. :35682 Max. :6082.0 Max. :15.0001
## median_house_value ocean_proximity
## Min. : 14999 <1H OCEAN :9034
## 1st Qu.:119500 INLAND :6496
## Median :179700 ISLAND : 5
## Mean :206864 NEAR BAY :2270
## 3rd Qu.:264700 NEAR OCEAN:2628
## Max. :500001
numeric_cols <- sapply(data_dropna, is.numeric)
for (col in names(data_dropna)[numeric_cols]) {
print(
ggplot(data_dropna, aes(x = .data[[col]])) +
geom_histogram(bins = 50, fill = "skyblue", color = "white") +
labs(title = paste("Histogram of", col), x = col)
)
}
Trong bộ dữ liệu này, ta có thể dễ dàng nhận thấy rằng 2 biến median_house_value và housing_median_age có giá trị max nhiều bất thườn. Kiểm tra thấy rằng có 2054 quan trắc, chiếm 10% bộ dữ liệu. 10% là con số tương đối đáng kể, nhưng vì bộ dữ liệu có hơn 20.000 quan trắc nên nhóm chọn cách loại bỏ các điểm dữ liệu này đi
max_val_1 <- max(data_dropna$median_house_value)
max_val_2 <- max(data_dropna$housing_median_age)
nrow(data_dropna[data_dropna$median_house_value == max_val_1 | data_dropna$housing_median_age == max_val_2, ])
## [1] 2054
data_dropoutlier <- data_dropna[data_dropna$median_house_value != max_val_1,]
data_dropoutlier <- data_dropoutlier[data_dropoutlier$housing_median_age != max_val_2,]
nrow(data_dropoutlier)
## [1] 18379
ggplot(data_dropoutlier, aes(x = median_house_value)) +
geom_histogram(bins = 50, fill = "skyblue", color = "white") +
labs(title = "Histogram of median_house_value")
# housing_median_age
ggplot(data_dropoutlier, aes(x = housing_median_age)) +
geom_histogram(bins = 50, fill = "lightgreen", color = "white") +
labs(title = "Histogram of housing_median_age")
set.seed(123)
n <- nrow(data_dropoutlier)
train_index <- sample(1:n, size = 0.8 * n)
train_data <- data_dropoutlier[train_index, ]
test_data <- data_dropoutlier[-train_index, ]
Vẽ biểu đồ tương quan giữa các biến, ta có thể dễ dàng nhận thấy rằng, biến longtitude và latitude có mối quan hệ tuyến tính khá mạnh với nhau. Bên cạnh đó các biến total_rooms, total_bedrooms, populations và households cũng có quan hệ tuyến tính với nhau.
num_vars <- sapply(train_data, is.numeric)
df_num <- train_data[, num_vars]
plot(df_num)
Vẽ ma trận tương quan để thấy các biến có quan hệ tuyến tính rõ hơn
numeric_vars <- train_data[sapply(train_data, is.numeric)]
corr_matrix <- cor(numeric_vars, use = "complete.obs")
corrplot(corr_matrix,
method = "color",
type = "upper",
tl.col = "black",
tl.srt = 45,
addCoef.col = "black", # <-- HIỂN THỊ GIÁ TRỊ SỐ
number.cex = 0.7, # kích thước chữ số
col = colorRampPalette(c("blue","white","red"))(200))
Hồi quy biến median_house_value theo 9 biến còn lại
full_model <- lm(median_house_value ~.,data = train_data)
summary(full_model)
##
## Call:
## lm(formula = median_house_value ~ ., data = train_data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -276615 -37328 -8751 26484 576611
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.931e+06 8.649e+04 -22.331 < 2e-16 ***
## longitude -2.256e+04 1.002e+03 -22.523 < 2e-16 ***
## latitude -2.044e+04 9.861e+02 -20.724 < 2e-16 ***
## housing_median_age 7.002e+02 4.902e+01 14.284 < 2e-16 ***
## total_rooms -7.123e+00 8.335e-01 -8.547 < 2e-16 ***
## total_bedrooms 8.766e+01 7.042e+00 12.448 < 2e-16 ***
## population -2.869e+01 1.055e+00 -27.196 < 2e-16 ***
## households 4.053e+01 7.483e+00 5.416 6.18e-08 ***
## median_income 3.844e+04 4.286e+02 89.691 < 2e-16 ***
## ocean_proximityINLAND -4.185e+04 1.737e+03 -24.091 < 2e-16 ***
## ocean_proximityISLAND 1.832e+05 4.189e+04 4.373 1.23e-05 ***
## ocean_proximityNEAR BAY -1.572e+04 2.069e+03 -7.598 3.19e-14 ***
## ocean_proximityNEAR OCEAN 5.741e+03 1.606e+03 3.574 0.000352 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 59200 on 14690 degrees of freedom
## Multiple R-squared: 0.6236, Adjusted R-squared: 0.6233
## F-statistic: 2028 on 12 and 14690 DF, p-value: < 2.2e-16
Kiểm tra vif của model
vif(full_model)
## GVIF Df GVIF^(1/(2*Df))
## longitude 16.516648 1 4.064068
## latitude 18.805012 1 4.336475
## housing_median_age 1.321904 1 1.149741
## total_rooms 14.386444 1 3.792947
## total_bedrooms 37.983306 1 6.163060
## population 6.380520 1 2.525969
## households 35.373700 1 5.947579
## median_income 1.908382 1 1.381442
## ocean_proximity 3.725273 4 1.178677
Vì trong ma trận tương quan, latitude có quan hệ tương quan với biến phụ thuộc cao hơn lontitude nên sẽ loại biến longitude. Biến total_bedrooms cũng có điểm VIF cao nhất nên cũng loại để xem model có được cải thiện không.Tiếp tục loại bỏ biến households vì VIF của biến này vẫn cao. Cuối cùng còn lại 6 biến có VIF khoảng từ 1-2, trong đó có total_rooms với population có VIF > 4.
model_1 <- lm(median_house_value ~. - longitude - total_bedrooms - households, data = train_data)
vif(model_1)
## GVIF Df GVIF^(1/(2*Df))
## latitude 1.467499 1 1.211404
## housing_median_age 1.291398 1 1.136397
## total_rooms 4.887500 1 2.210769
## population 4.525194 1 2.127250
## median_income 1.331484 1 1.153900
## ocean_proximity 1.734392 4 1.071256
summary(model_1)
##
## Call:
## lm(formula = median_house_value ~ . - longitude - total_bedrooms -
## households, data = train_data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -237820 -39587 -10152 27064 372672
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.301e+04 1.043e+04 4.123 3.76e-05 ***
## latitude 7.123e+02 2.883e+02 2.471 0.013490 *
## housing_median_age 6.722e+02 5.071e+01 13.256 < 2e-16 ***
## total_rooms 1.037e+01 5.084e-01 20.405 < 2e-16 ***
## population -1.641e+01 9.296e-01 -17.650 < 2e-16 ***
## median_income 3.377e+04 3.746e+02 90.145 < 2e-16 ***
## ocean_proximityINLAND -7.527e+04 1.399e+03 -53.788 < 2e-16 ***
## ocean_proximityISLAND 1.999e+05 4.382e+04 4.561 5.13e-06 ***
## ocean_proximityNEAR BAY -7.585e+03 2.142e+03 -3.541 0.000399 ***
## ocean_proximityNEAR OCEAN 1.326e+04 1.641e+03 8.083 6.82e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 61950 on 14693 degrees of freedom
## Multiple R-squared: 0.5877, Adjusted R-squared: 0.5874
## F-statistic: 2327 on 9 and 14693 DF, p-value: < 2.2e-16
Chọn mô hình bằng phương pháp stepwise với chuẩn BIC.Sau khi chạy thử, ta thấy rằng mô hình cho BIC nhỏ nhất khi loại bỏ biến latitude đi. V mô hình tốt nhất là mô hình gồm 5 biến, housing_median_age, total_rooms, population, median_income và ocean_proximity
modZero <- lm(formula = median_house_value ~ 1, data = train_data)
best_model_1 <- MASS::stepAIC(model_1, direction = "both", scope = list(lower = modZero, upper = model_1), k = log(nrow(train_data)), trace = 10)
## Start: AIC=324556.5
## median_house_value ~ (longitude + latitude + housing_median_age +
## total_rooms + total_bedrooms + population + households +
## median_income + ocean_proximity) - longitude - total_bedrooms -
## households
##
## Df Sum of Sq RSS AIC
## - latitude 1 2.3434e+10 5.6421e+13 324553
## <none> 5.6398e+13 324556
## - housing_median_age 1 6.7445e+11 5.7072e+13 324722
## - population 1 1.1957e+12 5.7593e+13 324855
## - total_rooms 1 1.5982e+12 5.7996e+13 324958
## - ocean_proximity 4 1.4702e+13 7.1100e+13 327924
## - median_income 1 3.1191e+13 8.7589e+13 331020
##
## Step: AIC=324553
## median_house_value ~ housing_median_age + total_rooms + population +
## median_income + ocean_proximity
##
## Df Sum of Sq RSS AIC
## <none> 5.6421e+13 324553
## + latitude 1 2.3434e+10 5.6398e+13 324556
## - housing_median_age 1 6.6929e+11 5.7090e+13 324717
## - population 1 1.2401e+12 5.7661e+13 324863
## - total_rooms 1 1.6256e+12 5.8047e+13 324961
## - ocean_proximity 4 1.6522e+13 7.2943e+13 328291
## - median_income 1 3.1169e+13 8.7590e+13 331010
###Kiểm định mô hình:
Kiểm định ANOVA: Với mức tin cậy 5%, ta có thể kết luận cả năm biến này đều có ý nghĩa thống kê, đóng góp vào mô hình
anova(best_model_1)
kiểm định xem phân phối của sai số có giống với phân phối chuẩn không. Vì shapiro test chỉ dùng cho dữ liệu < 5000, nên trong bài này nhóm dùng kiểm định Kolmogorov–Smirnov.Với p_value rất nhỏ, nên có thể kết luận rằng phân phối của biến phụ thuộc không giống với phân phối chuẩn.
ks.test(residuals(best_model_1),
"pnorm",
mean(residuals(best_model_1)),
sd(residuals(best_model_1)))
##
## Asymptotic one-sample Kolmogorov-Smirnov test
##
## data: residuals(best_model_1)
## D = 0.092895, p-value < 2.2e-16
## alternative hypothesis: two-sided
Kiểm định phương sai đồng nhất Breusch–Pagan. P-value nhỏ -> phương sai không đồng nhất
bptest(best_model_1)
##
## studentized Breusch-Pagan test
##
## data: best_model_1
## BP = 428.79, df = 8, p-value < 2.2e-16
plot(best_model_1,1)
plot(best_model_1,2)
Biến đổi Box-cox.
boxcox(best_model_1,lambda = seq(0.15, 0.25, length = 10))
Với \(\lambda\) ~ 0.19 mô hình biến đổi
thành
best_model_1 <- lm(median_house_value^0.19 ~ housing_median_age + total_rooms + population + median_income + ocean_proximity, data = train_data)
summary(best_model_1)
##
## Call:
## lm(formula = median_house_value^0.19 ~ housing_median_age + total_rooms +
## population + median_income + ocean_proximity, data = train_data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.3481 -0.4066 -0.0523 0.3469 2.8386
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 8.791e+00 2.533e-02 347.049 < 2e-16 ***
## housing_median_age 4.162e-03 4.902e-04 8.490 < 2e-16 ***
## total_rooms 9.721e-05 4.908e-06 19.806 < 2e-16 ***
## population -1.432e-04 8.946e-06 -16.011 < 2e-16 ***
## median_income 3.314e-01 3.621e-03 91.497 < 2e-16 ***
## ocean_proximityINLAND -9.076e-01 1.226e-02 -74.054 < 2e-16 ***
## ocean_proximityISLAND 1.648e+00 4.237e-01 3.889 0.000101 ***
## ocean_proximityNEAR BAY -7.947e-02 1.872e-02 -4.245 2.20e-05 ***
## ocean_proximityNEAR OCEAN 7.778e-02 1.586e-02 4.903 9.52e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.5991 on 14694 degrees of freedom
## Multiple R-squared: 0.6284, Adjusted R-squared: 0.6282
## F-statistic: 3107 on 8 and 14694 DF, p-value: < 2.2e-16
plot(best_model_1,1)
plot(best_model_1,2)
Loại outlier với cook distance:
cooks <- cooks.distance(best_model_1)
threshold <- 4/nrow(train_data)
influential_points <- which(cooks > threshold)
train_removed_cook <- train_data[-influential_points,]
best_model_1 <- lm(median_house_value^0.19 ~ housing_median_age + total_rooms + population + median_income + ocean_proximity, data = train_removed_cook)
summary(best_model_1)
##
## Call:
## lm(formula = median_house_value^0.19 ~ housing_median_age + total_rooms +
## population + median_income + ocean_proximity, data = train_removed_cook)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.68821 -0.37133 -0.03994 0.33253 2.10089
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 8.812e+00 2.355e-02 374.108 < 2e-16 ***
## housing_median_age 3.410e-03 4.392e-04 7.765 8.74e-15 ***
## total_rooms 1.411e-04 5.276e-06 26.736 < 2e-16 ***
## population -2.167e-04 9.808e-06 -22.089 < 2e-16 ***
## median_income 3.267e-01 3.368e-03 96.983 < 2e-16 ***
## ocean_proximityINLAND -9.487e-01 1.093e-02 -86.807 < 2e-16 ***
## ocean_proximityNEAR BAY -1.291e-01 1.698e-02 -7.601 3.13e-14 ***
## ocean_proximityNEAR OCEAN 9.171e-02 1.424e-02 6.442 1.21e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.5123 on 13969 degrees of freedom
## Multiple R-squared: 0.7088, Adjusted R-squared: 0.7087
## F-statistic: 4858 on 7 and 13969 DF, p-value: < 2.2e-16
plot(best_model_1,4)
kiểm định lại:
ks.test(residuals(best_model_1),
"pnorm",
mean(residuals(best_model_1)),
sd(residuals(best_model_1)))
##
## Asymptotic one-sample Kolmogorov-Smirnov test
##
## data: residuals(best_model_1)
## D = 0.032902, p-value = 1.441e-13
## alternative hypothesis: two-sided
bptest(best_model_1)
##
## studentized Breusch-Pagan test
##
## data: best_model_1
## BP = 324.52, df = 7, p-value < 2.2e-16
Tuy p_value vẫn rất nhỏ nhưng D-statistic đã giảm 1 nửa từ 0.09 xuống còn 0.03 (dữ liệu đã gần phân phối chuẩn hơn) và BP giảm từ 428 xuống còn 324, có cải thiện phương sai đồng nhất. Nhìn chung biến đổi box-cox đã cải thiện phân phối của sai số giống với phân phối chuẩn hơn và tính đồng nhất phương sai đã được cải thiện.
Vậy với hồi quy bội, mô hình được chọn có \(R^2\) = 07088, mô hình giải thích 70 biến phụ thuộc. mô hình có dạng \(y=x^\frac{1}{0.19}\) và các biến độc lập đều có ý nghĩa thống kê.
summary(best_model_1)
##
## Call:
## lm(formula = median_house_value^0.19 ~ housing_median_age + total_rooms +
## population + median_income + ocean_proximity, data = train_removed_cook)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.68821 -0.37133 -0.03994 0.33253 2.10089
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 8.812e+00 2.355e-02 374.108 < 2e-16 ***
## housing_median_age 3.410e-03 4.392e-04 7.765 8.74e-15 ***
## total_rooms 1.411e-04 5.276e-06 26.736 < 2e-16 ***
## population -2.167e-04 9.808e-06 -22.089 < 2e-16 ***
## median_income 3.267e-01 3.368e-03 96.983 < 2e-16 ***
## ocean_proximityINLAND -9.487e-01 1.093e-02 -86.807 < 2e-16 ***
## ocean_proximityNEAR BAY -1.291e-01 1.698e-02 -7.601 3.13e-14 ***
## ocean_proximityNEAR OCEAN 9.171e-02 1.424e-02 6.442 1.21e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.5123 on 13969 degrees of freedom
## Multiple R-squared: 0.7088, Adjusted R-squared: 0.7087
## F-statistic: 4858 on 7 and 13969 DF, p-value: < 2.2e-16
library(glmnet)
## Loading required package: Matrix
## Loaded glmnet 4.1-10
khai báo biến y và biến x
y <- train_removed_cook$median_house_value^0.19
x <- model.matrix(
median_house_value ~ .,
data = train_removed_cook)[, -1]
y_test <- test_data$median_house_value^0.19
x_test <- model.matrix(
median_house_value ~ .,
data = test_data
)[, colnames(x)]
Hồi quy Ridge, cross validation để chọn ra lambda tối ưu
cv_ridge <- cv.glmnet(x, y, alpha = 0)
plot(cv_ridge)
y_pred_ridge <- predict(cv_ridge, newx = x, s = "lambda.min")
r2_ridge <- 1- sum((y-y_pred_ridge)^2) / sum((y-mean(y))^2)
mse_ridge <- mean((y - y_pred_ridge)^2)
mse_ridge
## [1] 0.2460247
r2_ridge
## [1] 0.7269046
Hồi quy Lasso, cross validation để tìm ra lambda tối ưu.
cv_lasso <- cv.glmnet(x, y, alpha = 1)
plot(cv_lasso)
y_pred_lasso <- predict(cv_lasso, newx = x, s = "lambda.min")
r2_lasso <- 1- sum((y-y_pred_lasso)^2) / sum((y-mean(y))^2)
mse_lasso <- mean((y - y_pred_lasso)^2)
mse_lasso
## [1] 0.2363514
r2_lasso
## [1] 0.7376423
Hồi quy Elastic Net
cv_elastic <- cv.glmnet(x, y, alpha = 0.5)
plot(cv_elastic)
y_pred_elastic <- predict(cv_elastic, newx = x, s = "lambda.min")
r2_elastic <- 1- sum((y-y_pred_elastic)^2) / sum((y-mean(y))^2)
mse_elastic <- mean((y - y_pred_elastic)^2)
mse_elastic
## [1] 0.2363743
r2_elastic
## [1] 0.7376168
##Chạy thử mô hình trên tập test
# Predict trên tập test
y_test_pred_trans <- as.numeric(
predict(cv_ridge, newx = x_test, s = "lambda.1se")
)
#Tính RMSE
rmse_test_ridge_trans <- sqrt(mean((y_test - y_test_pred_trans)^2))
#Tính R^2
r2_test_ridge_trans <- 1 -
sum((y_test - y_test_pred_trans)^2) /
sum((y_test - mean(y_test))^2)
rmse_test_ridge_trans
## [1] 0.569859
r2_test_ridge_trans
## [1] 0.6440245
# Predict trên tập test
y_test_pred_trans <- as.numeric(
predict(cv_lasso, newx = x_test, s = "lambda.1se")
)
#Tính RMSE
rmse_test_lasso_trans <- sqrt(mean((y_test - y_test_pred_trans)^2))
#R^2
r2_test_lasso_trans <- 1 -
sum((y_test - y_test_pred_trans)^2) /
sum((y_test - mean(y_test))^2)
rmse_test_lasso_trans
## [1] 0.5616713
r2_test_lasso_trans
## [1] 0.6541803
# Predict trên scale transform
y_test_pred_trans <- as.numeric(
predict(cv_elastic, newx = x_test, s = "lambda.1se")
)
#Tính RMSE
rmse_test_elastic_trans <- sqrt(mean((y_test - y_test_pred_trans)^2))
#Tính R^2
r2_test_elastic_trans <- 1 -
sum((y_test - y_test_pred_trans)^2) /
sum((y_test - mean(y_test))^2)
rmse_test_elastic_trans
## [1] 0.5622278
r2_test_elastic_trans
## [1] 0.6534947
test_data$y_predicted <- predict(best_model_1, newdata = test_data)
test_data$y_predicted <- test_data$y_predicted
test_data
#Kiểm chứng với mô hình best_model
test_data$y_trains <- test_data$median_house_value^0.19
test_data
rmse_test_best_trans <- sqrt(mean((test_data$y_trains - test_data$y_predicted)^2))
r2_test_best_trans <- 1 -
sum((test_data$y_trains - test_data$y_predicted)^2) /
sum((test_data$y_trains - mean(test_data$y_predicted))^2)
rmse_test_best_trans
## [1] 0.591993
r2_test_best_trans
## [1] 0.615842
test_data$y_predicted <- predict(best_model_1, newdata = test_data)
test_data$median_house_value_predicted <- test_data$y_predicted^(1/0.19)
plot(test_data$median_house_value_predicted,test_data$median_house_value)
abline(0, 1, col = "red", lwd = 2)
Do mô hình được huấn luyện để dự đoán trên không gian đã biến đổi
y^(0.19) nên các giá trị dự đoán thu được từ mô hình chưa phải là giá
nhà thực tế. Để chuyển các dự đoán này trở về đơn vị gốc của giá nhà,
cần áp dụng phép biến đổi ngược, tức là y^(1/0.19). Bước này đảm bảo
rằng RMSE và R^2 được tính trên thang đo giá nhà thực. Sau khi áp dụng
phép biến đổi ngược và đánh giá trên tập test, Elastic Net regression
cho kết quả tốt nhất với RMSE thấp nhất và R^2 cao nhất. Điều này cho
thấy Elastic Net không chỉ dự đoán chính xác hơn mà còn tổng quát hóa
tốt, đặc biệt trong bối cảnh dữ liệu có đa cộng tuyến giữa các biến đầu
vào. Do đó, Elastic Net được lựa chọn là mô hình cuối cùng cho bài toán
dự đoán giá nhà.