Nhắc lại hồi quy tuyến tính
Chúng ta sử dụng phương pháp hồi quy tuyến tính để dự đoán một giá trị nào đó, dựa vào các biến số cho trước và mối quan hệ của các biến số này. Ví dụ, dựa vào hồi quy tuyến tính, chúng ta có thể biết được:
Dự báo giá trị của biến số dựa trên biến số khác của kết quả quan sát thu được
Kiểm định được dạng mô hình và đánh giá khả năng rủi ro dự đoán của mô hình
Ví dụ: Với bộ số liệu về giá trị nhà ở Boston, chúng ta có thể dự báo giá trị của ngôi nhà dựa vào diện tích của ngôi nhà, số phòng tắm, vị trí của ngôi nhà…dựa trên các quan sát của các biến này.
Mô hình chung
Phương trình hồi quy tổng quát
Phương trình hồi quy
Phương trình hồi quy tuyến tính dạng đơn biến có dạng tổng quát là: \[y^{i} = \beta_{0}+\beta_{1}x^{i}+\epsilon^{i}\] trong đó:
\(y^{i}\) là biến phụ thuộc (dependent variable), chính là giá trị của ngôi nhà chúng ta cần dự báo;
\(x^{i}\) là biến độc lập (independent variable), có thể lấy một biến làm ví dụ là diện tích của ngôi nhà;
\(\epsilon^{i}\) là phần dư, đây chính là giá trị chênh lệch của giá trị thực \(y^{i}\) và giá trị \(\hat{y}^{i}=\beta_{0}+\beta_{1}x^{i}\)
\(\beta_{0}\) là hệ số chặn của mô hình hồi quy (Intercept), hay chính là giá trị khi \(\beta_{1}=0\).
\(\beta_{1}\) là hệ số góc (Slope) của biến \(x^{i}\)
Các phương pháp xác định mô hình
Dựa vào MSE hoặc RMSE
Để đưa ra mô hình hồi quy có giá trị dự đoán tốt nhất, chúng ta thường xác định hai giá trị \(\beta_{0}\) và \(\beta_{1}\) sao cho tổng bình phương các giá trị phần dư là nhỏ nhất, phần dư càng nhỏ thì chứng tỏ mô hình dự báo càng chính xác, gần với giá trị thực. Phương pháp này có tên là \(SSE\) - Sum Squared Error, mục đích là tìm min \(SSE\). Tuy nhiên \(SSE\) là đại lượng bình phương nên rất khó diễn giải, do đó thường chúng ta thường đưa ra đại lượng khác là \(RMSE\) để diễn giải tốt hơn \(RMSE\) - Root Mean Squared Error \[RMSE = \sqrt{\frac{SSE}{N}}\] Do đó \(RMSE\) sẽ là đơn vị chuẩn hóa, lấy căn bậc hai sẽ có đơn vị là đơn vị của biến phụ thuộc, nó sẽ đại diện cho sai số trung bình của một quan sát trong mẫu.
Dựa vào R bình phương của mô hình
Trong phương pháp này chúng ta sẽ tính một giá trị \(R^{2}\) của mô hình hồi quy: \[R^{2}=1-\frac{SSE}{SST}\] Trong đó giá trị SST chính là giá trị bình phương độ lệch với mô hình cơ sở (baseline model). Mô hình cơ sở là mô hình không sử dụng bất kỳ biến số nào để dự báo. Ví dụ: Dự báo giá nhà chúng ta coi giá trị dự báo bằng trung bình của tất cả các giá trị khác đã quan sát, không sử dụng biến độc lập (diện tích nhà).
Ý nghĩa của \(R^{2}\) như sau:
| Giá trị R bình phương | Ý nghĩa |
|---|---|
| 0 | Không tốt hơn mô hình baseline |
| 1 | Mô hình dự báo chính xác hoàn toàn |
Baseline Model
Giá trị của \(R^{2}\) càng lớn, mô hình dự báo càng chính xác \(0<SSE<SST\)
Với mô hình đa biến
Hoàn toàn tương tự như mô hình đơn biến, tuy nhiên chúng ta cần có một chú ý, đó là Không nên sử dụng quá nhiều biến số cho mô hình đa biến, bởi vì:
Xảy ra hiện tượng overfitting:
- Khả năng dự báo trong mẫu quan sát là tốt, tuy nhiên khi dự báo với bộ số liệu mới cho kết quả rất kém, mô hình quá khớp so với bộ dữ liệu đã xây dựng
- Xảy ra khi sử dụng nhiều biến số cùng một lúc trong khi số lượng quan sát lại rất ít
Hai phương pháp để khắc phục hiện tượng overfitting:
- Sử dụng \(R^{2}\) hiệu chỉnh hay còn gọi là hiệu chỉnh R bình phương
- Tách mẫu quan sát thành hai phần: Một bộ mẫu xây dựng mô hình và một bộ mẫu kiểm định.
Chú ý các giá trị ngoại biên (outliers): Là các quan sát có giá trị bất thường so với các giá trị còn lại
Hiện tượng đa cộng tuyến: Khi các biến độc lập có tương quan lớn với nhau, khi đó có thể làm sai lệch mô hình. Ví dụ số phòng ngủ và số phòng tắm….
Áp dụng vào bộ dữ liệu Boston
Bộ dữ liệu Boston nằm trong thư viện MASS, chúng ta cần cài đặt thư viện này trước
Xem xét tương quan giữa các biến
## 'data.frame': 506 obs. of 14 variables:
## $ crim : num 0.00632 0.02731 0.02729 0.03237 0.06905 ...
## $ zn : num 18 0 0 0 0 0 12.5 12.5 12.5 12.5 ...
## $ indus : num 2.31 7.07 7.07 2.18 2.18 2.18 7.87 7.87 7.87 7.87 ...
## $ chas : int 0 0 0 0 0 0 0 0 0 0 ...
## $ nox : num 0.538 0.469 0.469 0.458 0.458 0.458 0.524 0.524 0.524 0.524 ...
## $ rm : num 6.58 6.42 7.18 7 7.15 ...
## $ age : num 65.2 78.9 61.1 45.8 54.2 58.7 66.6 96.1 100 85.9 ...
## $ dis : num 4.09 4.97 4.97 6.06 6.06 ...
## $ rad : int 1 2 2 3 3 3 5 5 5 5 ...
## $ tax : num 296 242 242 222 222 222 311 311 311 311 ...
## $ ptratio: num 15.3 17.8 17.8 18.7 18.7 18.7 15.2 15.2 15.2 15.2 ...
## $ black : num 397 397 393 395 397 ...
## $ lstat : num 4.98 9.14 4.03 2.94 5.33 ...
## $ medv : num 24 21.6 34.7 33.4 36.2 28.7 22.9 27.1 16.5 18.9 ...
## medv lstat age rm
## medv 1.0000000 -0.7376627 -0.3769546 0.6953599
## lstat -0.7376627 1.0000000 0.6023385 -0.6138083
## age -0.3769546 0.6023385 1.0000000 -0.2402649
## rm 0.6953599 -0.6138083 -0.2402649 1.0000000Xây dựng mô hình hồi quy
Chúng ta sử dụng thư viện caTools để phân tách bộ dữ liệu thành 2 bộ mẫu, bộ mẫu xây dựng và bộ mẫu kiểm định
Định hình, phân chia dữ liệu
Xây dựng mô hình
Mô hình đưa ra dự đoán giá trị trung vị của các ngôi nhà dựa vào biến tuổi của ngôi nhà, % dân cư có địa vị thấp trong khu vực đó và số phòng của ngôi nhà.
model1 <- lm(medv ~ lstat + age + rm, data = mauxaydung)
model1##
## Call:
## lm(formula = medv ~ lstat + age + rm, data = mauxaydung)
##
## Coefficients:
## (Intercept) lstat age rm
## -3.6123591 -0.5883268 0.0007888 5.3638120summary(model1)##
## Call:
## lm(formula = medv ~ lstat + age + rm, data = mauxaydung)
##
## Residuals:
## Min 1Q Median 3Q Max
## -18.535 -3.648 -1.258 2.149 28.793
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -3.6123591 4.1895781 -0.862 0.389
## lstat -0.5883268 0.0693345 -8.485 7.16e-16 ***
## age 0.0007888 0.0144048 0.055 0.956
## rm 5.3638120 0.5912642 9.072 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 5.931 on 333 degrees of freedom
## Multiple R-squared: 0.6208, Adjusted R-squared: 0.6174
## F-statistic: 181.7 on 3 and 333 DF, p-value: < 2.2e-16# Sử dụng mô hình để tính giá trị dự báo
dubao_xaydung <- predict(model1, mauxaydung)
# Hoặc sử dụng hàm fitted
dubao_xaydung1 <- fitted(model1)
# Tính tổng phần dư bình phương của bộ mẫu xây dựng của mô hình
SSE_xaydung <- sum((dubao_xaydung - mauxaydung$medv)^2)
# Tổng phần dư bình phương đối với mô hình cơ sở
SST_xaydung <- sum((mean(mauxaydung$medv) - mauxaydung$medv)^2)
# Tính giá trị R bình phương
R2_xaydung <- 1 - (SSE_xaydung/SST_xaydung)Thực hiện dự báo trên bộ mẫu kiểm định
# Sử dụng mô hình để tính giá trị dự báo
dubao_kiemdinh <- predict(model1, maukiemdinh)
# Tính tổng phần dư bình phương của bộ mẫu kiểm định của mô hình
SSE_kiemdinh <- sum((dubao_kiemdinh - maukiemdinh$medv)^2)
# Tổng phần dư bình phương đối với mô hình cơ sở, tuy nhiên sử dụng lại mean() của bộ mẫu xây dựng để tạo ra bộ mẫu hoàn toàn mới, chưa gặp bao giờ để tránh hiện tượng overfitting
SST_kiemdinh <- sum((mean(mauxaydung$medv) - maukiemdinh$medv)^2)
# Tính giá trị R bình phương
R2_kiemdinh <- 1 - (SSE_kiemdinh/SST_kiemdinh)Sử dụng mô hình với nhiều biến số, đa biến
mohinh2 <- lm(medv ~ .,data = myboston)
mohinh2##
## Call:
## lm(formula = medv ~ ., data = myboston)
##
## Coefficients:
## (Intercept) crim zn indus chas nox
## 36.461352 -0.108762 0.048031 0.019932 2.705245 -17.541602
## rm age dis rad tax ptratio
## 3.839225 -0.001938 -1.493304 0.324925 -0.011598 -0.947985
## black lstat ID
## 0.009357 -0.526184 -0.002526summary(mohinh2)##
## Call:
## lm(formula = medv ~ ., data = myboston)
##
## Residuals:
## Min 1Q Median 3Q Max
## -15.8948 -2.7585 -0.4663 1.7963 26.0911
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 36.461352 5.100994 7.148 3.21e-12 ***
## crim -0.108762 0.032855 -3.310 0.001000 **
## zn 0.048031 0.013785 3.484 0.000538 ***
## indus 0.019932 0.061468 0.324 0.745871
## chas 2.705245 0.861298 3.141 0.001786 **
## nox -17.541602 3.822390 -4.589 5.66e-06 ***
## rm 3.839225 0.418422 9.175 < 2e-16 ***
## age -0.001938 0.013380 -0.145 0.884866
## dis -1.493304 0.199892 -7.471 3.68e-13 ***
## rad 0.324925 0.068111 4.771 2.43e-06 ***
## tax -0.011598 0.003807 -3.046 0.002443 **
## ptratio -0.947985 0.130822 -7.246 1.67e-12 ***
## black 0.009357 0.002685 3.485 0.000536 ***
## lstat -0.526184 0.050704 -10.377 < 2e-16 ***
## ID -0.002526 0.002080 -1.215 0.225046
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.743 on 491 degrees of freedom
## Multiple R-squared: 0.7414, Adjusted R-squared: 0.734
## F-statistic: 100.6 on 14 and 491 DF, p-value: < 2.2e-16