Simple Linear Regression

Mô hình Hồi quy đơn giản

Như đã biết, để có thể xác định mối tương quan giữa 2 biến số liên tục x và y chúng ta có thể sử dụng hệ số tương quan:

• Hệ số tương quan Pearson (Có tham số)
• Hệ số tương quan Spearman (Phi tham số)

Hệ số tương quan r (được chuẩn hóa từ hiệp biến của 2 biến số x và y) phản ánh mối tương quan giữa x và y. Tuy nhiên để có thể mô mình hóa được hệ số tương quan này, chúng ta cần sử dụng tới hồi quy tuyến tính. Ví dụ được đưa ra trong bài viết đề cập tới 2 biến liên tục trong bộ dữ liệu nhỏ giải thích về thời gian mang thai ảnh hướng tới cân nặng của trẻ sơ sinh (dữ liệu thai kỳ).

Một chút về lịch sử:

Chúng ta sẽ nhắc lại một chút về lịch sử mô hình hồi quy. Mô hình Hồi quy tuyến tính có xuất phát điểm phát minh từ nhà khoa học lừng danh Francis Galton (em họ của Charles Darwin). Sir Francis Galton (16 tháng 2 năm 1822 - 17 tháng 1 năm 1911) là một nhà thông thái, nhà nhân chủng học, ưu sinh học, nhà thám hiểm, nhà địa lý, nhà phát minh, nhà khí tượng học và nhà thống kê người Anh. Ông được phong tước hiệp sĩ năm 1909.

Trong sự nghiệp của mình, Galton viết trên 340 bài báo và sách. Ông cũng là người sáng tạo ra khái niệm hệ số tương quan trong thống kê và phát triển hồi quy về điểm trung bình. Ông là người đầu tiên áp dụng các phương pháp thống kê đối với nghiên cứu về sự khác biệt của con người và tính di truyền của trí thông minh, ông giới thiệu việc sử dụng bảng hỏi và khảo sát nhằm thu thập dữ liệu với con người.

Francis Galton là nhà tiên phong của thuyết ưu sinh và là người đã đặt tên thuyết ưu sinh (tiếng Anh: eugenics). Ông cũng là người sáng lập ra ngành trắc nghiệm học tâm lý. Ông phát minh ra phương pháp phân loại dấu vân tay, đóng vai trò hữu dụng đối với khoa học pháp lý.

Mấu chốt nhất trong nghiên cứu thống kê của Galton đó là đã phát hiện ra mối liên quan giữa chiều cao của những người có tầng lớp xã hội cao, trí thức liên quan đến IQ của họ vì vào thời đó những người có chiều cao hơn với những người khác đều đại đa số là quan chức và tầng lớp trí thức còn những người có chiều cao thấp hơn đều là tầng lớp lao động bình thường của xã hội. Tuy nhiên, mặc dù thiên tài nhưng Galoton lại không thể đưa ra được một cách chứng minh và định lượng chính xác về vấn đề ông nghiên cứu mà sau này chính Karl Pearson (cha đẻ của hệ số tương quan Pearson) đã giúp ông định lượng vấn đề này.

Dữ liệu thai kỳ:

Bộ dữ liệu bao gồm 2 biến

• Biến thời gian(period): Miêu tả thời gian (tuần) của bà mẹ mang thai
• Biến cân nặng(weight): Miêu tả cân nặng (kg) của thai nhi

Chúng ta sẽ mô phỏng mối tương quan của 2 biến này và tính hệ số tương quan Pearson, miêu tả qua ma trận tương quan để trực quan dữ liệu.

## Warning: package 'psych' was built under R version 4.0.2

## [1] 0.8175185

Phương trình Hồi quy

Dạng biểu diễn Để giới thiệu mô hình hồi quy tuyến tính, chúng ta cần đến vài ký hiệu. Gọi \((x, y)\) là cặp giá trị \(x\) và \(y\) củ đối tượng \(i(i =1,2,3,....,n)\) trong đó \(x\) coi là period và \(y\) là weight. Trong ví dụ trên \(n =17\) cá nhân. Mô hình hồi quy tuyến tính phát biểu rằng giá trị \(y_{i}\) là hàm số của một giá trị ban đầu và một hệ số liên quan với x:

\(y_{i} = \alpha + \beta*x_{i}\)

Trong đó, \(\alpha\) là giá trị ban đầu của \(y\), còn gọi là intercept; \(\beta\) là hệ số liên quan tới \(x\), và thường được gọi là slope hay hệ số góc của mô hình hồi quy. Trong mô hình này, \(\alpha\) và \(\beta\) là 2 hằng số. Có thể xem mô hình trên là đại diện cho dữ liệu về mối liên quan giữa \(x\) và \(y\). Tuy nhiên chúng ta thấy rằng, đường thẳng này đều không thể nối tất cả các điểm có giá trị \((x, y)\) lại với nhau được; sẽ có một số giá trị lệch khỏi mô hinh. Do đó, chúng ta thêm một yếu tố khác của mô hình và ký hiệu là \(\epsilon_{i}\).

\(y_{i} = \alpha + \beta*x_{i} + \epsilon_{i}\)

Chúng ta sẽ làm đơn giản ký hiệu bằng cách bỏ ký hiệu \(i\) cho mỗi quan sát:

\(y = \alpha + \beta*x + \epsilon\)

Đây là mô hình cho quần thể (population), cần dùng 17 quan sát được để ước tính 2 tham số cho mô hình quần thể. Một cách tổng quát:

Giá trị quan sát của \(y\) = Giá trị tiên lượng + Phần dư (residual)

\(y =\hat{y} + e\)

Giá trị \(\hat{y}\) còn được gọi là expected value. Phần \((a+b*x)\) là phần mô hình hệ thống, và phần nhiễu(noise) \(e\) Ước tính tham số: Giả sử chúng ta có 2 điểm \((x_{1}, y_{1})\) và \((x_{2}, y_{2})\) trong một khoảng không gian \(x-y\), và mục tiêu là tìm một đường thẳng để nối 2 điểm này. Cách đơn giản nhất là dùng slope (hệ số góc) và gradient(độ dốc).

alt

\(slope = \frac{\Delta y}{\Delta x}\) = \(\frac{y_{2}-y_{1}}{x_{2}-x_{1}}\)

Mục tiêu đưa ra là tìm \(\alpha\) và \(\beta\) sao cho:

\(min\Sigma (y-\alpha-\beta*x)^2\) là bé nhất.

Ước tính tham số bằng R

Trong R, dùng hàm \(lm()\) dùng để ước tính tham số cho mô hình hồi quy tuyến tính bằng phương pháp bình phương tối thiểu. Cú pháp: \(lm( y\sim x, data)\)

Trong đó, \(y\) là biến phụ thuộc, \(x\) là biến tiên lượng hay còn gọi là biến độc lập. Theo bộ dữ liệu thaiky phía trên, chúng ta có thể mô tả như sau:

\(weight = \alpha + \beta(period) + \epsilon\)

m = lm(weight ~ period, data = thaiky)
summary(m)

## 
## Call:
## lm(formula = weight ~ period, data = thaiky)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -958.6 -157.2    7.7  201.0  667.8 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -4020.05    1263.05  -3.183  0.00618 ** 
## period        180.46      32.82   5.498 6.13e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 414.4 on 15 degrees of freedom
## Multiple R-squared:  0.6683, Adjusted R-squared:  0.6462 
## F-statistic: 30.23 on 1 and 15 DF,  p-value: 6.13e-05

Theo kết quả đưa ra, mô hình hồi quy sẽ là: \(Weight = -4020.5 + 180.46 * Period\)

Diễn giải tham số:

Chúng ta có thể thể hiện mô hình này bằng đường thẳng như sau. Đường biểu diễn màu đỏ chính là giá trị tiên lượng (-4020.5 +180.46 * period) của mô hình:

plot(weight ~ period, pch = 16, data = thaiky)
abline(lm(weight ~ period), col = 'red')

Ý nghĩa của slope

Tham số quan trọng nhất của mô hình hồi quy chính là \(\beta\) hay slope. Để hiểu rõ hơn ý nghĩa của tham số này, chúng ta dùng mô hình dự đoán cho 2 quan sát sau:

1. Quan sát thêm 1: Nếu thai kỳ là 40 tuần, trọng lượng trẻ sơ sinh sẽ là: \(weight = -4020.5 + 180.46 * 40 = 3197.9 (gram)\)
2. Nếu thai kỳ là 41 tuần, trọng lượng trẻ sơ sinh sẽ là: \(weight = -4020.5 + 180.46 * 41 = 3378.36 (gram)\) Do đó, sự khác biệt về thai kỳ là 3378.36-3197.9 = 180.46 (gram)

   Đây chính là slope. Mô hình này có ý nghĩa là khi thai kỳ được tăng lên 1 đơn vị (tuần) thì trọng lượng trẻ sơ sinh sẽ tăng lên xấp xỉ 180.46 gram.

Ý nghĩa cửa \(\alpha\)

Hằng số intercept = -4020.5 có nghĩa là khi period = 0 thì trọng lượng là -4020.5. Điều này cho thấy rằng hằng số \(\alpha\) trong mô hình hồi quy tuyến tính không có ý nghĩa. Không có ý nghĩa bởi vì trọng lượng không thể là số âm, và trong thực tế thì không thể nào có sản phụ nào có thai kỳ là 0 tuần. Do đó, đôi khi chúng ta xem \(\alpha\) là một “hằng số giả định nào đó”. Tuy nhiên, chúng ta có thể biến \(\alpha\) có ý nghĩa bằng cách hoán chuyển biến period sang đơn vị \(z\):

\(zperiod = \frac{period - trung bình(period)}{Độ lệch chuẩn của period}\)

Giá trị trung bình của period là 38.5 tuần, độ lệch chuẩn là 3.1 tuần. Điều này có nghĩa là khi trung bình period là 38.5 tuần thì zperiod bằng 0. Hoán chuyển lại thang đo: \(zperiod = scale(period)\) Theo thang đo mới thì \(zperiod\) sẽ có trung bình (\(\overline{zperiod}\)) là 0 và độ lệch chuẩn là 1. Áp dụng mô hình hồi quy cho \(zperiod\):

\(m = lm(weight \sim zperiod, data = thaiky)\) \(summary(m)\)

thaiky$zperiod <- scale(thaiky$period)
m=lm(weight ~ zperiod, data = thaiky)
summary(m)

## 
## Call:
## lm(formula = weight ~ zperiod, data = thaiky)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -958.6 -157.2    7.7  201.0  667.8 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   2902.0      100.5  28.872 1.46e-14 ***
## zperiod        569.6      103.6   5.498 6.13e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 414.4 on 15 degrees of freedom
## Multiple R-squared:  0.6683, Adjusted R-squared:  0.6462 
## F-statistic: 30.23 on 1 and 15 DF,  p-value: 6.13e-05

plot(density(thaiky$zperiod), xlab = "zperiod", main = "The density of zperiod")

Kiểm tra giả định:

Mô hình hồi quy tuyến tính, cũng tương tự như phương pháp phân tích tương quan dựa vào 4 giả định sau:

• Thứ nhất, mối liên hệ giữa \(x\) và \(y\) phải là mối quan hệ tuyến tính và là tham số tuyến tính
• Thứ hai, phân bố của phần dư \(e_{i}\) phải tuân theo (hay xấp xỉ) luật phân bố chuẩn với trung bình \(0\) và phương sai \(s^2\), theo ký hiệu: \(e \sim N(0,\sigma^2)\)
• Thứ ba, phương sai của \(y\) không thay đổi theo giá trị \(x\)
• Thứ tư, các giá trị của \(y_{i}\) và \(y_{i-1}\) độc lập với nhau. Chúng ta hoàn toàn có thể kiểm tra các giả định này thông qua autoplot() trong package ggfortify

#install.packages("ggfortify")
library(ggfortify)

## Warning: package 'ggfortify' was built under R version 4.0.3

## Loading required package: ggplot2

## Warning: package 'ggplot2' was built under R version 4.0.2

## 
## Attaching package: 'ggplot2'

## The following objects are masked from 'package:psych':
## 
##     %+%, alpha

autoplot(m)

## Warning: `arrange_()` is deprecated as of dplyr 0.7.0.
## Please use `arrange()` instead.
## See vignette('programming') for more help
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_warnings()` to see where this warning was generated.

• Biểu đồ 1 bên trái phần trên đưa ra mối liên hệ giữa giá trị dự đoán (prediction) \(\hat y_{i}\) và phần dư \(e_{i}\), cho thấy các giá trị phần dư xoay quanh 0, tức đúng với giả định một về giá trị trung bình của phần dư
• Biểu đồ 2 bên phải phần trên đưa ra mối liên hệ giữa giá trị lý thuyết và thực tế của phần dư. Nếu phân phối của phần dư tuân theo luật phân phối Chuẩn thì các giá trị nằm trên đường lý thuyết (lệch với trục \(O_{x}\)) góc 45 độ). Do đó theo biểu đồ trên, phù hợp với giả định thứ hai.
• Biểu đồ phần dưới bên trái trình bày mối tương quan giữa giá trị dự đoán và căn số bậc hai của phần dư. Biểu đồ này cho chúng ta biết phương sai của phần dư có hay không liên quan đến giá trị biến trên trục x. Biểu đồ cho thấy không có mối liên quan, và điều này cho chúng ta biết rằng giả định về đồng dạng phương sai (homogeneity) là có thể chấp nhận được.
• Biểu đồ thứ 4 dưới cùng bên phải hàng dưới trình bày giá trị “Leverage” và phần dư chuẩn hóa. Biểu đồ này cho chúng ta biết có những giá trị có ảnh hưởng cao hay không (influential observation). Ngoại trừ một vài quan sát (số 17), phần còn lại đều có giá trị phần dư chuẩn hóa nằm trong khoảng từ -2 đến +2. Do đó có thể chấp nhận được rằng không có giá trị ngoại lai (outliers) ảnh hưởng đến mô hình.