Import thư viện cần thiết
library(lmtest)## Warning: package 'lmtest' was built under R version 4.4.3## Loading required package: zoo## Warning: package 'zoo' was built under R version 4.4.3##
## Attaching package: 'zoo'## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numericlibrary(ggplot2)## Warning: package 'ggplot2' was built under R version 4.4.3library(car)## Warning: package 'car' was built under R version 4.4.3## Loading required package: carData## Warning: package 'carData' was built under R version 4.4.3library(dplyr)## Warning: package 'dplyr' was built under R version 4.4.3##
## Attaching package: 'dplyr'## The following object is masked from 'package:car':
##
## recode## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionlibrary(purrr)## Warning: package 'purrr' was built under R version 4.4.3##
## Attaching package: 'purrr'## The following object is masked from 'package:car':
##
## somelibrary(moments)
library(tidyr)## Warning: package 'tidyr' was built under R version 4.4.3library(corrplot)## Warning: package 'corrplot' was built under R version 4.4.3## corrplot 0.95 loadedlibrary(glmnet)## Warning: package 'glmnet' was built under R version 4.4.3## Loading required package: Matrix##
## Attaching package: 'Matrix'## The following objects are masked from 'package:tidyr':
##
## expand, pack, unpack## Loaded glmnet 4.1-10library(sandwich)## Warning: package 'sandwich' was built under R version 4.4.3library(nortest)
library(tibble)## Warning: package 'tibble' was built under R version 4.4.3df <- read.csv('bmw.csv')Tiền xử lý dữ liệu
# Tổng quan dữ liệu
str(df)## 'data.frame': 10781 obs. of 9 variables:
## $ model : chr " 5 Series" " 6 Series" " 5 Series" " 1 Series" ...
## $ year : int 2014 2018 2016 2017 2014 2016 2017 2018 2017 2016 ...
## $ price : int 11200 27000 16000 12750 14500 14900 16000 16250 14250 14250 ...
## $ transmission: chr "Automatic" "Automatic" "Automatic" "Automatic" ...
## $ mileage : int 67068 14827 62794 26676 39554 35309 38538 10401 42668 36099 ...
## $ fuelType : chr "Diesel" "Petrol" "Diesel" "Diesel" ...
## $ tax : int 125 145 160 145 160 125 125 145 30 20 ...
## $ mpg : num 57.6 42.8 51.4 72.4 50.4 60.1 60.1 52.3 62.8 68.9 ...
## $ engineSize : num 2 2 3 1.5 3 2 2 1.5 2 2 ...Dữ liệu có 10,781 quan trắc với 9 biến.
# Kiểm tra giá trị biến model
sapply(df[, sapply(df, is.character)], function(x) length(unique(x)))## model transmission fuelType
## 24 3 5Bộ dữ liệu gồm 24 mẫu xe khác nhau, đào sâu vào biến model xem ta có thể gom nhóm chúng hay không.
unique(df$model)## [1] " 5 Series" " 6 Series" " 1 Series" " 7 Series" " 2 Series" " 4 Series"
## [7] " X3" " 3 Series" " X5" " X4" " i3" " X1"
## [13] " M4" " X2" " X6" " 8 Series" " Z4" " X7"
## [19] " M5" " i8" " M2" " M3" " M6" " Z3"Có thể thấy có 5 dòng BMW chính ở đây là dòng Series, dòng X, dòng Z, dòng M, dòng i. với:
- Dòng Series: Dòng sedan, coupe, hoặc gran coupe truyền thống.
- Dòng X: Dòng xe thể thao đa dụng (SAV/SAC), gầm cao.
- Dòng Z: Dòng roadster 2 cửa, 2 chỗ ngồi, tập trung vào cảm giác lái thể thao mui trần.
- Dòng M: Phiên bản hiệu suất cao, được tinh chỉnh động cơ.
- Dòng i: Dòng xe điện (BEV) và hybrid cắm điện (PHEV), tập trung công nghệ bền vững.
# Gom nhóm biến model
df <- df %>%
mutate(model_group = case_when(
grepl("Series", model) ~ "Series",
grepl("X", model) ~ "SUV (X)",
grepl("Z", model) ~ "Roadster (Z)",
grepl("i", model) ~ "Electric (i)",
grepl("M", model) ~ "Performance (M)",
TRUE ~ "Other"
))
table(df$model_group)##
## Electric (i) Performance (M) Roadster (Z) Series SUV (X)
## 60 210 115 7945 2451# Kiểm tra giá trị biến transmission
table(df$transmission)##
## Automatic Manual Semi-Auto
## 3588 2527 4666# Kiểm tra giá trị biến fuelType
table(df$fuelType)##
## Diesel Electric Hybrid Other Petrol
## 7027 3 298 36 3417# Kiểm tra giá trị biến model
table(df$model)##
## 1 Series 2 Series 3 Series 4 Series 5 Series 6 Series 7 Series 8 Series
## 1969 1229 2443 995 1056 108 106 39
## i3 i8 M2 M3 M4 M5 M6 X1
## 43 17 21 27 125 29 8 804
## X2 X3 X4 X5 X6 X7 Z3 Z4
## 288 551 179 468 106 55 7 108# Giá trị Electric của biến FuelType chỉ có 3 quan trắc nên nhóm quyết định bỏ đi
df <- df[df$fuelType != "Electric", ]# Kiểm tra biến year
table(df$year)##
## 1996 1997 1998 1999 2000 2001 2002 2003 2004 2005 2006 2007 2008 2009 2010 2011
## 1 1 1 4 2 3 6 2 12 6 14 16 23 30 41 51
## 2012 2013 2014 2015 2016 2017 2018 2019 2020
## 119 357 501 921 1881 1720 848 3485 733Vì biến year có ít quan trắc ở các năm cũ và cũng nhiều giá trị, nên nhóm chọn xét biến year ở dạng numeric
# Đổi kiểu dữ liệu biến định tính thành factor
df <- df %>%
mutate(across(where(is.character), as.factor))
str(df)## 'data.frame': 10778 obs. of 10 variables:
## $ model : Factor w/ 24 levels " 1 Series"," 2 Series",..: 5 6 5 1 7 5 5 2 4 5 ...
## $ year : int 2014 2018 2016 2017 2014 2016 2017 2018 2017 2016 ...
## $ price : int 11200 27000 16000 12750 14500 14900 16000 16250 14250 14250 ...
## $ transmission: Factor w/ 3 levels "Automatic","Manual",..: 1 1 1 1 1 1 1 2 2 1 ...
## $ mileage : int 67068 14827 62794 26676 39554 35309 38538 10401 42668 36099 ...
## $ fuelType : Factor w/ 4 levels "Diesel","Hybrid",..: 1 4 1 1 1 1 1 4 1 1 ...
## $ tax : int 125 145 160 145 160 125 125 145 30 20 ...
## $ mpg : num 57.6 42.8 51.4 72.4 50.4 60.1 60.1 52.3 62.8 68.9 ...
## $ engineSize : num 2 2 3 1.5 3 2 2 1.5 2 2 ...
## $ model_group : Factor w/ 5 levels "Electric (i)",..: 4 4 4 4 4 4 4 4 4 4 ...# Kiểm tra dữ liệu khuyết
sum(is.na(df))## [1] 0Không có dữ liệu khuyết
# Tổng quan phân phối các biến định lượng
numeric_vars <- c("price", "year", "mileage", "tax", "mpg", "engineSize")
summary(df[numeric_vars])## price year mileage tax
## Min. : 1200 Min. :1996 Min. : 1 Min. : 0.0
## 1st Qu.: 14950 1st Qu.:2016 1st Qu.: 5528 1st Qu.:135.0
## Median : 20475 Median :2017 Median : 18345 Median :145.0
## Mean : 22735 Mean :2017 Mean : 25499 Mean :131.7
## 3rd Qu.: 27940 3rd Qu.:2019 3rd Qu.: 38208 3rd Qu.:145.0
## Max. :123456 Max. :2020 Max. :214000 Max. :580.0
## mpg engineSize
## Min. : 5.50 Min. :0.000
## 1st Qu.: 45.60 1st Qu.:2.000
## Median : 53.30 Median :2.000
## Mean : 56.28 Mean :2.168
## 3rd Qu.: 62.80 3rd Qu.:2.000
## Max. :470.80 Max. :6.600Ghi nhận thấy có engineSize = 0. EngineSize thực tế chỉ = 0 khi là xe điện vì không có động cơ đốt trong.
n_total <- nrow(df)
n_issue <- sum(df$engineSize == 0 & df$fuelType != "Electric", na.rm = TRUE)
percent_issue <- n_issue / n_total * 100
data.frame(
condition = "engineSize = 0 & fuelType != Electric",
n_rows = n_issue,
percent = percent_issue
)## condition n_rows percent
## 1 engineSize = 0 & fuelType != Electric 45 0.4175172Các quan trắc không phù hợp này chỉ chiếm 0.4% bộ dữ liệu, nhóm chọn loại bỏ khỏi bộ dữ liệu .
df <- df %>%
filter(!(engineSize == 0 & fuelType != "Electric"))# Kiểm tra dữ liệu trùng
sum(duplicated(df))## [1] 117sum(duplicated(df))/nrow(df)## [1] 0.01090096# Có 117 quang trắc bị trùng hoàn toàn, chiếm 1% bộ dữ liệu, nhóm chọn loại bỏ các quan trắc trùng lặp.
df <- df %>% distinct()# Kiểm tra các giá trị ngoại lai
count_outliers_iqr <- function(x) {
q1 <- quantile(x, 0.25, na.rm = TRUE)
q3 <- quantile(x, 0.75, na.rm = TRUE)
iqr <- q3 - q1
lower <- q1 - 1.5 * iqr
upper <- q3 + 1.5 * iqr
sum(x < lower | x > upper, na.rm = TRUE)
}
# Tạo bảng kết quả
outlier_table <- tibble(variable = numeric_vars) %>%
mutate(
n_outliers = map_int(variable, ~ count_outliers_iqr(df[[.x]])),
outlier_percent = n_outliers / nrow(df) * 100
) %>%
arrange(desc(outlier_percent))
outlier_table## # A tibble: 6 × 3
## variable n_outliers outlier_percent
## <chr> <int> <dbl>
## 1 engineSize 4108 38.7
## 2 tax 2692 25.4
## 3 price 470 4.43
## 4 mileage 292 2.75
## 5 mpg 238 2.24
## 6 year 211 1.99# Vẽ boxplot cho các biến định lượng
par(mfrow = c(2, 3)) # chỉnh layout nếu có 6 biến
for (v in numeric_vars) {
boxplot(df[[v]],
main = paste("Boxplot of", v),
ylab = v,
col = "lightgray")
}
par(mfrow = c(1, 1)) # reset layout# Hệ số xác định trước khi xử lý dữ liệu ngoại lai
model <- lm(price ~ . - model_group, data = df)
summary(model)$r.squared## [1] 0.8688227#Xử lý dữ liệu ngoại lai với phương pháp Winsorization cho những biến có % ngoại lại < 5%
handle_outliers <- function(df, column) {
Q1 <- quantile(df[[column]], 0.25, na.rm = TRUE)
Q3 <- quantile(df[[column]], 0.75, na.rm = TRUE)
IQR <- Q3- Q1
lower_bound <- Q1- 1.5 * IQR
upper_bound <- Q3 + 1.5 * IQR
df[[column]] <- ifelse(df[[column]] < lower_bound, lower_bound, df[[column]])
df[[column]] <- ifelse(df[[column]] > upper_bound, upper_bound, df[[column]])
return(df)
}
df <- handle_outliers(df, 'mileage')
df <- handle_outliers(df, 'mpg')
df <- handle_outliers(df, 'year')#Hệ số xác định sau khi xử lý dữ liệu ngoại lại
model <- lm(price ~ . - model_group, data = df)
summary(model)$r.squared## [1] 0.8902131Có thể thấy, R^2 của mô hình hồi quy tuyến tính OLS có cải thiện sau khi xử lý dữ liệu ngoại lai, tăng từ 0.868 lên 0.89
par(mfrow = c(2, 3)) # chỉnh layout nếu có 6 biến
for (v in numeric_vars) {
boxplot(df[[v]],
main = paste("Boxplot of", v),
ylab = v,
col = "lightgray")
}
par(mfrow = c(1, 1)) # reset layoutsummary(df)## model year price transmission
## 3 Series:2430 Min. :2012 Min. : 1200 Automatic:3499
## 1 Series:1957 1st Qu.:2016 1st Qu.: 14889 Manual :2480
## 2 Series:1186 Median :2017 Median : 20335 Semi-Auto:4637
## 5 Series:1051 Mean :2017 Mean : 22712
## 4 Series: 986 3rd Qu.:2019 3rd Qu.: 27950
## X1 : 797 Max. :2020 Max. :123456
## (Other) :2209
## mileage fuelType tax mpg engineSize
## Min. : 1 Diesel:6979 Min. : 0.0 Min. :19.80 Min. :0.60
## 1st Qu.: 5656 Hybrid: 264 1st Qu.:135.0 1st Qu.:45.60 1st Qu.:2.00
## Median :18673 Other : 36 Median :145.0 Median :53.30 Median :2.00
## Mean :25250 Petrol:3337 Mean :131.9 Mean :53.73 Mean :2.18
## 3rd Qu.:38481 3rd Qu.:145.0 3rd Qu.:62.80 3rd Qu.:2.00
## Max. :87719 Max. :580.0 Max. :88.60 Max. :6.60
##
## model_group
## Electric (i) : 24
## Performance (M): 206
## Roadster (Z) : 114
## Series :7863
## SUV (X) :2409
##
## Trực quan hoá dữ liệu
# scatter plot
vars <- c("mileage", "year", "engineSize", "mpg", "tax")
for (v in vars) {
p <- ggplot(df, aes_string(x = v, y = "price")) +
geom_point(alpha = 0.4) +
geom_smooth(method = "lm", se = FALSE, color = "red") +
theme_minimal() +
labs(
title = paste("Price vs", v, "with Regression Line"),
x = v,
y = "Price"
)
print(p)
}## Warning: `aes_string()` was deprecated in ggplot2 3.0.0.
## ℹ Please use tidy evaluation idioms with `aes()`.
## ℹ See also `vignette("ggplot2-in-packages")` for more information.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
# Tương quan giữa price với các biến định tính
cat_vars <- c("model_group", "fuelType", "transmission","model")
for (v in cat_vars) {
p <- ggplot(df, aes_string(x = v, y = "price")) +
geom_boxplot(fill = "lightgreen", outlier.alpha = 0.4) +
theme_minimal() +
labs(
title = paste("Price by", v),
x = v,
y = "Price"
) +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
print(p)
}
# Tương quan giữa tax với các biến định tính
cat_vars <- c("model_group", "fuelType", "transmission","model")
for (v in cat_vars) {
p <- ggplot(df, aes_string(x = v, y = "tax")) +
geom_boxplot(fill = "lightgreen", outlier.alpha = 0.4) +
theme_minimal() +
labs(
title = paste("Annual road tax by", v),
x = v,
y = "Annual road tax"
) +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
print(p)
}
Các dòng xe điện có mức thuế hằng năm thấp hơn so với các dòng còn
lại.
# Tương quan giữa mpg với các biến định tính
cat_vars <- c("model_group", "fuelType", "transmission","model")
for (v in cat_vars) {
p <- ggplot(df, aes_string(x = v, y = "mpg")) +
geom_boxplot(fill = "lightgreen", outlier.alpha = 0.4) +
theme_minimal() +
labs(
title = paste("Fuel efficiency by", v),
x = v,
y = "Fuel efficiency"
) +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
print(p)
}
Các dòng xe điện có có mpg lớn hơn (đồng nghĩa với hiệu suất nhiên liệu tốt hơn) so với các dòng xe còn lại. Tuy nhiên các dòng xe điện không được đo bằng mpg và được đo bằng MPGe (Miles Per Gallon equivalent) với 1 gallon xăng ≈ 33.7 kWh điện.
# Tương quan giữa egineSize với các biến định tính
cat_vars <- c("model_group", "fuelType", "transmission","model")
for (v in cat_vars) {
p <- ggplot(df, aes_string(x = v, y = "engineSize")) +
geom_boxplot(fill = "lightgreen", outlier.alpha = 0.4) +
theme_minimal() +
labs(
title = paste("EngineSize by", v),
x = v,
y = "EngineSize"
) +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
print(p)
}
Các dòng Hiệu suất cao M và SUV có enginesize cao hơn so với các dòng còn lại.Bên cạnh đó, động cơ của các xe lái tự động hoặc bán tự động cũng cao hơn động cơ xe số sàn.
#Histogram biến phụ thuộc Price
hist(df$price,
breaks = 30,
main = "Histogram of Price",
xlab = "Price",
col = "lightblue")
# Có thể thấy biến price tương đối lệch phải
skewness(df$price)## [1] 1.590796# Biến đổi log cho biến price
df <- df %>%
mutate(log_price = log(price))#Histogram biến phụ thuộc Log Price
hist(df$log_price,
breaks = 30,
main = "Histogram of Log Price",
xlab = "Price",
col = "lightblue")
# Biến Price đã đỡ lệch hơn
skewness(df$log_price)## [1] -0.2325722# Ma trận tương quan giữa các biến
numeric_vars <- df %>%
select(where(is.numeric))
cor_matrix <- cor(numeric_vars, use = "complete.obs")
corrplot(cor_matrix,
method = "color",
type = "upper",
addCoef.col = "black",
tl.col = "black",
tl.srt = 45,
number.cex = 0.7)
Dựa vào ma trận tương quan, có thể quan sát được các cặp biến có quan hệ tuyến tính mạnh như year-mileage, mileage-price.
Mô hình hồi quy
Hồi quy tuyến tính OLS
# Bỏ biến Price, sử dụng log_price làm biến phụ thuộc
df <- df[, !(names(df) %in% "price")]# Hồi quy OLS với biến gốc model
model_ols <- lm(log_price ~ .-model_group, data = df)
summary(model_ols)$r.squared## [1] 0.9173462# Hồi quy OLS với biến model_group
model_ols <- lm(log_price ~ .-model, data = df)
summary(model_ols)$r.squared## [1] 0.8838328Hồi quy OLS với biến model_group cho hệ số xác định là 0.88 thấp hơn 0.91 khi hồi quy với biến gốp model. Tuy nhiên, vì sau khi phân nhóm các giá trị của biến model thành các biến model_group, mô hình sẽ dễ giải thích hơn nên nhóm chọn hồi quy với biến model_group
# Bỏ biến model, sử dụng biến model_group
df <- df[, !(names(df) %in% "model")]# Tách bộ dữ liệu thành train và test với tỷ lệ 80-20
set.seed(1502) # để tái lập kết quả
n <- nrow(df)
train_index <- sample(seq_len(n), size = 0.8 * n)
train_df <- df[train_index, ]
test_df <- df[-train_index, ]# Hồi quy OLS với toàn bộ các biến
model_ols <- lm(log_price ~ ., data = train_df)
summary(model_ols)##
## Call:
## lm(formula = log_price ~ ., data = train_df)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.02646 -0.09403 -0.00059 0.09578 2.18909
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -2.168e+02 3.119e+00 -69.495 < 2e-16 ***
## year 1.129e-01 1.544e-03 73.121 < 2e-16 ***
## transmissionManual -1.360e-01 5.353e-03 -25.404 < 2e-16 ***
## transmissionSemi-Auto -8.615e-05 4.436e-03 -0.019 0.985
## mileage -5.828e-06 1.373e-07 -42.446 < 2e-16 ***
## fuelTypeHybrid 3.011e-01 1.327e-02 22.681 < 2e-16 ***
## fuelTypeOther 2.228e-01 3.491e-02 6.380 1.86e-10 ***
## fuelTypePetrol -1.397e-01 6.163e-03 -22.670 < 2e-16 ***
## tax -5.169e-04 3.860e-05 -13.390 < 2e-16 ***
## mpg -8.629e-03 3.066e-04 -28.143 < 2e-16 ***
## engineSize 2.365e-01 4.930e-03 47.974 < 2e-16 ***
## model_groupPerformance (M) -6.073e-01 4.888e-02 -12.423 < 2e-16 ***
## model_groupRoadster (Z) -8.854e-01 4.968e-02 -17.821 < 2e-16 ***
## model_groupSeries -9.169e-01 4.648e-02 -19.725 < 2e-16 ***
## model_groupSUV (X) -7.589e-01 4.673e-02 -16.241 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1674 on 8477 degrees of freedom
## Multiple R-squared: 0.8826, Adjusted R-squared: 0.8824
## F-statistic: 4553 on 14 and 8477 DF, p-value: < 2.2e-16Hồi quy OLS với toàn bộ các cho hệ số xác định đạt được 0.88 với đa phần các biến đều có ý nghĩa thống kê, ngoại trừtransmissionSemi-Auto và fuelTypeElectric. Kiểm tra Đa cộng tuyến của các biến.
vif(model_ols)## GVIF Df GVIF^(1/(2*Df))
## year 3.160063 1 1.777657
## transmission 1.436648 2 1.094807
## mileage 3.213909 1 1.792738
## fuelType 3.167628 3 1.211869
## tax 1.704613 1 1.305608
## mpg 4.433772 1 2.105652
## engineSize 2.088852 1 1.445286
## model_group 1.664750 4 1.065783mpg đang có GVIF^(1/(2*Df)) 2.08, không quá đáng lo nhưng nhóm loại xem thử kết quả đạt được tốt hơn không?
model_ols <- lm(log_price ~ .-mpg, data = train_df)
summary(model_ols)##
## Call:
## lm(formula = log_price ~ . - mpg, data = train_df)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.96258 -0.10742 -0.00268 0.10418 2.16099
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -2.234e+02 3.252e+00 -68.698 < 2e-16 ***
## year 1.159e-01 1.611e-03 71.910 < 2e-16 ***
## transmissionManual -1.397e-01 5.595e-03 -24.971 < 2e-16 ***
## transmissionSemi-Auto 6.328e-03 4.632e-03 1.366 0.172
## mileage -5.977e-06 1.435e-07 -41.665 < 2e-16 ***
## fuelTypeHybrid 1.568e-01 1.280e-02 12.245 < 2e-16 ***
## fuelTypeOther 3.093e-02 3.580e-02 0.864 0.388
## fuelTypePetrol -2.176e-02 4.724e-03 -4.606 4.16e-06 ***
## tax -1.479e-04 3.796e-05 -3.897 9.81e-05 ***
## engineSize 3.091e-01 4.394e-03 70.336 < 2e-16 ***
## model_groupPerformance (M) -5.903e-01 5.111e-02 -11.550 < 2e-16 ***
## model_groupRoadster (Z) -8.200e-01 5.189e-02 -15.803 < 2e-16 ***
## model_groupSeries -8.996e-01 4.860e-02 -18.511 < 2e-16 ***
## model_groupSUV (X) -6.888e-01 4.879e-02 -14.119 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1751 on 8478 degrees of freedom
## Multiple R-squared: 0.8716, Adjusted R-squared: 0.8714
## F-statistic: 4429 on 13 and 8478 DF, p-value: < 2.2e-16vif(model_ols)## GVIF Df GVIF^(1/(2*Df))
## year 3.145686 1 1.773608
## transmission 1.429126 2 1.093371
## mileage 3.209113 1 1.791400
## fuelType 1.483725 3 1.067970
## tax 1.508013 1 1.228012
## engineSize 1.517744 1 1.231968
## model_group 1.413823 4 1.044238Sau khi loại bỏ biến mpg, nhóm thấy vẫn tồn tại 2 biến không có ý nghĩa thống kê với mức ý nghĩa 5%, nhưng hệ số xác định lại giảm, nhóm chọn giữ lại biến mpg. Nhìn vào ma trận tương quan ở trên, nhóm cũng thấy 2 biến year và milage có quan hệ tuyến tính mạnh nhên nhóm sẽ thử loại 1 trong 2 biến.
model_ols <- lm(log_price ~ .-year, data = train_df)
summary(model_ols)##
## Call:
## lm(formula = log_price ~ . - year, data = train_df)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.18648 -0.12344 0.00494 0.13107 1.97061
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.128e+01 6.835e-02 164.959 < 2e-16 ***
## transmissionManual -1.804e-01 6.791e-03 -26.571 < 2e-16 ***
## transmissionSemi-Auto 2.409e-02 5.648e-03 4.265 2.02e-05 ***
## mileage -1.363e-05 1.104e-07 -123.424 < 2e-16 ***
## fuelTypeHybrid 3.625e-01 1.691e-02 21.433 < 2e-16 ***
## fuelTypeOther 3.479e-01 4.453e-02 7.814 6.21e-15 ***
## fuelTypePetrol -1.616e-01 7.860e-03 -20.561 < 2e-16 ***
## tax -7.990e-04 4.904e-05 -16.292 < 2e-16 ***
## mpg -1.014e-02 3.906e-04 -25.962 < 2e-16 ***
## engineSize 2.083e-01 6.276e-03 33.183 < 2e-16 ***
## model_groupPerformance (M) -4.781e-01 6.238e-02 -7.664 2.01e-14 ***
## model_groupRoadster (Z) -8.088e-01 6.343e-02 -12.752 < 2e-16 ***
## model_groupSeries -7.932e-01 5.932e-02 -13.373 < 2e-16 ***
## model_groupSUV (X) -6.287e-01 5.962e-02 -10.545 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2138 on 8478 degrees of freedom
## Multiple R-squared: 0.8086, Adjusted R-squared: 0.8083
## F-statistic: 2755 on 13 and 8478 DF, p-value: < 2.2e-16model_ols <- lm(log_price ~ .-mileage, data = train_df)
summary(model_ols)##
## Call:
## lm(formula = log_price ~ . - mileage, data = train_df)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.09443 -0.10432 0.00269 0.10866 2.26932
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -3.195e+02 2.165e+00 -147.563 < 2e-16 ***
## year 1.639e-01 1.071e-03 153.007 < 2e-16 ***
## transmissionManual -1.247e-01 5.886e-03 -21.180 < 2e-16 ***
## transmissionSemi-Auto 2.104e-02 4.853e-03 4.336 1.47e-05 ***
## fuelTypeHybrid 3.057e-01 1.461e-02 20.915 < 2e-16 ***
## fuelTypeOther 1.686e-01 3.842e-02 4.390 1.15e-05 ***
## fuelTypePetrol -1.214e-01 6.769e-03 -17.938 < 2e-16 ***
## tax -3.672e-04 4.232e-05 -8.675 < 2e-16 ***
## mpg -9.132e-03 3.373e-04 -27.069 < 2e-16 ***
## engineSize 2.283e-01 5.425e-03 42.090 < 2e-16 ***
## model_groupPerformance (M) -6.940e-01 5.378e-02 -12.905 < 2e-16 ***
## model_groupRoadster (Z) -9.632e-01 5.467e-02 -17.620 < 2e-16 ***
## model_groupSeries -1.007e+00 5.113e-02 -19.691 < 2e-16 ***
## model_groupSUV (X) -8.493e-01 5.140e-02 -16.525 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1843 on 8478 degrees of freedom
## Multiple R-squared: 0.8577, Adjusted R-squared: 0.8574
## F-statistic: 3930 on 13 and 8478 DF, p-value: < 2.2e-16vif(model_ols)## GVIF Df GVIF^(1/(2*Df))
## year 1.253148 1 1.119441
## transmission 1.418541 2 1.091341
## fuelType 3.147493 3 1.210582
## tax 1.690382 1 1.300147
## mpg 4.427156 1 2.104081
## engineSize 2.085640 1 1.444174
## model_group 1.660790 4 1.065465sau khi loại 1 trong 2 biến, nhóm nhận thấy hệ số xác định giảm đáng kể, từ 0.88 còn 0.85 và 0.80 nên nhóm quyết định giữ lại mô hình đầy đủ.
# Dùng BIC để chọn mô hình.
mod_ols <- lm(log_price ~ ., data = train_df)
modBIC <- MASS::stepAIC(mod_ols, k = log(nrow(train_df)), trace = 0)
summary(modBIC)##
## Call:
## lm(formula = log_price ~ year + transmission + mileage + fuelType +
## tax + mpg + engineSize + model_group, data = train_df)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.02646 -0.09403 -0.00059 0.09578 2.18909
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -2.168e+02 3.119e+00 -69.495 < 2e-16 ***
## year 1.129e-01 1.544e-03 73.121 < 2e-16 ***
## transmissionManual -1.360e-01 5.353e-03 -25.404 < 2e-16 ***
## transmissionSemi-Auto -8.615e-05 4.436e-03 -0.019 0.985
## mileage -5.828e-06 1.373e-07 -42.446 < 2e-16 ***
## fuelTypeHybrid 3.011e-01 1.327e-02 22.681 < 2e-16 ***
## fuelTypeOther 2.228e-01 3.491e-02 6.380 1.86e-10 ***
## fuelTypePetrol -1.397e-01 6.163e-03 -22.670 < 2e-16 ***
## tax -5.169e-04 3.860e-05 -13.390 < 2e-16 ***
## mpg -8.629e-03 3.066e-04 -28.143 < 2e-16 ***
## engineSize 2.365e-01 4.930e-03 47.974 < 2e-16 ***
## model_groupPerformance (M) -6.073e-01 4.888e-02 -12.423 < 2e-16 ***
## model_groupRoadster (Z) -8.854e-01 4.968e-02 -17.821 < 2e-16 ***
## model_groupSeries -9.169e-01 4.648e-02 -19.725 < 2e-16 ***
## model_groupSUV (X) -7.589e-01 4.673e-02 -16.241 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1674 on 8477 degrees of freedom
## Multiple R-squared: 0.8826, Adjusted R-squared: 0.8824
## F-statistic: 4553 on 14 and 8477 DF, p-value: < 2.2e-16Sử dụng phương pháp stepwise chuẩn BIC cũng cho ra kết quả là mô hình đầy đủ.
Phân thích thành phần chính
# Chọn các biến định lượng và loại bỏ biến phụ thuộc
num_vars <- train_df %>%
select(where(is.numeric)) %>%
select(-log_price)# Chuẩn hoá dữ liệu
df_PCA <- scale(num_vars, center = TRUE, scale = TRUE)
boxplot(df_PCA)
# Phân tích thành phần chính
res <- princomp(df_PCA, cor = FALSE)
summary(res, loadings = T)## Importance of components:
## Comp.1 Comp.2 Comp.3 Comp.4 Comp.5
## Standard deviation 1.5059035 1.2737670 0.7541976 0.60110674 0.42313132
## Proportion of Variance 0.4536025 0.3245347 0.1137762 0.07227437 0.03581224
## Cumulative Proportion 0.4536025 0.7781372 0.8919134 0.96418776 1.00000000
##
## Loadings:
## Comp.1 Comp.2 Comp.3 Comp.4 Comp.5
## year 0.401 0.573 0.157 0.693
## mileage -0.441 -0.535 0.716
## tax 0.474 -0.336 -0.576 0.571
## mpg -0.534 0.265 0.177 0.782
## engineSize 0.367 -0.450 0.782 0.224eig_vals <- res$sdev^2
df_scree <- data.frame(
PC = seq_along(eig_vals),
Eigenvalue = eig_vals
)
ggplot(df_scree, aes(x = PC, y = Eigenvalue)) +
geom_point() +
geom_line() +
theme_minimal() +
labs(
title = "Scree Plot (Elbow Method)",
x = "Principal Component",
y = "Eigenvalue"
)
Kết quả PCA cho thấy thành phần chính thứ nhất (PC1) giải thích 45,38% phương sai dữ liệu và PC2 giải thích thêm 32,38%, nâng tổng phương sai tích luỹ lên 77,81%; khi bổ sung PC3,PC4 tổng phương sai được giải thích đạt khoảng 96.4%. Thành phần chính thứ 5 đóng góp rất ít vào việc giải thích phương sai dữ liệu nên nhóm chọn giữ lại 4 thành phần chính đầu tiên.
# Hồi quy trên các thành phần chính
train_pca <- as.data.frame(res$scores)
train_pca$log_price <- train_df$log_price
train_pca$model_group <- train_df$model_group
train_pca$fuelType <- train_df$fuelType
train_pca$transmission <- train_df$transmission
pcr_model <- lm(log_price ~ Comp.1 + Comp.2 + Comp.3 + Comp.4
+ model_group + fuelType + transmission,
data = train_pca)
summary(pcr_model)##
## Call:
## lm(formula = log_price ~ Comp.1 + Comp.2 + Comp.3 + Comp.4 +
## model_group + fuelType + transmission, data = train_pca)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.02817 -0.09642 -0.00056 0.09834 2.15045
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 10.830906 0.047061 230.144 < 2e-16 ***
## Comp.1 0.244128 0.001770 137.889 < 2e-16 ***
## Comp.2 0.134872 0.001545 87.287 < 2e-16 ***
## Comp.3 0.138979 0.002916 47.653 < 2e-16 ***
## Comp.4 -0.046525 0.004017 -11.583 < 2e-16 ***
## model_groupPerformance (M) -0.580626 0.049422 -11.748 < 2e-16 ***
## model_groupRoadster (Z) -0.863281 0.050241 -17.183 < 2e-16 ***
## model_groupSeries -0.889286 0.046989 -18.925 < 2e-16 ***
## model_groupSUV (X) -0.729263 0.047228 -15.441 < 2e-16 ***
## fuelTypeHybrid 0.303459 0.013428 22.599 < 2e-16 ***
## fuelTypeOther 0.240870 0.035301 6.823 9.51e-12 ***
## fuelTypePetrol -0.142055 0.006233 -22.791 < 2e-16 ***
## transmissionManual -0.142316 0.005397 -26.371 < 2e-16 ***
## transmissionSemi-Auto -0.001338 0.004487 -0.298 0.766
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1694 on 8478 degrees of freedom
## Multiple R-squared: 0.8798, Adjusted R-squared: 0.8796
## F-statistic: 4774 on 13 and 8478 DF, p-value: < 2.2e-16Hồi quy Lasso, Ridge, ElasticNet
# Ma trận thiết kế (dummy hóa factor)
x_train <- model.matrix(
log_price ~ .,
data = train_df
)[, -1] # bỏ intercept
y_train <- train_df$log_price# Hồi quy Ridge
set.seed(1502)
cv_ridge <- cv.glmnet(
x_train,
y_train,
alpha = 0 # Ridge
)plot(cv_ridge)
# Lambda tối ưu
cv_ridge$lambda.min## [1] 0.03795819cv_ridge$lambda.1se## [1] 0.07280041# Fit model với lambda tối ưu
model_ridge <- glmnet(
x_train,
y_train,
alpha = 0,
lambda = cv_ridge$lambda.min
)
coef(model_ridge)## 15 x 1 sparse Matrix of class "dgCMatrix"
## s0
## (Intercept) -1.947244e+02
## year 1.015597e-01
## transmissionManual -1.415035e-01
## transmissionSemi-Auto 1.323016e-02
## mileage -6.075703e-06
## fuelTypeHybrid 2.963869e-01
## fuelTypeOther 2.591024e-01
## fuelTypePetrol -1.070170e-01
## tax -3.305090e-04
## mpg -6.831552e-03
## engineSize 2.215935e-01
## model_groupPerformance (M) 1.820550e-01
## model_groupRoadster (Z) -9.667674e-02
## model_groupSeries -1.320980e-01
## model_groupSUV (X) 3.978248e-02# Hồi quy Lasso
set.seed(123)
cv_lasso <- cv.glmnet(
x_train,
y_train,
alpha = 1 # Lasso
)
plot(cv_lasso)
-log(cv_lasso$lambda.min)## [1] 9.434755-log(cv_lasso$lambda.1se)## [1] 5.620372model_lasso <- glmnet(
x_train,
y_train,
alpha = 1,
lambda = cv_lasso$lambda.1se
)
coef(model_lasso)## 15 x 1 sparse Matrix of class "dgCMatrix"
## s0
## (Intercept) -2.171763e+02
## year 1.126790e-01
## transmissionManual -1.354343e-01
## transmissionSemi-Auto .
## mileage -5.749504e-06
## fuelTypeHybrid 2.869734e-01
## fuelTypeOther 2.014039e-01
## fuelTypePetrol -1.073908e-01
## tax -3.197133e-04
## mpg -6.577560e-03
## engineSize 2.409311e-01
## model_groupPerformance (M) 1.167760e-01
## model_groupRoadster (Z) -1.008612e-01
## model_groupSeries -1.692307e-01
## model_groupSUV (X) .Lasso loại được 2 biến khỏi mô hình là biến model_SUV và biến động cơ bán tự động.
set.seed(1502)
cv_enet <- cv.glmnet(
x_train,
y_train,
alpha = 0.4 # Elastic-Net
)
plot(cv_enet)
model_enet <- glmnet(
x_train,
y_train,
alpha = 1,
lambda = cv_enet$lambda.1se
)
coef(model_enet)## 15 x 1 sparse Matrix of class "dgCMatrix"
## s0
## (Intercept) -2.195170e+02
## year 1.137307e-01
## transmissionManual -1.328526e-01
## transmissionSemi-Auto .
## mileage -5.500037e-06
## fuelTypeHybrid 2.121249e-01
## fuelTypeOther 4.156418e-02
## fuelTypePetrol -6.354587e-02
## tax -4.110579e-05
## mpg -3.853825e-03
## engineSize 2.509637e-01
## model_groupPerformance (M) 6.759572e-02
## model_groupRoadster (Z) -5.388218e-02
## model_groupSeries -1.740446e-01
## model_groupSUV (X) .Chọn mô hình
# So sánh các model
eval_metrics <- function(y_true, y_pred) {
mse <- mean((y_true - y_pred)^2)
rmse <- sqrt(mse)
r2 <- 1 - sum((y_true - y_pred)^2) /
sum((y_true - mean(y_true))^2)
c(R2 = r2, RMSE = rmse, MSE = mse)
}# ols
pred_ols <- predict(modBIC, newdata = test_df)
ols_metrics <- eval_metrics(test_df$log_price, pred_ols)
ols_metrics## R2 RMSE MSE
## 0.88861446 0.15651122 0.02449576x_test <- model.matrix(
log_price ~ .,
data = test_df
)[, -1]
y_test <- test_df$log_price# Ridge
pred_ridge <- predict(
cv_ridge,
s = "lambda.min",
newx = x_test
)
ridge_metrics <- eval_metrics(test_df$log_price, pred_ridge)
ridge_metrics## R2 RMSE MSE
## 0.87982747 0.16256748 0.02642818# Lasso
pred_lasso <- predict(
cv_lasso,
s = "lambda.1se",
newx = x_test
)
lasso_metrics <- eval_metrics(test_df$log_price, pred_lasso)
lasso_metrics## R2 RMSE MSE
## 0.8811265 0.1616864 0.0261425# Elastic-Net
pred_enet <- predict(
cv_enet,
s = "lambda.1se",
newx = x_test
)
enet_metrics <- eval_metrics(test_df$log_price, pred_enet)
enet_metrics## R2 RMSE MSE
## 0.88106886 0.16172563 0.02615518# PCR
test_scaled <- scale(
test_df[, colnames(num_vars)], # đúng thứ tự biến numeric
center = attr(df_PCA, "scaled:center"),
scale = attr(df_PCA, "scaled:scale")
)
test_scores <- predict(res, newdata = test_scaled)
test_pca <- as.data.frame(test_scores)
test_pca$model_group <- test_df$model_group
test_pca$fuelType <- test_df$fuelType
test_pca$transmission <- test_df$transmission
pred_pca <- predict(pcr_model, newdata = test_pca)
pca_metrics <- eval_metrics(test_df$log_price, pred_pca)
pca_metrics## R2 RMSE MSE
## 0.88590714 0.15840186 0.02509115metrics_table <- rbind(
OLS = round(ols_metrics,4),
Ridge = round(ridge_metrics,4),
Lasso = round(lasso_metrics,4),
Elastic_Net = round(enet_metrics,4),
PCA = round(pca_metrics,4)
)
metrics_table## R2 RMSE MSE
## OLS 0.8886 0.1565 0.0245
## Ridge 0.8798 0.1626 0.0264
## Lasso 0.8811 0.1617 0.0261
## Elastic_Net 0.8811 0.1617 0.0262
## PCA 0.8859 0.1584 0.0251Dựa vào kết quả test ta có thể thấy, model hồi quy ols đang thể hiện tốt nhất. Tuy Lasso có loại được 2 biến so với các mô hình còn lại, nhưng nhóm vẫn chọn mô hình OLS với số biến nhìn chung không đáng kể.
Kiểm tra giả định mô hình
# Kiểm tra các điểm ảnh hưởng mạnh
plot(modBIC,4)
Có thể thấy có các điểm ảnh hưởng mạnh như quan trắc 9116, 9922, 10161
# Loại bỏ các điểm ảnh hưởng mạnh
cooks_d <- cooks.distance(modBIC)
n <- nrow(train_df)
threshold <- 4 / n
influential_points <- which(cooks_d > threshold)
train_df_clean <- train_df[-influential_points, ]# Chạy lại model ols với dữ liệu đã loại bỏ các điểm ảnh hưởng mạnh
model_ols_clean <- lm(log_price ~., data = train_df_clean)
summary(model_ols_clean)##
## Call:
## lm(formula = log_price ~ ., data = train_df_clean)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.43415 -0.08435 -0.00051 0.08502 0.54796
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -2.054e+02 2.530e+00 -81.169 < 2e-16 ***
## year 1.075e-01 1.251e-03 85.876 < 2e-16 ***
## transmissionManual -1.220e-01 4.278e-03 -28.530 < 2e-16 ***
## transmissionSemi-Auto 7.872e-03 3.530e-03 2.230 0.0258 *
## mileage -5.141e-06 1.106e-07 -46.467 < 2e-16 ***
## fuelTypeHybrid 3.974e-01 1.361e-02 29.203 < 2e-16 ***
## fuelTypeOther 3.126e-01 4.991e-02 6.263 3.96e-10 ***
## fuelTypePetrol -1.959e-01 5.300e-03 -36.958 < 2e-16 ***
## tax -3.132e-04 3.387e-05 -9.249 < 2e-16 ***
## mpg -1.263e-02 2.784e-04 -45.380 < 2e-16 ***
## engineSize 1.848e-01 4.230e-03 43.685 < 2e-16 ***
## model_groupPerformance (M) -6.844e-01 9.418e-02 -7.267 4.03e-13 ***
## model_groupRoadster (Z) -9.300e-01 9.476e-02 -9.814 < 2e-16 ***
## model_groupSeries -9.711e-01 9.317e-02 -10.423 < 2e-16 ***
## model_groupSUV (X) -8.530e-01 9.324e-02 -9.148 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1306 on 8099 degrees of freedom
## Multiple R-squared: 0.9148, Adjusted R-squared: 0.9146
## F-statistic: 6207 on 14 and 8099 DF, p-value: < 2.2e-16Sau khi loại bỏ các điểm ảnh hưởng, hệ số xác định đã tăng từ 0.88 lên .091
# Kiểm định phân phối chuẩn của sai số
qqnorm(residuals(model_ols_clean))
qqline(residuals(model_ols_clean), col = "red")
Các quan trắc đa phần đều nằm trên đường lý thuyết
# Do dữ liệu > 5000 dòng nên không dùng Shapiro được, dùng Kolmogorov–Smirnov hiệu chỉnh thay thế
lillie.test(residuals(model_ols_clean))##
## Lilliefors (Kolmogorov-Smirnov) normality test
##
## data: residuals(model_ols_clean)
## D = 0.01308, p-value = 0.003194Phần dư có phân phối chuẩn với độ tin cậy 0.001
# Kiểm định phương sai đồng nhất
bptest(model_ols_clean)##
## studentized Breusch-Pagan test
##
## data: model_ols_clean
## BP = 405.09, df = 14, p-value < 2.2e-16plot(model_ols_clean,1)
# Biến đổi box-cox để xem các kiểm định có được cải thiện không
L <- MASS::boxcox(model_ols_clean)
optimal_lambda_idx <-which.max(L$y)
optimal_lambda <-L$x[optimal_lambda_idx]
optimal_lambda## [1] -0.1818182train_df_clean$Y_boxcox<-(train_df_clean$log_price^optimal_lambda-1)/optimal_lambda
new_modBIC<-lm(Y_boxcox ~ . - log_price, data = train_df_clean)
summary(new_modBIC)##
## Call:
## lm(formula = Y_boxcox ~ . - log_price, data = train_df_clean)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.035011 -0.005464 0.000082 0.005646 0.034308
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.249e+01 1.671e-01 -74.712 < 2e-16 ***
## year 7.163e-03 8.265e-05 86.667 < 2e-16 ***
## transmissionManual -8.806e-03 2.826e-04 -31.163 < 2e-16 ***
## transmissionSemi-Auto 4.146e-04 2.332e-04 1.778 0.0754 .
## mileage -3.517e-07 7.309e-09 -48.124 < 2e-16 ***
## fuelTypeHybrid 2.535e-02 8.990e-04 28.201 < 2e-16 ***
## fuelTypeOther 2.014e-02 3.297e-03 6.109 1.05e-09 ***
## fuelTypePetrol -1.234e-02 3.501e-04 -35.262 < 2e-16 ***
## tax -1.267e-05 2.237e-06 -5.664 1.53e-08 ***
## mpg -7.775e-04 1.839e-05 -42.277 < 2e-16 ***
## engineSize 1.215e-02 2.794e-04 43.463 < 2e-16 ***
## model_groupPerformance (M) -4.129e-02 6.221e-03 -6.638 3.40e-11 ***
## model_groupRoadster (Z) -5.623e-02 6.260e-03 -8.984 < 2e-16 ***
## model_groupSeries -5.881e-02 6.154e-03 -9.557 < 2e-16 ***
## model_groupSUV (X) -5.108e-02 6.159e-03 -8.294 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.008628 on 8099 degrees of freedom
## Multiple R-squared: 0.9161, Adjusted R-squared: 0.916
## F-statistic: 6318 on 14 and 8099 DF, p-value: < 2.2e-16bptest(new_modBIC)##
## studentized Breusch-Pagan test
##
## data: new_modBIC
## BP = 372.28, df = 14, p-value < 2.2e-16plot(new_modBIC,1)
Biến đổi box cox vẫn không cải thiện được phương sai đồng nhất
# Kiểm định tính chuẩn sau khi biến đổi box-cox. Phần dư không còn tính chuẩn với độ tin cậy 0.0001.
lillie.test(residuals(new_modBIC))##
## Lilliefors (Kolmogorov-Smirnov) normality test
##
## data: residuals(new_modBIC)
## D = 0.015872, p-value = 8.736e-05Do không cải thiện được phương sai đồng nhất và phân dư không còn phân phối chuẩn với độ tin cậy 1/1000 nên nhóm chọn mô hình không có biến đổi box-cox
# Kiểm định tính độc lập của phần dư Durbin-Watson
dwtest(model_ols_clean)##
## Durbin-Watson test
##
## data: model_ols_clean
## DW = 2.0013, p-value = 0.5231
## alternative hypothesis: true autocorrelation is greater than 0pred_log_price <- predict(model_ols_clean, newdata = test_df)
pred_price <- exp(pred_log_price)
actual_price <- exp(test_df$log_price)
df_plot <- data.frame(
Actual_Price = actual_price,
Predicted_Price = pred_price
)
ggplot(df_plot, aes(x = Actual_Price, y = Predicted_Price)) +
geom_point(alpha = 0.4, color = "steelblue") +
geom_abline(slope = 1, intercept = 0,
linetype = "dashed", color = "red") +
theme_minimal() +
labs(
title = "Actual vs Predicted Price (OLS Clean Model)",
x = "Actual Price",
y = "Predicted Price"
)
Đề xuất sử dụng Roburst Standard Error hiệu chỉnh sai số chuẩn trong trường hợp phương sai sai số không đồng nhất. Việc này không làm thay đổi các hệ số ước lượng OLS, nhưng giúp các thống kê kiểm định (t-statistic, p-value và khoảng tin cậy) trở nên đáng tin cậy hơn.
coeftest(model_ols_clean, vcov = vcovHC(model_ols_clean, type = "HC1"))##
## t test of coefficients:
##
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -2.0536e+02 2.4925e+00 -82.3918 < 2e-16 ***
## year 1.0745e-01 1.2333e-03 87.1242 < 2e-16 ***
## transmissionManual -1.2205e-01 4.1982e-03 -29.0715 < 2e-16 ***
## transmissionSemi-Auto 7.8715e-03 3.6048e-03 2.1836 0.02902 *
## mileage -5.1413e-06 1.0364e-07 -49.6091 < 2e-16 ***
## fuelTypeHybrid 3.9745e-01 1.2061e-02 32.9526 < 2e-16 ***
## fuelTypeOther 3.1260e-01 1.6942e-02 18.4509 < 2e-16 ***
## fuelTypePetrol -1.9587e-01 5.9859e-03 -32.7217 < 2e-16 ***
## tax -3.1324e-04 3.0918e-05 -10.1313 < 2e-16 ***
## mpg -1.2634e-02 3.1595e-04 -39.9881 < 2e-16 ***
## engineSize 1.8480e-01 4.7555e-03 38.8605 < 2e-16 ***
## model_groupPerformance (M) -6.8439e-01 1.3792e-02 -49.6217 < 2e-16 ***
## model_groupRoadster (Z) -9.3001e-01 1.3225e-02 -70.3240 < 2e-16 ***
## model_groupSeries -9.7111e-01 9.9466e-03 -97.6323 < 2e-16 ***
## model_groupSUV (X) -8.5305e-01 1.0604e-02 -80.4448 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1