Câu 1 (5 điểm):

Sử dụng tập airquality, dự đoán Ozone dựa vào Wind và Temp bằng mô hình hồi quy tuyến tính.

Kiểm tra tính đa cộng tuyến bằng VIF.

1.1 Chuẩn bị dữ liệu

# Load dữ liệu
data(airquality)

# Loại bỏ NA
aq <- na.omit(airquality)

head(aq)
##   Ozone Solar.R Wind Temp Month Day
## 1    41     190  7.4   67     5   1
## 2    36     118  8.0   72     5   2
## 3    12     149 12.6   74     5   3
## 4    18     313 11.5   62     5   4
## 7    23     299  8.6   65     5   7
## 8    19      99 13.8   59     5   8

1.2 Xây dựng mô hình hồi quy tuyến tính

model_lm <- lm(Ozone ~ Wind + Temp, data = aq)

summary(model_lm)
## 
## Call:
## lm(formula = Ozone ~ Wind + Temp, data = aq)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -42.156 -13.216  -3.123  10.598  98.492 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -67.3220    23.6210  -2.850  0.00524 ** 
## Wind         -3.2948     0.6711  -4.909 3.26e-06 ***
## Temp          1.8276     0.2506   7.294 5.29e-11 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 21.73 on 108 degrees of freedom
## Multiple R-squared:  0.5814, Adjusted R-squared:  0.5736 
## F-statistic: 74.99 on 2 and 108 DF,  p-value: < 2.2e-16

1.3 Kiểm tra đa cộng tuyến bằng VIF

library(car)

vif(model_lm)
##    Wind    Temp 
## 1.32837 1.32837

Câu 2 (5 điểm)

Dùng tập Heart (từ thư viện mlbench), xây dựng mô hình hồi quy logistic dự đoán AHD.

Tính ROC-AUC score

2.1 Load dữ liệu

library(mlbench)
data(PimaIndiansDiabetes)

str(PimaIndiansDiabetes)
## 'data.frame':    768 obs. of  9 variables:
##  $ pregnant: num  6 1 8 1 0 5 3 10 2 8 ...
##  $ glucose : num  148 85 183 89 137 116 78 115 197 125 ...
##  $ pressure: num  72 66 64 66 40 74 50 0 70 96 ...
##  $ triceps : num  35 29 0 23 35 0 32 0 45 0 ...
##  $ insulin : num  0 0 0 94 168 0 88 0 543 0 ...
##  $ mass    : num  33.6 26.6 23.3 28.1 43.1 25.6 31 35.3 30.5 0 ...
##  $ pedigree: num  0.627 0.351 0.672 0.167 2.288 ...
##  $ age     : num  50 31 32 21 33 30 26 29 53 54 ...
##  $ diabetes: Factor w/ 2 levels "neg","pos": 2 1 2 1 2 1 2 1 2 2 ...

2.2 Chia train/test

set.seed(123)

n <- nrow(PimaIndiansDiabetes)
train_index <- sample(1:n, 0.7*n)

train <- PimaIndiansDiabetes[train_index, ]
test  <- PimaIndiansDiabetes[-train_index, ]

2.3 Logistic regression

model_logit <- glm(diabetes ~ ., 
                   data = train, 
                   family = binomial)

summary(model_logit)
## 
## Call:
## glm(formula = diabetes ~ ., family = binomial, data = train)
## 
## Coefficients:
##              Estimate Std. Error z value Pr(>|z|)    
## (Intercept) -8.405409   0.841872  -9.984  < 2e-16 ***
## pregnant     0.103471   0.037973   2.725  0.00643 ** 
## glucose      0.035730   0.004563   7.830 4.89e-15 ***
## pressure    -0.012707   0.006057  -2.098  0.03590 *  
## triceps      0.003563   0.008088   0.440  0.65959    
## insulin     -0.001710   0.001060  -1.613  0.10671    
## mass         0.088735   0.017954   4.942 7.72e-07 ***
## pedigree     0.696250   0.334761   2.080  0.03754 *  
## age          0.017015   0.011066   1.538  0.12415    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 694.17  on 536  degrees of freedom
## Residual deviance: 509.76  on 528  degrees of freedom
## AIC: 527.76
## 
## Number of Fisher Scoring iterations: 5
prob <- predict(model_logit, test, type = "response")
head(prob)
##          1          3          4          9         15         17 
## 0.72751772 0.77340088 0.04425949 0.71110216 0.61359933 0.37579888

2.4 Dự đoán xác suất

roc_obj <- roc(test$diabetes, prob)
## Setting levels: control = neg, case = pos
## Setting direction: controls < cases
auc(roc_obj)
## Area under the curve: 0.8445

2.5 ROC – AUC

plot(roc_obj, col="blue", main="ROC Curve")