Tiểu luận môn học PTDLĐT
Nguyễn Châu Quý
2025-07-13
CHƯƠNG 1: TỔNG QUAN NGHIÊN CỨU
1.1 Lý do chọn đề tài
Trong bối cảnh nền kinh tế toàn cầu ngày càng bộc lộ những bất bình đẳng sâu sắc về phân phối thu nhập và cơ hội, việc nghiên cứu các yếu tố quyết định đến mức thu nhập cá nhân không chỉ mang ý nghĩa học thuật, mà còn có giá trị thực tiễn sâu sắc trong hoạch định chính sách và phân tích cấu trúc xã hội. Các công trình kinh tế – xã hội đương đại (Atkinson, 2015; Piketty, 2014) đã nhấn mạnh rằng thu nhập không đơn thuần là kết quả của năng lực lao động, mà còn là biểu hiện của sự tương tác phức tạp giữa các yếu tố cá nhân (giới tính, học vấn, nghề nghiệp, tình trạng hôn nhân), yếu tố thể chế (chính sách thị trường lao động, thuế, bình đẳng giới), và đặc điểm xã hội – nhân khẩu học.
Trong bức tranh đó, trình độ học vấn và số giờ làm việc mỗi tuần nổi lên như hai yếu tố nền tảng, có vai trò quyết định đến khả năng đạt được mức thu nhập đủ sống và vươn lên. Tuy nhiên, phần lớn các nghiên cứu trước đây chủ yếu tập trung vào sự khác biệt về thu nhập trung bình, ít công trình xem xét khả năng đạt một ngưỡng thu nhập cụ thể – chẳng hạn như mức $50.000/năm – vốn thường được xem là chỉ dấu để nhận diện tầng lớp trung lưu trong xã hội hiện đại. Điều này tạo ra một khoảng trống đáng kể trong việc hiểu biết xác suất kinh tế, từ đó ảnh hưởng đến tính chính xác trong hoạch định chính sách giáo dục, lao động và phúc lợi.
Để góp phần lấp đầy khoảng trống đó, nghiên cứu này lựa chọn khai thác bộ dữ liệu Adult Census Income của Hoa Kỳ – một tập dữ liệu dân số đã được chuẩn hóa và ứng dụng rộng rãi trong nhiều nghiên cứu học thuật có giá trị (Kohavi & Becker, 1996; Ding và cộng sự, 2019). Dữ liệu này không chỉ cung cấp thông tin định lượng phong phú về thu nhập, học vấn, thời gian làm việc mà còn cho phép xây dựng mô hình dự báo xác suất thu nhập một cách hệ thống và có kiểm định.
Đề tài: “Phân tích các yếu tố ảnh hưởng đến khả năng đạt được thu nhập cao” được lựa chọn không chỉ vì tính cấp thiết và phù hợp với bối cảnh nghiên cứu hiện nay, mà còn vì khả năng đóng góp về mặt phương pháp luận, dữ liệu thực chứng và hàm ý chính sách. Kỳ vọng của nghiên cứu là cung cấp một cái nhìn rõ ràng, định lượng và có cơ sở khoa học về vai trò của học vấn và thời gian lao động đối với cơ hội kinh tế, từ đó góp phần thúc đẩy thảo luận học thuật cũng như hoạch định chính sách xã hội theo hướng công bằng và hiệu quả hơn.
1.2 Giới thiệu dữ liệu
Bài nghiên cứu sử dụng bộ dữ liệu Adult Census Income, được trích xuất từ cuộc điều tra dân số năm 1994 do Cục Điều tra Dân số Hoa Kỳ thực hiện. Tập dữ liệu bao gồm 30162 quan sát với 15 biến đặc trưng, phản ánh các khía cạnh đa dạng về nhân khẩu học, kinh tế và xã hội của các cá nhân trưởng thành tại Hoa Kỳ.
Trong số 15 biến, có 9 biến định tính và 6 biến định lượng, cụ thể như sau:
Các biến định lượng bao gồm:
Age: Độ tuổi của cá nhân tham gia khảo sát.
Final weight: Trọng số mẫu, phản ánh mức độ đại diện của mỗi cá nhân trong tổng thể dân số.
Education-num: Mức số hóa của trình độ học vấn, biểu thị thứ bậc học vấn.
Capital-gain: Thu nhập từ lợi nhuận vốn trong năm.
Capital-loss: Khoản lỗ từ các khoản đầu tư tài chính.
Hours-per-week: Số giờ làm việc trung bình mỗi tuần.
Các biến định tính bao gồm:
Workclass: loại hình đơn vị làm việc của cá nhân.
Private: Doanh nghiệp tư nhân.
Self-emp-not-inc: Tự làm chủ nhưng không đăng ký doanh nghiệp (ví dụ: tiểu thương, lao động tự do).
Self-emp-inc: Tự làm chủ có đăng ký doanh nghiệp (công ty TNHH hoặc tương đương).
Federal-gov:Làm việc cho chính phủ liên bang.
Local-gov: Làm việc cho chính quyền địa phương.
State-gov: Làm việc cho chính quyền bang.
Without-pay: Làm việc không lương (ví dụ: tình nguyện, phụ giúp gia đình).
Never-worked: Chưa từng làm việc.
Education: Trình độ học vấn cao nhất mà cá nhân đạt được.
Preschool: Mẫu giáo.
1st-4th: Lớp 1-4.
5th-6th: Lớp 5-6.
7th-8th:Lớp 7-8.
9th: Lớp 9.
10th: Lớp 10.
11th: Lớp 11.
12th: Lớp 12.
HS-grad: Tốt nghiệp THPT.
Some-college: Học một số tín chỉ đại học.
Assoc-voc: Bằng cao đẳng nghề.
Assoc-acdm: Cao đẳng học thuật.
Bachelors: Bằng cử nhân.
Masters: Bằng thạc sĩ.
Prof-school: Bằng chuyên nghiệp (cao hơn thạc sĩ).
Doctorate: Bằng tiến sĩ.
Marital status: Tình trạng hôn nhân của cá nhân.
Married-civ-spouse: Đã kết hôn và sống chung với vợ/chồng.
Married-spouse-absent: Đã kết hôn nhưng không sống cùng vợ/chồng (có thể vì ly thân, xa cách).
Divorced: Đã ly hôn.
Separated: Đã ly thân.
Widowed: Góa vợ/chồng.
Never-married:Chưa từng kết hôn.
Occupation: Nghề nghiệp hiện tại của cá nhân.
Tech-support: Hỗ trợ kỹ thuật.
Craft-repair:Thợ thủ công, sửa chữa (cơ khí, điện, mộc, v.v.).
Other-service: Dịch vụ khác (phục vụ, trông trẻ, vệ sinh,…).
Sales: Bán hàng và tiếp thị.
Exec-managerial: Quản lý điều hành (CEO, giám đốc, trưởng phòng).
Prof-specialty: Chuyên gia (bác sĩ, kỹ sư, luật sư, giảng viên đại học).
Handlers-cleaners: Công nhân bốc vác, dọn dẹp.
Machine-op-inspct: Vận hành và kiểm tra máy móc.
Adm-clerical: Hành chính – văn phòng.
Farming-fishing: Nông – ngư nghiệp.
Transport-moving: Vận chuyển (tài xế, điều phối).
Priv-house-serv: Phục vụ trong nhà riêng (giúp việc, bảo mẫu).
Protective-serv: Dịch vụ bảo vệ (cảnh sát, lính cứu hỏa).
Armed-Forces: Quân đội.
Relationship: Mô tả vai trò của cá nhân trong mối quan hệ với các thành viên trong hộ gia đình, chẳng hạn như vợ/chồng, con cái, hay chủ hộ.
Husband: Chồng.
Wife: Vợ.
Own-child: Con ruột của chủ hộ.
Not-in-family: Không có quan hệ huyết thống với chủ hộ (người ở trọ, bạn bè,…).
Other-relative: Họ hàng khác (cháu, chú bác,…).
Unmarried: Người sống độc thân (không phải là vợ/chồng).
Race: Biểu thị chủng tộc của cá nhân.
White: Người da trắng.
Black: Người da đen.
Asian-Pac-Islander: Người gốc Á hoặc đảo Thái Bình Dương.
Amer-Indian-Eskimo: Người da đỏ bản địa hoặc người Eskimo.
Other: Chủng tộc khác.
Sex: Phản ánh giới tính của cá nhân. Bao gồm Female là nữ và Male là nam.
Native country: Quốc gia nơi cá nhân sinh ra.
Income: Thể hiện mức thu nhập hàng năm của cá nhân, phân loại thành hai nhóm: “<=50K” và “>50K”.
1.3 Mục tiêu nghiên cứu
Mục tiêu của bài nghiên cứu là phân tích các yếu tố ảnh hưởng đến khả năng đạt được thu nhập cao trong dân số Hoa Kỳ, cụ thể là vượt ngưỡng 50.000 USD/năm.
1.4 Phương pháp nghiên cứu
Để đạt được các mục tiêu trên, nghiên cứu sử dụng phương pháp định lượng với quy trình phân tích cụ thể như sau:
Nguồn dữ liệu: Dữ liệu được lấy từ bộ Adult Census Income Dataset (Hoa Kỳ), một tập dữ liệu điều tra dân số đã qua xử lý và được sử dụng phổ biến trong các nghiên cứu về thu nhập.
Phân tích mô tả: Trình bày phân phối thu nhập theo học vấn và giờ làm việc, Relative Risk (RR), Odds Ratio (OR), kiểm định Chi-bình phương để xác định mối liên hệ thống kê.
Xây dựng mô hình: Ước lượng các mô hình hồi quy nhị phân gồm: Logistic Regression (Logit), Probit Regression, Complementary Log-log Regression (Cloglog). Từ đó so sánh hệ số, ý nghĩa thống kê và khả năng dự báo (dựa trên độ chính xác, độ nhạy, độ đặc hiệu, Brier Score và AIC).
Đánh giá mô hình và kết luận: So sánh hiệu năng dự báo giữa các mô hình để xác định mô hình phù hợp nhất; rút ra kết luận về tác động của học vấn và giờ làm việc đến thu nhập.
1.5 Ý nghĩa đề tài
Về mặt lý luận, nghiên cứu này góp phần bổ sung vào kho tàng tri thức về kinh tế lao động, đặc biệt trong lĩnh vực phân tích thu nhập cá nhân theo tiếp cận xác suất. Thay vì chỉ dừng lại ở việc so sánh mức thu nhập trung bình giữa các nhóm dân cư, đề tài tiếp cận vấn đề từ góc độ xác suất đạt được một ngưỡng thu nhập cụ thể, qua đó phản ánh tốt hơn mức độ bất bình đẳng cơ hội trong xã hội. Việc ứng dụng các mô hình hồi quy nhị phân như Logit, Probit và Complementary Log-log cho phép lượng hóa tác động của trình độ học vấn và thời lượng lao động đối với thu nhập một cách chính xác, đồng thời mở rộng khả năng ứng dụng các công cụ thống kê hiện đại vào nghiên cứu xã hội học định lượng.
Về mặt thực tiễn, đề tài mang lại những gợi ý chính sách hữu ích cho các nhà hoạch định chiến lược giáo dục, lao động và an sinh xã hội. Việc xác định rằng học vấn cao và số giờ làm việc đầy đủ có tương quan mạnh mẽ với khả năng đạt thu nhập cao sẽ giúp định hướng các chính sách đào tạo nghề, nâng cao kỹ năng cho người lao động, cũng như xây dựng khung pháp lý bảo vệ người làm việc bán thời gian hoặc thiếu việc làm. Trong bối cảnh thị trường lao động đang chịu nhiều biến động do tự động hóa và chuyển đổi số, kết quả nghiên cứu có thể đóng vai trò tham chiếu quan trọng cho việc thiết kế các chương trình tái đào tạo và chuyển đổi nghề nghiệp phù hợp với từng nhóm đối tượng dân cư.
Ngoài ra, nghiên cứu cũng có thể hỗ trợ các tổ chức xã hội, các doanh nghiệp và các cơ sở giáo dục trong việc đưa ra các quyết định liên quan đến tuyển dụng, đào tạo nguồn nhân lực và đánh giá hiệu quả chính sách nhân sự dựa trên bằng chứng định lượng rõ ràng.
1.6 Kết cấu bài nghiên cứu
Chương 1: Tổng quan nghiên cứu.
Chương 2: Phân tích các yếu tố ảnh hưởng đến khả năng đạt được thu nhập cao.
Chương 3: Kết luận.
CHƯƠNG 2: PHÂN TÍCH CÁC YẾU TỐ ẢNH HƯỞNG ĐẾN KHẢ NĂNG ĐẠT ĐƯỢC THU NHẬP CAO
2.1 Phân tích ảnh hưởng của trình độ học vấn đến khả năng đạt được mức thu nhập cao
Trong nền kinh tế hiện đại, học vấn được xem là chìa khóa mở ra cơ hội việc làm tốt hơn và mức thu nhập cao hơn. Người có trình độ học vấn cao thường có khả năng đảm nhận các vị trí đòi hỏi kiến thức chuyên môn sâu, kỹ năng giải quyết vấn đề, tư duy phản biện và khả năng sáng tạo. Ngược lại, những người có trình độ học vấn thấp thường bị giới hạn trong việc tiếp cận các vị trí có mức lương cao hoặc nghề nghiệp ổn định. Do đó, việc phân tích mối liên hệ giữa trình độ học vấn và khả năng đạt được thu nhập cao là bước đầu tiên để kiểm chứng luận điểm này.
2.1.1 Bảng tần số - tần suất cho biến thu nhập
Đầu tiên, nghiên cứu tiến hành lập bảng tần số cho biến income nhằm mục tiêu mô tả sơ bộ phân bố của thu nhập trong toàn bộ mẫu khảo sát. Việc này cho phép nhận diện tỷ trọng tương đối giữa nhóm cá nhân có thu nhập vượt quá 50000 USD/năm và nhóm có thu nhập thấp hơn ngưỡng này, từ đó cung cấp cái nhìn tổng quan ban đầu về sự bất bình đẳng thu nhập trong tập dữ liệu.
library(dplyr)
library(knitr)
# Tải các thư viện cần thiết
library(dplyr)
library(knitr)
library(kableExtra)
# Tạo bảng tần số và tỷ lệ phần trăm cho biến 'income'
bang_tanso_income <- b %>%
group_by(income) %>%
summarise(
`Tần số` = n(),
`Tần suất (%)` = round(100 * n() / nrow(b), 2)
) %>%
arrange(desc(`Tần số`))
# Hiển thị bảng với định dạng đẹp
kable(
bang_tanso_income,
caption = "Bảng 2.1: Phân bố tần số và tần suất của biến thu nhập",
col.names = c("Mức thu nhập", "Tần số", "Tần suất (%)"),
align = "lcc",
format = "html"
) %>%
kable_styling(
bootstrap_options = c("striped", "hover", "condensed", "responsive"),
full_width = FALSE,
position = "center"
) %>%
row_spec(0, bold = TRUE, color = "white", background = "#0073C2FF") %>% # tiêu đề xanh đậm
column_spec(2:3, width = "5em") # căn chỉnh cột| Mức thu nhập | Tần số | Tần suất (%) |
|---|---|---|
| <=50K | 22654 | 75.11 |
| >50K | 7508 | 24.89 |
Để có cái nhìn trực quan về cách phân bố mức thu nhập trong tập dữ liệu, tác giả xây dựng biểu đồ cột thể hiện số lượng cá nhân thuộc từng nhóm thu nhập (<=50K và >50K). Biểu đồ không chỉ minh họa sự chênh lệch rõ rệt giữa hai nhóm mà còn cung cấp thông tin về tỷ lệ phần trăm của từng nhóm trong tổng thể dữ liệu, từ đó giúp hiểu rõ hơn về phân phối thu nhập trong mẫu khảo sát.
# Tải thư viện cần thiết
library(ggplot2)
library(dplyr)
# Tính toán số lượng và nhãn
bang_plot <- b %>%
count(income) %>%
mutate(
ty_le = round(100 * n / sum(n), 1),
nhan = paste0(n, " (", ty_le, "%)")
)
# Vẽ biểu đồ
ggplot(bang_plot, aes(x = income, y = n, fill = income)) +
geom_bar(stat = "identity", width = 0.5, show.legend = FALSE) +
geom_text(aes(label = nhan), vjust = -0.5, size = 4.5, fontface = "bold") +
scale_fill_manual(values = c("<=50K" = "#00AFBB", ">50K" = "#E7B800")) +
labs(
title = "Biểu đồ 2.1: Phân bố mức thu nhập",
x = "Mức thu nhập",
y = "Số lượng cá nhân"
) +
scale_y_continuous(expand = expansion(mult = c(0, 0.1))) +
theme_minimal(base_size = 14) +
theme(
plot.title = element_text(face = "bold", hjust = 0.5, size = 16),
axis.title.x = element_text(face = "bold"),
axis.title.y = element_text(face = "bold"),
axis.text = element_text(color = "black")
)
Kết quả bảng tần số và biểu đồ trên cho thấy có sự chênh lệch rõ rệt giữa hai nhóm thu nhập. Cụ thể, có tới 22654 cá nhân (chiếm 75.11%) có mức thu nhập không vượt quá 50000 USD/năm, trong khi chỉ có 7508 cá nhân (tương ứng 24.89%) có mức thu nhập trên 50.000 USD/năm. Điều này cho thấy phần lớn người dân trong mẫu khảo sát có mức thu nhập ở mức thấp. Phân bố thu nhập này phản ánh sự bất cân xứng và có thể là dấu hiệu cho thấy sự phân hóa về thu nhập trong xã hội tại thời điểm khảo sát.
2.1.2 Bảng tần số - tần suất cho biến trình độ học vấn
Trong nghiên cứu này, tác giả tiến hành phân tích biến trình độ học vấn (education) nhằm tìm hiểu mức độ phổ biến của các cấp bậc học vấn trong tập dữ liệu, đồng thời chuẩn bị cho các phân tích tiếp theo. Để thuận tiện cho việc phân tích và so sánh, trình độ học vấn của các cá nhân được phân thành hai nhóm chính:
Nhóm học vấn cao (Higher Education): Bao gồm các cá nhân có bằng cấp từ cử nhân (Bachelors) trở lên, bao gồm cả thạc sĩ (Masters), trường chuyên nghiệp (Prof-school) và tiến sĩ (Doctorate). Vì đây là những trình độ sau trung học phổ thông, có tính chuyên môn và học thuật cao, thường yêu cầu thời gian đào tạo dài và được công nhận rộng rãi trên thị trường lao động.
Nhóm học vấn thấp (Lower Education): Bao gồm các cá nhân có trình độ học vấn còn lại như tốt nghiệp phổ thông, học nghề, cao đẳng hoặc chưa hoàn thành bậc đại học.
Việc phân nhóm này được thực hiện dựa trên sự phân biệt rõ ràng giữa các bậc học phổ thông và giáo dục bậc cao, nhằm làm nổi bật sự khác biệt trong tiềm năng thu nhập và nghề nghiệp của hai nhóm.
Bảng tần số cho biến trình độ học vấn sẽ được trình bày sau đây.
# Nhóm trình độ học vấn thành 2 nhóm
b$edu_group <- ifelse(b$education %in% c("Bachelors", "Masters", "Doctorate", "Prof-school"),
"Higher", "Lower")
# Thư viện hỗ trợ bảng đẹp
library(dplyr)
library(knitr)
library(kableExtra)
# Tạo bảng tần số và tỷ lệ phần trăm cho biến 'edu_group'
bang_tanso_edu_group <- b %>%
group_by(edu_group) %>%
summarise(
`Tần số` = n(),
`Tần suất (%)` = round(100 * n() / nrow(b), 2)
) %>%
arrange(desc(`Tần số`))
# Hiển thị bảng với định dạng đẹp
kable(
bang_tanso_edu_group,
caption = "Bảng 2.2: Phân bố tần số và tần suất theo nhóm trình độ học vấn",
col.names = c("Nhóm trình độ học vấn", "Tần số", "Tần suất (%)"),
align = "lcc",
format = "html"
) %>%
kable_styling(
bootstrap_options = c("striped", "hover", "condensed", "responsive"),
full_width = FALSE,
position = "center"
) %>%
row_spec(0, bold = TRUE, color = "white", background = "#0073C2FF") %>%
column_spec(2:3, width = "5em")| Nhóm trình độ học vấn | Tần số | Tần suất (%) |
|---|---|---|
| Lower | 22574 | 74.84 |
| Higher | 7588 | 25.16 |
Để trực quan hóa số liệu từ bảng phân bố tần số và tần suất theo nhóm trình độ học vấn, biểu đồ dưới đây thể hiện số lượng cá nhân ở mỗi nhóm học vấn (Higher Education và Lower Education). Việc sử dụng biểu đồ cột giúp người đọc dễ dàng nhận thấy sự chênh lệch giữa hai nhóm, từ đó có cái nhìn rõ hơn về cấu trúc học vấn của mẫu khảo sát. Biểu đồ không chỉ hỗ trợ việc so sánh trực quan mà còn làm nổi bật vai trò của trình độ học vấn trong các phân tích tiếp theo.
# Thư viện trực quan hóa
library(ggplot2)
# Chuẩn bị dữ liệu biểu đồ từ bảng đã tạo
bang_plot_edu_group <- b %>%
count(edu_group) %>%
mutate(
ty_le = round(100 * n / sum(n), 1),
nhan = paste0(n, " (", ty_le, "%)")
)
# Vẽ biểu đồ cột
ggplot(bang_plot_edu_group, aes(x = edu_group, y = n, fill = edu_group)) +
geom_bar(stat = "identity", width = 0.5, show.legend = FALSE) +
geom_text(aes(label = nhan), vjust = -0.5, size = 4.5, fontface = "bold") +
scale_fill_manual(values = c("Higher" = "#00AFBB", "Lower" = "#E7B800")) +
labs(
title = "Biểu đồ 2.2: Phân bố nhóm trình độ học vấn",
x = "Nhóm trình độ học vấn",
y = "Số lượng cá nhân"
) +
scale_y_continuous(expand = expansion(mult = c(0, 0.1))) +
theme_minimal(base_size = 14) +
theme(
plot.title = element_text(face = "bold", hjust = 0.5, size = 16),
axis.title = element_text(face = "bold"),
axis.text = element_text(color = "black")
)
Dựa vào bảng và đồ thị trên, có thể thấy rằng phần lớn cá nhân trong bộ dữ liệu thuộc nhóm trình độ học vấn thấp với 22574 người, chiếm 74.84% tổng mẫu. Trong khi đó, nhóm học vấn cao chỉ chiếm 25.16% tương ứng với 7588 người.
Điều này cho thấy sự phân bổ không đồng đều giữa hai nhóm học vấn, và gợi ý rằng nhóm có trình độ học vấn thấp đang chiếm ưu thế trong dữ liệu.
2.1.3 Thống kê mô tả cho hai biến thu nhập và trình độ học vấn
Để có cái nhìn tổng quan ban đầu về mối quan hệ giữa trình độ học vấn và mức thu nhập, bảng tần suất chéo giữa hai biến income (thu nhập) và edu_group (nhóm trình độ học vấn) đã được xây dựng.
Bảng dưới đây trình bày tần số và tần suất (tính theo %) của các nhóm kết hợp giữa học vấn và thu nhập.
# Tạo bảng tần suất chéo giữa 'income' và 'edu_group'
bang_cheo <- b %>%
group_by(edu_group, income) %>%
summarise(
`Tần số` = n(),
.groups = "drop"
) %>%
group_by(edu_group) %>%
mutate(`Tần suất (%)` = round(100 * `Tần số` / sum(`Tần số`), 2)) %>%
ungroup()
# Hiển thị bảng chéo với định dạng đẹp
library(knitr)
library(kableExtra)
kable(
bang_cheo,
caption = "Bảng 2.3: Bảng tần suất chéo giữa mức thu nhập và nhóm trình độ học vấn",
col.names = c("Nhóm trình độ học vấn", "Mức thu nhập", "Tần số", "Tần suất (%)"),
align = "lccc",
format = "html"
) %>%
kable_styling(
bootstrap_options = c("striped", "hover", "condensed", "responsive"),
full_width = FALSE,
position = "center"
) %>%
row_spec(0, bold = TRUE, color = "white", background = "#0073C2FF") %>%
column_spec(3:4, width = "6em")| Nhóm trình độ học vấn | Mức thu nhập | Tần số | Tần suất (%) |
|---|---|---|---|
| Higher | <=50K | 3858 | 50.84 |
| Higher | >50K | 3730 | 49.16 |
| Lower | <=50K | 18796 | 83.26 |
| Lower | >50K | 3778 | 16.74 |
Nhằm minh họa trực quan kết quả từ bảng tần suất chéo, biểu đồ dưới đây trình bày số lượng cá nhân ở từng nhóm trình độ học vấn.
# Thư viện vẽ biểu đồ
library(ggplot2)
# Vẽ biểu đồ cột cụm với số in đậm
ggplot(bang_cheo, aes(x = edu_group, y = `Tần số`, fill = income)) +
geom_bar(stat = "identity", position = "dodge", width = 0.6) +
geom_text(aes(label = `Tần số`),
position = position_dodge(width = 0.6),
vjust = -0.3, size = 3.8, fontface = "bold") + # Thêm fontface
scale_fill_manual(values = c("<=50K" = "#0073C2FF", ">50K" = "#EFC000FF")) +
labs(
title = "Biểu đồ 2.3: Phân bố mức thu nhập theo nhóm trình độ học vấn",
x = "Nhóm trình độ học vấn",
y = "Số lượng cá nhân",
fill = "Mức thu nhập"
) +
theme_minimal(base_size = 13) +
theme(
plot.title = element_text(face = "bold", hjust = 0.5, size = 14),
axis.title.x = element_text(face = "bold"),
axis.title.y = element_text(face = "bold")
)
Kết quả phân tích mối quan hệ giữa trình độ học vấn và mức thu nhập cho thấy sự khác biệt rõ rệt giữa hai nhóm học vấn. Cụ thể, trong nhóm có trình độ học vấn cao (Higher Education), tỷ lệ người có thu nhập trên 50K chiếm 49.16%, gần tương đương với tỷ lệ người có thu nhập dưới hoặc bằng 50K (50.84%). Điều này cho thấy phân bố thu nhập ở nhóm học vấn cao tương đối cân đối. Ngược lại, ở nhóm có trình độ học vấn thấp (Lower Education), phần lớn cá nhân (83.26%) có thu nhập dưới hoặc bằng 50K, trong khi chỉ có 16.74% đạt mức thu nhập cao hơn.
2.1.4 Relative Risk và khoảng ước lượng cho Relative Risk
Trong phần này, tác giả sẽ tiến hành tính toán Relative Risk (RR) nhằm đánh giá mức độ ảnh hưởng của từng nhóm trình độ học vấn đến xác suất đạt được mức thu nhập trên 50K.
library(epitools)
# Tạo bảng chéo giữa trình độ học vấn và thu nhập
table_edu_income <- table(b$edu_group, b$income)
# Đảm bảo biến theo đúng thứ tự: [Higher, Lower] x [<=50K, >50K]
# Tính Relative Risk
rr_result <- riskratio(table_edu_income)
# In kết quả
print(rr_result)## $data
##
## <=50K >50K Total
## Higher 3858 3730 7588
## Lower 18796 3778 22574
## Total 22654 7508 30162
##
## $measure
## risk ratio with 95% C.I.
## estimate lower upper
## Higher 1.0000000 NA NA
## Lower 0.3404646 0.328092 0.3533037
##
## $p.value
## two-sided
## midp.exact fisher.exact chi.square
## Higher NA NA NA
## Lower 0 0 0
##
## $correction
## [1] FALSE
##
## attr(,"method")
## [1] "Unconditional MLE & normal approximation (Wald) CI"Kết quả tính toán cho thấy Relative Risk (RR) của nhóm Lower (trình độ học vấn thấp) so với nhóm Higher (trình độ học vấn cao) là 0.3405, nhỏ hơn 1. Điều này có nghĩa là: những người có trình độ học vấn thấp chỉ có khả năng đạt mức thu nhập cao (>50K) bằng khoảng 34.1% so với những người có trình độ học vấn cao.
Nói cách khác, những người có trình độ học vấn cao có khả năng đạt thu nhập cao cao hơn gần 3 lần so với nhóm học vấn thấp.
Khoảng tin cậy 95% cho RR nằm trong khoảng từ 0.3281 đến 0.3533, nghĩa là với độ tin cậy 95%, tỉ lệ nguy cơ thực sự đạt thu nhập cao của nhóm học vấn thấp nằm trong khoảng này.
2.1.5 Odd Ratio và khoảng ước lượng cho Odds ratio
Để đánh giá mối liên hệ giữa trình độ học vấn và khả năng đạt được mức thu nhập cao (>50K), phân tích Odds Ratio được thực hiện.
library(epitools)
# Tính Odds Ratio
or_result <- oddsratio(table_edu_income)
# In kết quả
print(or_result)## $data
##
## <=50K >50K Total
## Higher 3858 3730 7588
## Lower 18796 3778 22574
## Total 22654 7508 30162
##
## $measure
## odds ratio with 95% C.I.
## estimate lower upper
## Higher 1.0000000 NA NA
## Lower 0.2078888 0.1963902 0.2200951
##
## $p.value
## two-sided
## midp.exact fisher.exact chi.square
## Higher NA NA NA
## Lower 0 0 0
##
## $correction
## [1] FALSE
##
## attr(,"method")
## [1] "median-unbiased estimate & mid-p exact CI"Giá trị Odds Ratio = 0.2079, nghĩa là odds đạt thu nhập cao (tức là tỷ lệ giữa số người có thu nhập >50K so với số người ≤50K) trong nhóm có trình độ học vấn thấp chỉ bằng khoảng 20.8% so với odds tương ứng ở nhóm có trình độ học vấn cao.
Khoảng tin cậy 95% của Odds Ratio nằm trong khoảng từ 0.1964 đến 0.2201, điều này cho thấy rằng, với độ tin cậy 95%, odds đạt thu nhập cao ở nhóm học vấn thấp chỉ dao động từ khoảng 19.6% đến 22.0% so với nhóm học vấn cao.
2.1.6 Thống kê suy diễn
2.1.6.1 Kiểm định tính độc lập
Để kiểm tra xem liệu trình độ học vấn của cá nhân có mối quan hệ thống kê với khả năng đạt được mức thu nhập cao hay không, nghiên cứu tiến hành áp dụng kiểm định Chi-bình phương về tính độc lập giữa hai biến phân loại: trình độ học vấn và thu nhập. Kiểm định này cho phép đánh giá giả thuyết liệu sự khác biệt về thu nhập giữa các nhóm trình độ học vấn có xảy ra một cách ngẫu nhiên hay thực sự phản ánh một mối liên hệ có ý nghĩa thống kê trong tổng thể.
Đặt giả thuyết:
\[ \left\{ \begin{array}{ll} H_0: \text{Trình độ học vấn của cá nhân không có mối quan hệ với mức thu nhập.} \\ H_1: \text{Trình độ học vấn của cá nhân có mối quan hệ với mức thu nhập.} \end{array} \right. \]
# Kiểm định Chi bình phương
chi_result <- chisq.test(table_edu_income)
print(chi_result)##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: table_edu_income
## X-squared = 3191, df = 1, p-value < 2.2e-16Vì p-value < 2.2e-16 < 0.05, ta bác bỏ giả thuyết \(H_0\). Điều này có nghĩa là trình độ học vấn của cá nhân có mối quan hệ với mức thu nhập.
2.1.6.2 So sánh tỷ lệ thu nhập >50K giữa 2 nhóm trình độ học vấn
Tiếp theo tác giả tiến hành so sánh cụ thể tỷ lệ thu nhập >50K giữa hai nhóm trình độ học vấn tiêu biểu.
- Giả thuyết \(H_0\): Tỷ lệ thu nhập cao của nhóm học vấn cao không lớn hơn nhóm học vấn thấp.
\[ H_0: p_{Higher} ≤ p_{Lower} \]
- Giả thuyết \(H_1\): Tỷ lệ thu nhập cao của nhóm học vấn cao lớn hơn nhóm học vấn thấp.
\[ H_1: p_{Higher} > p_{Lower} \]
Trong đó:
\(p_{Lower}\) là tỷ lệ người có thu nhập cao trong nhóm có trình độ học vấn thấp.
\(p_{Higher}\) là tỷ lệ người có thu nhập cao trong nhóm có trình độ học vấn cao.
# Tính số người có thu nhập cao theo từng nhóm
success <- table_edu_income[, ">50K"]
# Tổng số người trong từng nhóm
n_total <- rowSums(table_edu_income)
# Thực hiện kiểm định hiệu hai tỷ lệ (prop.test chỉ hỗ trợ hai phía, dùng binom.test nếu một phía)
prop.test(success, n_total, alternative = "greater", correct = TRUE)##
## 2-sample test for equality of proportions with continuity correction
##
## data: success out of n_total
## X-squared = 3191, df = 1, p-value < 2.2e-16
## alternative hypothesis: greater
## 95 percent confidence interval:
## 0.3138303 1.0000000
## sample estimates:
## prop 1 prop 2
## 0.4915656 0.1673607Kết quả kiểm định cho thấy tỷ lệ người có thu nhập cao ở nhóm học vấn cao là 49.16%, trong khi nhóm học vấn thấp chỉ là 16.74%. Hiệu tỷ lệ ước lượng là 32.42%, với khoảng tin cậy 95% từ 31.38% đến 100%. Đồng thời, giá trị p của kiểm định là nhỏ hơn 2.2e-16, thấp hơn nhiều so với mức ý nghĩa 5%. Do đó, tác giả bác bỏ giả thuyết \(H_0\) và kết luận rằng tỷ lệ người có thu nhập cao ở nhóm học vấn cao lớn hơn một cách có ý nghĩa thống kê so với nhóm học vấn thấp.
2.1.7 Mô hình hồi quy cho dữ liệu nhị phân
2.1.7.1 Mô hình logit
Sau khi thực hiện phân tích mô tả và kiểm định mối liên hệ giữa trình độ học vấn và mức thu nhập, để đánh giá sâu hơn ảnh hưởng của học vấn đến xác suất có thu nhập cao (>50K), tác giả tiến hành xây dựng mô hình hồi quy Logit.
Trong mô hình này, biến phụ thuộc là một biến nhị phân income_binary, nhận giá trị:
1 nếu thu nhập > 50K.
0 nếu thu nhập ≤ 50K.
# Tạo biến thu nhập nhị phân
b$income_binary <- ifelse(b$income == ">50K", 1, 0)
# Biến edu_group đã tạo từ trước với 2 nhóm: Higher và Lower
b$edu_group <- factor(b$edu_group, levels = c("Higher", "Lower")) # Đặt Higher làm nhóm tham chiếu
# Mô hình logit với biến edu_group
model_logit <- glm(income_binary ~ edu_group, data = b, family = binomial(link = "logit"))
summary(model_logit)##
## Call:
## glm(formula = income_binary ~ edu_group, family = binomial(link = "logit"),
## data = b)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.03374 0.02296 -1.469 0.142
## edu_groupLower -1.57071 0.02907 -54.028 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 33851 on 30161 degrees of freedom
## Residual deviance: 30909 on 30160 degrees of freedom
## AIC: 30913
##
## Number of Fisher Scoring iterations: 3Ta có kết từ hàm hồi quy thông qua mô hình logit như sau:
\[ \ln\left(\frac{\hat{\pi}}{1 - \hat{\pi}}\right) = -0.03374 - 1.57071 \cdot \text{edu_group}_{\text{Lower}} \]
Hệ số chặn (Intercept) \(\hat{\beta}_0\) là -0.03374, đại diện cho log-odds đạt mức thu nhập cao (tức thu nhập >50K) của nhóm có trình độ học vấn cao. Tuy nhiên, hệ số này không có ý nghĩa thống kê vì p_value = 0.142 > 0.05.
Đối với hệ số \(\hat{\beta}_1\) = -1.57071, thể hiện chênh lệch log-odds đạt mức thu nhập cao giữa nhóm có trình độ học vấn thấp so với nhóm có trình độ học vấn cao là −1.57071. Nói cách khác, nhóm học vấn thấp có log-odds đạt thu nhập cao thấp hơn 1.57071 đơn vị so với nhóm học vấn cao.
Nếu xét theo tỷ lệ odds (odds ratio), thì odds đạt thu nhập cao của nhóm có học vấn thấp chỉ bằng khoảng 20.8% so với nhóm có học vấn cao. Điều này có nghĩa là khả năng đạt thu nhập cao của nhóm học vấn thấp thấp hơn khoảng 4.8 lần so với nhóm học vấn cao.
Hệ số này có ý nghĩa thống kê rất mạnh vì p_value < 2e−16 < 0.05, cho thấy mối liên hệ giữa trình độ học vấn và khả năng đạt thu nhập cao là đáng tin cậy về mặt thống kê.
Từ mô hình logit ta dự báo xác suất có thu nhập cao tương ứng với mỗi nhóm học vấn như sau:
library(knitr)
library(kableExtra)
# Tạo dữ liệu mới theo nhóm học vấn
new_data <- data.frame(edu_group = factor(c("Higher", "Lower"), levels = c("Higher", "Lower")))
# Dự báo xác suất từ mô hình logit
new_data$predicted_prob <- predict(model_logit, newdata = new_data, type = "response")
# Hiển thị bảng đẹp với màu xanh cho tiêu đề và căn chỉnh cột
kable(new_data,
col.names = c("Nhóm trình độ học vấn", "Xác suất thu nhập cao (>50K)"),
digits = 4,
align = "c",
caption = "Bảng 2.4: Dự báo xác suất có thu nhập cao theo nhóm trình độ học vấn") %>%
kable_styling(full_width = FALSE, position = "center") %>%
row_spec(0, bold = TRUE, color = "white", background = "#0073C2FF") %>%
column_spec(2, width = "10em")| Nhóm trình độ học vấn | Xác suất thu nhập cao (>50K) |
|---|---|
| Higher | 0.4916 |
| Lower | 0.1674 |
Bảng dự báo xác suất có thu nhập cao theo nhóm trình độ học vấn được trực quan hóa bằng biểu đồ sau:
library(ggplot2)
# Tạo nhãn phần trăm làm chú thích
new_data$label <- paste0(new_data$edu_group, "\n(", round(new_data$predicted_prob * 100, 1), "%)")
# Vẽ biểu đồ tròn
ggplot(new_data, aes(x = "", y = predicted_prob, fill = edu_group)) +
geom_col(width = 1, color = "white") +
coord_polar("y") +
geom_text(aes(label = label),
position = position_stack(vjust = 0.5),
size = 5, color = "black", fontface = "bold") +
scale_fill_manual(values = c("#0099CC", "#FF9966")) +
theme_void() +
ggtitle("Biểu đồ 2.4: Xác suất thu nhập cao theo nhóm học vấn") +
theme(plot.title = element_text(hjust = 0.5, size = 16, face = "bold"),
legend.title = element_blank(),
legend.position = "bottom")
Kết quả cho thấy:
Đối với nhóm có trình độ học vấn cao (Higher), xác suất đạt thu nhập cao vào khoảng 49.16%. Điều này có nghĩa là, một cá nhân có trình độ học vấn cao có khoảng 49.16% khả năng đạt được mức thu nhập trung trở lên.
Trong khi đó, nhóm có trình độ học vấn thấp (Lower) chỉ có xác suất khoảng 16.74% đạt mức thu nhập cao. Như vậy, cơ hội đạt thu nhập cao của nhóm này thấp hơn so với nhóm học vấn cao.
Sự khác biệt rõ rệt về xác suất này cho thấy rằng trình độ học vấn đóng vai trò quan trọng trong việc quyết định khả năng đạt được mức thu nhập cao. Nhóm có trình độ cao gần như gấp ba lần cơ hội đạt thu nhập tốt so với nhóm học vấn thấp.
2.1.7.2 Mô hình probit
Trong phần này, tác giả sử dụng mô hình probit để đánh giá mối quan hệ giữa trình độ học vấn và xác suất đạt được mức thu nhập cao
# Mô hình Probit với biến edu_group
model_probit <- glm(income_binary ~ edu_group, data = b, family = binomial(link = "probit"))
summary(model_probit)##
## Call:
## glm(formula = income_binary ~ edu_group, family = binomial(link = "probit"),
## data = b)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.02114 0.01439 -1.469 0.142
## edu_groupLower -0.94350 0.01748 -53.989 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 33851 on 30161 degrees of freedom
## Residual deviance: 30909 on 30160 degrees of freedom
## AIC: 30913
##
## Number of Fisher Scoring iterations: 3Mô hình hồi quy Probit được ước lượng như sau:
\[ \hat{\pi} = \Phi(-0.02114 - 0.94350 \cdot \text{edu_group}_{\text{Lower}}) \]
Hệ số chặn \(\hat{\beta}_0\) = -0.02114, phản ánh z-score tương ứng với xác suất đạt thu nhập cao của nhóm có trình độ học vấn cao. Xác suất tương ứng được tính là \(\Phi(−0.02114)\) ≈ 0.4916, tức khoảng 49.2%. Tuy nhiên, hệ số này không có ý nghĩa thống kê (p_value = 0.142 > 0.05).
Hệ số \(\hat{\beta}_1\) = -0.94350 phản ánh rằng nhóm có trình độ học vấn thấp có z-score thấp hơn nhóm học vấn cao là 0.94350 đơn vị. Hệ số này có ý nghĩa thống kê rất cao (p_value < 2e−16).
Ta có thể tính xác suất đạt thu nhập cao ở nhóm học vấn thấp bằng: \(\Phi( -0.02114 − 0.94350\) ≈ 0.1670, tức là chỉ khoảng 16.7% người có trình độ học vấn thấp đạt được mức thu nhập cao.
Từ mô hình Probit, ta dự báo xác suất đạt mức thu nhập cao tương ứng với từng nhóm học vấn như sau:
library(knitr)
library(kableExtra)
# Tạo dữ liệu mới theo nhóm học vấn
new_data <- data.frame(edu_group = factor(c("Higher", "Lower"), levels = c("Higher", "Lower")))
# Dự báo xác suất từ mô hình logit
new_data$predicted_prob <- predict(model_logit, newdata = new_data, type = "response")
# Hiển thị bảng đẹp với tiêu đề màu xanh và cột thứ 2 rộng 10em
kable(new_data,
col.names = c("Nhóm trình độ học vấn", "Xác suất thu nhập cao (>50K)"),
digits = 4,
align = "c",
caption = "Bảng 2.5: Dự báo xác suất có thu nhập cao theo nhóm trình độ học vấn") %>%
kable_styling(full_width = FALSE, position = "center") %>%
row_spec(0, bold = TRUE, color = "white", background = "#0073C2FF") %>%
column_spec(2, width = "10em")| Nhóm trình độ học vấn | Xác suất thu nhập cao (>50K) |
|---|---|
| Higher | 0.4916 |
| Lower | 0.1674 |
Dự báo xác suất có thu nhập cao theo nhóm trình độ học vấn từ mô hình probit được trực quan hóa như sau:
# Cài đặt và nạp thư viện
library(ggplot2)
# Tạo bộ dữ liệu mới gồm 2 nhóm trình độ học vấn
new_data <- data.frame(
edu_group = factor(c("Higher", "Lower"), levels = c("Higher", "Lower"))
)
# Dự báo xác suất từ mô hình Probit đã huấn luyện trước đó
new_data$predicted_prob <- predict(model_probit, newdata = new_data, type = "response")
# Tạo nhãn phần trăm cho chú thích
new_data$label <- paste0(new_data$edu_group, "\n(", round(new_data$predicted_prob * 100, 1), "%)")
# Vẽ biểu đồ tròn
ggplot(new_data, aes(x = "", y = predicted_prob, fill = edu_group)) +
geom_col(width = 1, color = "white") +
coord_polar("y") +
geom_text(aes(label = label),
position = position_stack(vjust = 0.5),
size = 5, color = "black", fontface = "bold") +
scale_fill_manual(values = c("#0099CC", "#FF9966")) + # phối màu đẹp, rõ ràng
theme_void() +
ggtitle("Biểu đồ 2.5: Dự báo xác suất có thu nhập cao theo nhóm trình độ học vấn") +
theme(plot.title = element_text(hjust = 0.5, size = 16, face = "bold"),
legend.title = element_blank(),
legend.position = "bottom")
Kết quả cho thấy:
Nhóm có trình độ học vấn cao (Higher): Xác suất đạt mức thu nhập cao ước tính là khoảng 49.16%. Điều này có nghĩa là, một cá nhân thuộc nhóm học vấn cao có khoảng 49.16% khả năng đạt được mức thu nhập cao.
Nhóm có trình độ học vấn thấp (Lower): Xác suất đạt mức thu nhập cao chỉ vào khoảng 16.74%. Như vậy, khả năng đạt thu nhập cao của nhóm học vấn thấp là thấp hơn đáng kể so với nhóm học vấn cao.
Vậy, trình độ học vấn đóng vai trò quan trọng trong việc xác định khả năng đạt được mức thu nhập cao. Nhóm học vấn cao có khả năng đạt thu nhập cao gần gấp 3 lần so với nhóm học vấn thấp.
2.1.7.2 Mô hình cloglog
# Chạy mô hình cloglog
model_cloglog <- glm(income_binary ~ edu_group,
data = b,
family = binomial(link = "cloglog"))
# Tóm tắt kết quả mô hình
summary(model_cloglog)##
## Call:
## glm(formula = income_binary ~ edu_group, family = binomial(link = "cloglog"),
## data = b)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.39094 0.01669 -23.43 <2e-16 ***
## edu_groupLower -1.30648 0.02332 -56.02 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 33851 on 30161 degrees of freedom
## Residual deviance: 30909 on 30160 degrees of freedom
## AIC: 30913
##
## Number of Fisher Scoring iterations: 5Dựa trên kết quả mô hình hồi quy sử dụng liên kết complementary log-log (cloglog), phương trình ước lượng có dạng như sau:
\[ \log\left(-\log(1 - \hat{\pi})\right) = -0.39094 - 1.30648 \cdot \text{edu_group}_{\text{Lower}} \]
Hệ số chặn \(\hat{\beta}_0\) = -0.39094 đại diện cho giá trị log–log của xác suất đạt thu nhập cao đối với nhóm trình độ học vấn cao. Từ đó, ta tính được xác suất:
\(\hat{\pi}_{\text{Higher}} = 1 - \exp(-\exp(-0.39094)) \approx 0.7686\)
Nói cách khác, khoảng 76.9% những người thuộc nhóm có trình độ học vấn cao có khả năng đạt được mức thu nhập cao.
Hệ số \(\hat{\beta}_1\) = -1.30648 thể hiện sự chênh lệch về log(-log(1 - p)) giữa nhóm có trình độ học vấn thấp và nhóm học vấn cao. Cụ thể, giá trị âm cho thấy nhóm có học vấn thấp có khả năng đạt thu nhập cao thấp hơn so với nhóm học vấn cao.
Ta có thể ước tính xác suất đạt thu nhập cao cho nhóm học vấn thấp như sau:
\(\hat{\pi}_{\text{Lower}} = 1 - \exp(-\exp(-0.39094 - 1.30648)) \approx 0.1676\)
Điều này cho thấy chỉ khoảng 16.8% những người thuộc nhóm có học vấn thấp có khả năng đạt mức thu nhập cao – thấp hơn so với nhóm học vấn cao.
Với mức ý nghĩa 5%, cả hai hệ số đều có ý nghĩa thống kê rất cao vì p_value < 2e−16 < 0.05. Điều này chứng tỏ trình độ học vấn là yếu tố ảnh hưởng có ý nghĩa thống kê đến khả năng đạt thu nhập cao.
2.1.8 Đánh giá mô hình
2.1.8.1 AIC
Để đánh giá mức độ phù hợp của mô hình hồi quy Probit sử dụng trong nghiên cứu, chỉ số Akaike Information Criterion (AIC) được sử dụng như một tiêu chí quan trọng. AIC phản ánh sự đánh đổi giữa độ phù hợp của mô hình và độ phức tạp của nó. Cụ thể, mô hình có AIC càng thấp thì càng được xem là có hiệu quả hơn trong việc cân bằng giữa khả năng giải thích dữ liệu và tránh tình trạng overfitting.
model_logit <- glm(income_binary ~ edu_group, data = b, family = binomial(link = "logit"))
model_probit <- glm(income_binary ~ edu_group, data = b, family = binomial(link = "probit"))
model_cloglog <- glm(income_binary ~ edu_group,
data = b,
family = binomial(link = "cloglog"))library(knitr)
library(kableExtra)
# Tạo bảng so sánh AIC
aic_compare <- data.frame(
Model = c("Logit", "Probit", "Cloglog"),
AIC = c(
AIC(model_logit),
AIC(model_probit),
AIC(model_cloglog)
)
)
# Làm tròn giá trị AIC để dễ nhìn
aic_compare$AIC <- round(aic_compare$AIC, 2)
# Hiển thị bảng đẹp với màu xanh header
kable(aic_compare,
col.names = c("Mô hình", "Chỉ số AIC"),
align = "c",
caption = "Bảng 2.6. So sánh chỉ số AIC giữa các mô hình nhị phân") %>%
kable_styling(full_width = FALSE, position = "center", bootstrap_options = c("striped", "hover")) %>%
row_spec(0, bold = TRUE, color = "white", background = "#0073C2FF") %>% # màu xanh đậm cho dòng tiêu đề
column_spec(1, bold = TRUE) %>% # tô đậm cột mô hình
column_spec(2, width = "10em") # căn chỉnh chiều rộng cột AIC| Mô hình | Chỉ số AIC |
|---|---|
| Logit | 30913.33 |
| Probit | 30913.33 |
| Cloglog | 30913.33 |
Cả ba mô hình Logit, Probit và Cloglog đều cho cùng một giá trị AIC là 30913.33. Điều này cho thấy mức độ phù hợp mô hình là tương đương giữa các hàm liên kết, và không có mô hình nào thể hiện ưu thế rõ rệt về mặt thống kê trong bối cảnh dữ liệu hiện tại.
2.1.7.2 Hệ số Brier Score
Chỉ số Brier Score được sử dụng để đánh giá mức độ chính xác của các dự đoán xác suất, thông qua việc đo lường sai số bình phương giữa xác suất dự báo và kết quả thực tế. Giá trị Brier Score càng thấp cho thấy mô hình dự báo càng chính xác.
# Gọi thư viện
library(DescTools)
library(knitr)
library(kableExtra)
# Tính Brier Score cho từng mô hình
brier_logit <- BrierScore(model_logit)
brier_probit <- BrierScore(model_probit)
brier_cloglog <- BrierScore(model_cloglog)
# Tạo bảng kết quả
brier_table <- data.frame(
Mô_hình = c("Logit", "Probit", "Cloglog"),
Brier_Score = round(c(brier_logit, brier_probit, brier_cloglog), 5)
)
# Trình bày bảng đẹp
kable(brier_table,
col.names = c("Mô hình", "Chỉ số Brier Score"),
align = "c",
caption = "Bảng 2.7: So sánh chỉ số Brier Score giữa các mô hình nhị phân") %>%
kable_styling(full_width = FALSE, position = "center") %>%
row_spec(0, bold = TRUE, color = "white", background = "#0073C2FF") %>%
column_spec(2, width = "10em")| Mô hình | Chỉ số Brier Score |
|---|---|
| Logit | 0.16717 |
| Probit | 0.16717 |
| Cloglog | 0.16717 |
Trong phân tích này, cả ba mô hình Logit, Probit và Cloglog đều đạt cùng giá trị Brier Score là 0.16717. Điều này cho thấy độ chính xác trong dự báo xác suất của các mô hình là tương đương nhau. Không có mô hình nào nổi bật hơn về hiệu suất dự đoán, phản ánh rằng chúng đều hoạt động ở mức độ tương tự trong việc ước lượng khả năng một cá nhân có thu nhập cao.
2.1.7.3 Ma trận nhầm lẫn
Để đánh giá khả năng của mô hình trong việc phân biệt đúng những cá nhân có thu nhập cao, nghiên cứu sử dụng ma trận nhầm lẫn.
a. Ma trận nhầm lẫn cho mô hình logit
prob_pred <- predict(model_logit, type = "response")
#Tìm ngưỡng tối ưu theo chỉ số Youden
library(pROC)
roc_obj <- roc(b$income_binary, prob_pred)
opt_threshold <- coords(roc_obj, "best", ret = "threshold", best.method = "youden")
threshold_value <- as.numeric(opt_threshold[1])
cat("Ngưỡng tối ưu là:", round(threshold_value, 4), "\n")## Ngưỡng tối ưu là: 0.3295# Dự đoán phân loại theo ngưỡng tối ưu
pred_class <- ifelse(prob_pred > threshold_value, "TRUE", "FALSE")
actual_class <- ifelse(b$income_binary == 1, "TRUE", "FALSE")
# Chuyển về factor với đúng thứ tự levels
pred_class <- factor(pred_class, levels = c("FALSE", "TRUE"))
actual_class <- factor(actual_class, levels = c("FALSE", "TRUE"))
# Tạo confusion matrix, tập trung vào phân biệt người có thu nhập cao
library(caret)
confusionMatrix(pred_class, actual_class, positive = "TRUE")## Confusion Matrix and Statistics
##
## Reference
## Prediction FALSE TRUE
## FALSE 18796 3778
## TRUE 3858 3730
##
## Accuracy : 0.7468
## 95% CI : (0.7419, 0.7517)
## No Information Rate : 0.7511
## P-Value [Acc > NIR] : 0.9563
##
## Kappa : 0.3253
##
## Mcnemar's Test P-Value : 0.3660
##
## Sensitivity : 0.4968
## Specificity : 0.8297
## Pos Pred Value : 0.4916
## Neg Pred Value : 0.8326
## Prevalence : 0.2489
## Detection Rate : 0.1237
## Detection Prevalence : 0.2516
## Balanced Accuracy : 0.6633
##
## 'Positive' Class : TRUE
## b. Ma trận nhầm lẫn cho mô hình probit
# B2: Dự đoán xác suất từ mô hình
prob_pred <- predict(model_probit, type = "response")
# B3: Xác định ngưỡng tối ưu bằng chỉ số Youden
library(pROC)
roc_obj <- roc(b$income_bin, prob_pred)
opt_threshold <- coords(roc_obj, "best", ret = "threshold", best.method = "youden")
threshold_value <- as.numeric(opt_threshold[1])
cat("Ngưỡng tối ưu là:", round(threshold_value, 4), "\n")## Ngưỡng tối ưu là: 0.3295# B4: Dự đoán phân loại dựa vào ngưỡng tối ưu
pred_class <- ifelse(prob_pred > threshold_value, "TRUE", "FALSE")
actual_class <- ifelse(b$income_bin == 1, "TRUE", "FALSE")
# B5: Chuyển thành factor để tạo confusion matrix
pred_class <- factor(pred_class, levels = c("FALSE", "TRUE"))
actual_class <- factor(actual_class, levels = c("FALSE", "TRUE"))
# B6: Tạo ma trận nhầm lẫn
library(caret)
confusionMatrix(pred_class, actual_class, positive = "TRUE")## Confusion Matrix and Statistics
##
## Reference
## Prediction FALSE TRUE
## FALSE 18796 3778
## TRUE 3858 3730
##
## Accuracy : 0.7468
## 95% CI : (0.7419, 0.7517)
## No Information Rate : 0.7511
## P-Value [Acc > NIR] : 0.9563
##
## Kappa : 0.3253
##
## Mcnemar's Test P-Value : 0.3660
##
## Sensitivity : 0.4968
## Specificity : 0.8297
## Pos Pred Value : 0.4916
## Neg Pred Value : 0.8326
## Prevalence : 0.2489
## Detection Rate : 0.1237
## Detection Prevalence : 0.2516
## Balanced Accuracy : 0.6633
##
## 'Positive' Class : TRUE
## c. Ma trận nhầm lẫn cho mô hình cloglog
prob_pred <- predict(model_cloglog, type = "response")
# B2: Tìm ngưỡng tối ưu theo chỉ số Youden
library(pROC)
roc_obj <- roc(b$income_binary, prob_pred)
opt_threshold <- coords(roc_obj, "best", ret = "threshold", best.method = "youden")
threshold_value <- as.numeric(opt_threshold[1])
cat("Ngưỡng tối ưu là:", round(threshold_value, 4), "\n")## Ngưỡng tối ưu là: 0.3295# B3: Dự đoán phân loại theo ngưỡng tối ưu
pred_class <- ifelse(prob_pred > threshold_value, "TRUE", "FALSE")
actual_class <- ifelse(b$income_binary == 1, "TRUE", "FALSE")
# B4: Chuyển về factor với đúng thứ tự levels
pred_class <- factor(pred_class, levels = c("FALSE", "TRUE"))
actual_class <- factor(actual_class, levels = c("FALSE", "TRUE"))
# B5: Tạo confusion matrix, tập trung vào phân biệt người có thu nhập cao
library(caret)
confusionMatrix(pred_class, actual_class, positive = "TRUE")## Confusion Matrix and Statistics
##
## Reference
## Prediction FALSE TRUE
## FALSE 18796 3778
## TRUE 3858 3730
##
## Accuracy : 0.7468
## 95% CI : (0.7419, 0.7517)
## No Information Rate : 0.7511
## P-Value [Acc > NIR] : 0.9563
##
## Kappa : 0.3253
##
## Mcnemar's Test P-Value : 0.3660
##
## Sensitivity : 0.4968
## Specificity : 0.8297
## Pos Pred Value : 0.4916
## Neg Pred Value : 0.8326
## Prevalence : 0.2489
## Detection Rate : 0.1237
## Detection Prevalence : 0.2516
## Balanced Accuracy : 0.6633
##
## 'Positive' Class : TRUE
## Kết luận
Cả ba mô hình Logit, Probit và Cloglog đều cho cùng một ngưỡng phân loại tối ưu là 0.3295, tức là khi xác suất dự đoán lớn hơn 0.3295 thì cá nhân được phân loại là “có thu nhập cao” (TRUE), ngược lại là “không có thu nhập cao” (FALSE).
Ma trận nhầm lẫn cho thấy cả ba mô hình đều cho ra kết quả dự đoán giống nhau. Cụ thể, mô hình đã dự đoán đúng 3.730 cá nhân có thu nhập cao và 18.796 cá nhân không có thu nhập cao. Tuy nhiên, mô hình vẫn mắc sai lệch khi phân loại: 3.778 người có thu nhập cao bị dự đoán nhầm là không có, trong khi 3.858 người không có thu nhập cao lại bị dự đoán sai là có thu nhập cao.
Về các chỉ số đánh giá mô hình:
Độ chính xác tổng thể (Accuracy) là 74.68%, tức là tức là gần 3/4 số dự đoán là chính xác. Tuy nhiên, nó vẫn thấp hơn tỷ lệ nhóm chiếm đa số (No Information Rate = 75.11%), cho thấy mô hình chưa vượt trội so với việc luôn dự đoán theo nhóm phổ biến.
Độ nhạy (Sensitivity) đạt 49.68%, nghĩa là mô hình nhận diện được khoảng 1 nửa số người thực sự có thu nhập cao – cho thấy khả năng phát hiện nhóm mục tiêu còn hạn chế.
Độ đặc hiệu (Specificity) đạt 82.97%, cho thấy mô hình phân biệt khá tốt những người không có thu nhập cao.
Giá trị tiên đoán dương (PPV) là 49.16%, nghĩa là trong số các cá nhân được dự đoán là có thu nhập cao, chỉ khoảng một nửa thực sự đúng.
Giá trị tiên đoán âm (NPV) cao, đạt 83.26%, phản ánh rằng mô hình khá đáng tin cậy khi dự đoán ai đó không có thu nhập cao.
Balanced Accuracy là 66.33%, là trung bình giữa độ nhạy và độ đặc hiệu – cho thấy mô hình có mức độ chính xác tương đối ổn định giữa hai nhóm đối tượng.
2.2 Ảnh hưởng của số giờ làm việc trung bình mỗi tuần đến khả năng đạt được thu nhập cao
Trong nghiên cứu này, biến “số giờ làm việc mỗi tuần” (hours.per.week) được sử dụng để đánh giá mối quan hệ giữa thời lượng lao động và khả năng đạt được mức thu nhập cao. Tuy nhiên, do đặc điểm phân bố đa dạng của biến này trong dữ liệu gốc, tác giả tiến hành gộp nhóm để dễ dàng nhận diện và so sánh.
Cụ thể, biến hours.per.week được phân loại thành hai nhóm:
Nhóm làm việc nhiều giờ (“Highhour”): gồm những cá nhân làm việc từ 40 giờ trở lên mỗi tuần, phản ánh nhóm lao động toàn thời gian hoặc cường độ cao;
Nhóm làm việc ít giờ (“Lowhour”): gồm những cá nhân làm việc dưới 40 giờ mỗi tuần, thường thuộc nhóm bán thời gian hoặc làm việc bán thời gian không thường xuyên.
2.2.1 Bảng tần số - tần suất cho biến số giờ làm việc trung bình mỗi tuần
Bảng dưới đây thể hiện phân bố tần số và tần suất của các cá nhân theo từng nhóm giờ làm :
# Tạo biến nhóm làm việc: Lowhour (<40 giờ) và Highhour (>=40 giờ)
b$hours_group <- ifelse(b$hours.per.week < 40, "Lowhour", "Highhour")
b$hours_group <- factor(b$hours_group, levels = c("Lowhour", "Highhour"))
library(dplyr)
library(knitr)
library(kableExtra)
# Tạo bảng tần số và tần suất cho nhóm làm việc
bang_tanso_hours_group <- b %>%
group_by(hours_group) %>%
summarise(
`Tần số` = n(),
`Tần suất (%)` = round(100 * n() / nrow(b), 2)
)
kable(
bang_tanso_hours_group,
caption = "Bảng 2.8: Phân bố tần số và tần suất theo nhóm số giờ làm việc mỗi tuần",
col.names = c("Nhóm giờ làm việc", "Tần số", "Tần suất (%)"),
align = "lcc",
format = "html"
) %>%
kable_styling(
bootstrap_options = c("striped", "hover", "condensed", "responsive"),
full_width = FALSE,
position = "center"
) %>%
row_spec(0, bold = TRUE, color = "white", background = "#0073C2FF") %>%
column_spec(2:3, width = "5em")| Nhóm giờ làm việc | Tần số | Tần suất (%) |
|---|---|---|
| Lowhour | 6714 | 22.26 |
| Highhour | 23448 | 77.74 |
Mặc dù bảng tần số và tần suất đã cho thấy sự chênh lệch rõ rệt giữa hai nhóm số giờ làm việc, nhưng việc trình bày dưới dạng bảng có thể chưa trực quan. Do đó, tác giả sử dụng biểu đồ để trực quan hóa dữ liệu, giúp người đọc dễ dàng nhận diện và so sánh tỷ lệ giữa các nhóm một cách sinh động và rõ ràng hơn.
# Thư viện trực quan hóa
library(ggplot2)
library(dplyr)
# Chuẩn bị dữ liệu biểu đồ từ biến hours_group
bang_plot_hours_group <- b %>%
count(hours_group) %>%
mutate(
ty_le = round(100 * n / sum(n), 1),
nhan = paste0(n, " (", ty_le, "%)")
)
ggplot(bang_plot_hours_group, aes(x = hours_group, y = n, fill = hours_group)) +
geom_bar(stat = "identity", width = 0.5, show.legend = FALSE) +
geom_text(aes(label = nhan), vjust = -0.5, size = 4.5, fontface = "bold") +
scale_fill_manual(values = c("Lowhour" = "#F8766D", "Highhour" = "#00BFC4")) +
labs(
title = "Biểu đồ 2.6: Phân bố nhóm số giờ làm việc mỗi tuần",
x = "Nhóm giờ làm việc",
y = "Số lượng cá nhân"
) +
scale_y_continuous(expand = expansion(mult = c(0, 0.1))) +
theme_minimal(base_size = 14) +
theme(
plot.title = element_text(face = "bold", hjust = 0.5, size = 16),
axis.title = element_text(face = "bold"),
axis.text = element_text(color = "black")
)
Kết quả bảng tần số và đồ thị cho thấy phần lớn cá nhân trong mẫu khảo sát thuộc nhóm làm việc nhiều giờ (Highhour), chiếm tới 77.74%, tương ứng với 23448 người. Trong khi đó, nhóm làm việc ít giờ (Lowhour) chỉ chiếm 22.26%, tương ứng với 6714 người.
2.2.2 Thống kê mô tả cho hai biến thu nhập và số giờ làm việc
Tác giả đã lập bảng tần số giữa hai biến: thu nhập và số giờ làm việc. Việc trình bày bảng tần số giúp làm rõ cách phân bố thu nhập theo nhóm giờ làm việc, từ đó cung cấp cái nhìn trực quan ban đầu về mối quan hệ giữa thời lượng làm việc và khả năng đạt mức thu nhập cao.
# Thư viện cần thiết
library(dplyr)
library(knitr)
library(kableExtra)
# Tạo bảng tần suất chéo giữa income và hours_group
bang_tansuat_cheo <- b %>%
group_by(hours_group, income) %>%
summarise(`Tần số` = n(), .groups = "drop") %>%
group_by(hours_group) %>%
mutate(`Tần suất (%)` = round(100 * `Tần số` / sum(`Tần số`), 2)) %>%
ungroup()
# Hiển thị bảng đẹp
kable(
bang_tansuat_cheo,
caption = "Bảng 2.9: Phân bố tần suất chéo giữa nhóm giờ làm việc và mức thu nhập",
col.names = c("Nhóm giờ làm việc", "Mức thu nhập", "Tần số", "Tần suất (%)"),
align = "lccc",
format = "html"
) %>%
kable_styling(
bootstrap_options = c("striped", "hover", "condensed", "responsive"),
full_width = FALSE,
position = "center"
) %>%
row_spec(0, bold = TRUE, color = "white", background = "#0073C2FF") %>%
column_spec(3:4, width = "5em")| Nhóm giờ làm việc | Mức thu nhập | Tần số | Tần suất (%) |
|---|---|---|---|
| Lowhour | <=50K | 6057 | 90.21 |
| Lowhour | >50K | 657 | 9.79 |
| Highhour | <=50K | 16597 | 70.78 |
| Highhour | >50K | 6851 | 29.22 |
Bảng tần suất chéo giữa nhóm giờ làm việc và mức thu nhập cho thấy sự khác biệt đáng kể về tỷ lệ đạt thu nhập cao giữa hai nhóm. Tuy nhiên, để dễ hình dung và so sánh trực quan hơn mức chênh lệch này, tác giả tiếp tục minh họa dữ liệu bằng biểu đồ cột, giúp làm nổi bật mối liên hệ giữa số giờ làm việc và khả năng đạt thu nhập cao.
# Thư viện cần thiết
library(ggplot2)
library(dplyr)
# Tạo dữ liệu cho biểu đồ: tần số & tần suất chéo
plot_data <- b %>%
group_by(hours_group, income) %>%
summarise(count = n(), .groups = "drop") %>%
group_by(hours_group) %>%
mutate(percent = round(100 * count / sum(count), 1))
# Vẽ biểu đồ cột nhóm
ggplot(plot_data, aes(x = hours_group, y = count, fill = income)) +
geom_bar(stat = "identity", position = position_dodge(width = 0.7), width = 0.6) +
geom_text(aes(label = paste0(count, "\n(", percent, "%)")),
position = position_dodge(width = 0.7),
vjust = -0.6, size = 4, fontface = "bold") +
scale_fill_manual(values = c("<=50K" = "#E7B800", ">50K" = "#00AFBB"),
labels = c("≤50K", ">50K")) +
labs(
title = "Biểu đồ 2.7: Phân bố mức thu nhập theo nhóm giờ làm việc",
x = "Nhóm giờ làm việc mỗi tuần",
y = "Số lượng cá nhân",
fill = "Mức thu nhập"
) +
scale_y_continuous(expand = expansion(mult = c(0, 0.1))) +
theme_minimal(base_size = 14) +
theme(
plot.title = element_text(face = "bold", hjust = 0.5, size = 16),
axis.title = element_text(face = "bold"),
axis.text = element_text(color = "black"),
legend.title = element_text(face = "bold")
)
Kết quả từ bảng và đồ thị cho thấy:
Trong nhóm làm việc nhiều giờ (Highhour), có 29.22% người đạt mức thu nhập trên 50K, cao hơn so với nhóm làm việc ít giờ. Điều này cho thấy làm việc với thời lượng lớn hơn có liên quan đến khả năng đạt thu nhập cao hơn.
Ngược lại, ở nhóm làm việc ít giờ (Lowhour), chỉ 9.79% người có thu nhập trên 50K, trong khi 90.21% còn lại có thu nhập không quá 50K. Kết quả này phản ánh rằng làm việc ít giờ có thể gắn liền với khả năng đạt thu nhập thấp hơn.
2.2.3 Relative Risk và khoảng ước lượng cho Relative Risk
Trong phần này, tác giả tiến hành tính toán Relative Risk nhằm phân tích mức độ ảnh hưởng của số giờ làm việc trung bình mỗi tuần đến khả năng đạt được mức thu nhập cao của các cá nhân.
# Gọi thư viện cần thiết
library(epitools)
# Đảm bảo biến có thứ tự đúng: [Highhour, Lowhour] x [<=50K, >50K]
b$hours_group <- factor(b$hours_group, levels = c("Highhour", "Lowhour"))
# Tạo bảng chéo 2x2 giữa nhóm giờ làm việc và thu nhập
table_hours_income <- table(b$hours_group, b$income)
# Tính Relative Risk
rr_hours_result <- riskratio(table_hours_income)
# In kết quả
print(rr_hours_result)## $data
##
## <=50K >50K Total
## Highhour 16597 6851 23448
## Lowhour 6057 657 6714
## Total 22654 7508 30162
##
## $measure
## risk ratio with 95% C.I.
## estimate lower upper
## Highhour 1.000000 NA NA
## Lowhour 0.334916 0.3106196 0.3611128
##
## $p.value
## two-sided
## midp.exact fisher.exact chi.square
## Highhour NA NA NA
## Lowhour 0 6.382839e-267 2.958662e-231
##
## $correction
## [1] FALSE
##
## attr(,"method")
## [1] "Unconditional MLE & normal approximation (Wald) CI"Kết quả tính toán cho thấy Relative Risk (RR) của nhóm Lowhour so với nhóm Highhour là 0.3349, nhỏ hơn 1. Điều này có nghĩa là: Những người làm việc dưới 40 giờ mỗi tuần (Lowhour) có khả năng đạt mức thu nhập cao (>50K) chỉ bằng khoảng 33.5% so với những người làm việc từ 40 giờ trở lên mỗi tuần (Highhour). Nói cách khác, khả năng đạt thu nhập cao của nhóm làm việc nhiều giờ (Highhour) cao hơn gấp khoảng 3 lần so với nhóm làm việc ít giờ (Lowhour).Khoảng tin cậy 95% cho RR nằm trong khoảng từ 0.3106 đến 0.3611, nghĩa là với độ tin cậy 95%, nguy cơ thực sự nằm trong khoảng đó.
2.2.4 Odds Ratio và khoảng ước lượng cho Odds ratio
Để lượng hóa mối quan hệ giữa số giờ làm việc trung bình mỗi tuần và khả năng đạt mức thu nhập cao (trên 50K USD/năm), tác giả áp dụng phương pháp phân tích tỷ số odds (Odds Ratio)
# Gọi thư viện
library(epitools)
# Đảm bảo biến giờ làm việc được đặt đúng thứ tự: [Highhour, Lowhour]
b$hours_group <- factor(b$hours_group, levels = c("Highhour", "Lowhour"))
# Tạo bảng chéo 2x2 giữa nhóm giờ làm việc và thu nhập
table_hours_income <- table(b$hours_group, b$income)
# Tính Odds Ratio
or_hours_result <- oddsratio(table_hours_income)
# In kết quả
print(or_hours_result)## $data
##
## <=50K >50K Total
## Highhour 16597 6851 23448
## Lowhour 6057 657 6714
## Total 22654 7508 30162
##
## $measure
## odds ratio with 95% C.I.
## estimate lower upper
## Highhour 1.0000000 NA NA
## Lowhour 0.2628365 0.2411887 0.2860239
##
## $p.value
## two-sided
## midp.exact fisher.exact chi.square
## Highhour NA NA NA
## Lowhour 0 6.382839e-267 2.958662e-231
##
## $correction
## [1] FALSE
##
## attr(,"method")
## [1] "median-unbiased estimate & mid-p exact CI"Giá trị Odds Ratio = 0.2628, nghĩa là odds đạt thu nhập cao (tức là tỷ lệ giữa số người có thu nhập >50K so với số người ≤50K) trong nhóm làm việc dưới 40 giờ mỗi tuần (Lowhour) chỉ bằng khoảng 26.3% so với odds tương ứng ở nhóm làm việc từ 40 giờ trở lên mỗi tuần (Highhour).
Khoảng tin cậy 95% của Odds Ratio (OR) cho nhóm làm việc ít giờ (Lowhour) là từ 0.2412 đến 0.2860. Điều này có nghĩa là odds đạt thu nhập cao của nhóm làm việc ít giờ chỉ bằng khoảng 24.1% đến 28.6%
2.2.5 Thống kê suy diễn
2.2.5.1 Kiểm định tính độc lập
Để đánh giá mối quan hệ giữa số giờ làm việc trung bình mỗi tuần và khả năng đạt được mức thu nhập cao, nghiên cứu tiến hành kiểm định Chi-bình phương nhằm kiểm tra tính độc lập giữa hai biến phân loại: nhóm số giờ làm việc và mức thu nhập. Mục tiêu là xác định liệu sự khác biệt về thời lượng lao động hàng tuần có liên quan một cách có ý nghĩa thống kê đến khả năng đạt được thu nhập trên 50.000 USD mỗi năm hay không.
Đưa ra giả thuyết: \[ \left\{ \begin{array}{ll} H_0: \text{Số giờ làm việc trung bình mỗi tuần của cá nhân không có mối quan hệ với mức thu nhập.} \\ H_1: \text{Số giờ làm việc trung bình mỗi tuần của cá nhân có mối quan hệ với mức thu nhập.} \end{array} \right. \]
# Thực hiện kiểm định Chi bình phương
chi_result <- chisq.test(table_hours_income)
# In kết quả kiểm định
print(chi_result)##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: table_hours_income
## X-squared = 1053.2, df = 1, p-value < 2.2e-16Vì p_value < 2.2e-16 < 0.05, ta bác bỏ giả thuyết \(H_0\) tại mức ý nghĩa 5%. Vậy số giờ làm việc trung bình mỗi tuần của cá nhân có mối quan hệ với mức thu nhập.
2.2.5.2 So sánh tỷ lệ thu nhập >50K giữa 2 nhóm giờ làm việc
- Giả thuyết \(H_0\): Tỷ lệ người có thu nhập cao ở nhóm làm việc nhiều giờ không cao hơn so với nhóm làm việc ít giờ.
\[ H_0: p_{Highhour} ≤ p_{Lowhour} \]
- Giả thuyết \(H_1\):Tỷ lệ người có thu nhập cao ở nhóm làm việc nhiều giờ cao hơn so với nhóm làm việc ít giờ.
\[ H_1: p_{Highhour} > p_{Lowhour} \]
# Tạo bảng chéo giữa nhóm giờ làm việc và mức thu nhập
table_hours_income <- table(b$hours_group, b$income)
# Số người có thu nhập >50K trong mỗi nhóm giờ làm việc
success <- table_hours_income[, ">50K"]
# Tổng số người trong mỗi nhóm
n_total <- rowSums(table_hours_income)
# Kiểm định tỷ lệ: nhóm Highhour có tỷ lệ thu nhập cao hơn nhóm Lowhour
prop.test(success, n_total, alternative = "greater", correct = TRUE)##
## 2-sample test for equality of proportions with continuity correction
##
## data: success out of n_total
## X-squared = 1053.2, df = 1, p-value < 2.2e-16
## alternative hypothesis: greater
## 95 percent confidence interval:
## 0.1865179 1.0000000
## sample estimates:
## prop 1 prop 2
## 0.29217844 0.09785523Kết quả kiểm định tỷ lệ hai mẫu cho thấy tỷ lệ người có thu nhập >50K ở nhóm làm việc nhiều giờ (Highhour) là khoảng 29.22%, trong khi ở nhóm làm việc ít giờ (Lowhour) chỉ là khoảng 9.79%. Giá trị p_value < 2.2e-16 < 0.05. Vậy ta bác bỏ giả thuyết \(H_0\) tại mức ý nghĩa 5%. Điều này cho thấy có bằng chứng thống kê rõ ràng rằng tỷ lệ người có thu nhập cao ở nhóm làm việc nhiều giờ lớn hơn so với nhóm làm việc ít giờ. Khoảng tin cậy 95% cho hiệu số tỷ lệ nằm trong khoảng từ 18.65% đến 100%, càng củng cố cho kết luận về sự khác biệt rõ rệt này.
2.2.6 Mô hình hồi quy
2.2.6.1 Mô hình logit
# Chuyển biến highhour_bin thành factor với "Highhour" là nhóm tham chiếu
b$hours_group <- ifelse(b$hours.per.week >= 40, "Highhour", "Lowhour")
b$hours_group <- factor(b$hours_group, levels = c("Highhour", "Lowhour"))
# Tạo biến phụ thuộc income_bin
b$income_binary <- ifelse(b$income == ">50K", 1, 0)
# Ước lượng mô hình logit với hours_group là factor (Highhour là nhóm tham chiếu)
model_logit_income <- glm(income_binary ~ hours_group, data = b, family = binomial(link = "logit"))
summary(model_logit_income)##
## Call:
## glm(formula = income_binary ~ hours_group, family = binomial(link = "logit"),
## data = b)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.88483 0.01436 -61.62 <2e-16 ***
## hours_groupLowhour -1.33646 0.04351 -30.72 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 33851 on 30161 degrees of freedom
## Residual deviance: 32631 on 30160 degrees of freedom
## AIC: 32635
##
## Number of Fisher Scoring iterations: 4Ta có kết quả từ hàm hồi quy thông qua mô hình logit như sau:
\[ \ln\left(\frac{\hat{\pi}}{1 - \hat{\pi}}\right) = -0.88483 - 1.33646 \cdot \text{hours_group}_{\text{Lowhour}} \]
Hệ số chặn \(\hat{\beta}_0\) là -0.88483, điều này có nghĩa là khi không có sự tác động của nhóm làm việc ít giờ, thì log-odds đạt mức thu nhập cao của nhóm làm việc nhiều giờ là -0.88483.
Hệ số \(\hat{\beta}_1\) = -1.33646, thể hiện chênh lệch log-odds đạt thu nhập cao giữa nhóm làm việc ít giờ và nhóm làm việc nhiều giờ. Cụ thể, nhóm làm việc ít giờ có log-odds đạt thu nhập cao thấp hơn 1.33646 đơn vị so với nhóm làm việc nhiều giờ.
Nếu chuyển từ log-odds sang odds, ta thấy rằng odds đạt thu nhập cao của nhóm làm việc ít giờ chỉ bằng khoảng 26.3% so với nhóm làm việc nhiều giờ. Nói cách khác, khả năng đạt thu nhập cao ở nhóm làm ít giờ thấp hơn đáng kể.
Tất cả các hệ số đều có ý nghĩa thống kê rất mạnh với p_value nhỏ hơn 2e−16, tức là nhỏ hơn mức ý nghĩa 5%, cho thấy mối liên hệ giữa số giờ làm việc và khả năng đạt thu nhập cao là đáng tin cậy về mặt thống kê.
Từ mô hình logit tác giả dự báo xác suất có thu nhập cao tương ứng với mỗi nhóm học vấn như sau:
library(knitr)
library(kableExtra)
# Tạo dữ liệu mới để dự báo
new_data <- data.frame(hours_group = c("Highhour", "Lowhour"))
# Dự báo xác suất từ mô hình logit
new_data$predicted_prob <- predict(model_logit_income, newdata = new_data, type = "response")
# Làm đẹp bảng với tiêu đề màu xanh
kable(new_data,
col.names = c("Nhóm giờ làm việc", "Xác suất thu nhập cao (>50K)"),
digits = 4,
align = "c",
caption = "Bảng 2.10: Dự báo xác suất có thu nhập cao theo nhóm giờ làm việc") %>%
kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"),
full_width = FALSE, position = "center") %>%
row_spec(0, bold = TRUE, color = "white", background = "#0073C2FF")| Nhóm giờ làm việc | Xác suất thu nhập cao (>50K) |
|---|---|
| Highhour | 0.2922 |
| Lowhour | 0.0979 |
Để giúp trực quan hóa sự khác biệt về xác suất đạt thu nhập cao giữa các nhóm giờ làm việc, tác giả trình bày biểu đồ tròn sau. Biểu đồ này minh họa rõ ràng tỷ trọng dự báo thu nhập cao ở từng nhóm, qua đó hỗ trợ so sánh trực quan hơn so với bảng số liệu.
library(ggplot2)
# Tạo dữ liệu mới
new_data <- data.frame(
hours_group = factor(c("Highhour", "Lowhour"), levels = c("Highhour", "Lowhour"))
)
# Dự báo xác suất từ mô hình logit
new_data$predicted_prob <- predict(model_logit_income, newdata = new_data, type = "response")
# Tạo nhãn phần trăm
new_data$label <- paste0(ifelse(new_data$hours_group == "Highhour", "Nhiều giờ", "Ít giờ"),
"\n(", round(new_data$predicted_prob * 100, 1), "%)")
# Vẽ biểu đồ tròn
ggplot(new_data, aes(x = "", y = predicted_prob, fill = hours_group)) +
geom_col(width = 1, color = "white") +
coord_polar("y") +
geom_text(aes(label = label),
position = position_stack(vjust = 0.5),
size = 5, color = "black", fontface = "bold") +
scale_fill_manual(values = c("#66CC99", "#FF9999"),
labels = c("Nhiều giờ", "Ít giờ")) +
theme_void() +
ggtitle("Biểu đồ 2.8: Dự báo xác suất có thu nhập cao theo nhóm giờ làm việc") +
theme(plot.title = element_text(hjust = 0.5, size = 16, face = "bold"),
legend.title = element_blank(),
legend.position = "bottom")
Dựa trên kết quả dự báo từ mô hình hồi quy logistic, ta thấy rằng xác suất đạt mức thu nhập cao có sự khác biệt rõ rệt giữa hai nhóm giờ làm việc:
Đối với nhóm làm việc nhiều giờ (từ 40 giờ/tuần trở lên), xác suất đạt mức thu nhập cao là 0.2922, tức khoảng 29.22%.
Trong khi đó, ở nhóm làm việc ít giờ (dưới 40 giờ/tuần), xác suất này chỉ là 0.0979, tức khoảng 9.79%.
Kết quả này cho thấy những người làm việc nhiều giờ có khả năng đạt mức thu nhập cao gấp khoảng 3 lần so với những người làm việc ít giờ.
2.2.6.2 Mô hình probit
# Xây dựng mô hình Probit
model_probit_income <- glm(income_binary ~ hours_group,
family = binomial(link = "probit"),
data = b)
# Xem kết quả mô hình
summary(model_probit_income)##
## Call:
## glm(formula = income_binary ~ hours_group, family = binomial(link = "probit"),
## data = b)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.547032 0.008646 -63.27 <2e-16 ***
## hours_groupLowhour -0.746838 0.022702 -32.90 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 33851 on 30161 degrees of freedom
## Residual deviance: 32631 on 30160 degrees of freedom
## AIC: 32635
##
## Number of Fisher Scoring iterations: 4\[ \hat{\pi} = \Phi( -0.547032 -0.746838\cdot \text{hours_group}_{\text{Lowhour}}) \]
Hệ số chặn \(\hat{\beta}_0\) = -0.54703, tương ứng với xác suất đạt thu nhập cao của nhóm làm việc nhiều giờ. Xác suất tương ứng được tính là \(\Phi(−0.54704)\) ≈ 0.2922, tức khoảng tức khoảng 29.2%.
Hệ số \(\hat{\beta}_1\) = -0.74684 phản ánh rằng làm việc ít giờ có z-score thấp hơn nhóm làm việc nhiều giờ là 0.74684 đơn vị. Hai hệ số \(\hat{\beta}_0\) và \(\hat{\beta}_1\) có ý nghĩa thống kê rất cao (p_value < 2e−16).
Ta có thể tính xác suất đạt thu nhập cao ở nhóm học vấn thấp bằng: \(\Phi( -0.54703 − 0.74684\) ≈ 0.0983, tức là chỉ khoảng 9.8% người làm việc ít giờ đạt được mức thu nhập cao.
Từ mô hình Probit, ta dự báo xác suất đạt mức thu nhập cao tương ứng với từng nhóm số giờ làm việc như sau:
library(knitr)
library(kableExtra)
# Tạo dữ liệu mới theo nhóm số giờ làm việc
new_data <- data.frame(hours_group = factor(c("Highhour", "Lowhour"),
levels = c("Highhour", "Lowhour")))
# Dự báo xác suất từ mô hình Probit
new_data$predicted_prob <- predict(model_probit_income,
newdata = new_data,
type = "response")
# Hiển thị bảng kết quả
kable(new_data,
col.names = c("Nhóm số giờ làm việc", "Xác suất thu nhập cao (>50K)"),
digits = 4,
align = "c",
caption = "Bảng 2.11: Dự báo xác suất có thu nhập cao theo nhóm số giờ làm việc (Mô hình Probit)") %>%
kable_styling(full_width = FALSE, position = "center") %>%
row_spec(0, bold = TRUE, color = "white", background = "#0073C2FF") %>%
column_spec(2, width = "10em")| Nhóm số giờ làm việc | Xác suất thu nhập cao (>50K) |
|---|---|
| Highhour | 0.2922 |
| Lowhour | 0.0979 |
Tác giả thực hiện trực quan hóa dự báo xác suất có thu nhập cao theo nhóm số giờ làm việc (Mô hình Probit) bằng biểu đồ tròn.
# Cài đặt và nạp thư viện
library(ggplot2)
# Tạo bộ dữ liệu mới gồm 2 nhóm giờ làm việc
new_data <- data.frame(
hours_group = factor(c("Highhour", "Lowhour"), levels = c("Highhour", "Lowhour"))
)
# Dự báo xác suất từ mô hình Probit đã huấn luyện
new_data$predicted_prob <- predict(model_probit_income, newdata = new_data, type = "response")
# Đổi tên nhãn cho dễ hiểu
new_data$hours_label <- recode(new_data$hours_group,
"Highhour" = "Làm nhiều giờ",
"Lowhour" = "Làm ít giờ")
# Tạo nhãn phần trăm cho chú thích
new_data$label <- paste0(new_data$hours_label, "\n(", round(new_data$predicted_prob * 100, 1), "%)")
# Vẽ biểu đồ tròn
ggplot(new_data, aes(x = "", y = predicted_prob, fill = hours_label)) +
geom_col(width = 1, color = "white") +
coord_polar("y") +
geom_text(aes(label = label),
position = position_stack(vjust = 0.5),
size = 5, color = "black", fontface = "bold") +
scale_fill_manual(values = c("#0099CC", "#FF9966")) + # phối màu tương phản
theme_void() +
ggtitle("Biểu đồ 2.9: Dự báo xác suất có thu nhập cao theo nhóm giờ làm việc") +
theme(plot.title = element_text(hjust = 0.5, size = 16, face = "bold"),
legend.title = element_blank(),
legend.position = "bottom")
Kết quả cho thấy:
Nhóm làm việc từ 40 giờ trở lên mỗi tuần (Highhour): Xác suất đạt mức thu nhập cao được ước tính là khoảng 29.22%. Điều này cho thấy, một cá nhân làm việc toàn thời gian có gần 1/3 khả năng đạt được thu nhập cao.
Nhóm làm việc dưới 40 giờ mỗi tuần (Lowhour): Xác suất đạt thu nhập cao chỉ khoảng 9.79%, tức thấp hơn đáng kể so với nhóm làm việc nhiều giờ.
Như vậy, số giờ làm việc có ảnh hưởng rõ rệt đến khả năng đạt thu nhập cao. Cụ thể, những người làm việc từ 40 giờ trở lên có xác suất đạt thu nhập cao gấp gần 3 lần so với những người làm việc ít giờ hơn.
2.2.6.3 Mô hình clog-log
# Xây dựng mô hình clog-log
model_cloglog_income <- glm(income_binary ~ hours_group,
family = binomial(link = "cloglog"),
data = b)
# Xem kết quả mô hình
summary(model_cloglog_income)##
## Call:
## glm(formula = income_binary ~ hours_group, family = binomial(link = "cloglog"),
## data = b)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.06258 0.01214 -87.51 <2e-16 ***
## hours_groupLowhour -1.21064 0.04088 -29.62 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 33851 on 30161 degrees of freedom
## Residual deviance: 32631 on 30160 degrees of freedom
## AIC: 32635
##
## Number of Fisher Scoring iterations: 5Phương trình ước lượng của mô hình hồi quy với liên kết cloglog có dạng như sau:
\[ \log\left(-\log(1 - \hat{\pi})\right) = −1.06258 - 1.21064 \cdot \text{hours_group}_{\text{Lowhour}} \]
Hệ số chặn \(\hat{\beta}_0\) = -1.06258 đại diện cho giá trị log–log tương ứng với xác suất đạt thu nhập cao của nhóm làm việc nhiều giờ (từ 40 giờ trở lên mỗi tuần) .Xác suất này được tính như sau:
\(\hat{\pi}_{\text{Highhour}} = 1 - \exp(-\exp(-1.06258)) \approx 0.2916\)
Nói cách khác, khoảng 29.2% những người làm việc nhiều giờ có khả năng đạt mức thu nhập cao.
Hệ số \(\hat{\beta}_1 = -1.21064\) thể hiện sự chênh lệch về \(\log(-\log(1 - \hat{\pi}))\) giữa nhóm làm việc ít giờ và nhóm làm việc nhiều giờ. Dấu âm của hệ số cho thấy nhóm làm việc ít giờ có khả năng đạt thu nhập cao thấp hơn so với nhóm làm việc nhiều giờ.
Xác suất đạt thu nhập cao cho nhóm làm việc ít giờ được ước lượng như sau:
\(\hat{P}_{\text{lowhour}} = 1 - \exp(-\exp(−1.06258 - 1.21064)) \approx 0.0981\)
Điều này cho thấy chỉ khoảng 9.8% những người thuộc nhóm làm việc ít giờ có khả năng đạt được thu nhập cao – thấp hơn nhiều so với nhóm làm việc nhiều giờ. Tại mức ý nghĩa 5%, cả hai hệ số đều có ý nghĩa thống kê rất cao với giá trị p < 2e−16 < 0.05.
2.2.7 Đánh giá mô hình
2.2.7.1 AIC
# Nạp thư viện
library(knitr)
library(kableExtra)
# Tạo bảng dữ liệu AIC
aic_compare <- data.frame(
Model = c("Logit", "Probit", "Cloglog"),
AIC = c(
AIC(model_logit_income),
AIC(model_probit_income),
AIC(model_cloglog_income)
)
)
# Hiển thị bảng đẹp
kable(aic_compare,
col.names = c("Mô hình", "Giá trị AIC"),
digits = 2,
align = "c",
caption = "Bảng 2.12: So sánh giá trị AIC giữa các mô hình Logit, Probit và Cloglog") %>%
kable_styling(full_width = FALSE, position = "center") %>%
row_spec(0, bold = TRUE, color = "white", background = "#0073C2FF") %>%
column_spec(2, width = "10em")| Mô hình | Giá trị AIC |
|---|---|
| Logit | 32635.03 |
| Probit | 32635.03 |
| Cloglog | 32635.03 |
Cả ba mô hình Logit, Probit và Cloglog đều cho cùng một giá trị AIC là 32635.03. Điều này cho thấy mức độ phù hợp mô hình là tương đương giữa các hàm liên kết, và không có mô hình nào thể hiện ưu thế rõ rệt về mặt thống kê trong bối cảnh dữ liệu hiện tại.
2.2.7.2 Hệ số Brier Score
Chỉ số Brier Score được sử dụng để đánh giá mức độ chính xác của các dự đoán xác suất, thông qua việc đo lường sai số bình phương giữa xác suất dự báo và kết quả thực tế. Giá trị Brier Score càng thấp cho thấy mô hình dự báo càng chính xác.
# Gọi thư viện cần thiết
library(DescTools) # Để dùng hàm BrierScore()
library(knitr) # Để tạo bảng
library(kableExtra) # Để trình bày bảng đẹp
# Tính Brier Score bằng mô hình GLM trực tiếp
brier_logit <- BrierScore(model_logit_income)
brier_probit <- BrierScore(model_probit_income)
brier_cloglog <- BrierScore(model_cloglog_income)
# Tạo bảng kết quả
brier_table <- data.frame(
Mô_hình = c("Logit", "Probit", "Cloglog"),
Brier_Score = round(c(brier_logit, brier_probit, brier_cloglog), 5)
)
# Trình bày bảng đẹp
kable(brier_table,
col.names = c("Mô hình", "Chỉ số Brier Score"),
align = "c",
caption = "Bảng 2.13: So sánh chỉ số Brier Score giữa các mô hình nhị phân") %>%
kable_styling(full_width = FALSE, position = "center") %>%
row_spec(0, bold = TRUE, color = "white", background = "#0073C2FF") %>%
column_spec(2, width = "10em")| Mô hình | Chỉ số Brier Score |
|---|---|
| Logit | 0.18043 |
| Probit | 0.18043 |
| Cloglog | 0.18043 |
Trong phân tích này, cả ba mô hình Logit, Probit và Cloglog đều có cùng giá trị Brier Score là 0.18043. Điều này phản ánh rằng mức độ chính xác trong dự báo xác suất của ba mô hình là tương đương nhau, và không có mô hình nào thể hiện sự vượt trội rõ rệt về mặt hiệu suất dự báo.
2.2.7.3 Ma trận nhầm lẫn
Để đánh giá khả năng của mô hình trong việc phân biệt đúng những cá nhân có thu nhập cao, nghiên cứu sử dụng ma trận nhầm lẫn, với lớp dương được xác định là nhóm “TRUE” – tức là có thu nhập cao.
a. Ma trận nhầm lẫn cho mô hình logit
prob_pred <- predict(model_logit_income, type = "response")
#Tìm ngưỡng tối ưu theo chỉ số Youden
library(pROC)
roc_obj <- roc(b$income_binary, prob_pred)
opt_threshold <- coords(roc_obj, "best", ret = "threshold", best.method = "youden")
threshold_value <- as.numeric(opt_threshold[1])
cat("Ngưỡng tối ưu là:", round(threshold_value, 4), "\n")## Ngưỡng tối ưu là: 0.195# Dự đoán phân loại theo ngưỡng tối ưu
pred_class <- ifelse(prob_pred > threshold_value, "TRUE", "FALSE")
actual_class <- ifelse(b$income_binary == 1, "TRUE", "FALSE")
# Chuyển về factor với đúng thứ tự levels
pred_class <- factor(pred_class, levels = c("FALSE", "TRUE"))
actual_class <- factor(actual_class, levels = c("FALSE", "TRUE"))
# Tạo confusion matrix, tập trung vào phân biệt người có thu nhập cao
library(caret)
confusionMatrix(pred_class, actual_class, positive = "TRUE")## Confusion Matrix and Statistics
##
## Reference
## Prediction FALSE TRUE
## FALSE 6057 657
## TRUE 16597 6851
##
## Accuracy : 0.428
## 95% CI : (0.4224, 0.4336)
## No Information Rate : 0.7511
## P-Value [Acc > NIR] : 1
##
## Kappa : 0.1052
##
## Mcnemar's Test P-Value : <2e-16
##
## Sensitivity : 0.9125
## Specificity : 0.2674
## Pos Pred Value : 0.2922
## Neg Pred Value : 0.9021
## Prevalence : 0.2489
## Detection Rate : 0.2271
## Detection Prevalence : 0.7774
## Balanced Accuracy : 0.5899
##
## 'Positive' Class : TRUE
## b. Ma trận nhầm lẫn cho mô hình probit
# B1: Dự đoán xác suất từ mô hình Probit
prob_pred <- predict(model_probit_income, type = "response")
# B2: Tìm ngưỡng tối ưu theo chỉ số Youden
library(pROC)
roc_obj <- roc(b$income_binary, prob_pred)
opt_threshold <- coords(roc_obj, "best", ret = "threshold", best.method = "youden")
threshold_value <- as.numeric(opt_threshold[1])
cat("Ngưỡng tối ưu là:", round(threshold_value, 4), "\n")## Ngưỡng tối ưu là: 0.195# B3: Dự đoán phân loại theo ngưỡng tối ưu
pred_class <- ifelse(prob_pred > threshold_value, "TRUE", "FALSE")
actual_class <- ifelse(b$income_binary == 1, "TRUE", "FALSE")
# B4: Chuyển về factor với đúng thứ tự levels
pred_class <- factor(pred_class, levels = c("FALSE", "TRUE"))
actual_class <- factor(actual_class, levels = c("FALSE", "TRUE"))
# B5: Tạo confusion matrix, tập trung vào phân biệt người có thu nhập cao
library(caret)
confusionMatrix(pred_class, actual_class, positive = "TRUE")## Confusion Matrix and Statistics
##
## Reference
## Prediction FALSE TRUE
## FALSE 6057 657
## TRUE 16597 6851
##
## Accuracy : 0.428
## 95% CI : (0.4224, 0.4336)
## No Information Rate : 0.7511
## P-Value [Acc > NIR] : 1
##
## Kappa : 0.1052
##
## Mcnemar's Test P-Value : <2e-16
##
## Sensitivity : 0.9125
## Specificity : 0.2674
## Pos Pred Value : 0.2922
## Neg Pred Value : 0.9021
## Prevalence : 0.2489
## Detection Rate : 0.2271
## Detection Prevalence : 0.7774
## Balanced Accuracy : 0.5899
##
## 'Positive' Class : TRUE
## c. Ma trận nhầm lẫn cho mô hình cloglog
# B1: Dự đoán xác suất từ mô hình Probit
prob_pred <- predict(model_cloglog_income, type = "response")
# B2: Tìm ngưỡng tối ưu theo chỉ số Youden
library(pROC)
roc_obj <- roc(b$income_binary, prob_pred)
opt_threshold <- coords(roc_obj, "best", ret = "threshold", best.method = "youden")
threshold_value <- as.numeric(opt_threshold[1])
cat("Ngưỡng tối ưu là:", round(threshold_value, 4), "\n")## Ngưỡng tối ưu là: 0.195# B3: Dự đoán phân loại theo ngưỡng tối ưu
pred_class <- ifelse(prob_pred > threshold_value, "TRUE", "FALSE")
actual_class <- ifelse(b$income_binary == 1, "TRUE", "FALSE")
# B4: Chuyển về factor với đúng thứ tự levels
pred_class <- factor(pred_class, levels = c("FALSE", "TRUE"))
actual_class <- factor(actual_class, levels = c("FALSE", "TRUE"))
# B5: Tạo confusion matrix, tập trung vào phân biệt người có thu nhập cao
library(caret)
confusionMatrix(pred_class, actual_class, positive = "TRUE")## Confusion Matrix and Statistics
##
## Reference
## Prediction FALSE TRUE
## FALSE 6057 657
## TRUE 16597 6851
##
## Accuracy : 0.428
## 95% CI : (0.4224, 0.4336)
## No Information Rate : 0.7511
## P-Value [Acc > NIR] : 1
##
## Kappa : 0.1052
##
## Mcnemar's Test P-Value : <2e-16
##
## Sensitivity : 0.9125
## Specificity : 0.2674
## Pos Pred Value : 0.2922
## Neg Pred Value : 0.9021
## Prevalence : 0.2489
## Detection Rate : 0.2271
## Detection Prevalence : 0.7774
## Balanced Accuracy : 0.5899
##
## 'Positive' Class : TRUE
## Kết luận
Cả ba mô hình Logit, Probit và Cloglog đều cho cùng một ngưỡng phân loại tối ưu là 0.195 – nghĩa là nếu xác suất dự đoán một cá nhân có thu nhập cao vượt quá 0.195, thì người đó được phân loại là “có thu nhập cao” (TRUE); ngược lại là “không có thu nhập cao” (FALSE).
Ma trận nhầm lẫn của cả ba mô hình đều cho ra kết quả giống nhau: mô hình dự đoán đúng 6851 trường hợp có thu nhập cao và 6057 trường hợp không có thu nhập cao. Tuy nhiên, mô hình cũng phân loại sai khá nhiều, với 16597 cá nhân không có thu nhập cao bị dự đoán nhầm là có thu nhập cao.
Các chỉ số đánh giá mô hình chi tiết như sau:
Độ chính xác tổng thể (Accuracy) là 42.8%, thấp hơn tỷ lệ nhóm chiếm đa số (No Information Rate = 75.11%), phản ánh mô hình không vượt trội trong toàn bộ mẫu.
Độ nhạy (Sensitivity) rất cao, đạt 91.25%, cho thấy mô hình có khả năng phát hiện tốt nhóm mục tiêu là người có thu nhập cao.
Độ đặc hiệu (Specificity) thấp, chỉ 26.74%, phản ánh khả năng nhận diện đúng nhóm không có thu nhập cao còn hạn chế.
Giá trị tiên đoán dương (PPV) là 29.22%, nghĩa là trong số những người được dự đoán là có thu nhập cao, chỉ khoảng 29% là đúng.
Giá trị tiên đoán âm (NPV) cao, đạt 90.21%, tức mô hình khá tin cậy khi dự đoán ai đó không có thu nhập cao.
Balanced Accuracy đạt 58.99%, thể hiện hiệu suất trung bình giữa hai nhóm.
2.3 Phân tích tác động của trình độ học vấn, số giờ làm việc trung bình mỗi tuần đến khả năng một cá nhân có thu nhập cao
2.3.1 Mô hình hồi quy logit
Trong phần phân tích này, tác giả xây dựng một mô hình hồi quy logit nhị phân nhằm ước lượng xác suất một cá nhân có thu nhập >50K dựa trên các yếu tố bao gồm: trình độ học vấn, số giờ làm việc trung bình mỗi tuần.
# Mô hình Logit với 2 biến độc lập
model_logit <- glm(income_binary ~ edu_group + hours_group,
data = b,
family = binomial(link = "logit"))
summary(model_logit)##
## Call:
## glm(formula = income_binary ~ edu_group + hours_group, family = binomial(link = "logit"),
## data = b)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.14730 0.02411 6.11 9.99e-10 ***
## edu_groupLower -1.53209 0.02962 -51.72 < 2e-16 ***
## hours_groupLowhour -1.26396 0.04488 -28.16 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 33851 on 30161 degrees of freedom
## Residual deviance: 29926 on 30159 degrees of freedom
## AIC: 29932
##
## Number of Fisher Scoring iterations: 5Ta có kết quả từ hàm hồi quy thông qua mô hình logit như sau:
\[ \ln\left(\frac{\hat{\pi}}{1 - \hat{\pi}}\right) = 0.14730 - 1.53209 \cdot \text{edu_group}_{\text{Lower}} - 1.26396 \cdot \text{hours_group}_{\text{Lowhour}} \]
Kết quả từ mô hình hồi quy Logit cho thấy tất cả các biến Hệ số chặn (Intercept)\(\hat{\beta}_0\) là 0.14730, đại diện cho log-odds đạt mức thu nhập cao (tức thu nhập >50K) của nhóm có trình độ học vấn cao và làm việc từ 40 giờ trở lên mỗi tuần. Hệ số này có ý nghĩa thống kê rất cao vì p_value = 9.99e−10 < 0.05
Hệ số \(\hat{\beta}_1\) = -1.53209, thể hiện chênh lệch log-odds đạt mức thu nhập cao giữa nhóm có trình độ học vấn thấp so với nhóm có trình độ học vấn cao (trong cùng điều kiện về số giờ làm việc) là −1.53209. Nói cách khác, nhóm học vấn thấp có log-odds đạt thu nhập cao thấp hơn 1.53209 đơn vị so với nhóm học vấn cao.
Nếu xét theo tỷ lệ odds (odds ratio), thì odds đạt thu nhập cao của nhóm có học vấn thấp chỉ bằng khoảng 21.6% so với nhóm có học vấn cao. Điều này có nghĩa là khả năng đạt thu nhập cao của nhóm học vấn thấp thấp hơn khoảng 4.6 lần so với nhóm học vấn cao.
Đối với hệ số \(\hat{\beta}_2\) = -1.26396, thể hiện chênh lệch log-odds đạt mức thu nhập cao giữa nhóm làm việc dưới 40 giờ mỗi tuần so với nhóm làm việc từ 40 giờ trở lên (trong cùng điều kiện về trình độ học vấn) là −1.26396. Điều này có nghĩa là nhóm làm việc ít giờ có log-odds đạt thu nhập cao thấp hơn 1.26396 đơn vị so với nhóm làm việc đủ giờ. Hệ số này có ý nghĩa thống kê rất mạnh vì p_value < 2e−16 < 0.05, cho thấy mối liên hệ giữa trình độ học vấn và khả năng đạt thu nhập cao là đáng tin cậy về mặt thống kê.
Nếu xét theo tỷ lệ odds, thì odds đạt thu nhập cao của nhóm làm việc dưới 40 giờ chỉ bằng khoảng 28.3% so với nhóm làm việc đủ giờ. Nói cách khác, khả năng đạt thu nhập cao của nhóm làm việc ít giờ thấp hơn khoảng 3.5 lần so với nhóm làm việc đủ giờ. Hệ số này cũng có ý nghĩa thống kê rất mạnh vì p_value < 2e−16 < 0.05, cho thấy thời lượng làm việc là yếu tố ảnh hưởng rõ rệt đến khả năng đạt thu nhập cao.
library(knitr)
library(kableExtra)
# Tạo dữ liệu tổ hợp
new_data <- data.frame(
edu_group = factor(c("Higher", "Lower", "Higher", "Lower"), levels = c("Higher", "Lower")),
hours_group = factor(c("Highhour", "Highhour", "Lowhour", "Lowhour"), levels = c("Highhour", "Lowhour"))
)
# Dự báo xác suất từ mô hình logit
new_data$predicted_prob <- predict(model_logit, newdata = new_data, type = "response")
# Đặt nhãn cho hàng
new_data$Nhóm <- c(
"Học vấn cao – Làm ≥ 40 giờ",
"Học vấn thấp – Làm ≥ 40 giờ",
"Học vấn cao – Làm < 40 giờ",
"Học vấn thấp – Làm < 40 giờ"
)
# Chọn và sắp xếp lại cột để hiển thị
final_table <- new_data[, c("Nhóm", "predicted_prob")]
# Hiển thị bảng đẹp
kable(final_table,
col.names = c("Nhóm phân tích", "Xác suất thu nhập cao (>50K)"),
digits = 4,
align = "c",
caption = "Bảng 2.14: Dự báo xác suất có thu nhập cao theo trình độ học vấn và số giờ làm việc") %>%
kable_styling(full_width = FALSE, position = "center") %>%
row_spec(0, bold = TRUE, color = "white", background = "#0073C2FF") %>%
column_spec(2, width = "10em")| Nhóm phân tích | Xác suất thu nhập cao (>50K) |
|---|---|
| Học vấn cao – Làm ≥ 40 giờ | 0.5368 |
| Học vấn thấp – Làm ≥ 40 giờ | 0.2002 |
| Học vấn cao – Làm < 40 giờ | 0.2466 |
| Học vấn thấp – Làm < 40 giờ | 0.0661 |
Để giúp dễ dàng nhận diện sự khác biệt và so sánh giữa các nhóm, tác giả tiến hành trực quan hóa dữ liệu bằng biểu đồ cột.
# Gọi thư viện trực quan
library(ggplot2)
# Vẽ biểu đồ cột dự báo xác suất
ggplot(final_table, aes(x = Nhóm, y = predicted_prob, fill = Nhóm)) +
geom_bar(stat = "identity", width = 0.6) +
geom_text(aes(label = round(predicted_prob, 4)), vjust = -0.5, size = 4) +
scale_fill_brewer(palette = "Blues") +
labs(
title = "Biểu đồ 2.10: Dự báo xác suất thu nhập cao theo học vấn và giờ làm việc",
x = "Nhóm phân tích",
y = "Xác suất thu nhập cao (>50K)"
) +
theme_minimal(base_size = 13) +
theme(
legend.position = "none",
axis.text.x = element_text(angle = 25, hjust = 1),
plot.title = element_text(face = "bold", color = "#0073C2FF", size = 8) # giảm kích thước chữ ở đây
)
Kết quả từ mô hình logit cho thấy trình độ học vấn và số giờ làm việc đều ảnh hưởng rõ rệt đến xác suất có thu nhập cao. Người học vấn cao, làm ≥ 40 giờ/tuần có xác suất đạt thu nhập cao cao nhất (53.68%), trong khi nhóm học vấn thấp, làm < 40 giờ/tuần có xác suất thấp nhất (6.61%). Ngoài ra, với cùng một trình độ học vấn, người làm việc ≥ 40 giờ luôn có xác suất cao hơn so với người làm < 40 giờ.
2.3.2 Mô hình hồi quy probit
# Mô hình Logit với 2 biến độc lập
model_probit <- glm(income_binary ~ edu_group + hours_group,
data = b,
family = binomial(link = "probit"))
summary(model_probit)##
## Call:
## glm(formula = income_binary ~ edu_group + hours_group, family = binomial(link = "probit"),
## data = b)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.08389 0.01492 5.622 1.89e-08 ***
## edu_groupLower -0.91974 0.01772 -51.897 < 2e-16 ***
## hours_groupLowhour -0.70965 0.02379 -29.832 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 33851 on 30161 degrees of freedom
## Residual deviance: 29912 on 30159 degrees of freedom
## AIC: 29918
##
## Number of Fisher Scoring iterations: 5Kết quả hồi quy từ mô hình Probit nhằm đánh giá xác suất một cá nhân có thu nhập cao cho thấy rằng các yếu tố như trình độ học vấn, số giờ làm việc trung bình mỗi tuần đều có ảnh hưởng có ý nghĩa thống kê ở mức ý nghĩa 5% đến khả năng đạt thu nhập cao. Phương trình hồi quy probit ước lượng có dạng:
\[ \hat{\pi} = \Phi(0.08389 - 0.91974 \cdot \text{edu_group}_{\text{Lower}} - 0.70965 \cdot \text{hours_group}_{\text{Lowhour}}) \]
Hệ số chặn \(\hat{\beta}_0\) = 0.08389, tương ứng với xác suất đạt thu nhập cao của nhóm có trình độ học vấn cao và làm việc từ 40 giờ trở lên mỗi tuần. Xác suất tương ứng được tính là \(\Phi(0.08389)\) ≈ 0.5334, tức khoảng tức khoảng 53.3%.
Hệ số \(\hat{\beta}_1\) = -0.91974 phản ánh rằng nhóm có trình độ học vấn thấp có z-score thấp hơn nhóm học vấn cao là 0.91974 đơn vị (khi số giờ làm việc là như nhau).
Ta có thể tính xác suất đạt thu nhập cao ở nhóm học vấn thấp, làm việc ≥ 40 giờ bằng: \(\Phi( 0.08389 − 0.91974)\) ≈ 0.2016, tức là chỉ khoảng 20.2% người có trình độ học vấn thấp (nhưng vẫn làm đủ từ 40 giờ/tuần) đạt được mức thu nhập cao.
Hệ số \(\hat{\beta}_2\) cho thấy người làm dưới 40 giờ mỗi tuần có z-score thấp hơn 0.70965 đơn vị so với người làm ≥ 40 giờ (cùng trình độ học vấn). Xác suất đạt thu nhập cao ở nhóm có trình độ học vấn cao nhưng làm dưới 40 giờ/tuần là: \(\Phi(0.08389 − 0.70965)\) ≈ 0.266, tức khoảng 26.6%.
Cuối cùng, với nhóm vừa có trình độ học vấn thấp và làm việc dưới 40 giờ, z-score là: \(\Phi(0.08389 − 0.91974 − 0.70965)\) ≈ 0.0618, tức chỉ khoảng 6.2% có khả năng đạt thu nhập cao.
Tiếp theo, tác giả dự báo xác suất đạt mức thu nhập cao với mô hình probit
# Thư viện cần thiết
library(knitr)
library(kableExtra)
# Tạo dữ liệu tổ hợp
new_data <- data.frame(
edu_group = factor(c("Higher", "Lower", "Higher", "Lower"), levels = c("Higher", "Lower")),
hours_group = factor(c("Highhour", "Highhour", "Lowhour", "Lowhour"), levels = c("Highhour", "Lowhour"))
)
# Dự báo xác suất từ mô hình probit
new_data$predicted_prob <- predict(model_probit, newdata = new_data, type = "response")
# Gán nhãn nhóm
new_data$Nhóm <- c(
"Học vấn cao – Làm ≥ 40 giờ",
"Học vấn thấp – Làm ≥ 40 giờ",
"Học vấn cao – Làm < 40 giờ",
"Học vấn thấp – Làm < 40 giờ"
)
# Sắp xếp lại bảng
final_table <- new_data[, c("Nhóm", "predicted_prob")]
# Hiển thị bảng đẹp
kable(final_table,
col.names = c("Nhóm phân tích", "Xác suất thu nhập cao (>50K)"),
digits = 4,
align = "c",
caption = "Bảng 2.15: Dự báo xác suất có thu nhập cao theo trình độ học vấn và số giờ làm việc (Mô hình Probit)") %>%
kable_styling(full_width = FALSE, position = "center") %>%
row_spec(0, bold = TRUE, color = "white", background = "#0073C2FF") %>%
column_spec(2, width = "10em")| Nhóm phân tích | Xác suất thu nhập cao (>50K) |
|---|---|
| Học vấn cao – Làm ≥ 40 giờ | 0.5334 |
| Học vấn thấp – Làm ≥ 40 giờ | 0.2016 |
| Học vấn cao – Làm < 40 giờ | 0.2657 |
| Học vấn thấp – Làm < 40 giờ | 0.0611 |
Để giúp dễ dàng nhận diện sự khác biệt và so sánh giữa các nhóm, tác giả tiến hành trực quan hóa dữ liệu bằng biểu đồ cột.
# Thư viện trực quan hóa
library(ggplot2)
# Vẽ biểu đồ cột thể hiện xác suất thu nhập cao theo từng nhóm
ggplot(new_data, aes(x = Nhóm, y = predicted_prob, fill = Nhóm)) +
geom_bar(stat = "identity", width = 0.6) +
geom_text(aes(label = round(predicted_prob, 3)),
vjust = -0.5, size = 4, fontface = "bold") +
scale_fill_brewer(palette = "Set2") +
labs(title = "Biểu đồ 2.11: Xác suất thu nhập cao theo học vấn và số giờ làm việc (Probit)",
x = "Nhóm phân tích",
y = "Xác suất dự báo (>50K)") +
ylim(0, 1) +
theme_minimal() +
theme(
axis.text.x = element_text(angle = 20, hjust = 1),
plot.title = element_text(face = "bold", hjust = 0.5, size = 10),
legend.position = "none"
)
Dựa trên kết quả dự báo từ mô hình hồi quy Probit, ta thấy rằng xác suất đạt mức thu nhập cao có sự chênh lệch rõ rệt giữa các nhóm phân tích, đặc biệt theo trình độ học vấn và số giờ làm việc hàng tuần.
Cụ thể, những người có học vấn cao và làm việc từ 40 giờ/tuần trở lên có xác suất đạt thu nhập cao lên đến 0.5334 (tức 53.34%), trong khi nhóm có học vấn thấp nhưng làm cùng số giờ chỉ đạt 0.2016 (20.16%). Như vậy, trong cùng điều kiện về thời gian làm việc, người có học vấn cao có khả năng đạt thu nhập cao gấp hơn 2.6 lần so với người có học vấn thấp.
Tương tự, khi so sánh giữa các mức giờ làm việc trong cùng một nhóm học vấn, những người có học vấn cao làm việc dưới 40 giờ/tuần có xác suất đạt thu nhập cao là 0.2657, cao hơn khoảng 4.3 lần so với người có học vấn thấp làm việc dưới 40 giờ (chỉ đạt 0.0611).
Tổng thể, mô hình cho thấy xác suất đạt mức thu nhập cao tăng đáng kể khi cá nhân có trình độ học vấn cao hơn và làm việc nhiều giờ hơn mỗi tuần.
2.3.3 Mô hình hồi quy cloglog
# Xây dựng mô hình cloglog
model_cloglog <- glm(income_binary ~ edu_group + hours_group,
data = b,
family = binomial(link = "cloglog"))
# Xem kết quả mô hình
summary(model_cloglog)##
## Call:
## glm(formula = income_binary ~ edu_group + hours_group, family = binomial(link = "cloglog"),
## data = b)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.25838 0.01707 -15.14 <2e-16 ***
## edu_groupLower -1.24812 0.02342 -53.29 <2e-16 ***
## hours_groupLowhour -1.09673 0.04081 -26.87 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 33851 on 30161 degrees of freedom
## Residual deviance: 29950 on 30159 degrees of freedom
## AIC: 29956
##
## Number of Fisher Scoring iterations: 5Mô hình hồi quy với hàm liên kết cloglog được ước lượng như sau:
\[ \log\left(-\log(1 - \hat{\pi})\right) = -0.25838 - 1.24812 \cdot \text{edu_group}_{\text{Lower}} - 1.09673 \cdot \text{hours_group}_{\text{Lowhour}} \]
Hệ số chặn \(\hat{\beta}_0\) = -0.25838 đại diện cho giá trị log–log tương ứng với xác suất đạt thu nhập cao của nhóm có trình độ học vấn cao và làm việc từ 40 giờ trở lên mỗi tuần.
Xác suất tương ứng được tính theo công thức nghịch đảo của liên kết cloglog như sau:
\[ \hat{\pi}_{\text{HighEdu, HighHour}} = 1 - \exp\left(-\exp(-0.25838)\right) \approx 0.4181, \]
tức khoảng 41.8%. Hệ số này có ý nghĩa thống kê rất cao (p-value < 2e−16).
Hệ số \(\hat{\beta}_1 = -1.24812\) phản ánh rằng nhóm có trình độ học vấn thấp có log–log của xác suất đạt thu nhập cao thấp hơn nhóm học vấn cao là 1.24812 đơn vị (giữ nguyên số giờ làm việc). Hệ số này có ý nghĩa thống kê rất cao (p-value < 2e−16).
Ta có thể tính xác suất đạt thu nhập cao ở nhóm học vấn thấp (nhưng vẫn làm ≥ 40 giờ/tuần) bằng:
\[ \hat{\pi}_{\text{LowEdu, HighHour}} = 1 - \exp\left(-\exp(-0.25838 - 1.24812)\right) \approx 0.1440, \]
tức chỉ khoảng 14.4% có khả năng đạt thu nhập cao.
Tương tự, hệ số \(\hat{\beta}_2 = -1.09673\) cho thấy người làm việc dưới 40 giờ có log–log xác suất thấp hơn 1.09673 đơn vị so với người làm ≥ 40 giờ (giữ nguyên trình độ học vấn). Hệ số này có ý nghĩa thống kê rất cao (p-value < 2e−16).
Ta có thể tính xác suất đạt thu nhập cao của nhóm học vấn cao nhưng làm việc ít giờ là:
\[ \hat{\pi}_{\text{HighEdu, LowHour}} = 1 - \exp\left(-\exp(-0.25838 - 1.09673)\right) \approx 0.1655, \]
tức khoảng 16.6%.
Cuối cùng, với nhóm vừa có trình độ học vấn thấp và làm việc dưới 40 giờ mỗi tuần:
\[ \hat{\pi}_{\text{LowEdu, LowHour}} = 1 - \exp\left(-\exp(-0.25838 - 1.24812 - 1.09673)\right) \approx 0.0563, \]
tức chỉ khoảng 5.6% người thuộc nhóm này có khả năng đạt thu nhập cao.
2.3.4 Đánh giá mô hình
2.3.4.1 AIC
Để lựa chọn mô hình hồi quy nhị phân phù hợp nhất cho việc phân tích xác suất đạt mức thu nhập cao, ta tiến hành so sánh mức độ phù hợp của ba mô hình phổ biến gồm: Logit, Probit và Cloglog. Việc so sánh được thực hiện dựa trên chỉ số AIC (Akaike Information Criterion), một tiêu chí thống kê thường được sử dụng nhằm đánh giá độ phù hợp của mô hình trong khi vẫn kiểm soát mức độ phức tạp. Phần dưới đây trình bày kết quả so sánh AIC giữa các mô hình.
library(knitr)
library(kableExtra)
# Tạo bảng so sánh AIC
aic_compare <- data.frame(
Model = c("Logit", "Probit", "Cloglog"),
AIC = c(
AIC(model_logit),
AIC(model_probit),
AIC(model_cloglog)
)
)
# Làm tròn giá trị AIC để dễ nhìn
aic_compare$AIC <- round(aic_compare$AIC, 2)
# Trình bày bảng đẹp
kable(aic_compare,
col.names = c("Mô hình", "Chỉ số AIC"),
align = "c",
caption = "Bảng 2.16: So sánh chỉ số AIC giữa các mô hình nhị phân") %>%
kable_styling(full_width = FALSE, position = "center") %>%
row_spec(0, bold = TRUE, color = "white", background = "#0073C2FF") %>%
column_spec(2, width = "10em")| Mô hình | Chỉ số AIC |
|---|---|
| Logit | 29931.97 |
| Probit | 29917.91 |
| Cloglog | 29956.00 |
Dựa trên kết quả so sánh chỉ số AIC giữa ba mô hình nhị phân gồm Logit, Probit và Cloglog, ta thấy rằng mô hình Probit có chỉ số AIC thấp nhất (29917.91), cho thấy đây là mô hình phù hợp nhất trong ba mô hình được xem xét.
2.3.4.2 Hệ số Brier Score
Chỉ số Brier Score được sử dụng để đánh giá mức độ chính xác của các dự đoán xác suất, thông qua việc đo lường sai số bình phương giữa xác suất dự báo và kết quả thực tế. Giá trị Brier Score càng thấp cho thấy mô hình dự báo càng chính xác..
# Cài đặt nếu chưa có
# install.packages("DescTools")
library(DescTools)
library(knitr)
library(kableExtra)
# Tính Brier Score cho từng mô hình
brier_logit <- BrierScore(model_logit)
brier_probit <- BrierScore(model_probit)
brier_cloglog <- BrierScore(model_cloglog)
# Tạo bảng so sánh Brier Score
brier_compare <- data.frame(
Mô_hình = c("Logit", "Probit", "Cloglog"),
Brier_Score = c(brier_logit, brier_probit, brier_cloglog)
)
# Trình bày bảng đẹp
kable(brier_compare,
col.names = c("Mô hình", "Chỉ số Brier Score"),
digits = 4,
align = "c",
caption = "Bảng 2.17: So sánh chỉ số Brier Score giữa các mô hình") %>%
kable_styling(full_width = FALSE, position = "center") %>%
row_spec(0, bold = TRUE, color = "white", background = "#0073C2FF") %>%
column_spec(2, width = "10em")| Mô hình | Chỉ số Brier Score |
|---|---|
| Logit | 0.1624 |
| Probit | 0.1623 |
| Cloglog | 0.1625 |
Dựa vào chỉ số Brier Score, mô hình Probit cho kết quả dự báo tốt nhất với giá trị thấp nhất (0.1623), nhỉnh hơn so với Logit (0.1624) và Cloglog (0.1625). Điều này nhất quán với kết quả AIC và cho thấy Probit là lựa chọn phù hợp nhất trong phân tích xác suất thu nhập cao.
2.3.4.3 Ma trận nhầm lẫn
Để đánh giá khả năng của mô hình trong việc phân biệt đúng những cá nhân có thu nhập cao, nghiên cứu sử dụng ma trận nhầm lẫn, với lớp dương được xác định là nhóm “TRUE” – tức là có thu nhập cao.
a. Ma trận nhầm lẫn cho mô hình logit
prob_pred <- predict(model_logit, type = "response")
# B2: Tìm ngưỡng tối ưu theo chỉ số Youden
library(pROC)
roc_obj <- roc(b$income_binary, prob_pred)
opt_threshold <- coords(roc_obj, "best", ret = "threshold", best.method = "youden")
threshold_value <- as.numeric(opt_threshold[1])
cat("Ngưỡng tối ưu là:", round(threshold_value, 4), "\n")## Ngưỡng tối ưu là: 0.2234# B3: Dự đoán phân loại theo ngưỡng tối ưu
pred_class <- ifelse(prob_pred > threshold_value, "TRUE", "FALSE")
actual_class <- ifelse(b$income_binary == 1, "TRUE", "FALSE")
# B4: Chuyển về factor với đúng thứ tự levels
pred_class <- factor(pred_class, levels = c("FALSE", "TRUE"))
actual_class <- factor(actual_class, levels = c("FALSE", "TRUE"))
# B5: Tạo confusion matrix, tập trung vào phân biệt người có thu nhập cao
library(caret)
confusionMatrix(pred_class, actual_class, positive = "TRUE")## Confusion Matrix and Statistics
##
## Reference
## Prediction FALSE TRUE
## FALSE 18796 3778
## TRUE 3858 3730
##
## Accuracy : 0.7468
## 95% CI : (0.7419, 0.7517)
## No Information Rate : 0.7511
## P-Value [Acc > NIR] : 0.9563
##
## Kappa : 0.3253
##
## Mcnemar's Test P-Value : 0.3660
##
## Sensitivity : 0.4968
## Specificity : 0.8297
## Pos Pred Value : 0.4916
## Neg Pred Value : 0.8326
## Prevalence : 0.2489
## Detection Rate : 0.1237
## Detection Prevalence : 0.2516
## Balanced Accuracy : 0.6633
##
## 'Positive' Class : TRUE
## Kết quả cho thấy ngưỡng phân loại tối ưu là 0.2234, tức là nếu xác suất dự đoán lớn hơn 0.2234 thì cá nhân đó được phân loại là “có thu nhập cao” (TRUE), ngược lại là “không có thu nhập cao” (FALSE).
Ma trận nhầm lẫn tương ứng cho thấy mô hình dự đoán đúng 3730 trường hợp có thu nhập cao và 18796 trường hợp không có thu nhập cao. Tuy nhiên, mô hình cũng phân loại sai khá nhiều, với 3778 người có thu nhập cao bị dự đoán nhầm là không có, và 3858 người không có thu nhập cao bị dự đoán nhầm là có thu nhập cao.
Các chỉ số đánh giá mô hình chi tiết như sau:
Độ chính xác tổng thể (Accuracy) là 74.68%, gần tương đương với tỷ lệ nhóm chiếm đa số (No Information Rate = 75.11%), phản ánh mô hình không cải thiện đáng kể so với việc chỉ dự đoán nhóm phổ biến.
Độ nhạy (Sensitivity) đạt 49.68%, tức mô hình chỉ phát hiện đúng khoảng một nửa số người có thu nhập cao.
Độ đặc hiệu (Specificity) là 82.97%, cho thấy mô hình có khả năng phân loại khá tốt nhóm không có thu nhập cao.
Giá trị tiên đoán dương (PPV) là 49.16%, nghĩa là trong số những cá nhân được dự đoán là có thu nhập cao, chỉ khoảng 49% thực sự đúng.
Giá trị tiên đoán âm (NPV) đạt 83.26%, cho thấy mô hình khá đáng tin cậy khi dự đoán một cá nhân không thuộc nhóm thu nhập cao.
Balanced Accuracy đạt 66.33%, thể hiện hiệu suất trung bình giữa hai nhóm là chấp nhận được.
b. Ma trận nhầm lẫn cho mô hình probit
# B1: Dự đoán xác suất từ mô hình Probit
prob_pred <- predict(model_probit, type = "response")
# B2: Tìm ngưỡng tối ưu theo chỉ số Youden
library(pROC)
roc_obj <- roc(b$income_binary, prob_pred)
opt_threshold <- coords(roc_obj, "best", ret = "threshold", best.method = "youden")
threshold_value <- as.numeric(opt_threshold[1])
cat("Ngưỡng tối ưu là:", round(threshold_value, 4), "\n")## Ngưỡng tối ưu là: 0.2337# B3: Dự đoán phân loại theo ngưỡng tối ưu
pred_class <- ifelse(prob_pred > threshold_value, "TRUE", "FALSE")
actual_class <- ifelse(b$income_binary == 1, "TRUE", "FALSE")
# B4: Chuyển về factor với đúng thứ tự levels
pred_class <- factor(pred_class, levels = c("FALSE", "TRUE"))
actual_class <- factor(actual_class, levels = c("FALSE", "TRUE"))
# B5: Tạo confusion matrix, tập trung vào phân biệt người có thu nhập cao
library(caret)
confusionMatrix(pred_class, actual_class, positive = "TRUE")## Confusion Matrix and Statistics
##
## Reference
## Prediction FALSE TRUE
## FALSE 18796 3778
## TRUE 3858 3730
##
## Accuracy : 0.7468
## 95% CI : (0.7419, 0.7517)
## No Information Rate : 0.7511
## P-Value [Acc > NIR] : 0.9563
##
## Kappa : 0.3253
##
## Mcnemar's Test P-Value : 0.3660
##
## Sensitivity : 0.4968
## Specificity : 0.8297
## Pos Pred Value : 0.4916
## Neg Pred Value : 0.8326
## Prevalence : 0.2489
## Detection Rate : 0.1237
## Detection Prevalence : 0.2516
## Balanced Accuracy : 0.6633
##
## 'Positive' Class : TRUE
## Kết quả phân tích cho thấy ngưỡng phân loại tối ưu là 0.2337, tức là nếu xác suất dự đoán vượt quá 23.37% thì cá nhân đó được phân loại là “có thu nhập cao” (TRUE), ngược lại là “không có thu nhập cao” (FALSE).
Ma trận nhầm lẫn tương ứng cho thấy mô hình dự đoán đúng 6851 trường hợp có thu nhập cao và 6057 trường hợp không có thu nhập cao. Tuy nhiên, mô hình cũng phân loại sai khá nhiều, với 16597 cá nhân không có thu nhập cao bị dự đoán nhầm là có thu nhập cao.
Các chỉ số đánh giá mô hình chi tiết như sau:
Độ chính xác tổng thể (Accuracy) là 74.68%, tức mô hình dự đoán đúng khoảng 3/4 trường hợp. Tuy nhiên, độ chính xác này không cao hơn đáng kể so với No Information Rate (NIR = 75.11%), tức là chỉ cần đoán “không có thu nhập cao” cho mọi người cũng đạt độ chính xác tương đương.
Độ nhạy (Sensitivity) là 49.68%, nghĩa là mô hình chỉ phát hiện đúng được khoảng 50% người có thu nhập cao, còn lại 50% bị bỏ sót.
Độ đặc hiệu (Specificity) đạt 82.97%, tức mô hình nhận diện khá tốt những người không có thu nhập cao, với hơn 82% được phân loại đúng.
Giá trị tiên đoán dương (PPV) là 49.16%, nghĩa là trong số những người được dự đoán là có thu nhập cao, chỉ khoảng một nửa thực sự đúng.
-Giá trị tiên đoán âm (NPV) cao, đạt 83.26%, cho thấy mô hình khá tin cậy khi dự đoán ai đó không có thu nhập cao.
- Balanced Accuracy đạt 66.33%, phản ánh hiệu suất trung bình giữa hai nhóm, và cho thấy mô hình có khả năng phân biệt tốt hơn so với ngẫu nhiên nhưng vẫn còn hạn chế.
Tóm lại, mô hình Probit thể hiện khả năng nhận diện khá tốt nhóm không có thu nhập cao (nhờ độ đặc hiệu cao), nhưng lại bỏ sót gần một nửa số người có thu nhập cao. Mặc dù độ chính xác tổng thể đạt khoảng 74.7%, nhưng không cao hơn đáng kể so với đoán theo nhóm chiếm đa số, cho thấy hiệu suất mô hình vẫn còn hạn chế trong việc phân biệt rõ ràng hai nhóm thu nhập.
c. Ma trận nhầm lẫn cho mô hình cloglog
# B1: Dự đoán xác suất từ mô hình Probit
prob_pred <- predict(model_cloglog, type = "response")
# B2: Tìm ngưỡng tối ưu theo chỉ số Youden
library(pROC)
roc_obj <- roc(b$income_binary, prob_pred)
opt_threshold <- coords(roc_obj, "best", ret = "threshold", best.method = "youden")
threshold_value <- as.numeric(opt_threshold[1])
cat("Ngưỡng tối ưu là:", round(threshold_value, 4), "\n")## Ngưỡng tối ưu là: 0.2131# B3: Dự đoán phân loại theo ngưỡng tối ưu
pred_class <- ifelse(prob_pred > threshold_value, "TRUE", "FALSE")
actual_class <- ifelse(b$income_binary == 1, "TRUE", "FALSE")
# B4: Chuyển về factor với đúng thứ tự levels
pred_class <- factor(pred_class, levels = c("FALSE", "TRUE"))
actual_class <- factor(actual_class, levels = c("FALSE", "TRUE"))
# B5: Tạo confusion matrix, tập trung vào phân biệt người có thu nhập cao
library(caret)
confusionMatrix(pred_class, actual_class, positive = "TRUE")## Confusion Matrix and Statistics
##
## Reference
## Prediction FALSE TRUE
## FALSE 18796 3778
## TRUE 3858 3730
##
## Accuracy : 0.7468
## 95% CI : (0.7419, 0.7517)
## No Information Rate : 0.7511
## P-Value [Acc > NIR] : 0.9563
##
## Kappa : 0.3253
##
## Mcnemar's Test P-Value : 0.3660
##
## Sensitivity : 0.4968
## Specificity : 0.8297
## Pos Pred Value : 0.4916
## Neg Pred Value : 0.8326
## Prevalence : 0.2489
## Detection Rate : 0.1237
## Detection Prevalence : 0.2516
## Balanced Accuracy : 0.6633
##
## 'Positive' Class : TRUE
## Kết quả cho thấy ngưỡng phân loại tối ưu là 0.2131, tức là nếu xác suất dự đoán lớn hơn 0.2131 thì cá nhân đó được phân loại là “có thu nhập cao” (TRUE), ngược lại là “không có thu nhập cao” (FALSE).
Ma trận nhầm lẫn tương ứng cho thấy mô hình dự đoán đúng 3730 trường hợp có thu nhập cao và 18796 trường hợp không có thu nhập cao. Tuy nhiên, mô hình cũng có một tỷ lệ phân loại sai đáng kể, với 3778 người có thu nhập cao bị dự đoán nhầm là không có, và 3858 người không có thu nhập cao bị dự đoán nhầm là có.
Các chỉ số đánh giá mô hình chi tiết như sau:
Độ chính xác tổng thể (Accuracy) là 74.68%, rất gần với tỷ lệ nhóm phổ biến (No Information Rate = 75.11%), cho thấy mô hình không cải thiện đáng kể so với việc luôn dự đoán nhóm chiếm đa số.
Độ nhạy (Sensitivity) đạt 49.68%, tức là chỉ khoảng một nửa số người có thu nhập cao được mô hình phát hiện đúng.
Độ đặc hiệu (Specificity) đạt 82.97%, phản ánh khả năng nhận diện nhóm không có thu nhập cao là khá tốt.
Giá trị tiên đoán dương (PPV) là 49.16%, nghĩa là trong số những người được mô hình dự đoán là có thu nhập cao, chỉ có khoảng 49% thực sự đúng.
Giá trị tiên đoán âm (NPV) đạt 83.26%, cho thấy mô hình khá đáng tin cậy khi dự đoán một cá nhân không thuộc nhóm thu nhập cao.
Balanced Accuracy đạt 66.33%, phản ánh hiệu suất trung bình giữa hai nhóm.
Kết luận
Trong ba mô hình, Probit có độ nhạy cao nhất (91.25%), tức là phát hiện đúng phần lớn người có thu nhập cao, dù độ chính xác tổng thể còn thấp (42.8%) do phân loại sai nhóm còn lại. Ngược lại, Logit và Cloglog có độ chính xác tổng thể cao hơn (74.68%) nhưng độ nhạy chỉ khoảng 49.7%, bỏ sót nhiều trường hợp thu nhập cao. Vậy mô hình Probit là lựa chọn phù hợp nhất khi nhận diện người có thu nhập cao.
CHƯƠNG 3: KẾT LUẬN
3.1 Kết luận
Nghiên cứu này đã phân tích mối quan hệ giữa thu nhập và hai yếu tố then chốt: trình độ học vấn và số giờ làm việc mỗi tuần. Các phân tích thống kê mô tả, kiểm định giả thuyết và hồi quy nhị phân đều cho thấy cả hai yếu tố này đều có ảnh hưởng rõ rệt đến khả năng đạt mức thu nhập cao (>50K).
Cụ thể, nhóm có trình độ học vấn cao có tỷ lệ thu nhập cao là 49.16%, trong khi nhóm học vấn thấp chỉ đạt 16.74%. Kiểm định Chi-bình phương cho p-value < 2.2e-16, khẳng định mối liên hệ có ý nghĩa thống kê giữa học vấn và thu nhập. Chỉ số RR = 0.3405 cho thấy người học vấn thấp chỉ có 34.1% khả năng đạt thu nhập cao so với nhóm học vấn cao. Đồng thời, OR = 0.2079 cho thấy odds đạt mức thu nhập >50K ở nhóm học vấn thấp chỉ bằng 20.8% so với nhóm học vấn cao. Các mô hình hồi quy nhị phân (Logit, Probit, Cloglog) đều nhất quán, với hệ số âm có ý nghĩa thống kê rất cao (p < 2e-16). Trong mô hình Logit, odds đạt thu nhập cao ở nhóm học vấn thấp thấp hơn khoảng 4.8 lần. Kết quả dự đoán xác suất từ mô hình Probit cũng sát với thống kê mô tả: 49.16% (Higher) và 16.74% (Lower).
Tương tự, số giờ làm việc mỗi tuần cũng cho thấy ảnh hưởng đáng kể đến thu nhập. Nhóm làm việc từ 40 giờ trở lên có tỷ lệ thu nhập cao là 29.22%, trong khi nhóm làm dưới 40 giờ chỉ là 9.79%. Kiểm định Chi-bình phương tiếp tục cho p < 2.2e-16. Chỉ số RR = 0.3349 cho thấy người làm ít giờ chỉ có 33.5% khả năng đạt thu nhập cao so với người làm đủ giờ, còn OR = 0.2628 thể hiện odds chỉ bằng 26.3%. Mô hình Logit cho thấy odds đạt thu nhập cao của nhóm làm việc ít giờ thấp hơn khoảng 3.8 lần. Mô hình Probit dự đoán xác suất là 29.22% (Highhour) và 9.79% (Lowhour), sát với số liệu thực tế.
Khi kết hợp cả hai yếu tố, sự chênh lệch rõ rệt hơn: nhóm có học vấn cao và làm việc đủ giờ (Higher & Highhour) có xác suất đạt thu nhập cao cao nhất (53.68%), trong khi nhóm học vấn thấp và làm việc ít giờ (Lower & Lowhour) có xác suất thấp nhất (6.61%). Cả ba mô hình hồi quy đều cho thấy các hệ số tương tác có ý nghĩa thống kê mạnh mẽ, khẳng định ảnh hưởng đồng thời của học vấn và thời gian lao động đến thu nhập.
Về mặt phân loại, cả ba mô hình Logit, Probit và Cloglog đều cho kết quả tương tự nhau với độ chính xác 74.68%, độ nhạy 49.68% và độ đặc hiệu 82.97%. Điều này cho thấy mô hình dự đoán tốt nhóm không có thu nhập cao, nhưng bỏ sót gần một nửa số người có thu nhập cao.
Tóm lại, tuy hiệu quả phân loại của ba mô hình là tương đương, mô hình Probit được đánh giá là phù hợp hơn cả khi xét đến tiêu chí tối ưu hóa xác suất dự đoán (AIC và Brier Score), dù khả năng phát hiện người có thu nhập cao vẫn còn hạn chế.
3.2 Hàm ý chính sách
Kết quả nghiên cứu cho thấy rằng trình độ học vấn và số giờ làm việc mỗi tuần đều có ảnh hưởng đáng kể đến xác suất đạt được mức thu nhập cao. Do đó, một số hàm ý chính sách có thể được rút ra như sau:
Thứ nhất, cần tiếp tục mở rộng cơ hội tiếp cận giáo dục đại học và sau đại học, đặc biệt đối với các nhóm dân cư có hoàn cảnh khó khăn. Việc nâng cao trình độ học vấn không chỉ giúp cải thiện kỹ năng lao động mà còn làm gia tăng xác suất đạt được thu nhập cao một cách rõ rệt (từ 16.74% ở nhóm học vấn thấp lên 49.16% ở nhóm học vấn cao).
Thứ hai, chính sách lao động cần khuyến khích việc làm toàn thời gian và ổn định, nhất là trong các lĩnh vực có năng suất cao. Nhóm làm việc ≥40 giờ/tuần có xác suất thu nhập cao gần gấp 3 lần so với nhóm làm việc ít giờ (29.22% so với 9.79%), cho thấy vai trò của cường độ lao động trong tích lũy thu nhập.
Thứ ba, cần có các chương trình đào tạo lại và nâng cao kỹ năng gắn liền với nhu cầu thị trường, nhằm hỗ trợ những người có học vấn thấp hoặc làm việc bán thời gian có thể tiếp cận các vị trí có thu nhập tốt hơn.
Cuối cùng, việc sử dụng mô hình dự báo xác suất trong hoạch định chính sách lao động – giáo dục nên được thúc đẩy nhằm xác định chính xác nhóm có nguy cơ thu nhập thấp, từ đó can thiệp sớm và hiệu quả hơn bằng các chính sách phúc lợi hoặc hỗ trợ nghề nghiệp phù hợp.
3.2 Hạn chế và hướng nghiên cứu tiếp theo
Mặc dù nghiên cứu đã cung cấp các bằng chứng thực nghiệmvề vai trò của trình độ học vấn và số giờ làm việc đối với xác suất đạt thu nhập cao, vẫn tồn tại một số hạn chế cần được xem xét.
Thứ nhất, mô hình hồi quy nhị phân (Logit, Probit, Cloglog) chỉ cho phép phân tích mối quan hệ định hướng xác suất mà chưa giải thích được tác động cận biên của từng biến giải thích trong các ngữ cảnh khác nhau. Việc chưa đưa vào các biến tương tác hoặc kiểm soát thêm các yếu tố khác như kinh nghiệm làm việc, ngành nghề cụ thể, vị trí địa lý hoặc tình trạng hôn nhân có thể dẫn đến sai lệch do biến bị thiếu.
Thứ hai, nghiên cứu chỉ sử dụng dữ liệu cắt ngang từ một thời điểm cố định, do đó không thể rút ra kết luận về quan hệ nhân quả theo thời gian. Thu nhập là một biến kinh tế có thể thay đổi mạnh mẽ theo chu kỳ đời sống và điều kiện kinh tế vĩ mô, điều mà dữ liệu cắt ngang không phản ánh được.
Thứ ba, độ chính xác dự báo của các mô hình vẫn còn hạn chế, đặc biệt ở mô hình Probit tuy có độ nhạy cao nhưng độ chính xác tổng thể còn thấp (42.8%), cho thấy khả năng phân loại nhóm thu nhập thấp và cao chưa tối ưu.
Do đó, hướng nghiên cứu tiếp theo nên mở rộng theo các hướng sau:
Mở rộng bộ biến: Bổ sung thêm các yếu tố như thâm niên làm việc, vị trí công tác, lĩnh vực nghề nghiệp, kỹ năng mềm và kỹ năng số, hoặc biến định lượng về chất lượng giáo dục để mô hình có khả năng giải thích mạnh hơn.
Áp dụng mô hình nâng cao: Ngoài mô hình nhị phân truyền thống, có thể xem xét các mô hình phi tuyến nâng cao như mô hình hồi quy bán tham số (semiparametric), mô hình cây quyết định (decision tree), hoặc các kỹ thuật học máy nhằm cải thiện độ chính xác phân loại.
Phân tích dị biệt theo nhóm: Nghiên cứu sâu hơn sự khác biệt theo giới tính, chủng tộc, khu vực cư trú (nông thôn – thành thị), hoặc kết hợp phân tầng theo nhóm ngành nghề để đánh giá tác động không đồng nhất của các yếu tố giải thích.
Những hướng phát triển này không chỉ giúp cải thiện mô hình dự báo mà còn góp phần nâng cao hiểu biết về cơ chế hình thành thu nhập trong thị trường lao động hiện đại.