Mô hình hồi quy tuyến tính đơn giản

Tên: Mai Huy

MSSV: 43.01.104.065

Số thứ tự: 08

1) Nhập thư viện và tìm hiểu tập dữ liệu

# Load thư viện Mass chứa tập dữ liệu Boston
library(MASS)
# Liệt kê tên các biến thuộc tính của dữ liệu Boston
names(Boston)

 [1] "crim"    "zn"      "indus"   "chas"    "nox"     "rm"      "age"     "dis"     "rad"     "tax"    
[11] "ptratio" "black"   "lstat"   "medv"   

# Tìm hiểu thêm về thông tin tập dữ liệu Boston
?Boston

Sử dụng câu lệnh ?Boston giúp chúng ta tìm hiểu rõ hơn về các biến thuộc tính trong dữ liệu Boston (có 506 dòng và 14 cột), chẳng hạn như là:

• medv là giá trị nhà cửa trung bình

• rm là số phòng trung bình mỗi nhà

• age là trung bình độ tuổi của các căn nhà

• lstat là tỉ lệ hộ gia đình có địa vị kinh tế xã hội thấp

2) Phân tích, dự đoán, cách thức hoạt động một mô hình hồi quy đơn giản

# lm(medv~lstat) là hàm dùng để fit một mô hình hồi quy tuyến tính đơn giản với medv là giá trị đầu ra cần được dự đoán  và lstat là giá trị đầu vào trong tập dữ liệu Boston
# Tuy nhiên hệ thống sẽ báo lỗi do R không nhận dạng được 2 biến medv và lstat 
lm.fit=lm(medv~lstat)

Error in eval(predvars, data, env) : object 'medv' not found

# Do đó chúng ta cần truyền vào tập dữ liệu Boston để R có nhận dạng được 2 biến medv và lstat
# R đã nhận dạng được 2 biến medv và lstat trong dữ liệu Boston
lm.fit=lm(medv~lstat, data=Boston)

Dưới đây là cách không cần truyền vào tập dữ liệu Boston nhưng R vẫn có thể nhận dạng được các biến của dữ liệu

# attach dùng để khiến cho những biến feature trong dữ liệu có sẵn trong Rstudio theo tên
attach(Boston)
# Giờ đây hàm lm() vẫn có thể sử dụng được khi không cần truyền vào tập dữ liệu nào cả 
lm.fit=lm(medv~lstat)

# Biến lm.fit dùng để hiển thị một số thông tin cơ bản về mô hình đầu ra như thông tin về hàm lm() đã sử dũng trước đó cũng như là những hệ số B0 = 34.55 và B1= -0.95 với medv = 34.55-0.95 lstat
lm.fit

Call:
lm(formula = medv ~ lstat)

Coefficients:
(Intercept)        lstat  
      34.55        -0.95  

# Và để hiển thị thông tin chi tiết hơn , chúng ta sử dụng hàm summary(lm.fit)
# Hàm sẽ cho chúng ta biết thông tin về p-values, thống kê t(t-statistic), độ lệch chuẩn cho những hệ số của biến lstat , cũng như là thống kê R^2, sai số chuẩn (standard errors) và thống kê F (F-statistic) cho mô hình
summary(lm.fit)

Call:
lm(formula = medv ~ lstat)

Residuals:
    Min      1Q  Median      3Q     Max 
-15.168  -3.990  -1.318   2.034  24.500 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 34.55384    0.56263   61.41   <2e-16 ***
lstat       -0.95005    0.03873  -24.53   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 6.216 on 504 degrees of freedom
Multiple R-squared:  0.5441,    Adjusted R-squared:  0.5432 
F-statistic: 601.6 on 1 and 504 DF,  p-value: < 2.2e-16

# Xem các biến được lưu trữ trong lm.fit
names(lm.fit)

 [1] "coefficients"  "residuals"     "effects"       "rank"          "fitted.values" "assign"       
 [7] "qr"            "df.residual"   "xlevels"       "call"          "terms"         "model"        

# Xem giá trị của biến hệ số (coefficients) cho biết thông tin về B0 và B1
lm.fit$coefficients

(Intercept)       lstat 
 34.5538409  -0.9500494 

# Cũng được dùng để Xem giá trị của biến hệ số (coefficients) cho biết thông tin về B0 và B1
coef(lm.fit)

(Intercept)       lstat 
 34.5538409  -0.9500494 

# Hàm confint cho biết khoảng tin cậy của các hệ số B0 và B1 đang ước lượng
confint(lm.fit)

                2.5 %     97.5 %
(Intercept) 33.448457 35.6592247
lstat       -1.026148 -0.8739505

Hàm predict tạo ra khoảng tin cậy

# Cột fit là giá trị medv được dự đoán khi đưa vào dữ liệu (data) mới là biến lstat được gán cho 3 giá trị lần lượt là 5,10,15, lwr và upr là khoảng trên và khoảng dưới của khoảng tin cậy.
predict(lm.fit, data.frame(lstat=c(5,10,15)), interval="confidence")

       fit      lwr      upr
1 29.80359 29.00741 30.59978
2 25.05335 24.47413 25.63256
3 20.30310 19.73159 20.87461

Hàm predict tạo ra khoảng dự đoán

# Cột fit là giá trị medv được dự đoán khi đưa vào dữ liệu (data) mới là biến lstat được gán cho 3 giá trị lần lượt là 5,10,15, lwr và upr là khoảng trên và khoảng dưới của khoảng dự đoán.
predict(lm.fit, data.frame(lstat=c(5,10,15)),
interval="prediction")

       fit       lwr      upr
1 29.80359 17.565675 42.04151
2 25.05335 12.827626 37.27907
3 20.30310  8.077742 32.52846

3) Biểu đồ phân tán minh hoạ mô hình hồi quy đơn giản

# attach dùng để khiến cho những biến feature trong dữ liệu có sẵn trong Rstudio theo tên
attach(Boston)

The following objects are masked from Boston (pos = 3):

    age, black, chas, crim, dis, indus, lstat, medv, nox, ptratio, rad, rm, tax, zn

The following objects are masked from Boston (pos = 4):

    age, black, chas, crim, dis, indus, lstat, medv, nox, ptratio, rad, rm, tax, zn

The following objects are masked from Boston (pos = 5):

    age, black, chas, crim, dis, indus, lstat, medv, nox, ptratio, rad, rm, tax, zn

The following objects are masked from Boston (pos = 6):

    age, black, chas, crim, dis, indus, lstat, medv, nox, ptratio, rad, rm, tax, zn

The following objects are masked from Boston (pos = 7):

    age, black, chas, crim, dis, indus, lstat, medv, nox, ptratio, rad, rm, tax, zn

# Hàm plot để tạo một đồ thị phân tán với giá trị đầu vào x= lstat và đầu ra y= medv
plot(lstat,medv)
# abline tạo ra 1 đường tuyến tính fit với mô hình hồi quy hiện tại với giá trị đầu vào là lstat và đầu ra là medv
abline(lm.fit)

# Điều chỉnh độ rộng hay độ dày của đường tuyến tính để chúng ta nhìn thấy đường thẳng rõ nét hơn
abline(lm.fit,lwd=3)

# Điều chỉnh độ rộng hay độ dày của đường tuyến tính để chúng ta nhìn thấy đường thẳng rõ nét hơn
plot(lstat,medv)
abline(lm.fit,lwd=3)

# Thay đổi màu sắc đường tuyến tính thành màu đỏ
plot(lstat,medv)
abline(lm.fit,lwd=3,col="red")

# Thay đổi các điểm dữ liệu của đồ thị phân tán thành màu xanh với đầu vào là lstat và đầu ra là medv 
plot(lstat,medv,col="blue")

# Thay đổi hình dạng những điểm dữ liệu dưới dạng những viên đạn và có vòng tròn nhỏ hơn
plot(lstat,medv,pch=20)

# Thay đổi hình dạng những điểm dữ liệu những dấu + 
plot(lstat,medv,pch="+")

# Biểu diễn những điểm dữ liệu với x và y từ 1->20 và hình dạng những điểm dữ liệu cũng thay đổi liên tục từ dạng 1->20
plot(1:20,1:20,pch=1:20)

#Biểu đồ phân tán của các biến  dữ liệu với các dữ liệu cột được chia thành cửa sổ 2x2 bằng cách sử dụng hàm mfrow=c(2,2)
par(mfrow=c(2,2))

plot(lm.fit)

• Biểu đồ 1 (Residuals vs Fitted) giúp quan sát những điểm giá trị được dự đoán(fit) và giá trị sai lệch dư thừa của những điểm dự đoán này so với những điểm thật

• Biểu đồ 2 (Normal Q-Q) giúp so sánh 2 phân phối xác suất bằng cách vẽ các lượng tử của chúng với nhau

• Biểu đồ 3 (Scale- Location) tương tự như biểu đồ 1 nhưng khác là lấy căn bậc 2 giá trị dư thừa chuẩn hoá

• Biểu đồ 4 (Residuals vs Leverage) giúp phát hiện những ngoại lệ trong mô hình hồi quy đơn giản, cho phép quan sát sự lan rộng của những khoảng dư thừa thay đổi theo giá trị leverage( sự nhạy cảm của giá trị dự doán đối với sự thay đổi của những giá trị thật)

# Mô hình phân tán điểm dữ liệu được dự đoán và giá trị sai lệch dư thừa được tính từ mô hình hồi quy sử dụng hàm residuals()
plot(predict(lm.fit),residuals(lm.fit))

# Mô hình phân tán điểm dữ liệu được dự đoán và giá trị sai lệch dư thừa được tính từ phân phối student sử dụng hàm rstudent()
plot(predict(lm.fit),rstudent(lm.fit))

#Trên cơ sở các mô hình plot các giá trị dư thừa, có một số bằng chứng về phi tuyến tính. Thống kê đòn bẩy(leverage statistics) có thể được tính cho bất kỳ giá trị dự đoán nào bằng cách sử dụng Hàm hatvalues ()
plot(hatvalues(lm.fit))

# Chỉ mục cho chúng ta dòng nào có giá trị thống kê đòn đòn bẩy (leverage statistics) là lớn nhất
which.max(hatvalues(lm.fit))

375 
375