Pearson

Pearson’s product-moment correlation (hay Pearson’s correlation) là một trong những phương pháp để ước lượng mối quan hệ giữa 2 biến liên tục.

Điều kiện áp dụng:

• 2 biến cần tính tương quan là 2 biến liên tục, có phân phối chuẩn.
• Mỗi giá trị của biến này sẽ liên quan với 1 giá trị của biến kia theo từng cặp.
• Phải có mối quan hệ tuyến tính giữa 2 biến này.
• Không có outliers đáng kể.

Hệ số tương quan Pearson (r) có giá trị từ -1 cho đến 1, được đánh giá như sau:

• r = 0, không có mối tương quan
• r = 1, tương quan tuyệt đối
• r < 0, tương quan nghịch
• r > 0, tương quan thuận

Và

• 0 < |r| < 0.2, tương quan rất yếu
• 0.2 < |r| < 0.4, tương quan yếu
• 0.4 < |r| < 0.6, tương quan trung bình
• 0.6 < |r| < 0.8, tương quan mạnh
• 0.8 < |r| < 1, tương quan rất mạnh

Thực hành trên R

Số liệu mẫu

# Ta lấy ngẫu nhiên mẫu 2000 trong 50000 trường hợp đầu tiên của dataset diamonds làm ví dụ
df = ggplot2::diamonds[sample(1:50000, 2000),]
head(df)

## # A tibble: 6 × 10
##   carat cut       color clarity depth table price     x     y     z
##   <dbl> <ord>     <ord> <ord>   <dbl> <dbl> <int> <dbl> <dbl> <dbl>
## 1  1.5  Very Good J     I1       63.3    58  5141  7.27  7.23  4.59
## 2  0.74 Fair      E     SI1      66      56  2069  5.6   5.67  3.72
## 3  0.38 Ideal     E     SI1      61.9    55   700  4.64  4.7   2.89
## 4  0.52 Ideal     D     SI2      62.5    56  1289  5.08  5.12  3.19
## 5  0.33 Very Good D     VS1      63.2    56  1109  4.45  4.44  2.81
## 6  2.05 Good      G     SI2      63.7    61 16148  7.99  8.03  5.1

(Lưu ý là hàm sample() là một hàm ngẫu nhiên, do đó mỗi lần chạy sẽ cho các số liệu khác nhau, dẫn tới kết quả sẽ khác nhau giữa từng lần chạy nhé)

Ta sẽ thử áp dụng Pearson’s correlation cho 2 biến liên tục là depth và table

Bước 1. Kiểm tra phân bố chuẩn của 2 biến này bằng biểu đồ boxplot

(Có nhiều phương pháp để kiểm tra phân bố chuẩn, tuy nhiên nếu dùng boxplot ta sẽ tiện kiểm tra outliers luôn)

# Boxplot của biến depth
p1 = ggplot(df) +
  geom_boxplot(aes(y = depth))

# Boxplot của biến table
p2 = ggplot(df) +
  geom_boxplot(aes(y = table))

# Print 2 boxplot cùng nhau
ggarrange(p1, p2, 
          labels = c("Biến depth", "Biến table"),
          ncol = 2, nrow = 1)

B2. Kiểm tra mối quan hệ tuyến tính giữa 2 biến

# Sử dụng biểu đồ phân tán để kiểm tra
df %>% 
  ggplot(aes(x = depth, y = table)) +
  geom_point()

Ta thấy có xu hướng khi depth tăng thì table giảm, theo 1 đường thẳng (tuyến tính)

B3. Không có outlier đáng kể.

Từ 2 biểu đồ trên cho thấy có outliers (biến ngoại lai). Để áp dụng Pearson’s correlation ta phải loại bỏ các outliers này.

# Tìm ra các giá trị outliers
out_depth = boxplot(df$depth, plot = F)$out
out_table = boxplot(df$table, plot = F)$out

# Loại bỏ các trường hợp outliers này ra khỏi data frame
df = df[-which(df$depth %in% out_depth), ]
df = df[-which(df$table %in% out_table), ]

Sau đó kiểm tra lại

df %>% 
  ggplot(aes(x = depth, y = table)) +
  geom_point()

B4. Tính hệ số tương quan Pearson

cor(df$depth, df$table, method = "pearson")

## [1] -0.219185

Sau khi có hệ số tương quan. Ta so sánh giá trị này với phần đánh giá ở trên nhé!

Kiểm định hệ số tương quan (r = 0)

cor.test(df$depth, df$table, method = "pearson")

## 
##  Pearson's product-moment correlation
## 
## data:  df$depth and df$table
## t = -9.7664, df = 1890, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.2616653 -0.1758572
## sample estimates:
##       cor 
## -0.219185

Giả sử ta thu được kết quả như sau:

r = -0.256, 95%CI: [-0.298; -0.213], p < 0.001

Kết luận đầy đủ là:

depth và table có mối tương quan nghịch, mức độ tương quan là yếu và mối tương quan này có ý nghĩa thống kê.

Spearman

Spearman’s rank-order correlation (hay Spearman’s correlation) là phương pháp nonparametric để ước lượng mối quan hệ của 2 biến thứ bậc (ordinal). Có thể dùng thay thế cho Pearson’s correlation khi 2 biến liên tục không có phân bố chuẩn.

Điều kiện áp dụng:

• 2 biến cần tính tương quan là 2 biến ordinal hoặc continuous.
• Mỗi giá trị của biến này sẽ liên quan với 1 giá trị của biến kia theo từng cặp.
• Phải có mối quan hệ đơn điệu (monotonic) giữa 2 biến này.

Quan hệ đơn điệu tức là nếu biến này tăng thì biến kia một là tăng, hai là giảm, không thể vừa tăng và vừa giảm.

Xem ví dụ bên dưới để biết thế nào là không đơn điệu:

data <- data.frame(hours=c(6, 9, 12, 14, 30, 35, 40, 47, 51, 55, 60),
                   happiness=c(14, 28, 50, 70, 89, 94, 90, 75, 59, 44, 27))

ggplot(data, aes(x = hours, y = happiness)) +
  geom_point() +
  stat_smooth(method = "lm", formula = y ~ x + I(x^2), size = 1)

Bạn có thể thấy, ban đầu x tăng thì y tăng, nhưng lúc sau x tăng thì y giảm. Đó là ví dụ về tính không đơn điệu.

Áp dụng Spearman’s correlation

Ta có thể dùng lại ví dụ của Pearson’s correlation để thử áp dụng, bằng cách thay đổi method sang spearman.

cor(df$depth, df$table, method = "spearman")

## [1] -0.1865611

Hệ số tương quan của phương pháp Spearman gọi là rho. Vì rho = -0.21, ta có thể nói có mối tương quan nghịch mức độ yếu giữa 2 biến depth và table.

Kiểm định hệ số tương quan (rho = 0)

cor.test(df$depth, df$table, method = "spearman")

## Warning in cor.test.default(df$depth, df$table, method = "spearman"): Cannot
## compute exact p-value with ties

## 
##  Spearman's rank correlation rho
## 
## data:  df$depth and df$table
## S = 1339374826, p-value = 2.805e-16
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##        rho 
## -0.1865611

Giả sử ta thu được rho = -0.21, p <0.001, ta kết luận như sau:

Vì p < 0.001, ta có thể mạnh dạn bác bỏ H0 (rho = 0). Ta kết luận hệ số tương quan giữa depth và table có ý nghĩa thống kê (rho = -0.21, p < 0.001)

Kendall

Ở ví dụ chạy tương quan Spearman, bạn sẽ thấy cảnh báo “Cannot compute exact p-value with ties”. Vậy ties là gì?

Tied data là khi số liệu bị trùng lặp.

Ví dụ: Ta có số liệu 1,2,3,4,4,4,5,5,6,7,8,9. Khi đó 4 và 5 bị lặp lại. Ta gọi là ties.

Ties ảnh hưởng đến phân bậc thứ hạng (rank) trong các kiểm định phi tham số (nonparametric).

Kendall’s correlation giúp kiểm tra hệ số tương quan khi số liệu có nhiều ties, cỡ mẫu nhỏ.

Áp dụng tính hệ số tương quan Kendall (tau):

cor(df$depth, df$table, method = "kendall")

## [1] -0.1372866

Kiểm định hệ số tương quan (tau = 0)

cor.test(df$depth, df$table, method = "kendall")

## 
##  Kendall's rank correlation tau
## 
## data:  df$depth and df$table
## z = -8.3118, p-value < 2.2e-16
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
##        tau 
## -0.1372866

Nhận xét tương tự 2 ví dụ về Pearson’s và Spearman’s correlation nhé!

Tuyến tính (Linear)

Là phương trình bậc 1, mô tả mối quan hệ giữa biến độc lập và biến phụ thuộc thông qua 2 giá trị b0 và b.

• b0 còn gọi là hằng số, hệ số chặn, constant, intercept,…
• b còn gọi là hệ số, hệ số góc, độ nghiên, coefficient, slobe,…

Đồ thị của phương trình tuyến tính là một đường thẳng.

Khi xem xét mối quan hệ giữa 2 biến, người ta chú ý đến hệ số b hơn vì nó quy định mức độ ảnh hưởng của biến độc lập lên biến phụ thuộc:

• b = 0, biến độc lập không có tác động gì tới biến phụ thuộc.
• b < 0, giá trị của biến độc lập tăng sẽ làm giảm giá trị của biến phụ thuộc và ngược lại (tương quan nghịch).
• b > 0, giá trị của biến độc lập tăng làm cho giá trị của biến phụ thuộc cũng tăng theo (tương quan thuận).
• b càng lớn, ảnh hưởng của biến độc lập lên biến phụ thuộc càng nhiều.
• b càng bé, ảnh hưởng của biến độc lập lên biến phụ thuộc càng ít.

data = read_csv("income.data.csv")
head(data)

## # A tibble: 6 × 3
##    ...1 income happiness
##   <dbl>  <dbl>     <dbl>
## 1     1   3.86      2.31
## 2     2   4.98      3.43
## 3     3   4.92      4.60
## 4     4   3.21      2.79
## 5     5   7.20      5.60
## 6     6   3.73      2.46

p1 = data %>% ggplot(aes(y = income)) + geom_boxplot()
p2 = data %>% ggplot(aes(y = happiness)) + geom_boxplot()

ggarrange(p1, p2, 
          labels = c("Biến income", "Biến happiness"),
          ncol = 2, nrow = 1)

data %>% 
  ggplot(aes(x = income, y = happiness)) +
  geom_point()

mo_hinh = lm(formula = happiness~income, data = data)
summary(mo_hinh)

## 
## Call:
## lm(formula = happiness ~ income, data = data)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2.02479 -0.48526  0.04078  0.45898  2.37805 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.20427    0.08884   2.299   0.0219 *  
## income       0.71383    0.01854  38.505   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.7181 on 496 degrees of freedom
## Multiple R-squared:  0.7493, Adjusted R-squared:  0.7488 
## F-statistic:  1483 on 1 and 496 DF,  p-value: < 2.2e-16

par(mfrow=c(2,2)) # Chia output ra thành bảng 2x2 (chứa 4 biểu đồ)
plot(mo_hinh)

Logistics (Logit)

Tên đầy đủ: Hồi quy logistic nhị phân (Binomial logistics regression, Binomial logit model, logit model, logistics regression).

Khác với hồi quy tuyến tính, vì biến phụ thuộc trong hồi quy logistic chỉ nhận 2 giá trị (biến nhị phân) nên ta không thể áp dụng phương trình y = b0 + bx được.

Để áp dụng, người ta biến đổi y thành phương trình logarit của OR bằng giá trị y’, lúc này y’ của chúng ta không còn là biến nhị phân nữa mà đã là biến liên tục. Lúc này, ta có thể áp dụng hồi quy tuyến tính. y’ = b0 +bx, với y’ là logarit tự nhiên của p/(1-p) (OR của y). Khi đó, ta đánh giá hệ số b giống như mô hình hồi quy tuyến tính.

Hệ số b thể hiện mối quan hệ giữa x và y’. Để xem xét trực tiếp mối quan hệ giữa biến độc lập và biến phụ thuộc (x và y) mà không phải gián tiếp qua y’, ta chỉ cần lấy e mũ cơ số b (exp(b)). Giá trị exp(b) lúc này gọi là OR của biến độc lập. Việc đánh giá exp(b) (hay OR) của biến độc lập có khác so với hệ số b trong hồi quy tuyến tính một chút:

• OR luôn là số dương.
• OR = 1, biến độc lập không có ảnh hưởng gì đến biến phụ thuộc.
• OR < 1, biến độc lập tăng làm giảm (nguy cơ) giá trị của biến phụ thuộc.

Vì biến phụ thuộc lúc này chỉ có 2 giá trị, nếu ta chọn 1 giá trị làm giá trị tham chiếu (gọi là A), thì OR < 1 có nghĩa là biến độc lập tăng làm giảm nguy cơ biến phụ thuộc có giá trị còn lại kia (giá trị B). OR > 1, biến độc lập tăng làm tăng nguy cơ biến phụ thuộc có giá trị là B.

^^ Ta có thể diễn tả theo cách ngược lại như sau: OR < 1, biến độc lập tăng làm tăng nguy cơ biến phụ thuộc có giá trị A. OR > 1, biến độc lập tăng làm giảm nguy cơ biến phụ thuộc có giá trị A. Dùng cách nào cũng được. Tuy nhiên để thống nhất với nhau, ta sẽ dùng cách đầu tiên.

head(mtcars)

##                    mpg cyl disp  hp drat    wt  qsec vs am gear carb
## Mazda RX4         21.0   6  160 110 3.90 2.620 16.46  0  1    4    4
## Mazda RX4 Wag     21.0   6  160 110 3.90 2.875 17.02  0  1    4    4
## Datsun 710        22.8   4  108  93 3.85 2.320 18.61  1  1    4    1
## Hornet 4 Drive    21.4   6  258 110 3.08 3.215 19.44  1  0    3    1
## Hornet Sportabout 18.7   8  360 175 3.15 3.440 17.02  0  0    3    2
## Valiant           18.1   6  225 105 2.76 3.460 20.22  1  0    3    1

mo_hinh = glm(data = mtcars, 
              formula = vs~wt,
              family = "binomial")
summary(mo_hinh)

## 
## Call:
## glm(formula = vs ~ wt, family = "binomial", data = mtcars)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -1.9003  -0.7641  -0.1559   0.7223   1.5736  
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)   
## (Intercept)   5.7147     2.3014   2.483  0.01302 * 
## wt           -1.9105     0.7279  -2.625  0.00867 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 43.860  on 31  degrees of freedom
## Residual deviance: 31.367  on 30  degrees of freedom
## AIC: 35.367
## 
## Number of Fisher Scoring iterations: 5

# Lấy các giá trị của log OR biến vs
or_vs = predict(mo_hinh, type = "response")

# Hàm biến đổi
logit = log(or_vs/(1-or_vs))

# Kiểm tra mối quan hệ của hàm biến đổi này và biến wt
ggplot() +geom_point(aes(x = mtcars$wt, y = logit))

exp(-1.911)

## [1] 0.1479324