Tái chọn mẫu ngẫu nhiên là một chủ đề không mới trong môn Thống kê, tuy nhiên gần đây nó có tính « thời sự » hơn khi Machine learning trở thành một công cụ suy diễn thống kê mới. Tái chọn mẫu giữ cho các mô hình thống kê được an toàn trước nguy cơ « overfiting » và được khuyến cáo nên áp dụng thường quy cho mọi phân tích (suy đến cùng, thống kê cổ điển cũng chính là các mô hình).
Một quy trình dựng mô hình hiện đại bao gồm 2 công đoạn : Tái chọn mẫu trên tập dữ liệu riêng (Trainset) và Kiểm chứng độc lập trên tập dữ liệu khác (Testset).
Trước kia, những phương pháp tái chọn mẫu như Bootstrap, K-folds Crossvalidation, Monte Carlo … ít được thực hiện đại trà, một phần do chủ quan, nhưng nguyên nhân chính là do sự phức tạp của chúng. Chỉ có người làm Stats chuyên nghiệp mới có khả năng viết R function để làm các quy trình này. Chưa có một chuẩn phổ quát trong R cho phép người dùng thiết lập, tùy chỉnh và kiểm soát toàn bộ các quy trình tái chọn mẫu.
Một lần nữa, hai tác giả Max Kuhn và Hadley Wickham lại cùng nhau giải quyết vấn đề này. Họ đã tạo ra một package rất tiện dụng có tên là rsample mà Nhi sẽ giới thiệu với các bạn trong bài hôm nay. Thực ra, bản chất của mọi quy trình tái chọn mẫu rất đơn giản : mục tiêu của chúng ta là tạo ra nhiều phiên bản dữ liệu bằng cách rút ngẫu nhiên các phần tử từ tập dữ liệu gốc. Mô hình/phân tích thống kê sẽ được dựng trên mỗi phiên bản này, từ đó một chuỗi kết quả sẽ được khảo sát để lượng giá kích thước hiệu ứng, phẩm chất/độ tin cậy của mô hình…
Bước đầu tiên của tái chọn mẫu là lập ra một danh sách mô tả chi tiết về các phiên bản dữ liệu này, bao gồm : số lượt lấy mẫu (bao nhiêu phiên bản sẽ được tạo ra), mỗi phần tử có thể xuất hiện nhiều lần hay không (bootstrap= chọn mẫu có lặp lại), có dành một phần để kiểm định hay không ? (và tỉ lệ bao nhiêu ? ), cấu trúc khối (kiểm chứng chéo K khối), có bảo toàn tỉ lệ phân bố của một biến nào đó hay không ? (thí dụ phân vị của 1 biến kết quả liên tục, tỉ lệ 2 nhãn giá trị của 1 biến kết quả nhị phân). Chính những đòi hỏi kỹ thuật chi tiết nhỏ nhặt này làm khó người dùng R, vì họ phải viết function để đáp ứng đồng thời yêu cầu về tính ngẫu nhiên và các điều kiện kèm theo.
Package rsample được tạo ra để giải quyết công đoạn khó nhất này, nó tạo ra một nền tảng phổ quát, cho phép người dùng kiểm soát đến từng chi tiết của tất cả quy trình tái chọn mẫu bằng những hàm có cú pháp tương tự. Chỉ 1 dòng code, ta có thể tạo danh sách tái chọn mẫu. Danh sách dữ liệu này lại có cấu trúc đồng nhất cho mọi quy trình, gắn kết dễ dàng với những hàm tiện ích khác từ broom, purrr, vv. và có thể dùng để lưu trữ trực tiếp kết quả thống kê/kiểm định mô hình từ một hàm map().
Trước khi bắt đầu, Nhi sẽ tải về máy dữ liệu Hypothyroid (bệnh Suy nhược tuyến giáp) của viện Garvan (Úc), gồm hơn 2000 bệnh nhân.
NHi sẽ dùng caret để cắt ngẫu nhiên dữ liệu gốc thành 2 phần với tỉ lệ 80% và 20% : Phần lớn (trainset) sẽ dùng cho phân tích của chúng ta :
library(tidyverse)
df=read.csv("https://www.openml.org/data/get_csv/53534/hypothyroid.csv",na.strings = "?")
df=df%>%dplyr::select(TSH,T3,TT4,T4U,FTI,binaryClass)%>%na.omit()
# Spliting
idx=caret::createDataPartition(y=df$binaryClass, p=0.8,list=FALSE)
trainset=df[idx,]
testset=df[-idx,]
trainset%>%head(5)%>%knitr::kable()
TSH | T3 | TT4 | T4U | FTI | binaryClass | |
---|---|---|---|---|---|---|
1 | 1.30 | 2.5 | 125 | 1.14 | 109 | P |
5 | 0.72 | 1.2 | 61 | 0.87 | 70 | P |
9 | 0.60 | 2.2 | 123 | 0.93 | 132 | P |
10 | 2.40 | 1.6 | 83 | 0.89 | 93 | P |
11 | 1.10 | 2.2 | 115 | 0.95 | 121 | P |
Đầu tiên, Nhi giả định ta muốn tạo ra một danh sách tái chọn mẫu sử dụng phương pháp Bootstrap. Package rsample cho phép làm điều này cực kì đơn giản, với 1 dòng code duy nhất bằng hàm bootstraps
library(rsample)
# A rsample object
set.seed(2105)
boot_obj<-bootstraps(trainset,times = 10,strata="binaryClass")
Nhi vừa áp dụng một quy trình bootstrap ngắn, với 10 lượt tái chọn mẫu, trên tập trainset, với điều kiện phụ là bảo toàn tỉ lệ biến binaryClass, và lưu danh sách phiên bản dữ liệu vào 1 object có tên là boot_obj
Ta biết rằng tập trainset có kích thước 177 Kb
library(mlbench)
library(pryr)
##
## Attaching package: 'pryr'
## The following objects are masked from 'package:purrr':
##
## compose, partial
object_size(trainset)
## 177 kB
object danh sách dữ liệu tái chọn mẫu có kích thước là 277 Kb (kích thước này sẽ tỉ lệ với số lượt lấy mẫu)
object_size(boot_obj)
## 277 kB
Như vậy mỗi phiên bản dữ liệu trong 1 lượt lấy mẫu chỉ có kích thước là 27.7 kb, rất nhỏ.
object_size(boot_obj)/nrow(boot_obj)
## 27.7 kB
Với 10 lượt lấy mẫu, kích thước của danh sách so với dữ liệu gốc chỉ tăng gấp rưỡi
(object_size(boot_obj)/object_size(trainset))%>%as.numeric()
## [1] 1.560318
Bây giờ, ta sẽ tìm hiểu cấu trúc bên trong object danh sách này nhé:
Đây là một cấu trúc dữ liệu 2 chiều, mỗi hàng tương ứng với 1 lượt bootstrap, nó có 2 cột, cột splits là một list, chứa tập dữ liệu được rút ngẫu nhiên từ trainset, cột id có vai trò định danh.
boot_obj
## # Bootstrap sampling using stratification
## # A tibble: 10 x 2
## splits id
## <list> <chr>
## 1 <S3: rsplit> Bootstrap01
## 2 <S3: rsplit> Bootstrap02
## 3 <S3: rsplit> Bootstrap03
## 4 <S3: rsplit> Bootstrap04
## 5 <S3: rsplit> Bootstrap05
## 6 <S3: rsplit> Bootstrap06
## 7 <S3: rsplit> Bootstrap07
## 8 <S3: rsplit> Bootstrap08
## 9 <S3: rsplit> Bootstrap09
## 10 <S3: rsplit> Bootstrap10
Ta thử trích 3 hàng đầu tiên của splits, và nhận thông tin như sau:
boot_obj$splits[c(1:3)]
## $`1`
## <2203/800/2203>
##
## $`2`
## <2203/785/2203>
##
## $`3`
## <2203/819/2203>
3 con số này lần lượt có tên là: analysis = tập dữ liệu chính, dùng để dựng mô hình, assessment = tập dữ liệu dùng kiểm định mô hình (chúng có nội dung hoàn toàn độc lập), và cuối cùng là tập dữ liệu gốc = trainset.
Như ta biết, quy trình bootstrap tạo ra một phiên bản dữ liệu có cùng kích thước với dữ liệu gốc, nhưng chứa những phần tử được rút ngẫu nhiên và có lặp lại từ dữ liệu gốc này. Đó là ý nghĩa con số 2203 = tập analysis; package rsample còn đi xa hơn, khi không chỉ làm bootstrap cho phân tích chính/mô hình mà còn bootstrap lần thứ 2 để tạo ra 1 tập kiểm định nhỏ độc lập với tập analysis, tập kiểm định assessment này chứa khoảng 800 phần tử.
Tiếp theo mới là điều thú vị: khi áp dụng hàm analysis trên 1 hàng của object danh sách, nó sẽ trích xuất phiên bản dữ liệu (tập analysis)
boot_obj$splits[[1]]%>%analysis()%>%head(5)%>%knitr::kable()
TSH | T3 | TT4 | T4U | FTI | binaryClass | |
---|---|---|---|---|---|---|
1 | 1.3 | 2.5 | 125 | 1.14 | 109 | P |
9 | 0.6 | 2.2 | 123 | 0.93 | 132 | P |
9.1 | 0.6 | 2.2 | 123 | 0.93 | 132 | P |
10 | 2.4 | 1.6 | 83 | 0.89 | 93 | P |
10.1 | 2.4 | 1.6 | 83 | 0.89 | 93 | P |
Tương tự, hàm assessment sẽ trích xuất phiên bản dữ liệu kiểm định:
boot_obj$splits[[2]]%>%assessment()%>%head(5)%>%knitr::kable()
TSH | T3 | TT4 | T4U | FTI | binaryClass | |
---|---|---|---|---|---|---|
1 | 1.300 | 2.5 | 125 | 1.14 | 109 | P |
15 | 3.300 | 1.8 | 109 | 0.91 | 119 | P |
29 | 1.900 | 1.5 | 113 | 1.06 | 106 | P |
31 | 0.035 | 2.5 | 119 | 1.55 | 76 | P |
34 | 1.700 | 1.9 | 95 | 1.05 | 90 | P |
Như vậy, bạn cũng có thể hình dung ta sẽ khai thác object danh sách này như thế nào khi dựng và kiểm định mô hình hàng loạt, thí dụ bằng hàm map hay vòng lặp for loop.
Ta kiểm tra thêm 1 chi tiết nữa, đó là điều kiện bảo toàn tỉ lệ biến kết quả (binaryClass): tỉ lệ này là như nhau ở mỗi phiên bản dữ liệu, và bằng với tỉ lệ trong tập trainset
boot_obj$splits[[5]]%>%as.data.frame()%>%.$binaryClass%>%table()
## .
## N P
## 188 2015
table(trainset$binaryClass)
##
## N P
## 179 2024
Thật không tệ ! Ta có thể yên tâm khi áp dụng quy trình này cho những bài toán với data imbalance.
Đầu tiên, Nhi sẽ minh họa một phân tích so sánh bằng bootstrap.
Nhi viết 1 function để tính khác biệt trung vị của biến TT4 giữa 2 phân nhóm P và N, áp dụng cho mỗi lượt bootstrap
Thử function này trên trainset, nó chạy tốt:
# A simple comparison
med_diff <- function(splits,target="TT4") {
temp<- analysis(splits)
res<-temp%>%group_by(binaryClass)%>%summarise_at(target,mean,na.rm=T)
median(res[2,2]%>%as.numeric())-median((res[1,2]%>%as.numeric()))
}
med_diff(trainset,"TT4")
## [1] 42.97813
Sau đó, chỉ với 2 dòng code, ta có thể làm một so sánh trung vị bằng bootstrap 1000 lần:
set.seed(2105)
bigboot <- bootstraps(trainset,times=1000)
bigboot$TT4_meddif<- map_dbl(bigboot$splits,med_diff)
Như bạn thấy, Nhi tạo danh sách dữ liệu bigboot trước, sau đó dùng chính danh sách này để lưu kết quả khác biệt trung vị bằng 1 hàm map của package purrr. Rất đơn giản, nhanh chóng
Tiếp theo ta có thể làm gì tùy thích, thí dụ tính 97.5%CI của khác biệt trung vị
quantile(bigboot$TT4_meddif,
probs = c(0.025, 0.05, 0.500, 0.95,0.975))
## 2.5% 5% 50% 95% 97.5%
## 37.62887 38.47300 43.00951 47.44953 48.21423
Hay làm một phản nghiệm với ngưỡng null hypothesis là +35
threshold = +35
WVPlots::ShadedDensity(bigboot,"TT4_meddif", threshold,
title="threshold= +35")+theme_bw()
Hay vẽ density plot cho kết quả bootstrap
bigboot%>%mutate(group=rep(c(1:10),100))%>%
ggplot(aes(x=TT4_meddif,col=factor(group)))+
geom_density(alpha=0.02,fill="red",show.legend = F)+
scale_color_brewer(palette="Reds")+
geom_vline(xintercept=35,linetype=2,col="black")+
theme_bw()
## Warning in RColorBrewer::brewer.pal(n, pal): n too large, allowed maximum for palette Reds is 9
## Returning the palette you asked for with that many colors
Bây giờ ta sẽ phức tạp hóa vấn đề lên một chút, cũng là so sánh TT4 giữa 2 phân nhóm P/N nhưng qua 1 mô hình hồi quy tuyến tính bằng hàm glm()
Mô hình so sánh có dạng: TT4 ~ binaryClass -1
Ta viết 1 hàm glm_coefs cho phép trích xuất coefficients của mô hình này, cho mỗi lượt bootstrap.
Sau đó gắn hàm này vào danh sách bigboot bằng 1 hàm map
# Model
mod_form=as.formula(TT4~binaryClass-1)
glm_coefs <- function(splits, ...) {
mod <- glm(..., data = analysis(splits),family="gaussian")
as.data.frame(t(coef(mod)))
}
coefs <- map(.x = bigboot$splits,
.f = glm_coefs,
mod_form)
object coefs là 1 list, nên hàm lapply cho phép trích xuất kết quả của 2 vectors phân nhóm P,N thật đơn giản:
coefdat=data_frame(Group=rep(c(1:10),100),
ClassN=lapply(coefs,`[[`, "binaryClassN")%>%as.numeric(),
ClassP=lapply(coefs,`[[`, "binaryClassP")%>%as.numeric())
library(ggridges)
coefdat%>%gather(ClassN,ClassP,key="Coef",value="Estimated")%>%
ggplot(aes(y=Coef,x=Estimated,col=factor(Group)))+
geom_density_ridges(alpha=0.02,fill="red",show.legend = F)+
scale_color_brewer(palette="Reds")+
theme_bw()
## Picking joint bandwidth of 0.553
## Warning in RColorBrewer::brewer.pal(n, pal): n too large, allowed maximum for palette Reds is 9
## Returning the palette you asked for with that many colors
Ta vừa mới khảo sát trung bình TT4 cho mỗi phân nhóm P/N trên 1000 phiên bản dữ liệu.
Trong 2 thí dụ nêu trên, ta chỉ mới dùng đến tập analysis trong danh sách bootstrap. Thí dụ sau đây Nhi sẽ thử một phương pháp tái chọn mẫu khác là K_folds cross validation, với cấu trúc 10x10
Giả định ta muốn dựng 1 mô hình logistic để phân loại 2 nhãn P/N dựa vào2 predictors là TSH và FTI
Ta có thể dùng toàn bộ tập trainset
# cross validation a logistic
trainset$Class=as.numeric(trainset$binaryClass)-1
lmod0=glm(data=trainset,binaryClass~TSH+FTI,family="binomial")
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
summary(lmod0)
##
## Call:
## glm(formula = binaryClass ~ TSH + FTI, family = "binomial", data = trainset)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.8503 0.1514 0.2208 0.2866 7.6591
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.542680 0.523932 1.036 0.3
## TSH -0.211228 0.019724 -10.709 < 2e-16 ***
## FTI 0.031448 0.005226 6.018 1.77e-09 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1241.69 on 2202 degrees of freedom
## Residual deviance: 658.02 on 2200 degrees of freedom
## AIC: 664.02
##
## Number of Fisher Scoring iterations: 8
Bây giờ, ta sẽ áp dụng một quy trình kiểm định chéo 10x10, tức chia trainset thành 10 blocks, dùng 9 blocks dựng mô hình, kiểm định trên block còn lại, và lặp lại 10 lần.
Hàm vfold_cv() cho phép tạo danh sách dữ liệu cho quy trình kiểm chứng chéo này
# KFCV
mod_form <- as.formula(Class~TSH+FTI)
rs_obj <- vfold_cv(trainset, V=10, repeats = 10)
rs_obj$splits[1]
## $`1`
## <1982/221/2203>
rs_obj
## # 10-fold cross-validation repeated 10 times
## # A tibble: 100 x 3
## splits id id2
## <list> <chr> <chr>
## 1 <S3: rsplit> Repeat01 Fold01
## 2 <S3: rsplit> Repeat01 Fold02
## 3 <S3: rsplit> Repeat01 Fold03
## 4 <S3: rsplit> Repeat01 Fold04
## 5 <S3: rsplit> Repeat01 Fold05
## 6 <S3: rsplit> Repeat01 Fold06
## 7 <S3: rsplit> Repeat01 Fold07
## 8 <S3: rsplit> Repeat01 Fold08
## 9 <S3: rsplit> Repeat01 Fold09
## 10 <S3: rsplit> Repeat01 Fold10
## # ... with 90 more rows
Như ta thấy, cấu trúc của danh sách có thay đổi chút ít nhưng tương tự như bootstrap: cột splits chứa những phiên bản dữ liệu, tập analysis gốm 9 blocks, kích thước khoảng 1980/2203, tập assessment là 1 block còn lại, khoảng 220/2203. Ta không áp dụng điều kiện bảo tồn tỉ lệ P/N
Tiếp theo, Nhi viết 1 hàm holdout_results với nội dung : dựng mô hình logistic cho mỗi tập analysis, sau đó dùng hàm augment của broom package để phân loại cho tập assessment tương ứng, so sánh kết quả phân loạibởi mô hình và giá trị thực của biến binaryClass, rồi đếm tần suất phân loại chính xác P/N
Như trên, ta dùng hàm map để áp dụng hàm holdout này cho danh sách tái chọn mẫu
holdout_results <- function(splits, ...) {
mod <- glm(..., data = analysis(splits), family = binomial)
holdout <- assessment(splits)
res <- broom::augment(mod, newdata = holdout)
lvls <- levels(holdout$binaryClass)
predictions <- factor(ifelse(res$.fitted > 0, lvls[2], lvls[1]),
levels = lvls)
res$correct <- predictions == holdout$binaryClass
res
}
rs_obj$results <- map(rs_obj$splits,
holdout_results,
mod_form)
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
rs_obj$accuracy <- map_dbl(rs_obj$results, function(x) mean(x$correct))
quantile(rs_obj$accuracy,probs = c(0.05, 0.500, 0.95))
## 5% 50% 95%
## 0.9272727 0.9502262 0.9730718
Kết quả cho thấy 95%CI của accuracy từ 92.73% đến 96.38%, trung vị accuracy=95%
Lưu ý rằng bản thân dữ liệu này bị bất xứng (imbalance), cụ thể tỉ lệ P/N là 92%, do đó mô hình logistic nêu trên chưa thể gọi là tối ưu vì Accuracy chỉ cao hơn 92% một chút.
table(trainset$binaryClass)
##
## N P
## 179 2024
2024/2203
## [1] 0.9187472
Ta có thể thay Accuracy bằng một tiêu chí kiểm định bất kì tùy thích, thí dụ Fscore, Recall, Precision, vv, chỉ cần đặt tiêu chí này vào hàm holdout
Kết quả kiểm chứng chéo 10x10 có thể vẽ thành biểu đồ như sau:
library(viridis)
## Loading required package: viridisLite
rs_obj%>%ggplot(aes(x=id2,y=id))+geom_tile(aes(fill=accuracy))+
scale_fill_viridis(option="D")+
scale_y_discrete("Repeat")+
scale_x_discrete("Folds")
Ta có thể thay kiểm chứng chéo bằng tái chọn mẫu Monte Carlo, kết quả Accuracy tương tự
# Monte Carlo cross-validation
mc_obj <- mc_cv(trainset, prop=0.8, times =100)
mc_obj$splits[1]
## [[1]]
## <1763/440/2203>
mc_obj$results <- map(mc_obj$splits,
holdout_results,
mod_form)
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
mc_obj$accuracy <- map_dbl(mc_obj$results, function(x) mean(x$correct))
mc_obj%>%mutate(Iter=as.numeric(row.names(.)))%>%
ggplot(aes(x=Iter,y=accuracy))+
geom_path(col="red")+
geom_point(shape=21,fill="red",col="red4")+
geom_hline(linetype=2,col="blue",size=1,yintercept = median(mc_obj$accuracy))+
theme_bw()
Bài thực hành đến đây là hết. Các bạn vừa làm quen với package rsample của Max Kuhn và Wickham. Như những package khác của họ, rsample không chỉ là một công cụ chuyên biệt mà là một nền tảng phổ quát có thể tái sử dụng nhiều lần và trở thành một bộ phận của ngữ pháp R, cho phép liên kết với những công cụ khác để tạo thành quy trình khép kín, đơn giản và phổ quát.
Việc tách rời danh sách phiên bản dữ liệu, hàm phân tích lõi bên trong và quy trình mô tả/đồ họa cho phép người dùng tùy chỉnh và sáng tạo tùy theo mục tiêu của mình. Cách làm việc được khuyến khích đó là sử dụng hàm map của package purrr, hàm tidy và augmente của package broom mà Nhi đã giới thiệu trước đây.
Ưu điểm của package rsample đó là nó đơn giản tối đa quá trình tạo dữ liệu tái chọn mẫu, vừa đảm bảo tính ngẫu nhiên, vừa cho phép tùy chỉnh và quan trọng nhất là lưu trữ các phiên bản dữ liệu này một cách tường minh, cho phép tái lập kết quả chứ không chạy ngầm như những package khác.