Bệnh tiểu đường, hay còn gọi là đái tháo đường, đang trở thành một trong những thách thức sức khỏe cộng đồng lớn nhất của thế kỷ 21. Theo Liên đoàn Đái tháo đường Quốc tế (IDF), số người mắc bệnh tiểu đường trên toàn cầu đã vượt qua con số 537 triệu vào năm 2021 và dự kiến sẽ tăng lên 783 triệu vào năm 2045. Sự gia tăng này không chỉ gây ra gánh nặng lớn về kinh tế, xã hội và y tế mà còn làm suy giảm nghiêm trọng chất lượng cuộc sống của người bệnh. Bệnh tiểu đường là nguyên nhân hàng đầu dẫn đến các biến chứng nghiêm trọng như bệnh tim mạch, suy thận, mù lòa, tổn thương thần kinh, và cắt cụt chi.
Tại Việt Nam, tỷ lệ mắc bệnh tiểu đường đang tăng nhanh, đặc biệt ở các khu vực đô thị hóa, nơi lối sống hiện đại và chế độ ăn uống không lành mạnh ngày càng phổ biến. Theo thống kê của Bộ Y tế Việt Nam, khoảng 7% dân số trưởng thành mắc bệnh tiểu đường, và đáng lo ngại hơn, gần 60% trong số đó không được chẩn đoán cho đến khi xuất hiện biến chứng. Điều này phản ánh sự thiếu hụt trong nhận thức cộng đồng, các chương trình tầm soát, và các chiến lược phòng ngừa hiệu quả.
Việc xác định sớm và chính xác các yếu tố nguy cơ liên quan đến bệnh tiểu đường là một bước quan trọng để giảm thiểu tác động của bệnh. Các yếu tố như di truyền, lối sống, chỉ số khối cơ thể (BMI), huyết áp, và các yếu tố môi trường đóng vai trò then chốt trong sự phát triển của bệnh. Bằng cách phân tích các yếu tố này, chúng ta có thể xây dựng các mô hình dự báo rủi ro, từ đó đề xuất các chiến lược can thiệp sớm như thay đổi lối sống, cải thiện chế độ ăn uống, và tăng cường hoạt động thể chất. Những can thiệp này không chỉ giúp giảm nguy cơ mắc bệnh mà còn hỗ trợ phát hiện bệnh ở giai đoạn đầu, khi việc điều trị còn hiệu quả và ít tốn kém hơn. Với những lý do cấp thiết nêu trên, đề tài “Phân tích một số yếu tố ảnh hưởng đến rủi ro mắc bệnh tiểu đường” được lựa chọn để nghiên cứu, nhằm cung cấp cơ sở khoa học cho các chiến lược phòng ngừa và kiểm soát bệnh.
Bệnh tiểu đường (Diabetes Mellitus) là một bệnh rối loạn chuyển hóa mạn tính, đặc trưng bởi tình trạng tăng đường huyết (mức đường trong máu cao hơn bình thường). Nguyên nhân là do cơ thể thiếu hụt insulin, kháng insulin, hoặc cả hai. Insulin là một hormone do tuyến tụy sản xuất, có vai trò giúp glucose từ máu đi vào tế bào để chuyển hóa thành năng lượng.
Có ba loại tiểu đường chính:
Đường huyết cao kéo dài có thể gây tổn thương nghiêm trọng cho nhiều cơ quan trong cơ thể, dẫn đến các biến chứng:
Cơ sở lý thuyết: HbA1c (Hemoglobin A1c) phản ánh mức đường huyết trung bình trong khoảng 2–3 tháng gần nhất, do đó nó được xem là chỉ số “vàng” trong chẩn đoán và theo dõi hiệu quả điều trị bệnh tiểu đường. Theo hướng dẫn của Hiệp hội Tiểu đường Hoa Kỳ (ADA, 2020), mức HbA1c \(\geq\) 6.5% là một trong những tiêu chuẩn để chẩn đoán bệnh tiểu đường, trong khi mức từ 5.7% đến 6.4% cho thấy tình trạng tiền tiểu đường. Do đó, giá trị HbA1c càng cao, mức độ kiểm soát đường huyết càng kém, dẫn đến nguy cơ mắc bệnh và các biến chứng liên quan càng lớn.
Nghiên cứu liên quan: Một nghiên cứu quan trọng của Selvin và cộng sự (2010) được công bố trên New England Journal of Medicine đã khẳng định HbA1c là một yếu tố dự báo mạnh mẽ không chỉ cho bệnh tiểu đường tuýp 2 mà còn cho các biến chứng tim mạch. Nghiên cứu này nhấn mạnh độ nhạy và độ đặc hiệu cao của chỉ số này trong chẩn đoán lâm sàng.
Cơ sở lý thuyết: Đường huyết lúc đói (Fasting Blood Sugar, FBS) là chỉ số đo nồng độ glucose trong máu sau ít nhất 8 giờ nhịn ăn. Đây là một xét nghiệm cơ bản và là tiêu chuẩn chẩn đoán quan trọng (mức FBS \(\geq\) 126 mg/dL theo ADA). Mức FBS tăng cao cho thấy cơ thể đang gặp vấn đề trong việc điều hòa glucose, một dấu hiệu đặc trưng của tình trạng kháng insulin hoặc suy giảm chức năng tế bào beta của tuyến tụy, thường thấy ở bệnh nhân tiền tiểu đường và tiểu đường tuýp 2.
Nghiên cứu liên quan: Nghiên cứu của Tirosh và cộng sự (2005) trên New England Journal of Medicine đã cung cấp bằng chứng cho thấy ngay cả khi mức FBS tăng dần trong phạm vi được coi là “bình thường cao” (90–125 mg/dL), nó vẫn liên quan chặt chẽ đến việc tăng nguy cơ phát triển tiểu đường tuýp 2 trong tương lai, đặc biệt ở nam giới.
Cơ sở lý thuyết: Huyết áp tâm thu (Systolic Blood Pressure) cao, đặc biệt là tình trạng tăng huyết áp (thường được định nghĩa là \(\geq\) 140 mmHg), có mối liên hệ sinh lý bệnh chặt chẽ với bệnh tiểu đường. Cả hai bệnh lý này thường cùng tồn tại và chia sẻ các cơ chế chung như kháng insulin, rối loạn chức năng nội mô và tình trạng viêm mãn tính. Tăng huyết áp không chỉ là yếu tố nguy cơ của bệnh tiểu đường mà còn làm gia tăng đáng kể nguy cơ biến chứng tim mạch ở những người đã mắc bệnh.
Nghiên cứu liên quan: Một bài tổng quan của Ferrannini & Cushman (2012) trên tạp chí Diabetes Care đã chỉ ra rằng tăng huyết áp tâm thu là một yếu tố nguy cơ độc lập đối với bệnh tiểu đường tuýp 2. Cụ thể, nghiên cứu ước tính rằng mỗi 10 mmHg tăng lên ở huyết áp tâm thu có thể làm tăng nguy cơ mắc bệnh tiểu đường khoảng 10–15%.
Cơ sở lý thuyết: Hút thuốc lá được công nhận là một yếu tố nguy cơ có thể thay đổi đối với nhiều bệnh mạn tính, bao gồm cả tiểu đường tuýp 2. Các chất độc trong khói thuốc có thể gây ra viêm hệ thống, tăng stress oxy hóa và trực tiếp làm suy giảm độ nhạy của tế bào với insulin. Do đó, người hút thuốc có nguy cơ mắc bệnh tiểu đường cao hơn đáng kể so với người không hút thuốc.
Nghiên cứu liên quan: Một phân tích tổng hợp quy mô lớn của Willi và cộng sự (2007) đăng trên JAMA đã kết luận rằng việc hút thuốc làm tăng nguy cơ mắc bệnh tiểu đường tuýp 2 lên đến 44% so với người không bao giờ hút thuốc. Nghiên cứu này cũng chỉ ra một mối quan hệ liều-phản ứng, nghĩa là nguy cơ càng tăng khi số lượng thuốc hút càng nhiều.
Cơ sở lý thuyết: Đi tiểu thường xuyên (hay đa niệu) là một trong những triệu chứng lâm sàng kinh điển của bệnh tiểu đường. Cơ chế của hiện tượng này là do nồng độ glucose trong máu tăng cao vượt ngưỡng tái hấp thu của thận. Lượng đường dư thừa sẽ được thải ra ngoài qua nước tiểu, kéo theo nước và gây ra hiện tượng bài niệu thẩm thấu, khiến người bệnh phải đi tiểu nhiều lần với lượng nước tiểu lớn. Triệu chứng này thường là dấu hiệu cho thấy bệnh đã tiến triển.
Nghiên cứu liên quan: Nghiên cứu của Gupta và cộng sự (2018) trên Journal of Diabetes Research đã xác nhận rằng đi tiểu nhiều là một triệu chứng có giá trị trong việc phát hiện các trường hợp tiểu đường chưa được chẩn đoán, đặc biệt khi nó xuất hiện cùng với các triệu chứng khác như khát nước nhiều (polydipsia) và mệt mỏi, giúp sàng lọc sớm các đối tượng nguy cơ cao.
Trong nghiên cứu này, biến phụ thuộc là tình trạng mắc bệnh tiểu đường, là một biến nhị phân (có/không, 1/0). Do đó, các mô hình hồi quy cho biến nhị phân là lựa chọn phù hợp.
Dựa trên những nghiên cứu trước về các yếu tố ảnh hưởng đến bệnh tiểu đường, em đề xuất mô hình nghiên cứu như sau:
library(DiagrammeR)
grViz("
digraph bankruptcy_model {
graph [layout = dot, rankdir = LR]
node [shape = box, style = filled, fillcolor = darkslateblue, fontcolor = white]
X1 [label = 'HbA1c']
X2 [label = 'FastingBloodSugar']
X3 [label = 'SystolicBP']
X4 [label = 'Smoking']
X5 [label = 'FrequentUrination']
Y [label = 'Diagnosis', shape = box, fillcolor = darkslateblue]
X1 -> Y
X2 -> Y
X3 -> Y
X4 -> Y
X5 -> Y
}
")Đây là mô hình được sử dụng phổ biến nhất để phân tích mối quan hệ giữa một biến phụ thuộc nhị phân và một hoặc nhiều biến độc lập (định lượng hoặc định tính).
Mô hình Logit không dự đoán trực tiếp giá trị 0 hay 1, mà dự đoán xác suất để biến phụ thuộc nhận giá trị 1, ký hiệu là \(P_i = P(Y_i=1)\).
Công thức của mô hình có dạng: \[P_i = E(Y=1|X_i) = \frac{1}{1 + e^{-(\beta_0 + \beta_1 X_{1i} + ... + \beta_k X_{ki})}} = \frac{e^{Z_i}}{1 + e^{Z_i}}\] Trong đó: * \(P_i\) là xác suất đối tượng \(i\) mắc bệnh tiểu đường. * \(X_{1i}, ..., X_{ki}\) là giá trị của các biến độc lập (yếu tố nguy cơ) của đối tượng \(i\). * \(\beta_0, \beta_1, ..., \beta_k\) là các tham số hồi quy cần ước lượng. * \(Z_i = \beta_0 + \beta_1 X_{1i} + ... + \beta_k X_{ki}\). * \(e\) là cơ số của logarit tự nhiên.
Để tuyến tính hóa mô hình, người ta sử dụng phép biến đổi logit: \[L_i = \ln\left(\frac{P_i}{1-P_i}\right) = \beta_0 + \beta_1 X_{1i} + ... + \beta_k X_{ki}\] Đại lượng \(\frac{P_i}{1-P_i}\) được gọi là tỷ số chênh (odds), cho biết khả năng xảy ra sự kiện (mắc bệnh) lớn gấp bao nhiêu lần khả năng không xảy ra sự kiện. \(L_i\) là logarit của tỷ số chênh.
Hệ số \(\beta_k\) được diễn giải thông qua tỷ số chênh lệch (Odds Ratio - OR), tính bằng \(e^{\beta_k}\). Nó cho biết khi biến độc lập \(X_k\) tăng một đơn vị (và các biến khác không đổi), tỷ số chênh (odds) mắc bệnh sẽ thay đổi (tăng hoặc giảm) một lượng bằng \(e^{\beta_k}\) lần.
Mô hình Probit là một lựa chọn thay thế cho mô hình Logit. Thay vì sử dụng hàm phân phối logistic, mô hình Probit sử dụng hàm phân phối tích lũy chuẩn (cumulative standard normal distribution function), ký hiệu là \(\Phi(Z)\).
\[P_i = P(Y=1|X_i) = \Phi(Z_i) = \int_{-\infty}^{Z_i} \frac{1}{\sqrt{2\pi}} e^{-t^2/2} dt\] Trong đó \(Z_i = \beta_0 + \beta_1 X_{1i} + ... + \beta_k X_{ki}\).
Mô hình Cloglog (Complementary log-log) là một dạng mở rộng của mô hình hồi quy nhị phân, tương tự như mô hình Logit và Probit, nhưng sử dụng một hàm liên kết khác. Mô hình này thường được sử dụng trong các trường hợp mà xác suất xảy ra sự kiện cần được mô tả bất đối xứng, đặc biệt là khi sự kiện xảy ra có xác suất rất nhỏ hoặc rất lớn.
Mô hình cloglog sử dụng hàm liên kết:
\[ g(p) = \log(-\log(1 - p)) \]
Trong đó: - \(p = \mathbb{P}(Y = 1 | X)\) là xác suất của biến phản hồi nhị phân bằng 1. - Hàm liên kết cloglog là bất đối xứng, khác với logit và probit vốn có tính đối xứng.
Khi đó, hàm nghịch đảo của cloglog là:
\[ p = 1 - e^{-e^{\eta}}, \quad \text{với } \eta = X\beta \]
Điều này dẫn đến công thức xác suất trong mô hình cloglog:
\[ \mathbb{P}(Y = 1 | X) = 1 - \exp\left(-\exp(X\beta)\right) \]
Do các mô hình Logit và Probit là phi tuyến tính, phương pháp Bình phương Tối thiểu Thông thường (OLS) không còn phù hợp. Thay vào đó, các tham số \(\beta\) được ước lượng bằng Phương pháp Hợp lý Cực đại (MLE).
Nguyên tắc của MLE là tìm ra một bộ giá trị của các tham số \((\hat{\beta}_0, \hat{\beta}_1, ..., \hat{\beta}_k)\) sao cho xác suất xảy ra của mẫu dữ liệu quan sát được là lớn nhất.
Hàm hợp lý (Likelihood Function) cho một mẫu gồm \(n\) quan sát độc lập được xây dựng như sau: \[L(\beta_0, ..., \beta_k) = \prod_{i=1}^{n} P_i^{Y_i} (1-P_i)^{1-Y_i}\] Trong đó: * \(Y_i = 1\) nếu đối tượng \(i\) mắc bệnh, và \(Y_i = 0\) nếu không. * \(P_i\) là xác suất đối tượng \(i\) mắc bệnh, được tính từ mô hình Logit hoặc Probit.
Vì hàm \(L\) là một tích phức tạp, người ta thường làm việc với logarit của nó, gọi là hàm Log-Likelihood: \[\ln L = \sum_{i=1}^{n} [Y_i \ln(P_i) + (1-Y_i) \ln(1-P_i)]\] Phương pháp MLE sẽ tìm các giá trị \(\hat{\beta}\) để tối đa hóa hàm \(\ln L\) này. Quá trình này thường được thực hiện bằng các thuật toán tối ưu hóa trên máy tính.
Phương pháp Random Forest là một thuật toán học máy (machine learning) được Breiman giới thiệu vào năm 2001. Cụ thể, trong bối cảnh hồi quy, Random Forest sẽ là nơi một biến phản hồi \(Y \in \mathbb{R}\) được giải thích bởi một tập hợp các biến đồng biến \(X = (X^{(1)}, \dots, X^{(k)})\). Một tập hợp gồm \(N\) tập dữ liệu, bao gồm biến phản hồi và các biến đồng biến liên quan có thể được sử dụng để huấn luyện một Random Forest.
Một RF có kích thước \(B\) được tạo thành từ \(B\) cây hồi quy. Một cây là một cấu trúc được tạo từ các nút nhị phân, được xây dựng dần dần từ trên xuống dưới cho đến khi đạt điều kiện dừng. Có hai loại nút trong cây: nút bên trong (internal) và nút cuối (terminal), nút cuối còn được gọi là lá (leaves). Tại một nút bên trong, một quy tắc nhị phân dưới dạng \(X^{(i)} \leq s\) hoặc \(X^{(i)} > s\) so sánh một biến đồng biến \(X^{(i)}\) với một ngưỡng \(s\). Kết quả của bài kiểm tra này chia không gian dự đoán và tập dữ liệu huấn luyện thành hai phần mới, mỗi phần nằm trong một nhánh khác nhau.
Khi xây dựng cây dựa trên một tập dữ liệu huấn luyện, chỉ số của biến đồng biến \(i\) và giá trị ngưỡng \(s\) được xác định sao cho:
\[ L^2 = \sum_{i=1}^{n}(y_i - \hat{y}_i)^2 \rightarrow \min \]
có nghĩa rằng với bất kỳ cặp \((i, s)\) nào làm cho giá trị dự báo gần với giá trị thực tế nhất sẽ được chọn. Một biến đồng biến có thể được sử dụng nhiều lần cho các lựa chọn khác nhau của \(i\) tại các cấp khác nhau trong quá trình xây dựng cây. Quá trình chia tiếp tục cho đến khi tất cả các quan sát trong một nút cụ thể của cây có cùng giá trị biến đồng biến, khi đó nút sẽ trở thành lá. Ngoài ra, khi nút có ít hơn \(N_{\text{min}}\) (thường \(N_{\text{min}} = 5\) trong bối cảnh hồi quy) quan sát, nút đó cũng trở thành lá.
Khi cây được xây dựng xong, một giá trị của biến phản hồi được gán cho mỗi lá, tương ứng với giá trị trung bình của các giá trị phản hồi của các tập dữ liệu có mặt tại lá đó. Với một tập dữ liệu đã quan sát và một tập biến đồng biến mới \(X\), việc dự đoán giá trị phản hồi \(Y\) chính là việc đi qua các nút của cây dựa trên giá trị của biến đồng biến \(X\) mới. Kết quả của quá trình dự đoán là giá trị được gán cho lá mà tập dữ liệu đó dừng lại sau khi đi theo con đường này.
Điểm chính tạo nên sự “Random” trong phương pháp này đó chính là việc tổng hợp (hoặc lấy mẫu) các cây hồi quy ngẫu nhiên. Một số lượng lớn các cây được huấn luyện trên các tập mẫu bootstrap (có kích thước bằng tập dữ liệu gốc, được tạo bằng cách lấy mẫu ngẫu nhiên có hoàn lại từ tập dữ liệu huấn luyện ban đầu), điều này có nghĩa là mỗi cây hồi quy sẽ được huấn luyện trên một tập bootstrap khác nhau.
Thêm vào đó, ở mỗi nút, một tập con \(m_{\text{try}}\) các đồng biến trong tổng \(k\) biến đồng biến có sẵn sẽ được chọn một cách hoàn toàn ngẫu nhiên, từ đó sẽ lựa chọn cặp \((i, s)\) tối ưu dựa trên tập \(m_{\text{try}}\), và tập này sẽ khác nhau một cách ngẫu nhiên qua từng nút.
Giá trị dự đoán cuối cùng của một Random Forest được xác định bằng
cách lấy trung bình các giá trị dự đoán của \(B\) cây. Đối với bài toán phân
loại (như dự báo là phá sản hay an toàn), kết quả dự đoán được
xác định bằng nguyên tắc bỏ phiếu đa số. Mỗi cây quyết
định sẽ đưa ra dự đoán về lớp (ví dụ: "phá sản" hoặc
"không phá sản"). Nhãn nào nhận được nhiều phiếu bầu nhất
từ các cây sẽ là dự đoán cuối cùng của mô hình.
Trong quá trình xây dựng các cây trong Random Forest, một khía cạnh quan trọng là việc xác định mức độ đóng góp của từng biến đồng biến đối với quá trình phân chia dữ liệu và dự đoán kết quả. Để định lượng điều này, chỉ số tầm quan trọng của biến (variable importance) được sử dụng, trong đó một phương pháp phổ biến là Gini Importance, còn gọi là Mean Decrease in Impurity (MDI).
Trong bối cảnh bài toán phân loại (ví dụ như phân loại doanh nghiệp có nguy cơ phá sản hay không), mỗi nút bên trong (internal node) của cây quyết định được xây dựng dựa trên việc chọn một biến đồng biến và một ngưỡng chia sao cho độ hỗn tạp (impurity) tại node đó được giảm nhiều nhất sau khi chia dữ liệu.
Một chỉ số thường dùng để đo độ hỗn tạp là Gini index, được tính như sau tại một nút \(t\):
\[ \text{Gini}(t) = 1 - \sum_{c=1}^{C} p_c^2 \]
Trong đó: - \(C\) là số lượng nhãn (ví dụ với biến phụ thuộc là nhị phân thì \(C = 2\)), - \(p_c\) là tỷ lệ mẫu của lớp \(c\) trong nút đó.
Mỗi lần chia dữ liệu, Random Forest sẽ chia sao cho độ hỗn tạp sẽ giảm đi tại nút đó. Khi một biến \(X^{(i)}\) được chọn để chia tại một nút, mức độ giảm độ hỗn tạp do biến đó tạo ra được tính theo công thức:
\[ \Delta i(t, X^{(i)}) = i(t) - \left[ \frac{N_L}{N_t} \cdot i(t_L) + \frac{N_R}{N_t} \cdot i(t_R) \right] \]
Trong đó: - \(i(t)\) là độ hỗn tạp tại nút hiện tại, - \(t_L, t_R\) là hai nút con sau khi chia, - \(N_t\) là số lượng mẫu tại nút hiện tại và \(N_L, N_R\) là số mẫu tại hai nút con.
Sau khi xây dựng toàn bộ cây, tổng độ giảm độ hỗn tạp mà một biến gây ra tại tất cả các nút mà nó được chọn làm tiêu chí chia sẽ được cộng dồn lại. Khi huấn luyện Random Forest với \(B\) cây, chỉ số tầm quan trọng của một biến sẽ được tính trung bình qua tất cả các cây. Cuối cùng, tầm quan trọng của biến được chuẩn hóa thành tỷ lệ phần trăm đóng góp để có thể so sánh giữa các biến.
Phương pháp được sử dụng để đánh giá hiệu suất của mô hình phân loại nhị phân dựa trên các thông tin thu được từ ma trận nhầm lẫn (confusion matrix). Một ma trận nhầm lẫn có thể được trình bày như sau:
| Thực tế: 1 | Thực tế: 0 | |
|---|---|---|
| Dự báo: 1 | TP | FP |
| Dự báo: 0 | FN | TN |
Trong đó:
- TP (True Positive): số lượng doanh nghiệp phá sản
trong thực tế và được dự báo đúng là phá sản.
- TN (True Negative): số lượng doanh nghiệp không phá
sản trong thực tế và được dự báo đúng là không phá sản.
- FP (False Positive): số lượng doanh nghiệp không phá
sản nhưng bị mô hình dự báo sai là phá sản.
- FN (False Negative): số lượng doanh nghiệp phá sản
nhưng bị mô hình dự báo sai là không phá sản.
Từ ma trận nhầm lẫn, các chỉ số đánh giá mô hình được tính như sau:
Accuracy = \((TP + TN) / (TP + TN + FP + FN)\): đo lường tỷ lệ dự báo chính xác tổng thể của mô hình. Theo nghiên cứu của Mu-Yen Chen (2011), Accuracy là một trong những chỉ số phổ biến nhất được sử dụng trong đánh giá mô hình.
Precision = \(TP / (TP + FP)\): đo lường tỷ lệ doanh nghiệp thực sự phá sản trong số các doanh nghiệp được mô hình dự báo là phá sản.
Recall (Sensitivity) = \(TP / (TP + FN)\): phản ánh khả năng phát hiện đúng các doanh nghiệp phá sản trong số toàn bộ các doanh nghiệp phá sản thực tế.
F1-score = \(2 \times \frac{Precision \times Recall}{Precision + Recall}\): là chỉ số trung bình điều hòa giữa Precision và Recall, đặc biệt hữu ích trong trường hợp dữ liệu mất cân bằng, giúp đánh giá độ chính xác tổng thể của mô hình một cách toàn diện.
Ngoài ra, một phương pháp đánh giá phổ biến khác là dựa trên ROC (Receiver Operating Characteristic) là một đường cong thể hiện hiệu suất của mô hình phân loại nhị phân ở các ngưỡng phân loại khác nhau.
Đường cong ROC được vẽ bằng cách biểu diễn: - Trục hoành (x-axis): Tỷ lệ dương giả (False Positive Rate – FPR) - Trục tung (y-axis): Tỷ lệ dương thật (True Positive Rate – TPR), còn gọi là Recall hoặc Sensitivity
Công thức: - \(\text{TPR} = \frac{TP}{TP + FN}\) - \(\text{FPR} = \frac{FP}{FP + TN}\)
Khi đó, chỉ số AUC (Area Under the Curve) là diện tích dưới đường cong ROC, được dùng để đo tổng thể khả năng phân biệt của mô hình.
Diễn giải: > AUC biểu thị xác suất mà mô hình sẽ xếp hạng một quan sát dương tính cao hơn một quan sát âm tính một cách ngẫu nhiên.
Tập dữ liệu gồm tổng cộng 1.879 quan sát đã được chia ngẫu nhiên thành hai phần: 70% (tương đương 1.314 quan sát) được sử dụng làm tập huấn luyện (train set) và 30% còn lại (565 quan sát) làm tập kiểm tra (test set). Việc phân chia được thực hiện theo phương pháp lấy mẫu ngẫu nhiên có phân tầng nhằm đảm bảo giữ nguyên tỷ lệ phân bố của biến phụ thuộc giữa hai tập, từ đó giúp mô hình tránh sai lệch do mất cân bằng dữ liệu. Các biến của mô hình bao gồm:
Biến phụ thuộc:
0 = Không mắc bệnh,
1 = Mắc bệnh tiểu đường.Biến độc lập:
SystolicBP: Huyết áp tâm thu, nằm trong khoảng từ 90 đến 180 mmHg.
FastingBloodSugar: Đường huyết lúc đói, nằm trong khoảng từ 70 đến 200 mg/dL.
HbA1c: Mức Hemoglobin A1c, nằm trong khoảng từ 4.0% đến 10.0%.
Smoking: Tình trạng hút thuốc, trong đó:
0: Không hút thuốc1: Có hút thuốcFrequentUrination: Triệu chứng đi tiểu thường xuyên, trong đó:
0: Không có triệu chứng1: Có triệu chứnglibrary(xlsx)
library(tidyverse)
data <- read.xlsx(file = "C:/Users/Khanh Hoa/OneDrive/UFM/R/Tuần 4 - PTDLDT/diabetes_data_org.xlsx", header = T, sheetIndex = 1) #load data
data1 <- select(data, c("Diagnosis","HbA1c","FastingBloodSugar","SystolicBP","Smoking","FrequentUrination"))
data1$Diagnosis <- as.factor(data1$Diagnosis)
data1$Smoking <- as.factor(data1$Smoking)
data1$FrequentUrination <- as.factor(data1$FrequentUrination)
library(rsample)
set.seed(111)
split <- initial_split(data1, prop = 0.7, strata = "Diagnosis") # Chia tập dữ liệu giữ tỷ lệ ban đầu của biến chẩn đoán
train_data <- training(split)
test_data <- testing(split)
library(kableExtra)
head(data1,5) %>%
kable(
caption = "<b style='display:block; text-align:center;'>Bảng 1: Năm dòng đầu tiên của dữ liệu</b>",
escape = FALSE, align = 'c'
) %>%
kable_styling(
position = "center", font_size = 13,
full_width = FALSE,
bootstrap_options = c("striped", "hover") #striped table
)| Diagnosis | HbA1c | FastingBloodSugar | SystolicBP | Smoking | FrequentUrination |
|---|---|---|---|---|---|
| 1 | 9.283631 | 163.68716 | 93 | 1 | 0 |
| 1 | 7.326871 | 188.34707 | 165 | 0 | 0 |
| 0 | 4.083426 | 127.70365 | 119 | 0 | 0 |
| 0 | 6.516645 | 82.68842 | 169 | 1 | 0 |
| 0 | 5.607221 | 90.74339 | 165 | 0 | 0 |
thong_ke_dinh_tinh <- function(data, var_name) {
#Convert to factor if needed
variable <- as.factor(data[[var_name]])
#freq
freq <- table(variable)
percent <- prop.table(freq) * 100
#result
result <- data.frame(
Variable = var_name,
Label = names(freq),
Frequency = as.vector(freq),
Percent = round(as.vector(percent), 2)
)
return(result)
}
# trình bày kết quả
tb1 <- thong_ke_dinh_tinh(train_data, "Diagnosis")
tb2 <- thong_ke_dinh_tinh(train_data, "Smoking")
tb3 <- thong_ke_dinh_tinh(train_data, "FrequentUrination")
# combine
rbind(tb1, tb2, tb3) %>%
kable(
caption = "<b style='display:block; text-align:center;'>Bảng 2: Thống kê mô tả các biến định tính tập train</b>",
escape = FALSE, align = 'c'
) %>%
kable_styling(
position = "center", font_size = 13,
full_width = FALSE,
bootstrap_options = c("striped", "hover") #striped table
)| Variable | Label | Frequency | Percent |
|---|---|---|---|
| Diagnosis | 0 | 788 | 59.97 |
| Diagnosis | 1 | 526 | 40.03 |
| Smoking | 0 | 955 | 72.68 |
| Smoking | 1 | 359 | 27.32 |
| FrequentUrination | 0 | 1060 | 80.67 |
| FrequentUrination | 1 | 254 | 19.33 |
library(ggplot2)
library(patchwork)
# Biểu đồ 1: Diagnosis
p1 <- ggplot(tb1, aes(x = as.factor(Label), y = Frequency, fill = as.factor(Label))) +
geom_bar(stat = "identity") +
coord_flip() +
labs(title = "Diagnosis", x = "Label", y = "Frequency") +
geom_text(aes(label = Frequency),
position = position_dodge(width = 0.9),
hjust = -0.1, colour = "#0099CC") +
theme_minimal() +
theme(legend.position = "none",
plot.title = element_text(hjust = 0.5, face = "bold", size = 14))
# Biểu đồ 2: Smoking
p2 <- ggplot(tb2, aes(x = as.factor(Label), y = Frequency, fill = as.factor(Label))) +
geom_bar(stat = "identity") +
labs(title = "Smoking", x = "Label", y = "Frequency") +
geom_text(aes(label = Frequency),
position = position_dodge(width = 0.9),
vjust = 1, colour = "#0099CC") +
theme_minimal() +
theme(legend.position = "none",
plot.title = element_text(hjust = 0.5, face = "bold", size = 14)) +
scale_fill_manual(values = c("0" = "#CCCC99", "1" = "#FFCCCC"))
# Biểu đồ 3: FrequentUrination
tb3$prop <- tb3$Frequency / sum(tb3$Frequency)
tb3$label_text <- paste0(tb3$Label, ": ", round(tb3$prop * 100, 1), "%")
p3 <- ggplot(tb3, aes(x = "", y = prop, fill = as.factor(Label))) +
geom_bar(stat = "identity", width = 1) +
coord_polar(theta = "y") +
labs(title = "Frequent Urination") +
theme_void() +
geom_text(aes(label = label_text), position = position_stack(vjust = 0.7), size = 3, color = "white") +
theme(
legend.position = "none",
plot.title = element_text(hjust = 0.5, face = "bold", size = 14)
) +
scale_fill_manual(values = c("0" = "#99CCFF", "1" = "#99CC33"))
(p1) / (p2 | p3)Thống kê mô tả ban đầu cho tập huấn luyện cho thấy sự phân bố của ba biến nhị phân liên quan đến sức khỏe. Trong số 1.314 quan sát, biến mục tiêu Diagnosis (chẩn đoán bệnh) có 59,97% số mẫu không mắc bệnh (label = 0) và 40,03% số mẫu mắc bệnh (label = 1), cho thấy dữ liệu không quá mất cân bằng. Biến Smoking (hút thuốc) ghi nhận 72,68% người không hút thuốc, trong khi 27,32% có hút thuốc, phản ánh tỷ lệ người hút thuốc tương đối thấp trong tập mẫu. Cuối cùng, FrequentUrination (tiểu nhiều) có 80,67% người không gặp triệu chứng, và chỉ 19,33% có dấu hiệu tiểu nhiều, cho thấy phần lớn mẫu không biểu hiện triệu chứng này
# trình bày kết quả
tb4 <- thong_ke_dinh_tinh(test_data, "Diagnosis")
tb5 <- thong_ke_dinh_tinh(test_data, "Smoking")
tb6 <- thong_ke_dinh_tinh(test_data, "FrequentUrination")
# combine
rbind(tb4, tb5, tb6) %>%
kable(
caption = "<b style='display:block; text-align:center;'>Bảng 3: Thống kê mô tả các biến định tính tập test</b>",
escape = FALSE, align = 'c'
) %>%
kable_styling(
position = "center", font_size = 13,
full_width = FALSE,
bootstrap_options = c("striped", "hover") #striped table
)| Variable | Label | Frequency | Percent |
|---|---|---|---|
| Diagnosis | 0 | 339 | 60.00 |
| Diagnosis | 1 | 226 | 40.00 |
| Smoking | 0 | 395 | 69.91 |
| Smoking | 1 | 170 | 30.09 |
| FrequentUrination | 0 | 448 | 79.29 |
| FrequentUrination | 1 | 117 | 20.71 |
# Biểu đồ 1: Diagnosis
p4 <- ggplot(tb4, aes(x = as.factor(Label), y = Frequency, fill = as.factor(Label))) +
geom_bar(stat = "identity") +
coord_flip() +
labs(title = "Diagnosis", x = "Label", y = "Frequency") +
geom_text(aes(label = Frequency),
position = position_dodge(width = 0.9),
hjust = -0.1, colour = "#0099CC") +
theme_minimal() +
theme(legend.position = "none",
plot.title = element_text(hjust = 0.5, face = "bold", size = 14))
# Biểu đồ 2: Smoking
p5 <- ggplot(tb5, aes(x = as.factor(Label), y = Frequency, fill = as.factor(Label))) +
geom_bar(stat = "identity") +
labs(title = "Smoking", x = "Label", y = "Frequency") +
geom_text(aes(label = Frequency),
position = position_dodge(width = 0.9),
vjust = 1, colour = "#0099CC") +
theme_minimal() +
theme(legend.position = "none",
plot.title = element_text(hjust = 0.5, face = "bold", size = 14)) +
scale_fill_manual(values = c("0" = "#CCCC99", "1" = "#FFCCCC"))
# Biểu đồ 3: FrequentUrination
tb6$prop <- tb6$Frequency / sum(tb6$Frequency)
tb6$label_text <- paste0(tb6$Label, ": ", round(tb6$prop * 100, 1), "%")
p6 <- ggplot(tb6, aes(x = "", y = prop, fill = as.factor(Label))) +
geom_bar(stat = "identity", width = 1) +
coord_polar(theta = "y") +
labs(title = "Frequent Urination") +
theme_void() +
geom_text(aes(label = label_text), position = position_stack(vjust = 0.7), size = 3, color = "white") +
theme(
legend.position = "none",
plot.title = element_text(hjust = 0.5, face = "bold", size = 14)
) +
scale_fill_manual(values = c("0" = "#99CCFF", "1" = "#99CC33"))
(p4) / (p5 | p6)Trong tập dữ liệu kiểm tra gồm 565 quan sát, tỷ lệ bệnh nhân được chẩn đoán mắc bệnh tiểu đường chiếm 40%, trong khi 60% còn lại không mắc bệnh, cho thấy sự phân bổ tương đối cân bằng giữa hai nhóm. Đối với biến Smoking, phần lớn người tham gia là không hút thuốc (chiếm khoảng 69.91%), trong khi chỉ có 30.09% là người hút thuốc. Về tần suất tiểu tiện, có đến 79.29% số người khảo sát không thường xuyên đi tiểu, trong khi 20.71% cho biết họ có triệu chứng tiểu tiện thường xuyên. Những thống kê này cho thấy sự tương đồng nhất định với tập huấn luyện, đảm bảo tính nhất quán và đại diện khi tiến hành xây dựng và đánh giá mô hình.
phan_tich_cap_bien <- function(data, var1, var2) {
tbl <- table(data[[var1]], data[[var2]])
percent_total <- round(prop.table(tbl) * 100, 2)
percent_row <- round(prop.table(tbl, margin = 1) * 100, 2)
percent_col <- round(prop.table(tbl, margin = 2) * 100, 2)
result <- list(
Bang_tan_suat = tbl,
Ti_le_theo_tong = percent_total,
Ti_le_theo_hang = percent_row,
Ti_le_theo_cot = percent_col
)
return(result)
}
result1 <- phan_tich_cap_bien(train_data, "Smoking", "Diagnosis")
#convert short into long
dia_sm1 <- as.data.frame(
as.table(result1$Bang_tan_suat)) %>%
rename(Smoking = Var1, Diagnosis = Var2, Frequency = Freq) %>%
mutate(
Percent_Total = as.vector(result1$Ti_le_theo_tong),
Percent_of_Smoking = as.vector(result1$Ti_le_theo_hang),
Percent_Diagnosis = as.vector(result1$Ti_le_theo_cot)
)
dia_sm1 <- data.frame(
Dataset = "Train set",
dia_sm1
)
result2 <- phan_tich_cap_bien(test_data, "Smoking", "Diagnosis")
#convert short into long
dia_sm2 <- as.data.frame(
as.table(result2$Bang_tan_suat)) %>%
rename(Smoking = Var1, Diagnosis = Var2, Frequency = Freq) %>%
mutate(
Percent_Total = as.vector(result2$Ti_le_theo_tong),
Percent_of_Smoking = as.vector(result2$Ti_le_theo_hang),
Percent_Diagnosis = as.vector(result2$Ti_le_theo_cot)
)
dia_sm2 <- data.frame(
Dataset = "Test set",
dia_sm2
)
rbind(dia_sm1, dia_sm2) %>%
kable(
caption = "<b style='display:block; text-align:center;'>Bảng 4: Thống kê mô tả cặp biến Diagnosis x Smoking</b>",
escape = FALSE, align = 'c'
) %>%
kable_styling(
position = "center", font_size = 13,
full_width = FALSE,
bootstrap_options = c("striped", "hover") #striped table
)| Dataset | Smoking | Diagnosis | Frequency | Percent_Total | Percent_of_Smoking | Percent_Diagnosis |
|---|---|---|---|---|---|---|
| Train set | 0 | 0 | 589 | 44.82 | 61.68 | 74.75 |
| Train set | 1 | 0 | 199 | 15.14 | 55.43 | 25.25 |
| Train set | 0 | 1 | 366 | 27.85 | 38.32 | 69.58 |
| Train set | 1 | 1 | 160 | 12.18 | 44.57 | 30.42 |
| Test set | 0 | 0 | 243 | 43.01 | 61.52 | 71.68 |
| Test set | 1 | 0 | 96 | 16.99 | 56.47 | 28.32 |
| Test set | 0 | 1 | 152 | 26.90 | 38.48 | 67.26 |
| Test set | 1 | 1 | 74 | 13.10 | 43.53 | 32.74 |
p7 <- ggplot(train_data, aes(x = Diagnosis, fill = Smoking)) +
geom_bar(position = "dodge") +
geom_text(stat = "count", aes(label = ..count..),
position = position_dodge(width = 0.9),
vjust = -0.5, colour = "#0099CC") +
labs(title = "Tập train", x = "Diagnosis", y = "Tần số") +
theme_minimal() +
theme(
plot.title = element_text(hjust = 0.5, face = "bold", size = 14),
axis.text.x = element_text(hjust = 1))
p8 <- ggplot(dia_sm2, aes(x = Diagnosis, y = Smoking, fill = Frequency)) +
geom_tile(color = "white") +
scale_fill_gradient(low = "white", high = "#0099CC") +
geom_text(aes(label = Frequency),
hjust = 0.5, colour = "#FFCCCC") +
labs(
title = "Tập test",
x = "Diagnosis",
y = "Smoking") +
theme_minimal() +
theme(
plot.title = element_text(hjust = 0.5, face = "bold", size = 14),
axis.text.x = element_text(angle = 0, hjust = 1)
)
p7 | p8Phân tích mối quan hệ giữa hành vi hút thuốc (Smoking) và khả năng được chẩn đoán mắc tiểu đường (Diagnosis) cho thấy xu hướng nhất quán trên cả tập huấn luyện và tập kiểm tra.
Trong tập huấn luyện, ước tính có khoảng 38.32% người không hút thuốc bị chẩn đoán mắc bệnh tiểu đường, trong khi tỷ lệ này ở nhóm có hút thuốc là 44.57%. Điều này cho thấy nhóm hút thuốc có nguy cơ được chẩn đoán mắc bệnh cao hơn so với nhóm không hút thuốc. Đồng thời, trong số những người được chẩn đoán mắc bệnh (526 người), phần lớn vẫn thuộc nhóm không hút thuốc (69.58%), điều này phản ánh rằng mặc dù tỷ lệ nguy cơ cao hơn trong nhóm hút thuốc, nhưng do số lượng người không hút thuốc nhiều hơn, nên họ chiếm tỷ lệ lớn hơn trong nhóm bệnh nhân.
Kết quả ở tập kiểm tra thể hiện xu hướng tương tự. Tỷ lệ người mắc bệnh trong nhóm không hút thuốc là 38.48%, thấp hơn so với 43.53% ở nhóm hút thuốc. Trong nhóm bệnh nhân được chẩn đoán, có đến 67.26% là người không hút thuốc, và 32.74% là người hút thuốc.
Nhìn chung, kết quả mô tả cho thấy hút thuốc có thể là yếu tố liên quan đến nguy cơ mắc bệnh tiểu đường, khi tỷ lệ người mắc bệnh trong nhóm hút thuốc luôn cao hơn nhóm không hút thuốc. Tuy nhiên, do số lượng người không hút thuốc chiếm đa số trong mẫu khảo sát, nên họ vẫn chiếm tỷ trọng lớn hơn trong tổng số ca được chẩn đoán mắc bệnh. Tiếp tục, để cung cấp thêm các đánh giá khách quan, chỉ số Relative Risk và Odds Ratio được tính toán bên dưới.
library(epitools)
tanso1 <- table(train_data$Smoking, train_data$Diagnosis)
rr1 <- epitab(tanso1, method = 'riskratio')
rr1 <- data.frame(rr1$tab)
or1 <- epitab(tanso1, method = 'oddsratio')
or1 <- data.frame(or1$tab)
tanso2 <- table(test_data$Smoking, test_data$Diagnosis)
rr2 <- epitab(tanso2, method = 'riskratio')
rr2 <- data.frame(rr2$tab)
or2 <- epitab(tanso2, method = 'oddsratio')
or2 <- data.frame(or2$tab)
data.frame(
Criteria = c("RR (P(Y=1;X=1)/P(Y=1;X=0))", "OR (Odds(X=1)/Odds(X=0))"),
Train_Set = c(rr1$riskratio[2],or1$oddsratio[2]),
Test_set = c(rr2$riskratio[2], or2$oddsratio[2])
) %>%
kable(
caption = "<b style='display:block; text-align:center;'>Bảng 5: Relative Risk và Odds Ratio của Diagnosis x Smoking</b>",
escape = FALSE, align = 'c'
) %>%
kable_styling(
position = "center", font_size = 13,
full_width = FALSE,
bootstrap_options = c("striped", "hover") #striped table
)## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
| Criteria | Train_Set | Test_set |
|---|---|---|
| RR (P(Y=1;X=1)/P(Y=1;X=0)) | 1.162915 | 1.131192 |
| OR (Odds(X=1)/Odds(X=0)) | 1.293901 | 1.232319 |
RR ~ 1.16 (train) và 1.13 (test), Có thể hiểu rằng người hút thuốc có nguy cơ được chẩn đoán mắc bệnh tiểu đường cao hơn khoảng 13–16% so với người không hút thuốc. Tương tự, OR ~ 1.29 (train) và 1.23 (test), cho thấy người hút thuốc có odds (tỷ lệ chọi) mắc bệnh cao hơn 23–29% so với người không hút thuốc.
result3 <- phan_tich_cap_bien(train_data, "FrequentUrination", "Diagnosis")
#convert short into long
dia_u1 <- as.data.frame(
as.table(result3$Bang_tan_suat)) %>%
rename(FrequentUrination = Var1, Diagnosis = Var2, Frequency = Freq) %>%
mutate(
Percent_Total = as.vector(result3$Ti_le_theo_tong),
Percent_of_FrequentUrination = as.vector(result3$Ti_le_theo_hang),
Percent_Diagnosis = as.vector(result3$Ti_le_theo_cot)
)
dia_u1 <- data.frame(
Dataset = "Train set",
dia_u1
)
result4 <- phan_tich_cap_bien(test_data, "FrequentUrination", "Diagnosis")
#convert short into long
dia_u2 <- as.data.frame(
as.table(result4$Bang_tan_suat)) %>%
rename(FrequentUrination = Var1, Diagnosis = Var2, Frequency = Freq) %>%
mutate(
Percent_Total = as.vector(result4$Ti_le_theo_tong),
Percent_of_FrequentUrination = as.vector(result4$Ti_le_theo_hang),
Percent_Diagnosis = as.vector(result4$Ti_le_theo_cot)
)
dia_u2 <- data.frame(
Dataset = "Test set",
dia_u2
)
rbind(dia_u1, dia_u2) %>%
kable(
caption = "<b style='display:block; text-align:center;'>Bảng 6: Thống kê mô tả cặp biến Diagnosis x FrequentUrination</b>",
escape = FALSE, align = 'c'
) %>%
kable_styling(
position = "center", font_size = 13,
full_width = FALSE,
bootstrap_options = c("striped", "hover") #striped table
)| Dataset | FrequentUrination | Diagnosis | Frequency | Percent_Total | Percent_of_FrequentUrination | Percent_Diagnosis |
|---|---|---|---|---|---|---|
| Train set | 0 | 0 | 669 | 50.91 | 63.11 | 84.90 |
| Train set | 1 | 0 | 119 | 9.06 | 46.85 | 15.10 |
| Train set | 0 | 1 | 391 | 29.76 | 36.89 | 74.33 |
| Train set | 1 | 1 | 135 | 10.27 | 53.15 | 25.67 |
| Test set | 0 | 0 | 291 | 51.50 | 64.96 | 85.84 |
| Test set | 1 | 0 | 48 | 8.50 | 41.03 | 14.16 |
| Test set | 0 | 1 | 157 | 27.79 | 35.04 | 69.47 |
| Test set | 1 | 1 | 69 | 12.21 | 58.97 | 30.53 |
p9 <- ggplot(train_data, aes(x = Diagnosis, fill = FrequentUrination)) +
geom_bar(position = "dodge") +
geom_text(stat = "count", aes(label = ..count..),
position = position_dodge(width = 0.9),
vjust = -0.5, colour = "#0099CC") +
labs(title = "Tập train", x = "Diagnosis", y = "Tần số") +
theme_minimal() +
theme(
plot.title = element_text(hjust = 0.5, face = "bold", size = 14),
axis.text.x = element_text(hjust = 1))
p10 <- ggplot(dia_u2, aes(x = Diagnosis, y = FrequentUrination, fill = Frequency)) +
geom_tile(color = "white") +
scale_fill_gradient(low = "white", high = "#0099CC") +
geom_text(aes(label = Frequency),
hjust = 0.5, colour = "#FFCCCC") +
labs(
title = "Tập test",
x = "Diagnosis",
y = "FrequentUrination") +
theme_minimal() +
theme(
plot.title = element_text(hjust = 0.5, face = "bold", size = 14),
axis.text.x = element_text(angle = 0, hjust = 1)
)
p9 | p10Dữ liệu mô tả cho thấy có mối liên hệ rõ rệt giữa triệu chứng tiểu tiện thường xuyên (Frequent Urination) và tình trạng được chẩn đoán mắc bệnh tiểu đường (Diagnosis). Trong cả hai tập dữ liệu huấn luyện và kiểm tra, phần lớn các trường hợp không có tiểu tiện thường xuyên thuộc về nhóm không mắc bệnh (chiếm khoảng 85% trong tập train và 86% trong tập test). Ngược lại, tỷ lệ người mắc bệnh trong nhóm có triệu chứng tiểu tiện thường xuyên cao hơn đáng kể so với nhóm không có triệu chứng này: cụ thể trong tập huấn luyện, 53.15% người có triệu chứng bị chẩn đoán mắc bệnh, so với chỉ 36.89% ở nhóm không triệu chứng. Tương tự, trong tập kiểm tra, tỷ lệ này là 58.97% so với 35.04%.
Ngoài ra, nếu xét theo tỷ lệ trong từng nhóm Diagnosis, có đến 25–30% người mắc bệnh thuộc nhóm có triệu chứng tiểu tiện thường xuyên, trong khi tỷ lệ này ở nhóm không mắc bệnh chỉ khoảng 15%. Những con số này cho thấy Frequent Urination là một dấu hiệu cảnh báo tiềm năng, có giá trị phân biệt tốt hơn so với biến Smoking trước đó. Kết quả này hỗ trợ giả định rằng tiểu tiện thường xuyên có thể là một trong những yếu tố dự báo quan trọng cho bệnh tiểu đường.
tanso1 <- table(train_data$FrequentUrination, train_data$Diagnosis)
rr1 <- epitab(tanso1, method = 'riskratio')
rr1 <- data.frame(rr1$tab)
or1 <- epitab(tanso1, method = 'oddsratio')
or1 <- data.frame(or1$tab)
tanso2 <- table(test_data$FrequentUrination, test_data$Diagnosis)
rr2 <- epitab(tanso2, method = 'riskratio')
rr2 <- data.frame(rr2$tab)
or2 <- epitab(tanso2, method = 'oddsratio')
or2 <- data.frame(or2$tab)
data.frame(
Criteria = c("RR (P(Y=1;X=1)/P(Y=1;X=0))", "OR (Odds(X=1)/Odds(X=0))"),
Train_Set = c(rr1$riskratio[2],or1$oddsratio[2]),
Test_set = c(rr2$riskratio[2], or2$oddsratio[2])
) %>%
kable(
caption = "<b style='display:block; text-align:center;'>Bảng 5: Relative Risk và Odds Ratio của Diagnosis x FrequentUrination</b>",
escape = FALSE, align = 'c'
) %>%
kable_styling(
position = "center", font_size = 13,
full_width = FALSE,
bootstrap_options = c("striped", "hover") #striped table
)## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
| Criteria | Train_Set | Test_set |
|---|---|---|
| RR (P(Y=1;X=1)/P(Y=1;X=0)) | 1.440884 | 1.682835 |
| OR (Odds(X=1)/Odds(X=0)) | 1.941048 | 2.664411 |
Trong tập huấn luyện, tỷ lệ mắc bệnh ở nhóm có triệu chứng tiểu tiện thường xuyên cao hơn nhóm không có triệu chứng khoảng 44% (RR = 1.44), đồng nghĩa với việc nhóm này có nguy cơ mắc bệnh cao hơn đáng kể. Xét theo Odds Ratio, Odds mắc bệnh của nhóm này cao hơn gần 1.94 lần so với nhóm không có triệu chứng.
Trong tập kiểm tra, kết quả thậm chí còn rõ rệt hơn với RR = 1.68 và OR = 2.66, tức là nhóm có triệu chứng tiểu tiện thường xuyên có Odds mắc bệnh cao hơn tới gần 2.7 lần theo tỷ lệ chênh odds.
train_quant_data <- select(train_data, c("HbA1c","FastingBloodSugar","SystolicBP"))
test_quant_data <- select(test_data, c("HbA1c","FastingBloodSugar","SystolicBP"))
rbind(
data.frame(
Dataset = "Train Set",
skimr::skim(train_quant_data)),
data.frame(
Dataset = "Test Set",
skimr::skim(test_quant_data))) %>%
kable(
caption = "<b style='display:block; text-align:center;'>Bảng 5: Relative Risk và Odds Ratio của Diagnosis x FrequentUrination</b>",
escape = FALSE, align = 'c'
) %>%
kable_styling(
position = "center", font_size = 13,
full_width = FALSE,
bootstrap_options = c("striped", "hover") #striped table
)## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
| Dataset | skim_type | skim_variable | n_missing | complete_rate | numeric.mean | numeric.sd | numeric.p0 | numeric.p25 | numeric.p50 | numeric.p75 | numeric.p100 | numeric.hist |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Train Set | numeric | HbA1c | 0 | 1 | 7.003278 | 1.745091 | 4.003088 | 5.458524 | 7.088102 | 8.478017 | 9.991193 | ▇▇▇▇▇ |
| Train Set | numeric | FastingBloodSugar | 0 | 1 | 135.561918 | 37.274788 | 70.074649 | 102.555304 | 137.408488 | 167.922448 | 199.935506 | ▇▇▇▇▇ |
| Train Set | numeric | SystolicBP | 0 | 1 | 134.370624 | 25.699454 | 90.000000 | 113.000000 | 134.000000 | 156.000000 | 179.000000 | ▇▇▇▇▇ |
| Test Set | numeric | HbA1c | 0 | 1 | 6.913003 | 1.725863 | 4.008058 | 5.384626 | 7.112198 | 8.302510 | 9.970911 | ▇▆▇▇▆ |
| Test Set | numeric | FastingBloodSugar | 0 | 1 | 134.373233 | 38.090575 | 70.171130 | 101.164192 | 137.031202 | 166.162554 | 199.758457 | ▇▆▇▇▇ |
| Test Set | numeric | SystolicBP | 0 | 1 | 133.306195 | 25.420647 | 90.000000 | 111.000000 | 134.000000 | 154.000000 | 179.000000 | ▇▇▇▇▆ |
Ba biến định lượng quan trọng gồm HbA1c, đường huyết lúc đói (Fasting Blood Sugar) và huyết áp tâm thu (Systolic Blood Pressure) đều không có giá trị thiếu và có sự phân bố tương đối đồng đều giữa tập huấn luyện và tập kiểm tra, cho thấy dữ liệu được chia tách hợp lý và đáng tin cậy để phân tích. Biến HbA1c có giá trị trung bình khoảng 7.0 (Train) và 6.91 (Test), cao hơn ngưỡng bình thường (dưới 5.7%), cho thấy nhiều bệnh nhân trong mẫu thuộc nhóm tiền đái tháo đường hoặc đái tháo đường. Độ lệch chuẩn khoảng 1.7 phản ánh mức độ phân tán trung bình, trong khi phân vị 25–75% dao động từ khoảng 5.4 đến 8.4. Biến Fasting Blood Sugar có giá trị trung bình khoảng 135 mg/dL ở cả hai tập, cũng cao hơn ngưỡng bình thường (dưới 100 mg/dL), với độ lệch chuẩn lớn hơn (~37–38), phản ánh sự phân tán rộng hơn giữa các bệnh nhân; khoảng phân vị giữa 25–75% nằm trong khoảng 101–167 mg/dL, cho thấy phần lớn bệnh nhân có đường huyết cao hơn mức an toàn. Đối với Systolic Blood Pressure, giá trị trung bình khoảng 134 mmHg, cao hơn ngưỡng lý tưởng (dưới 120 mmHg), cho thấy nhiều bệnh nhân trong mẫu đang ở trạng thái tiền tăng huyết áp hoặc tăng huyết áp nhẹ. Tổng thể, cả ba biến định lượng đều có giá trị trung bình vượt ngưỡng bình thường, đồng thời có sự phân bố ổn định giữa hai tập dữ liệu, khẳng định vai trò quan trọng của chúng trong việc mô tả đặc điểm sức khỏe và hỗ trợ phân loại nguy cơ bệnh lý trong nghiên cứu này.
model <- glm(Diagnosis ~ HbA1c + FastingBloodSugar + SystolicBP + Smoking + FrequentUrination, data = train_data, family = "binomial")
# Calculate VIF
vif_results <- car::vif(model)
# Convert VIF results to a data frame for plotting
vif_df <- data.frame(Variable = names(vif_results), VIF = vif_results)
# Set a threshold to indicate high VIF
high_vif_threshold <- 5
# Create a ggplot bar plot to visualize VIF values
ggplot(vif_df, aes(x = Variable, y = VIF)) +
geom_bar(stat = "identity", fill = "#99CCCC") +
geom_text(aes(label = round(VIF, 2)), vjust = -0.5, color = "#99CC00") +
geom_hline(yintercept = high_vif_threshold, linetype = "dashed", color = "red") +
scale_y_continuous(limits = c(0, max(vif_df$VIF) + 1)) +
labs(title = "Variance Inflation Factor (VIF)",
y = "VIF",
x = "Variable") +
theme_minimal() +
theme(
plot.title = element_text(hjust = 0.5, face = "bold", size = 14),
axis.text.x = element_text(angle = 45, hjust = 1)) +
ylim(c(0,10))Kết quả kiểm định đa cộng tuyến thông qua hệ số phóng đại phương sai (VIF) cho thấy tất cả các biến độc lập trong mô hình đều có VIF nhỏ hơn 2, với giá trị dao động từ 1.0 đến 1.4. Cụ thể, các biến như FastingBloodSugar, FrequentUrination, HbA1c, Smoking, và SystolicBP đều có VIF thấp, trong đó biến cao nhất là HbA1c với VIF = 1.38. Mức ngưỡng cảnh báo thường được sử dụng trong thực hành là 5 hoặc 10, và tất cả các biến trong mô hình đều nằm dưới ngưỡng này một cách đáng kể. Do đó, có thể kết luận rằng mô hình không có dấu hiệu đa cộng tuyến nghiêm trọng, các biến giải thích không bị phụ thuộc tuyến tính mạnh với nhau, giúp đảm bảo tính ổn định và độ tin cậy cho ước lượng hồi quy.
#logit
logit_model <- glm(Diagnosis ~ HbA1c + FastingBloodSugar + SystolicBP + Smoking + FrequentUrination, data = train_data, family = binomial(link = "logit"))
#probit
probit_model <- glm(Diagnosis ~ HbA1c + FastingBloodSugar + SystolicBP + Smoking + FrequentUrination, data = train_data, family = binomial(link = "probit"))
#cloglog
cloglog_model <- glm(Diagnosis ~ HbA1c + FastingBloodSugar + SystolicBP + Smoking + FrequentUrination, data = train_data, family = binomial(link = "cloglog"))
# show and compare
jtools::export_summs(list("Logistic" = logit_model,
"Probit" = probit_model,
"Cloglog" = cloglog_model),scale = TRUE)| Logistic | Probit | Cloglog | |
|---|---|---|---|
| (Intercept) | -1.08 *** | -0.61 *** | -1.24 *** |
| (0.10) | (0.06) | (0.08) | |
| HbA1c | 1.48 *** | 0.83 *** | 0.99 *** |
| (0.10) | (0.05) | (0.06) | |
| FastingBloodSugar | 1.55 *** | 0.88 *** | 1.02 *** |
| (0.10) | (0.05) | (0.06) | |
| SystolicBP | -0.10 | -0.06 | -0.01 |
| (0.08) | (0.04) | (0.05) | |
| Smoking1 | 0.25 | 0.16 | 0.20 |
| (0.17) | (0.10) | (0.11) | |
| FrequentUrination1 | 1.40 *** | 0.80 *** | 0.93 *** |
| (0.20) | (0.11) | (0.13) | |
| N | 1314 | 1314 | 1314 |
| AIC | 1101.03 | 1106.61 | 1114.85 |
| BIC | 1132.12 | 1137.69 | 1145.94 |
| Pseudo R2 | 0.55 | 0.54 | 0.54 |
| All continuous predictors are mean-centered and scaled by 1 standard deviation. The outcome variable is in its original units. *** p < 0.001; ** p < 0.01; * p < 0.05. | |||
pred_logit <- predict(logit_model, newdata = test_data, type = "response")
pred_probit <- predict(probit_model, newdata = test_data, type = "response")
pred_cloglog <- predict(cloglog_model, newdata = test_data, type = "response")
pred_logit1 <- predict(logit_model, newdata = train_data, type = "response")
pred_probit1 <- predict(probit_model, newdata = train_data, type = "response")
pred_cloglog1 <- predict(cloglog_model, newdata = train_data, type = "response")
library(pROC)
roc_logit <- roc(test_data$Diagnosis, pred_logit)
roc_probit <- roc(test_data$Diagnosis, pred_probit)
roc_cloglog <- roc(test_data$Diagnosis, pred_cloglog)
roc_logit1 <- roc(train_data$Diagnosis, pred_logit1)
roc_probit1 <- roc(train_data$Diagnosis, pred_probit1)
roc_cloglog1 <- roc(train_data$Diagnosis, pred_cloglog1)
ggroc(list(
Train_Logistic = roc_logit1,
Train_Probit = roc_probit1,
Train_Cloglog = roc_cloglog1
)) +
geom_segment(aes(x = 1, y = 0, xend = 0, yend = 1),
color = "gray", linetype = "dashed") +
labs(
title = "Đường cong ROC",
subtitle = paste0(
"AUC Train Logistic = ", round(auc(roc_logit1), 5), " | ",
"Probit = ", round(auc(roc_probit1), 5), " | ",
"Cloglog = ", round(auc(roc_cloglog1), 5)
)
) +
theme_minimal() +
theme(
plot.title = element_text(hjust = 0.5, face = "bold", size = 14),
plot.subtitle = element_text(hjust = 0.5)) +
scale_color_manual(values = c("#330000","#00DD00","#FF0000"))Trong mô hình Logistic, HbA1c (hệ số = 0.850, p < 0.001) có \(e^{0.850} \approx 2.34\), cho thấy mỗi đơn vị tăng (% HbA1c) làm tăng odds mắc bệnh lên 2.34 lần, khẳng định vai trò quan trọng trong chẩn đoán. FastingBloodSugar (hệ số = 0.042, p < 0.001) với \(e^{0.042} \approx 1.043\) làm tăng odds mắc bệnh 4.3% mỗi mg/dL tăng thêm. FrequentUrination1 (hệ số = 1.399, p < 0.001) có \(e^{1.399} \approx 4.05\), cho thấy người có triệu chứng đái nhiều có odds mắc bệnh cao gấp 4.05 lần, phản ánh tính đặc hiệu của triệu chứng (Gupta et al., 2018). Smoking1 (hệ số = 0.255, p = 0.131) với \(e^{0.255} \approx 1.29\) làm tăng odds 29%, nhưng không đạt ý nghĩa thống kê (p > 0.05). SystolicBP (hệ số = -0.004, p = 0.176) với \(e^{-0.004} \approx 0.996\) có tác động âm không đáng kể, cho thấy huyết áp tâm thu không phải yếu tố dự đoán mạnh trong mẫu này .
Trong mô hình Probit và Cloglog, các hệ số nhỏ hơn (ví dụ: HbA1c: 0.478 và 0.565; FastingBloodSugar: 0.023 và 0.027; FrequentUrination: 0.802 và 0.931) do đặc điểm hàm liên kết khác nhau, nhưng HbA1c, FastingBloodSugar, và FrequentUrination vẫn có ý nghĩa thống kê (p < 0.001), trong khi SystolicBP (p = 0.196 và 0.882) và Smoking (p = 0.093 và 0.071) không đạt ý nghĩa ở mức 0.05. Vấn đề này có thể được giải thích rằng tỷ lệ người hút thuốc (27,32%) và phạm vi SystolicBP (90–179 mmHg, trung bình 134,37 mmHg) có thể không đủ biến thiên để phát hiện hiệu ứng trong mẫu này. (4) Bối cảnh nghiên cứu: Các nghiên cứu trước thường phân tích trên dân số lớn hơn hoặc dài hạn, trong khi mẫu này có thể tập trung vào các yếu tố lâm sàng ngắn hạn. Do đó, SystolicBP và Smoking có thể vẫn quan trọng trong các bối cảnh khác, nhưng trong nghiên cứu này, chúng bị lấn át bởi các yếu tố mạnh hơn như HbA1c, FastingBloodSugar, và FrequentUrination.
So sánh ba mô hình dựa trên các chỉ số đánh giá, mô hình Logistic vượt trội với AIC thấp nhất (1101.0), BIC thấp nhất (1132.1), và Log-Likelihood cao nhất (-544.516), cho thấy nó phù hợp hơn với dữ liệu so với Probit (AIC = 1106.6; BIC = 1137.7; Log-Lik = -547.304) và Cloglog (AIC = 1114.9; BIC = 1145.9; Log-Lik = -551.427). Cả ba mô hình có RMSE giống nhau (0.36), cho thấy độ chính xác dự đoán tương đương, nhưng Logistic có lợi thế về khả năng giải thích và tính đơn giản.
Thêm vào đó, Biểu đồ ROC cho thấy cả ba mô hình hồi quy nhị phân — Logistic, Probit và Cloglog — đều đạt hiệu suất dự đoán rất tốt, với giá trị AUC trên tập huấn luyện là 0.8862, 0.8862 và 0.8859. Sự khác biệt giữa AUC của tập huấn luyện ở cả ba mô hình đều rất nhỏ, cho thấy các mô hình không bị overfitting và có tính ổn định cao. Đường cong ROC của các mô hình gần như trùng nhau, thể hiện hiệu suất phân loại tương đương và đều vượt trội so với đường phân loại ngẫu nhiên. Trong đó, Logistic và Probit gần như giống hệt nhau về kết quả, còn Cloglog chỉ thấp hơn rất nhẹ. Do đó, mô hình Logistic là lựa chọn phù hợp nhất vì không chỉ có hiệu suất tương đương với các mô hình khác mà còn dễ diễn giải hơn.
Với mô hình Logistic được xem là phù hợp nhất, để dễ dàng diễn giải, ta có thể tính toán hiệu ứng biên bằng cách lấy đạo hàm riêng của \(p\) theo biến độc lập \(X_i\) để đo lường sự thay đổi của xác suất \(p\) khi \(X_i\) thay đổi một đơn vị. Công thức tính hiệu ứng biên cho mô hình Logistic được trình bày như sau:
\[ \frac{\partial p}{\partial X_i} = p \times (1 - p) \times \beta_i, \]
trong đó, \(p \times (1 - p)\) là xác suất biên, phụ thuộc vào giá trị cụ thể của \(X_i\). Ở nghiên cứu này, hiệu ứng biên sẽ được tính cho từng quan sát của \(X_i\) và rồi lấy trung bình các giá trị đó ((Average Marginal Effect - AME).
library(margins)
summary(margins(logit_model)) %>%
kable(
caption = "<b style='display:block; text-align:center;'>Bảng 6: Hiệu ứng biên của mô hình Logistic</b>",
escape = FALSE, align = 'c'
) %>%
kable_styling(
position = "center", font_size = 13,
full_width = FALSE,
bootstrap_options = c("striped", "hover") #striped table
)| factor | AME | SE | z | p | lower | upper |
|---|---|---|---|---|---|---|
| FastingBloodSugar | 0.0055557 | 0.0002053 | 27.060891 | 0.0000000 | 0.0051533 | 0.0059581 |
| FrequentUrination1 | 0.1916521 | 0.0257056 | 7.455658 | 0.0000000 | 0.1412701 | 0.2420342 |
| HbA1c | 0.1134641 | 0.0043312 | 26.196764 | 0.0000000 | 0.1049751 | 0.1219532 |
| Smoking1 | 0.0342498 | 0.0227636 | 1.504584 | 0.1324311 | -0.0103661 | 0.0788657 |
| SystolicBP | -0.0005293 | 0.0003900 | -1.357237 | 0.1747059 | -0.0012937 | 0.0002351 |
Kết quả AME từ mô hình logistic cho thấy mức độ thay đổi trung bình trong xác suất mắc bệnh tiểu đường (Diagnosis = 1) khi mỗi biến độc lập tăng một đơn vị, giữ các biến khác cố định. HbA1c có AME = 0,1135 (SE = 0,0043, p < 0,0001), nghĩa là mỗi đơn vị tăng (% HbA1c) làm tăng xác suất mắc bệnh thêm 11,35 điểm phần trăm, với ý nghĩa thống kê cao (z = 26,1968), xác nhận vai trò quan trọng của HbA1c trong chẩn đoán tiểu đường. FastingBloodSugar với AME = 0,0056 (SE = 0,0002, p < 0,0001) cho thấy mỗi mg/dL tăng làm tăng xác suất mắc bệnh 0,56 điểm phần trăm, tuy hiệu ứng nhỏ hơn HbA1c nhưng vẫn rất đáng kể (z = 27,0609). FrequentUrination1 có AME = 0,1917 (SE = 0,0257, p < 0,0001), tức là những người có triệu chứng đái nhiều có xác suất mắc bệnh cao hơn 19,17 điểm phần trăm so với người không có triệu chứng, nhấn mạnh tính đặc hiệu của triệu chứng này. Hai biến còn lại không có ý nghĩa thống kê.
Để mở rộng phân tích dự đoán chẩn đoán bệnh tiểu đường, thuật toán Random Forest được áp dụng như một phương pháp học máy bổ sung nhằm đánh giá khả năng dự đoán của các biến độc lập (HbA1c, FastingBloodSugar, SystolicBP, Smoking, FrequentUrination). Không giống các mô hình hồi quy như Logistic, Probit, hay Cloglog, vốn cung cấp hệ số và odds ratio để diễn giải trực tiếp mối quan hệ giữa các biến độc lập và biến phụ thuộc Diagnosis, Random Forest là một mô hình phi tham số dựa trên tập hợp cây quyết định, khó diễn giải mối quan hệ tuyến tính hay xác suất cụ thể. Thay vào đó, nó cung cấp giá trị tầm quan trọng của biến (feature importance), đo lường mức độ đóng góp của mỗi biến trong việc cải thiện độ chính xác dự đoán của mô hình. Do đó, đánh giá Random Forest sẽ tập trung vào việc xác định tầm quan trọng của các biến, giúp nhận diện những yếu tố có ảnh hưởng lớn nhất đến chẩn đoán bệnh tiểu đường, đồng thời bổ sung góc nhìn phi tuyến tính so với các mô hình hồi quy trước đó.
library(randomForest)
rf_model <- randomForest(Diagnosis ~ HbA1c + FastingBloodSugar + SystolicBP + Smoking + FrequentUrination, data = train_data)
print(rf_model)Call: randomForest(formula = Diagnosis ~ HbA1c + FastingBloodSugar + SystolicBP + Smoking + FrequentUrination, data = train_data) Type of random forest: classification Number of trees: 500 No. of variables tried at each split: 2
OOB estimate of error rate: 10.12%
Confusion matrix: 0 1 class.error 0 745 43 0.05456853 1 90 436 0.17110266
# show importance of vars
library(waffle)
importance <- c(
HbA1c = rf_model$importance[1],
FastingBloodSugar = rf_model$importance[2],
SystolicBP = rf_model$importance[3],
Smoking = rf_model$importance[4],
FrequentUrination = rf_model$importance[5]
)
importance_pct <- round(importance / sum(importance) * 100)
waffle(importance_pct, rows = 10, size = 0.5,
colors = c("#99FFFF", "#f1c40f", "#FFCCCC", "#CC66FF", "#95a5a6"),
title = "Tầm quan trọng của các biến trong khả năng dự báo",
xlab = "1 square = 1%") +
labs(subtitle = paste0(
"Mean Decrease Gini: HbA1c = ", round(rf_model$importance[1],3), " | ",
"FastingBloodSugar = ", round(rf_model$importance[2],3), "\n",
"SystolicBP = ", round(rf_model$importance[3],3), " | ",
"Smoking = ", round(rf_model$importance[4],3), " | ",
"FrequentUrination = ", round(rf_model$importance[5],3))) +
theme_minimal() +
theme(
plot.title = element_text(hjust = 0.5, face = "bold", size = 14),
plot.subtitle = element_text(hjust = 0.5)) Kết quả phân tích tầm quan trọng biến trong mô hình Random Forest cho thấy hai yếu tố Fasting Blood Sugar và HbA1c là những biến quan trọng nhất, với giá trị MeanDecreaseGini lần lượt là 261.73 và 256.01, vượt trội so với các biến còn lại như SystolicBP, Smoking, hay FrequentUrination. Điều này hoàn toàn tương đồng với kết quả từ mô hình hồi quy logistic (và cả Bayes logistic trước đó), trong đó hai biến này cũng có hệ số ước lượng cao và có ý nghĩa thống kê mạnh, cho thấy ảnh hưởng lớn đến khả năng mắc bệnh tiểu đường. Sự thống nhất giữa hai mô hình — một mô hình tuyến tính (logistic) và một mô hình phi tuyến mạnh mẽ hơn (Random Forest) — củng cố độ tin cậy trong việc khẳng định HbA1c và Fasting Blood Sugar là hai chỉ báo sinh học then chốt trong việc chẩn đoán bệnh tiểu đường.
Để đánh giá hiệu suất dự báo một cách tổng quan, một thuật toán học máy phổ biến là Random Forest sẽ được cung cấp và so sánh với mô hình Logistic.
pred_logit <- predict(logit_model,newdata = test_data, type = "response")
pred_rf <- predict(rf_model, newdata = test_data, type = "prob")[, 2]
roc_logit <- roc(test_data$Diagnosis, pred_logit)
roc_rf <- roc(test_data$Diagnosis, pred_rf)
ggroc(list(
Logistic = roc_logit,
Random_Forest = roc_rf
)) +
geom_segment(aes(x = 1, y = 0, xend = 0, yend = 1),
color = "gray", linetype = "dashed") +
labs(
title = "Đường cong ROC",
subtitle = paste0(
"AUC Logistic = ", round(auc(roc_logit), 5), " | ",
"Random Forest = ", round(auc(roc_rf), 5)
)
) +
theme_minimal() +
theme(
plot.title = element_text(hjust = 0.5, face = "bold", size = 14),
plot.subtitle = element_text(hjust = 0.5)) +
scale_color_manual(values = c("#00DD00","#FF0000"))Biểu đồ ROC trên cho thấy cả hai mô hình Logistic Regression và Random Forest đều đạt hiệu quả phân loại tốt, với diện tích dưới đường cong (AUC) lần lượt là 0.90303 và 0.92997. Mô hình Random Forest vượt trội hơn nhẹ so với Logistic Regression, thể hiện qua đường cong nằm phía trên và gần góc trên bên trái hơn — điều này phản ánh khả năng phân biệt giữa hai lớp cao hơn. Sự chênh lệch AUC tuy không lớn nhưng cho thấy Random Forest tận dụng tốt các mối quan hệ phi tuyến và tương tác giữa các biến, phù hợp với đánh giá tầm quan trọng biến trước đó, trong khi Logistic Regression vẫn là một mô hình tuyến tính hiệu quả và dễ diễn giải. Tiếp tục, các chỉ tiêu đánh giá hiệu suất dự báo sẽ được cung cấp.
# gắn nhãn
pred_logit_label <- ifelse(pred_logit >= 0.5, 1, 0)
pred_rf_label <- ifelse(pred_rf >= 0.5, 1, 0)
# ma trận nhầm lẫn
library(caret)
# Logistic
cm_logit <- confusionMatrix(factor(pred_logit_label), factor(test_data$Diagnosis), positive = "1")
# Random Forest
cm_rf <- confusionMatrix(factor(pred_rf_label), factor(test_data$Diagnosis), positive = "1")
cm_logit_tb <- as.data.frame(cm_logit$table)
cm_rf_tb <- as.data.frame(cm_rf$table)
#visualize
cm_logit_v <- ggplot(cm_logit_tb, aes(x = Reference, y = Prediction, fill = Freq)) +
geom_tile(color = "white") +
scale_fill_gradient(low = "white", high = "#0099CC") +
geom_text(aes(label = Freq),
hjust = 0.5, colour = "#FFCCCC", size = 5) +
labs(
title = "Logistic",
x = "Reference",
y = "Prediction") +
theme_minimal() +
theme(
plot.title = element_text(hjust = 0.5, face = "bold", size = 14),
axis.text.x = element_text(angle = 0, hjust = 1)
)
cm_rf_v <- ggplot(cm_rf_tb, aes(x = Reference, y = Prediction, fill = Freq)) +
geom_tile(color = "white") +
scale_fill_gradient(low = "white", high = "#0099CC") +
geom_text(aes(label = Freq),
hjust = 0.5, colour = "#FFCCCC", size = 5) +
labs(
title = "Random Forest",
x = "Reference",
y = "Prediction") +
theme_minimal() +
theme(
plot.title = element_text(hjust = 0.5, face = "bold", size = 14),
axis.text.x = element_text(angle = 0, hjust = 1)
)
cm_logit_v | cm_rf_v data.frame(
Model = c("Logistic", "Random Forest"),
Accuracy = c(round(cm_logit[["overall"]][["Accuracy"]],3),round(cm_rf[["overall"]][["Accuracy"]],3)),
Precision = c(round(cm_logit[["byClass"]][["Precision"]],3),round(cm_rf[["byClass"]][["Precision"]],3)),
Recall = c(round(cm_logit[["byClass"]][["Recall"]],3),round(cm_rf[["byClass"]][["Recall"]],3)),
F1_Score = c(round(cm_logit[["byClass"]][["F1"]],3),round(cm_rf[["byClass"]][["F1"]],3))
) %>%
kable(
caption = "<b style='display:block; text-align:center;'>Bảng 7: Các chỉ tiêu đánh giá kết quả dự báo</b>",
escape = FALSE, align = 'c'
) %>%
kable_styling(
position = "center", font_size = 13,
full_width = FALSE,
bootstrap_options = c("striped", "hover") #striped table
)| Model | Accuracy | Precision | Recall | F1_Score |
|---|---|---|---|---|
| Logistic | 0.825 | 0.801 | 0.748 | 0.773 |
| Random Forest | 0.904 | 0.922 | 0.832 | 0.874 |
Kết quả dự báo cho thấy mô hình Random Forest vượt trội so với mô hình Logistic trên tất cả các chỉ số đánh giá. Với Accuracy đạt 0,904, Random Forest dự đoán đúng 90,4% các trường hợp, cao hơn đáng kể so với Logistic (0,825). Về Precision, Random Forest (0,922) cho thấy 92,2% các trường hợp được dự đoán mắc bệnh tiểu đường (Diagnosis = 1) là đúng, so với 80,1% của Logistic, thể hiện khả năng giảm thiểu dự đoán sai tích cực tốt hơn. Recall của Random Forest (0,832) cao hơn Logistic (0,748), nghĩa là nó phát hiện được 83,2% các trường hợp thực sự mắc bệnh, đảm bảo ít bỏ sót hơn. F1-Score, kết hợp giữa Precision và Recall, của Random Forest (0,874) cũng vượt trội so với Logistic (0,773), cho thấy sự cân bằng tốt hơn giữa độ chính xác và khả năng phát hiện. Mặc dù mô hình Logistic có ưu điểm về khả năng diễn giải (qua odds ratio), Random Forest chứng minh hiệu quả dự đoán cao hơn nhờ khả năng xử lý các mối quan hệ phi tuyến tính và tương tác phức tạp giữa các biến như HbA1c, FastingBloodSugar, và FrequentUrination.
Trong nghiên cứu này, ba mô hình hồi quy nhị phân gồm logistic, probit và complementary log-log (cloglog) cùng với một mô hình học máy là Random Forest đã được sử dụng để phân loại nguy cơ mắc bệnh tiểu đường dựa trên các yếu tố: HbA1c, đường huyết lúc đói (FastingBloodSugar), huyết áp tâm thu (SystolicBP), hút thuốc (Smoking) và tiểu nhiều (FrequentUrination).
Kết quả cho thấy cả ba mô hình hồi quy tuyến tính tổng quát đều cho kết quả khá tương đồng về hệ số và chiều hướng tác động của các biến. Tuy nhiên, về mặt hiệu năng, logistic regression có chỉ số AIC thấp nhất, điều này cho thấy logistic phù hợp hơn trong ngữ cảnh dữ liệu hiện tại và là lựa chọn cân bằng giữa khả năng giải thích và hiệu quả dự báo. So với các mô hình này, Random Forest nổi bật hơn với AUC lên đến 0.93, phản ánh khả năng phân biệt tốt hơn nhờ khai thác được các mối quan hệ phi tuyến và tương tác phức tạp giữa các biến.
Xét về tác động của các biến đến khả năng mắc bệnh, kết quả hồi quy logistic cho thấy HbA1c có hệ số dương và ý nghĩa cao, xác nhận rằng chỉ số HbA1c càng cao thì xác suất mắc tiểu đường càng lớn. Tương tự, đường huyết lúc đói cũng cho kết quả tương tự, là một trong những biến quan trọng nhất trong cả hồi quy lẫn Random Forest. Tiểu nhiều cũng có mối liên hệ dương với xác suất mắc bệnh, nhưng mức độ ảnh hưởng thấp hơn. Ngược lại, huyết áp tâm thu và vấn đề hút thuốc không cho thấy ý nghĩa thống kê rõ ràng trong một số mô hình, gợi ý rằng biến này có thể không đóng vai trò quyết định trong dữ liệu hiện tại.
Tổng thể, mô hình Random Forest cho hiệu suất dự báo cao nhất, trong khi hồi quy logistic vẫn giữ ưu điểm về khả năng giải thích và tính minh bạch. Kết hợp cả hai phương pháp này trong thực hành lâm sàng và nghiên cứu dịch tễ học có thể giúp nâng cao khả năng phát hiện sớm và kiểm soát bệnh tiểu đường, đặc biệt trong các hệ thống y tế với nguồn lực hạn chế.
Dựa trên kết quả phân tích, các khuyến nghị sau được đề xuất:
Ứng dụng lâm sàng: Tích hợp HbA1c, FastingBloodSugar, và FrequentUrination vào các công cụ sàng lọc sớm bệnh tiểu đường trong thực hành y tế, do chúng là các yếu tố dự đoán mạnh. Mô hình Random Forest nên được ưu tiên trong các hệ thống hỗ trợ chẩn đoán tự động nhờ độ chính xác cao (Accuracy = 0,904).
Cải thiện mô hình: Xem xét bổ sung các biến khác như BMI hoặc FamilyHistoryDiabetes, vốn có liên quan mạnh đến bệnh tiểu đường trong các nghiên cứu trước (WHO, 2021), để tăng cường khả năng dự đoán, đặc biệt khi SystolicBP và Smoking cho thấy tác động yếu.
Phân tích sâu hơn: Thực hiện kiểm tra tương tác giữa các biến (ví dụ: HbA1c và FastingBloodSugar) và áp dụng kỹ thuật cân bằng mẫu để xử lý tỷ lệ FrequentUrination thấp (19,33%) nhằm cải thiện độ nhạy của mô hình.
Nghiên cứu tiếp theo: Mở rộng nghiên cứu với dữ liệu lớn hơn và đa dạng hơn về dân tộc hoặc khu vực địa lý để đảm bảo tính tổng quát hóa của mô hình. Đồng thời, đánh giá thêm các thuật toán học máy khác như Gradient Boosting để so sánh hiệu quả với Random Forest.
Giáo dục và phòng ngừa: Tăng cường giáo dục sức khỏe về kiểm soát đường huyết và nhận biết sớm triệu chứng đái nhiều, đồng thời khuyến khích bỏ hút thuốc để giảm nguy cơ, dù tác động của Smoking trong mẫu này không mạnh.
Trần Mạnh Tường (2025), Giáo trình Phân tích Dữ liệu Định tính.
Selvin, E., Steffes, M. W., Zhu, H., Matsushita, K., Wagenknecht, L., Pankow, J., … & Brancati, F. L. (2010). Glycated hemoglobin, diabetes, and cardiovascular risk in nondiabetic adults. New England Journal of Medicine, 362(9), 800–811.
Tirosh, A., Shai, I., Tekes-Manova, D., Israeli, E., Pereg, D., Shochat, T., … & Rudich, A. (2005). Normal fasting plasma glucose levels and type 2 diabetes in young men. New England Journal of Medicine, 353(14), 1454–1462.
Ferrannini, E., & Cushman, W. C. (2012). Diabetes and hypertension: The bad companions. Diabetes Care, 35(8), 1613–1620.
Willi, C., Bodenmann, P., Ghali, W. A., Faris, P. D., & Cornuz, J. (2007). Active smoking and the risk of type 2 diabetes: A systematic review and meta-analysis. JAMA, 298(22), 2654–2664.
Gupta, S., Gupta, V., & Gupta, A. (2018). Prevalence of undiagnosed diabetes in an urban population: Role of symptoms in early detection. Journal of Diabetes Research, 2018, 1626257.
American Diabetes Association. (2020). Standards of Medical Care in Diabetes—2020. Diabetes Care, 43(Supplement 1), S14–S31.