Đánh giá mô hình phân loại

1 Giới thiệu

Trong thời gian gần đây, Machine learning (ML) bắt đầu được ứng dụng phổ biến trong nghiên cứu y học lâm sàng. Các bác sĩ bắt đầu sử dụng những algorithm mới lạ, thí dụ SVM, Decision Tree, Random Forest,… ML cũng mang đến những quy trình mới, thí dụ cross validation, kiểm định độc lập và nhiều tiêu chí đánh giá, so sánh, lựa chọn mô hình mới. Tuy nhiên các bác sĩ vẫn gặp nhiều khó khăn khi tiếp thu dòng kiến thức cực lớn này. Một mặt, giới bác sĩ vẫn còn giữ quán tính rất lớn khi đi theo cách làm việc truyền thống và bị chi phối bởi cách suy nghĩ cũ theo trường phái thống kê y sinh. Thí dụ, khi nhắc đến việc dựng Mô hình, đa số bác sĩ thường chỉ nghĩ đến mục tiêu Diễn dịch (interpretive), họ chưa quen với mục tiêu tiên lượng. Trong các nghiên cứu y học, khái niệm mô hình gần như đồng nghĩa với Hồi quy tuyến tính, vì đây từng là công cụ phổ biến nhất mà môn thống kê hỗ trợ. Ngoài ra, việc kiểm tra đánh giá hiệu năng/phẩm chất của mô hình; so sánh mô hình, chọn lọc mô hình tối ưu… chưa được quan tâm đúng mức, và nếu có thực hiện cũng chỉ giới hạn trên dữ liệu hiện hành chứ không dùng dữ liệu độc lập. Bác sĩ cũng chưa quen với việc nhìn classifier đa biến như một test chẩn đoán để liên hệ chúng với những khái niệm đặc thù của ngành Y, thí dụ Likelihood ratios. Mặt khác vấn đề càng phức tạp hơn khi giữa ML và Thống kê có quá nhiều dị biệt về thuật ngữ (thậm chí có những điều cùng nguyên lý nhưng thuật ngữ cũng hoàn toàn khác nhau), ngăn cản bác sĩ học kiến thức từ Machine Learning và diễn đạt những tiêu chí đặc thù của ngành Y cho các data scientist bên ngoài.

Việc đánh giá mô hình chắc chắn cần thiết, ngay cả khi mục tiêu của nghiên cứu không phải là tiên lượng. Khi đánh giá mô hình, ta mới có thể kiểm soát, điều chỉnh những algorithm đạt hiệu quả tối ưu, so sánh giữa nhiều mô hình và chọn ra mô hình chính xác nhất, và xác nhận điều này trên một dữ liệu độc lập. Tuy nhiên câu hỏi cần giải đáp: Kiểm tra bằng những tiêu chí nào ?

Do đó, Nhi viết bài thực hành này với hy vọng giải quyết càng nhiều càng tốt các vấn đề nêu trên, với 3 mục tiêu chính:

1. Liệt kê tất cả những chỉ số / tiêu chí có thể dùng để đánh giá hiệu năng, phẩm chất của một mô hình nhị phân, thí dụ Random Forest.

2. Phân tích các dị biệt về thuật ngữ giữa ML và Thống kê y học, và chuẩn hóa chúng cho ứng dụng lâm sàng.

3. Tính toán các chỉ số này trong R, theo cách hoàn toàn thủ công hoặc sử dụng các hàm viết sẵn từ các thư viện (R packages).

2 Thí dụ minh họa: Random Forest model

Để minh họa, ta phải có một mô hình. Nhi sẽ dùng lại bộ số liệu Heart disease. Đây là 1 tập hợp dữ liệu của hơn 600 bệnh nhân từ 4 bệnh viện (Cleveland,Budapest,Long Beach and Zurich). Ta sẽ dựng một mô hình tiên lượng cho bệnh Tim mạch dựa vào 14 biến số bao gồm Tuổi, Giới tính, Triệu chứng Đau ngực, cholesterol, fasting blood sugar test, và stress test gồm nhịp tim và đoạn ST của ECG. Algorithm được chọn là Random Forest. Công đoạn training sẽ được thực hiện bằng caret. Mô hình RF sẽ được huấn luyện trên 50% dữ liệu, 50% còn lại dùng để kiểm định độc lập. Quy trình huấn luyện do caret thực hiện tự động bằng cách tinh chỉnh tham số mtry, nhằm tối ưu hóa tiêu chí Accuracy.

library(tidyverse)

va=read.table("https://archive.ics.uci.edu/ml/machine-learning-databases/heart-disease/processed.va.data", sep =",",na.strings="?",strip.white=TRUE, fill = TRUE)%>%as_tibble()
hu=read.table("https://archive.ics.uci.edu/ml/machine-learning-databases/heart-disease/processed.hungarian.data", sep =",",na.strings="?",strip.white=TRUE, fill = TRUE)%>%as_tibble()
sw=read.table("https://archive.ics.uci.edu/ml/machine-learning-databases/heart-disease/processed.switzerland.data", sep =",",na.strings="?",strip.white=TRUE, fill = TRUE)%>%as_tibble()
cl=read.table("https://archive.ics.uci.edu/ml/machine-learning-databases/heart-disease/processed.cleveland.data", sep =",",na.strings="?",strip.white=TRUE, fill = TRUE)%>%as_tibble()

df=rbind(va,hu,sw,cl)

names(df)=c("Age","Sex","ChestPain",
            "RestBP","Chol","FBS",
            "RestECG","MaxHR","CPETAgina",
            "Oldpeak","Slope","CA","Thal","Class")

data<-df[,-c(11,12,13)]%>%filter(.,Chol!=0)%>%na.omit()

data$Sex%<>%as.factor()%>%
  recode_factor(.,`0` = "Female", 
                  `1` = "Male")

data$ChestPain%<>%as.factor()%>%
  recode_factor(.,`1` = "Typical", 
                `2` = "Atypical",
                `3` = "Non_aginal", 
                `4` = "asymptomatic" )

data$FBS%<>%as.factor()%>%
  recode_factor(.,`0` = "No", 
                `1` = "Yes")

data$RestECG%<>%as.factor()%>%
  recode_factor(.,`0` = "Normal", 
                `1` = "Abnormal_ST",
                `2` = "LVHypertrophy")

data$CPETAgina%<>%as.factor()%>%
  recode_factor(.,`0` = "No",
                `1` = "Yes")

data$Class%<>%as.factor()%>%
  recode_factor(.,`0` = "Negative",
                `1` = "Positive",
                `2` = "Positive", 
                `3` = "Positive",
                `4` = "Positive")

rm(cl,df,hu,sw,va)

library(caret)

set.seed(1234)
idTrain=createDataPartition(y=data$Class, p=0.5,list=FALSE)

trainset=data[idTrain,]
testset=data[-idTrain,]

#KRF

Control=trainControl(method= "repeatedcv",
                     number=5,
                     repeats=5,
                     classProbs=TRUE,
                     summaryFunction=multiClassSummary)

rfmod=caret::train(Class~.,
            data=trainset,
            method = "rf",
            trControl=Control,
            tuneLength=5)

Đầu tiên, ta có thể thấy việc đánh giá mô hình đã được caret áp dụng ngay từ công đoạn huấn luyện, caret lần lượt thử mtry=2,4,7,10… , mỗi lần như vậy mô hình được dựng trên 4 blocks dữ liệu và kiểm định trên block còn lại, lặp lại 5 lượt. Cuối cùng, tiêu chí Accuracy được chọn để quyết định mtry tối ưu=2.

Trong kết quả crossvalidation, ta cũng có thể thấy những chỉ số khác như AUC, Kappa, F1, Sensitivity, Specificity, Balanced accuracy. Chúng ta sẽ nói về chúng sâu hơn trong phần sau, các bạn kiên nhẫn.

rfmod

## Random Forest 
## 
## 331 samples
##  10 predictor
##   2 classes: 'Negative', 'Positive' 
## 
## No pre-processing
## Resampling: Cross-Validated (5 fold, repeated 5 times) 
## Summary of sample sizes: 265, 264, 265, 265, 265, 265, ... 
## Resampling results across tuning parameters:
## 
##   mtry  logLoss    AUC        prAUC      Accuracy   Kappa      F1       
##    2    0.4318672  0.8853035  0.8539108  0.7784629  0.5546227  0.7944399
##    4    0.4466792  0.8752756  0.8445773  0.7753411  0.5488683  0.7895396
##    7    0.4597419  0.8651575  0.8335033  0.7747536  0.5479423  0.7877628
##   10    0.4725115  0.8603050  0.8257138  0.7638525  0.5260140  0.7774332
##   13    0.4744599  0.8598974  0.8149111  0.7596005  0.5169086  0.7758893
##   Sensitivity  Specificity  Pos_Pred_Value  Neg_Pred_Value  Precision
##   0.8175462    0.7357258    0.7773326       0.7892797       0.7773326
##   0.8072605    0.7406452    0.7768561       0.7819277       0.7768561
##   0.8003697    0.7469355    0.7793826       0.7767778       0.7793826
##   0.7887731    0.7367339    0.7710711       0.7640531       0.7710711
##   0.7980504    0.7175000    0.7588405       0.7687769       0.7588405
##   Recall     Detection_Rate  Balanced_Accuracy
##   0.8175462  0.4296036       0.7766360        
##   0.8072605  0.4241397       0.7739528        
##   0.8003697  0.4205393       0.7736526        
##   0.7887731  0.4144597       0.7627535        
##   0.7980504  0.4192994       0.7577752        
## 
## Accuracy was used to select the optimal model using the largest value.
## The final value used for the model was mtry = 2.

plot(rfmod)

3 Dữ liệu kết quả tiên lượng

Việc kiểm định mô hình sẽ được thực hiện độc lập trên dữ liệu testset (n=330), đầu tiên, chúng ta tạo ra một dataframe có nội dung như sau:

dfpred=data_frame(Truth=testset$Class,
                  BinTruth=if_else(testset$Class=="Positive",1,0),
                  ClassNeg=predict(rfmod,newdata=testset,type="prob")%>%.[[1]],
                  ClassPos=predict(rfmod,newdata=testset,type="prob")%>%.[[2]],
                  ClassLab=predict(rfmod,newdata=testset))

dfpred%>%head()

## # A tibble: 6 x 5
##   Truth    BinTruth ClassNeg ClassPos ClassLab
##   <fctr>      <dbl>    <dbl>    <dbl> <fctr>  
## 1 Negative     0      0.560     0.440 Negative
## 2 Positive     1.00   0.0480    0.952 Positive
## 3 Positive     1.00   0.100     0.900 Positive
## 4 Negative     0      0.446     0.554 Positive
## 5 Positive     1.00   0.592     0.408 Negative
## 6 Negative     0      0.854     0.146 Negative

1. Truth: Giá trị thực tế quan sát được,dưới dạng factor.

2. BinTruth: Giá trị thực tế dưới dạng xác suất tuyệt đối: Positive=1, Negative=0.

3. ClassNeg: Xác suất tiên lượng cho nhãn Negative.

4. ClassPos: Xác suất tiên lượng cho nhãn Positive.

5. ClassLab: Kết quả phân loại của mô hình.

Từ dữ liệu này, chúng ta có thể xác định được tất cả các tiêu chí về phẩm chất mô hình.

4 Confusion matrix

Đầu tiên, ta bàn về một hình thức tóm tắt kết quả về hiệu năng mô hình rất lý thú: Confusion matrix hay error matrix. Confusion matrix là thuật ngữ của giới Machine learning, tương đương với khái niệm contingency table hay cross-table bên giới thống kê. Nó có bản chất là một bảng phân phối tần số 2 chiều (bảng chéo) cho phép trình bày tỉ lệ tương hợp và bất xứng giữa Thực tế và kết quả phân loại của quy luật cần kiểm tra (mô hình).Trong bài này, ta chỉ xét trường hợp đơn giản nhất của confusion matrix áp dụng cho bài toán nhị phân (Binary classification).Khi đó, Confusion matrix trình bày tần suất của 4 tổ hợp: TP (True Positive), TN (True Positive), FP (False positive) và FN (False negetive). True (Đúng) chỉ sự tương hợp, khi kết quả phân loại phù hợp với giá trị thực tế; False (Sai) chỉ sự bất xứng hay nhầm lẫn, khi mô hình phân loại nhầm so với thực tế.

Ý nghĩa quan trọng của confusion matrix là ở chỗ TP,TN,FP,FN là 4 chỉ số cơ bản, từ đó cho phép suy ra kết quả của hầu hết những chỉ số khác.

Confusion matrix có thể tính tự động từ caret:

cf=caret::confusionMatrix(reference=dfpred$Truth,
                       data=dfpred$ClassLab,
                       positive="Positive",
                       mode="everything")

cf$table

##           Reference
## Prediction Negative Positive
##   Negative      139       32
##   Positive       34      125

Tuy nhiên, ta có thể tự dựng confusion matrix thủ công như sau:

confmat=table(Classification=dfpred$ClassLab,Truth=dfpred$Truth)

confmat

##               Truth
## Classification Negative Positive
##       Negative      139       32
##       Positive       34      125

Ta cũng có thể tính thủ công trực tiếp từng chỉ số một từ dữ liệu kiểm định:

TP=with(dfpred,sum(ClassLab==Truth & ClassLab=="Positive"))
TN=with(dfpred,sum(ClassLab==Truth & ClassLab!="Positive"))
FN=with(dfpred,sum(ClassLab!=Truth & Truth=="Positive"))
FP=with(dfpred,sum(ClassLab!=Truth & Truth!="Positive"))

cbind(TP,TN,FN,FP)

##       TP  TN FN FP
## [1,] 125 139 32 34

Lưu ý: Positive/negative chỉ là tên gọi của 2 nhãn giá trị, mang tính quy ước và có ý nghĩa tương đối tùy theo mục tiêu của người dùng (Dương tính = điều ta quan tâm tìm kiếm và Âm tính = loại trừ (tất cả mọi) thứ còn lại).

Trong y học lâm sàng còn dịch là “Thật” và “Giả”, nhưng Nhi cho rằng phát biểu “Mô hình phân loại sai, mô hình chẩn đoán lầm, loại trừ nhầm” thì dễ hiểu hơn so với “kết quả âm tính giả/dương tính giả”. Khi áp dụng cho mục tiêu chẩn đoán, TP có thể dịch là “chẩn đoán/phát hiện đúng/trúng”, TN là “loại trừ đúng”, FP là “chẩn đoán nhầm/sai”, FN là “bỏ sót, loại trừ nhầm”.

5 Tính hữu dụng:

Đầu tiên, ta sẽ kiểm tra xem mô hình hữu dụng đến mức nào, thông qua các chỉ số như sau:

1. Sensitivity - Độ nhạy

True Positive rate (TPR) : tỉ lệ phân loại Positive đúng trên tổng số các trường hợp Positive:

\[TPR = Sensitivity = Recall = \frac{TP}{TP+FN} = \frac{TP}{\sum Positive}\]

Tỉ lệ TPR còn có tên gọi khác là Sensitivity (độ nhạy), hit rate (tỉ lệ trúng đích), và Recall. Do trên lâm sàng ta quen với khái niệm Sensitivity, Nhi sẽ sử dụng nó như thuật ngữ duy nhất. Recall là thuật ngữ ít thông dụng hơn nhưng có nhiều tài liệu sử dụng nó nên các bạn cần ghi nhớ.

2. Specificity: Độ đặc hiệu

True negative rate (TNR) : tỉ lệ loại trừ đúng trên tổng số các trường hợp Negative hay còn gọi là Specificity (độ đặc hiệu) trong y học lâm sàng.

\[TNR = Specificity = \frac{TN}{TN+FP} = \frac{TP}{\sum Negative}\] Nhận xét:

Như đã nói ở trên, ý nghĩa của Positive/negative là tương đối, tương tự cho Sensitivity/Specificity. Nếu người bác sĩ đảo vị trí của mục tiêu phát hiện/loại trừ thì Sens/specs sẽ đổi vị trí cho nhau, nhưng chúng đều đo lường tính hữu dụng của mô hình (test chẩn đoán), hay khả năng cho phép đáp ứng mục tiêu phát hiện/loại trừ bệnh lý của người bác sĩ. Một mô hình có sensitivity cao cũng giống như một bác sĩ “can đảm” nói “Có” trước bệnh nhân, và specificity cao là sự can đảm nói “Không” khi loại trừ đối tượng không mắc bệnh.

Sensitivity và Specificity luôn phải đi cặp với nhau, vì chúng bổ sung cho nhau. Mỗi chỉ số này chỉ mới kiểm tra 1 hàng trong Confusion matrix (1 ứng dụng : phát hiện/loại trừ) và do đó, độc lập với ứng dụng còn lại. Một vài mô hình có thể rất tốt cho ứng dung phát hiện nhưng kém cho ứng dụng loại trừ và ngược lại. Mặt khác, một mô hình “cực đoan” (null model) - luôn luôn trả lời “Có” hoặc “Không” sẽ cho ra giá trị zero cho Sensitivity hoặc Specificity. Khi một tác giả nào đó chỉ trình bày Sensitivity hoặc Specificity nhưng cố tình lờ đi chỉ số còn lại, anh ta/chị ta đang muốn che giấu sự thật, và một nửa sự thật không phải là sự thật.

6 Tính chính xác

Khái niệm “chính xác” rất khó để diễn đạt, vì ở đây ta có đến 2 nhãn giá trị, do đó ta có thể bàn về tính chính xác chung chung, tổng quát với chỉ số Accuracy:

1. Accuracy : (ACC)

Có thể tạm dịch thuật ngữ Acuracy như “độ chính xác tổng quát”, vì nó đơn giản là tỉ lệ của tất cả trường hợp phân loại Đúng (không phân biệt negative/positive) trên toàn bộ trường hợp trong mẫu kiểm định.

\[Accuracy = \frac{TP+TN}{TP+TN+FP+FN}\]

Đây là tiêu chí phổ biến nhất (thường được nghĩ đến đầu tiên) khi kiểm định hiệu năng của mô hình phân loại, tuy nhiên giá trị thực dụng của nó thường kém vì nó không đặc hiệu cho một mục tiêu nào cả.

2. Balanced accuracy:

Khi 2 nhãn Positive/Negative bị mất cân đối, một tiêu chí khác phù hợp hơn là Balanced accuracy (BAC) hay độ chính xác sau cân bằng:

\[BAC = \frac{TPR+TNR}{2}\]

Cho mục tiêu chuyên biệt hơn, ta có thể sử dụng:

3. Precision (độ chính xác) hay Positive predictive value (PPV):

Là tỉ lệ thực sự positive trên tổng số các trường hợp được mô hình dán nhãn “Positive”. Precision là một thuật ngữ bên Data science và trong y học nó tương đương với khái niệm « PPV ». Nó đo lường tính “xác định”, hay khả năng phân loại Positive chính xác của mô hình. Nhi thích dùng Precision hơn PPV, vì thuật ngữ “value” tuy chuyên biệt cho Positive, nhưng không nói lên ý nghĩa “khả năng” (performance), còn Positive hay Negative chỉ là tên gọi quy ước.

\[Precision = PPV = \frac{TP}{TP+FP}\] 4) Tương tự: Negative predictive value (NPV) đo lường khả năng loại trừ chính xác:

\[NPV = \frac{TN}{TN+FN}\]

Lưu ý: Cẩn thận phân biệt giữa PPV/sensitivity và Recall/Precision. Precision đo lường tính “chuẩn xác”, Recall đo lường tính “hữu dụng”. Mẹo để nhớ là “Dám quyết định (Recall cao) và quyết định đúng (Precision cao)”. Cũng như cẩn thận khi dùng thuật ngữ “Độ chính xác” vì gây nhầm lẫn giữa Accuracy và Precision (Accuracy không phân biệt Positive/negative). Có thể dịch Precision là “Khả năng xác định”.

Precision không khảo sát mô hình môt cách độc lập, mà đặt mô hình vào một bối cảnh (dữ liệu). Do đó ta không nói chung chung:”mô hình chính xác”, mà là :” Mô hình chính xác đối với mẫu/dữ liệu hiện thời”.

5. F score

Do Recall (Sensitivity) và Precision là hai khái niệm khác nhau, ta có chỉ số F1 score cho phép đánh giá cân bằng giữa 2 phẩm chất này.

F1 score: Được định nghĩa như trung bình điều hòa (harmonic mean) giữa Precision và Recall (PPV và Sens) . \[F1 = 2\times \frac{1}{\frac{1}{Recall}+\frac{1}{Precision}}\] \[F1 = 2\times \frac{Precision \times Recall}{Precision + Recall}\] Hoặc:

\[F1 = \frac{2\times TP}{2\times TP+FP+FN}\]

F1 là trường hợp đặc biệt của công thức tổng quát Fbeta:

\[F_{\beta } = (1+\beta ^{2})\times \frac{Precision\times Recall}{(\beta ^{2}\times Precision)+Recall}\]

Hay:

\[F_{\beta } = \frac{(1+\beta ^{2})\times TP}{(1+\beta ^{2})\times TP + (\beta ^{2}\times FN) + FP}\]

Từ công thức này ta còn tính được chỉ số F2 (đặt trọng số cao hơn cho Recall so với Precision) và F0.5 (trọng số cao hơn cho Precision so với Recall).

Như vậy F1 được dùng khi ta quan tâm đồng đều vai trò của cả Precision và Recall, nói cách khác ta muốn Mô hình (quy luật chẩn đoán) vừa Nhạy, vừa chính xác. Việc lựa chọn giữa Recall (Sensitivity) và Precision (PPV) tùy thuộc vào mục tiêu ứng dụng của mô hình: người bác sĩ muốn Tầm soát bệnh hay muốn Xác định bệnh ? Nếu xem cả 2 đều quan trọng như nhau thì F1 là tiêu chí thích hợp nhất khi kiểm định mô hình phân loại vì khi F1 đạt tối ưu thì cả Precision và Recall đều phải tối ưu, ngược lại chỉ cần một trong 2 có giá trị thấp thì F1 sẽ thấp.

6. G index:

Cũng liên quan đến sự kết hợp Precision và Recall, ta có G measure hay còn gọi là Fowlkes–Mallows index (1983). G sử dụng trung bình nhân:

\[G = \sqrt{Recall \times Precision} = \sqrt{PPV \times TPR}\]

7 Sai sót và nhầm lẫn

Một cách tiếp cận khác để đánh giá mô hình, đó là ta quan tâm đến nguy cơ nhầm lẫn, sai sót. Bạn có thể nghĩ đơn giản về sự bù trừ giữa chính xác/nhầm lẫn, nhưng quan hệ này phức tạp hơn cho những mục tiêu chuyên biệt, và bài toán phân loại nhiều nhãn giá trị, khi đó tính chính xác tổng quát không đảm bảovề tỉ lệ nhầm lẫn tối thiểu.Trong y học,nguy cơ phân loại nhầm là tiêu chí quan trọng cần kiểm tra .

Ngay từ Confusion matrix, ta đã có tần suất nhầm lẫn tuyệt đối : FN và FP. Từ đó, ta có thể tính tỉ lệ sai sót

1. False negative rate: Tỉ lệ loại trừ nhầm (FNR) và

2. False positive rate = Tỉ lệ phát hiện nhầm (FPR).

\[FNR = \frac{FN}{\sum Positive} = \frac{FN}{FN+TP} = 1-TPR\]

\[FPR = \frac{FP}{\sum Negative} = \frac{FP}{FP+TN} = 1-TNR\]

Một mô hình tốt cần có FNR và FPR thấp. Lưu ý là tầm quan trọng của FNR và FPR tùy thuộc vào bản chất của vấn đề. Trong bài toán xếp loại tín dụng, việc phân loại nhầm một hồ sơ vay vốn từ tốt thành xấu sẽ gây hậu quả ít nghiêm trọng hơn trường hợp ngược lại. Trong Y học thì khác, việc loại trừ nhầm lẫn một căn bệnh đôi khi sẽ gây hậu quả nghiêm trọng vì bệnh nhân có thể tử vong khi không được điều trị kịp thời, trong khi việc chẩn đoán nhầm người bình thường thành bệnh nhân thì hậu quả có thể ít nghiêm trọng hơn (do vẫn còn cơ hội kiểm chứng lại với những xét nghiệm bổ sung).

Hai chỉ số khác hiếm gặp hơn là:

3. False discovery rate (FDR), tỉ lệ phát hiện nhầm, có ý nghĩa nghịch với Precision (PPV).

\[FDR = \frac{FP}{FP+TP} = 1-PPV = 1 – Precision\] Và

4. False omission rate (FOR): tỉ lệ loại trừ nhầm

\[FDR = \frac{FN}{FN+TN} = 1-NPV\]

Lưu ý: tạm thời chúng ta tách khái niệm FPR, FNR, FDR ra khỏi khái niệm Sai lầm type I và Type II trong kiểm định giả thuyết thống kê, dù chúng có liên hệ với nhau. Thuật ngữ Sai lầm type I và II không giúp ích cho việc kiểm định mô hình phân loại.

5. Tương tự BAC, ta cũng có thể tính Balanced Error rate (tỉ lệ sai lầm sau cân bằng):

\[BER = \frac{FPR+FNR}{2}\]

6. Brier score:

Điểm số Brier (BS) được giới thiệu bởi Glenn W. Brier vào năm 1950. BS được xác định như giá trị trung bình của các bình phương sai số giữa xác suất dự báo bởi mô hình và xác suất thực tế. Lưu ý: BS chuyên biệt cho từng nhãn giá trị, thí dụ nếu đặt Positive làm mục tiêu, thì xác suất thực tế của một trường hợp Positive sẽ = 1, và ngược lại trường hợp negative có xác suất thực tế =0. Nếu trên một cá thể Positive mà mô hình ước tính ra xác suất là 0.79 chẳng hạn, thì bình phương khác biệt cho case đó là (0.79-1)^2 = 0.0441.

Một mô hình tốt sẽ ước tính xác suất với sai biệt nhỏ nhất có thể, cho tất cả nhãn giá trị và tất cả trường hợp trong mẫu kiểm định, khi đó Brier score sẽ thấp. Giá trị tối ưu của BS là 0, và giá trị tồi tệ nhất là 1. Do đó Brier score được xem như 1 tiêu chí để đánh giá khả năng Sai lầm của mô hình (tương tự FNR và FPR).

Ta có thể tính BS hoàn toàn thủ công hoặc sử dụng hàm bierscore trong package scoring

7. Logarithmic loss

Một chỉ số khác là logloss, được sử dụng nhiều trong giới Machine learning, logloss được định nghĩa như trung bình của -log(xác suất dự báo).

Cho bài toán nhị phân, ta mặc định tính logloss theo xác suất cho nhãn Positive.

\[logloss = \frac{\sum -log(p_{i})}{N}\]

Logloss nhắm đến mục tiêu đo lường sai biệt của mô hình, giá trị tối thiểu (mô hình chính xác nhất) của logloss=0, nhưng không có ngưỡng cao nhất.

8. Tương tự, ta có chỉ số Logarithmic scoring rule (LSR)

\[LSR = mean(log(p_{i}))\] LSR trái dấu với logloss nhưng ý nghĩa tương đương, mô hình càng chính xác thì LSR càng gần 0. Mô hình càng kém chính xác thì LSR tiến về -Inf

9. Một chỉ số khác có thể tính là SSR : Spherical Scoring Rule,

SSR được xác định bằng trung bình của tích xác suất tiên lượng của một nhãn cho mỗi trường hợp i và tổng xác suất tiên lượng cho (các) nhãn j còn lại cho trường hợp i này :

\[SSR = Mean(P_{i}(\sum P_{ij}))\]

Giá trị SSR có thể được tính bởi hàm sphscore của package scoring.

10. Mean misclassification error (MMCE)

Cuối cùng: Nếu xem kết quả kiểm định mô hình là một biến logical: Y=1 nếu phân loại nhầm, Y=0 nếu phân loại đúng, ta sẽ thống kê được sai số trung bình: Mean misclassification error (MMCE) :

\[mmce = \frac{\sum Prediction \neq Truth}{N}\]

8 Mô hình như test chẩn đoán:

Một khái niệm khác, rất phổ biến trong thống kê y học nhưng không được quan tâm bên giới Machine learning, đó là Likelihood ratios (Positive và Negative), tức tỉ số khả dĩ Dương và Âm.

LR+ được xác định bằng độ nhạy (của test chẩn đoán) chia cho (1-độ đặc hiệu). LR- được tính bằng (1-độ nhạy) chia cho độ đặc hiệu:

\[LR+ = \frac{SEN}{1-SPEC}\] \[LR- = \frac{1-SEN}{SPEC}\]

Một cách tổng quát, LR đo lường mức độ liên hệ giữa quy luật chẩn đoán và trạng thái bệnh lý. Thí dụ, LR+ cao (>1) cho thấy kết quả xét nghiệm dương tính có liên quan đến sự hiện diện của bệnh lý. Cạnh đó, LR còn là phương tiện trung gian cho phép biện luận, so sánh hiệu năng của một xét nghiệm chẩn đoán này với một xét nghiệm khác. Việc so sánh này dựa vào Pretest probability và Posttest probability.

Posttest probability của một xét nghiệm có thể được xác định theo 3 bước:

\[ pretest Odds = Pretest Prob / (1-Pretest Prob)\]

\[ Posttest Odds = Pretest odds * LR \]

\[Post test Prob = Posttest odds/(1+Post test odds)\]

Như vậy, một xét nghiệm mới nếu muốn được công nhận là tốt hơn so với một xét nghiệm cũ, xét nghiệm mới này phải có Posttest Prob cao hơn so với Pretest Prob, tức là LR phải > 1 và càng cao càng tốt. Người ta ước tính rằng để Post test probability cho chẩn đoán Có bệnh tăng hơn 15% thì LR+ phải lớn hơn 2, tương tự để post test prob cho phép loại trừ bệnh thì LR- phải nhỏ hơn 0,2. LR càng gần 1 thì hiệu năng của xét nghiệm càng thấp (so với Pretest prob). Post test Prob của xét nghiệm mới này lại có thể được dùng như Pretest Prob cho một bước tiếp theo trong quy trình chẩn đoán.

Trong thí dụ này, giả sử ta có prevalence của bệnh lý Tim mạch được ước tính từ trainset là 47.6% và dùng nó như PretestProb, từ LR+ = 4.05 ta có thể ước tính được Posttest Prob của quy luật chẩn đoán mới (mô hình Random Forest) là 0.79, tức là với quy luật mới ta đã làm thay đổi Pretest Prob lên tới +65.5 %, cho thấy mô hình mới là một công cụ thực sự có ích so với suy đoán ngẫu nhiên..

9 Tính tương hợp Tiên lượng/Thực tế:

1. Trị số Kappa của Cohen:

Kappa là một trị số thống kê do Jacob Cohen đề xuất năm 1960 nhằm đo lường mức độ tương hợp giữa 2 hay nhiều kết quả đo lường (inter-rater agreement). Kappa được định nghĩa như:

\[\kappa = 1- \frac{1-p_{o}}{1-p_{e}}\]

Trong đó po là tỉ lệ tương hợp quan sát được trên thực tế, và pe là tỉ lệ tương hợp ước đoán theo phân bố ngẫu nhiên. Nếu sự tương hợp là tuyệt đối, Kappa=1, Nếu sự tương hợp kém (so với phân bố ngẫu nhiên) thì kappa giảm về 0.

Cho mục tiêu kiểm định phẩm chất mô hình, ta có thể áp dụng hàm cohen.kappa của package psych để tính giá trị Kappa và Weighted Kappa cho cả kết quả là Label lẫn Xác suất dự báo.

2. Hệ số tương quan Matthews (MCC)

MCC là một phương pháp riêng của giới Machine learning, dùng để đánh giá phẩm chất của mô hình phân loại nhị phân, được giới thiệu bởi Brian W. Matthews vào năm 1975. Mục tiêu của MCC là khắc phục vấn đề dữ liệu bị mất cân bằng. MCC có bản chất là một hệ số tương quan giữa giá trị Thực tế và Kết quả dự báo của mô hình, nó được xác định như sau:

\[MCC = \frac{(TP*TN)-(FP*FN)}{\sqrt{(TP+FP)*(TP+FN)*(TN+FP)*(TN+FN)}}\]

MCC có liên hệ trực tiếp với trị số thống kê Chi bình phương của bảng 2x2 :

\[\left | MCC \right | = \sqrt{\frac{\chi ^2}{N}}\]

Giá trị của MCC dao động trong khoảng -1 đến +1, MCC=+1 biểu thị cho “kết quả phân loại hoàn hảo”, MCC=0 cho thấy mô hình vô dụng (không hơn gì sự phán đoán ngẫu nhiên) còn MCC=-1 cho thấy mô hình không những vô dụng hoàn toàn mà còn tuyệt đối sai, vì kết quả phân loại hoàn toàn trái nghịch với quan sát thực tế.

MCC có thể được phân tích thành 2 nhân tố:

3. Chỉ số Markedness (MK) được định nghĩa như:

\[MK = PPV + NPV – 1\]

4. Chỉ số Informedness hay còn gọi là Bookmaker Informedness (BM):

\[BM = TPR + TNR – 1\]

Nếu tinh ý ta sẽ thấy BM chính là trị số Youden’s J trong phân tích ROC, vì:

\[Youden’s J = Sensitivity + Specificity – 1\]

MCC chính là trung bình nhân của Informedness (BM) và Markedness (MK)

\[\left | MCC \right | = \sqrt{BM *MK }\] # Phân tích ROC và AUC

Phương pháp cuối cùng mà chúng ta sẽ làm đó là phân tích ROC, cụ thể là diện tích dưới ROC (AUC). Khi nhắc đến đường cong ROC (Receiver operating characteristic), các bạn bác sĩ dễ dàng liên hệ tới phương pháp để đánh giá hiệu quả của một test chẩn đoán đơn biến. Mục tiêu là thiết lập một biểu đồ thể hiện sự dao động của TPR tùy theo FPR khi dịch chuyển qua các ngưỡng cắt khác nhau trên thang điểm chẩn đoán (Score).

Một cách tổng quát, phân tích ROC có thể áp dụng cho mọi mô hình nhị phân với Score chính là kết quả xác suất mà mô hình này tiên lượng.

Theo mặc định, ngưỡng cắt của xác suất phân loại là 0.5 nhưng trong một số hoàn cảnh, ngưỡng này có thể được điều chỉnh (gọi là cân bằng, calibration). Thí dụ nếu ta muốn tạo ra một quy luật chẩn đoán nhạy hơn so với mặc định, ta có thể giảm ngưỡng cắt còn 0.4 hay 0.3 cho nhãn Positive, như vậy sẽ làm tăng đồng thời FPR và TPR, và giảm đồng thời FNR, TNR theo luật bù trừ, làm giảm thiểu tỉ lệ bỏ sót chẩn đoán. Đường cong ROC cho phép khảo sát tất cả khả năng về ngưỡng cắt khác nhau và ảnh hưởng lên TPR, FPR. ROC cho phép đánh giá hiệu năng tổng quát của mô hình. Mô hình tối ưu sẽ có TPR cao và FPR thấp (gần góc trên bên trái ROC) và AUC gần 1 nhất có thể.

10 Cách làm gian lận

Đầu tiên, Nhi sẽ chia sẻ với các bạn cách tính nhanh nhất các tiêu chí đánh giá mô hình mà không cần bận tâm đến công thức của chúng. Nhi gọi đây là cách làm “đểu”, “cheating” hay “Múa rìu qua mắt thợ”, nhờ vào package mlr. Trong mlr có tích hợp hầu hết các chỉ số đánh giá mô hình, ngay cả những chỉ số hiếm gặp nhất. Nhưng vấn đề đó là các hàm này chỉ áp dụng được cho model do chính mlr huấn luyện và chỉ nhận object prediction của chính mlr. Do đó, ta phảican thiệp để đánh tráo nội dung prediction của một model bất kì, bằng kết quả prediction của mô hình mà ta đang muốn kiểm tra.

Chỉ với 5 bước, ta có kết quả như sau:

library(mlr)

# Step 1:
task=makeClassifTask(data=testset,target="Class")

# Step 2
dummylrner=makeLearner("classif.rpart",predict.type ="prob")
dummymod=train(dummylrner,task)
dummypred=predict(dummymod,task)

# Step 3

dummypred$data$truth<-dfpred$Truth
dummypred$data$prob.Negative<-dfpred$ClassNeg
dummypred$data$prob.Positive<-dfpred$ClassPos
dummypred$data$response<-dfpred$ClassLab

# Step 4

measure=list(mlr::acc,
             mlr::auc,
             mlr::bac,
             mlr::ber,
             mlr::brier,
             mlr::f1,
             mlr::gmean,
             mlr::gpr,
             mlr::kappa,
             mlr::logloss,
             mlr::lsr,
             mlr::mcc,
             mlr::mmce,
             mlr::qsr,
             mlr::ssr)

# Step 5

performance(dummypred,measure)%>%as.matrix()

##               [,1]
## acc      0.8000000
## auc      0.8711204
## bac      0.7998233
## ber      0.2001767
## brier    0.1461624
## f1       0.8081395
## gmean    0.7998150
## gpr      0.8081532
## kappa    0.5992935
## logloss  0.4498465
## lsr     -0.4498465
## mcc      0.5993377
## mmce     0.2000000
## qsr      0.7076752
## ssr      0.8365076

Oh yes ! Rất đơn giản phải không nè ?

Bây giờ ta sẽ làm “thật”, tức là tính thủ công lần lượt từng chỉ số :

11 Cách làm thật

11.1 Tính hữu dụng

TPR=TP/(TP+FN)
TNR=TN/(TN+FP)
SEN=TPR
SPEC=TNR
Recall=SEN

rbind(TPR,TNR,SEN,SPEC,Recall)

##             [,1]
## TPR    0.7961783
## TNR    0.8034682
## SEN    0.7961783
## SPEC   0.8034682
## Recall 0.7961783

11.2 Tính chính xác

ACC=(TP+TN)/(TP+FP+TN+FN)
BAC=(TPR+TNR)/2
PPV=TP/(TP+FP)
NPV=TN/(TN+FN)
Precision=PPV
G=sqrt(Recall*Precision)
F2=(1+2^2)*(Precision*Recall)/((2^2)*Precision+Recall)
F0.5=(1+(0.5)^2)*(Precision*Recall)/((0.5^2)*Precision+Recall)
F1=2*TP/(TP+TP+FP+FN)
Fscore=(1+1)*(Precision*Recall)/(1*Precision+Recall)

rbind(ACC,BAC,PPV,NPV,Precision,G,F1,F2,F0.5,Fscore)

##                [,1]
## ACC       0.8000000
## BAC       0.7998233
## PPV       0.7861635
## NPV       0.8128655
## Precision 0.7861635
## G         0.7911551
## F1        0.7911392
## F2        0.7941550
## F0.5      0.7881463
## Fscore    0.7911392

11.3 Nguy cơ sai lầm

FNR=FN/(FN+TP)
FPR=FP/(FP+TN)
FDR=FP/(TP+FP)
FOR=1-NPV
BER=with(dfpred,mean((ClassPos>0.5 & BinTruth==0) | (ClassPos<0.5 & BinTruth==1)))
logloss=with(dfpred,-mean(log10(ClassPos)))
LSR=with(dfpred,mean(log10(ClassPos)))
mmce=with(dfpred,mean(Truth!=ClassLab))

library(scoring)
brier= dfpred%>%scoring::brierscore(data=.,BinTruth~ClassPos)%>%mean()

SSR=dfpred%>%
  scoring::sphscore(data=.,BinTruth~ClassPos)%>%mean()

rbind(FNR,FPR,FDR,FOR,BER,logloss,LSR,mmce,brier,SSR)

##               [,1]
## FNR      0.2038217
## FPR      0.1965318
## FDR      0.2138365
## FOR      0.1871345
## BER      0.2000000
## logloss  0.4643810
## LSR     -0.4643810
## mmce     0.2000000
## brier    0.1461624
## SSR      0.1634924

11.4 Tính tương hợp

BM=TPR+TNR-1
MK=PPV+NPV-1
MCC=(TP*TN-FP*FN)/(sqrt((TP+FP)*(TP+FN)*(TN+FP)*(TN+FN)))
Kappa=psych::cohen.kappa(as.matrix(dfpred[,c(1,5)]))%>%.$weighted.kappa

rbind(BM,MK,MCC,Kappa)

##            [,1]
## BM    0.5996466
## MK    0.5990290
## MCC   0.5993377
## Kappa 0.5992935

11.5 Likelihood ratios

LRpos=SEN/(1-SPEC)
LRneg=(1-SEN)/SPEC

rbind(LRpos,LRneg)

##            [,1]
## LRpos 4.0511428
## LRneg 0.2536773

11.6 ROC-AUC

library(pROC)
myroc <- roc(Truth~ClassPos, data=dfpred)

AUC=myroc$auc

data_frame(fpr=1-myroc$specificities,
           tpr=myroc$sensitivities)%>%
  ggplot(aes(x=fpr,ymin=0,ymax=tpr))+
  geom_polygon(aes(y=tpr),fill="red",alpha=0.3)+
  geom_path(aes(y=tpr),col="red3",size=1)+
  geom_abline(linetype='dashed')+
  theme_bw()+
  coord_equal()+
  labs(x="False positive rate",y="True Positive rate")+
  ggtitle(paste0("AUC=",round(AUC,3)))

12 Tổng kết

Bài thực hành đến đây là hết. Chắc chắn nội dung của bài có nhiều thiếu sót, nhưng Nhi đã cố gắng hệ thống tất cả những gì Nhi học đươc từ giới Machine learning và diễn đạt nó theo hướng ứng dụng vào y học lâm sàng. Như các bạn đã thấy, việc tính thủ công những giá trị này hoàn toàn khả thi.Hy vọng bài viết giúp ích cho các bạn bác sĩ đang làm quen với ML. Chúc các bạn thành công.