Nhiệm vụ tuần 4:
Nhiệm vụ tuần 5:
Bộ dữ liệu Adult (hay Census Income 1994) ghi lại thông tin nhân khẩu học và việc làm của 48842 cá nhân ở Mỹ, với 14 đặc trưng bao gồm các biến phân loại (như workclass, education, marital-status, occupation…) và các biến số (như age, education-num, hours-per-week…)
Nhiệm vụ phân loại chính là dự đoán xem thu nhập hàng năm của một người có vượt trên $50.000 hay không không dựa trên dữ liệu điều tra dân số. Đây còn được gọi là tập dữ liệu “Thu nhập điều tra dân số”.
Theo thông tin từ UCI Machine Learning Repository, các bản ghi được trích xuất từ dữ liệu điều tra dân số năm 1994 bởi Barry Becker (với điều kiện lọc: tuổi >16, thu nhập cá nhân >100, trọng số dân số >1, giờ làm >0).
Đặc điểm của tập dữ liệu: Đa biến
Lĩnh vực chủ đề: Khoa học xã hội
Loại tính năng: Phân loại, Số nguyên
Các trường hợp: 48842
Các biến: Gồm 14 biến:
Đọc file: Đầu tiên ta đọc dữ liệu adult.csv vào R, xác định cột và loại bỏ các giá trị thiếu được mã hoá bằng ký tự “?”.
adult <- read.csv("C:/Users/HP/Downloads/adult.csv", header = TRUE, stringsAsFactors = FALSE)
colnames(adult) <- c("age","workclass","fnlwgt","education","education_num",
"marital_status","occupation","relationship","race","sex",
"capital_gain","capital_loss","hours_per_week",
"native_country","income")
Cấu trúc dữ liệu: Sử dụng str() để xem thông tin cấu trúc. Ta thấy tổng số hàng và kiểu dữ liệu của mỗi cột.
str(adult)
## 'data.frame': 32560 obs. of 15 variables:
## $ age : int 50 38 53 28 37 49 52 31 42 37 ...
## $ workclass : chr " Self-emp-not-inc" " Private" " Private" " Private" ...
## $ fnlwgt : int 83311 215646 234721 338409 284582 160187 209642 45781 159449 280464 ...
## $ education : chr " Bachelors" " HS-grad" " 11th" " Bachelors" ...
## $ education_num : int 13 9 7 13 14 5 9 14 13 10 ...
## $ marital_status: chr " Married-civ-spouse" " Divorced" " Married-civ-spouse" " Married-civ-spouse" ...
## $ occupation : chr " Exec-managerial" " Handlers-cleaners" " Handlers-cleaners" " Prof-specialty" ...
## $ relationship : chr " Husband" " Not-in-family" " Husband" " Wife" ...
## $ race : chr " White" " White" " Black" " Black" ...
## $ sex : chr " Male" " Male" " Male" " Female" ...
## $ capital_gain : int 0 0 0 0 0 0 0 14084 5178 0 ...
## $ capital_loss : int 0 0 0 0 0 0 0 0 0 0 ...
## $ hours_per_week: int 13 40 40 40 40 16 45 50 40 80 ...
## $ native_country: chr " United-States" " United-States" " United-States" " Cuba" ...
## $ income : chr " <=50K" " <=50K" " <=50K" " <=50K" ...
Thăm dò nhanh: Hiển thị một vài dòng đầu/cuối để hiểu sơ bộ dữ liệu.
head(adult)
## age workclass fnlwgt education education_num marital_status
## 1 50 Self-emp-not-inc 83311 Bachelors 13 Married-civ-spouse
## 2 38 Private 215646 HS-grad 9 Divorced
## 3 53 Private 234721 11th 7 Married-civ-spouse
## 4 28 Private 338409 Bachelors 13 Married-civ-spouse
## 5 37 Private 284582 Masters 14 Married-civ-spouse
## 6 49 Private 160187 9th 5 Married-spouse-absent
## occupation relationship race sex capital_gain capital_loss
## 1 Exec-managerial Husband White Male 0 0
## 2 Handlers-cleaners Not-in-family White Male 0 0
## 3 Handlers-cleaners Husband Black Male 0 0
## 4 Prof-specialty Wife Black Female 0 0
## 5 Exec-managerial Wife White Female 0 0
## 6 Other-service Not-in-family Black Female 0 0
## hours_per_week native_country income
## 1 13 United-States <=50K
## 2 40 United-States <=50K
## 3 40 United-States <=50K
## 4 40 Cuba <=50K
## 5 40 United-States <=50K
## 6 16 Jamaica <=50K
tail(adult)
## age workclass fnlwgt education education_num marital_status
## 32555 22 Private 310152 Some-college 10 Never-married
## 32556 27 Private 257302 Assoc-acdm 12 Married-civ-spouse
## 32557 40 Private 154374 HS-grad 9 Married-civ-spouse
## 32558 58 Private 151910 HS-grad 9 Widowed
## 32559 22 Private 201490 HS-grad 9 Never-married
## 32560 52 Self-emp-inc 287927 HS-grad 9 Married-civ-spouse
## occupation relationship race sex capital_gain
## 32555 Protective-serv Not-in-family White Male 0
## 32556 Tech-support Wife White Female 0
## 32557 Machine-op-inspct Husband White Male 0
## 32558 Adm-clerical Unmarried White Female 0
## 32559 Adm-clerical Own-child White Male 0
## 32560 Exec-managerial Wife White Female 15024
## capital_loss hours_per_week native_country income
## 32555 0 40 United-States <=50K
## 32556 0 38 United-States <=50K
## 32557 0 40 United-States >50K
## 32558 0 40 United-States <=50K
## 32559 0 20 United-States <=50K
## 32560 0 40 United-States >50K
Giá trị thiếu: Kiểm tra các biến định tính (categorical) xem có NA không. Ta phát hiện có NA ở các cột workclass, occupation, native_country (giá trị “?”). Ví dụ:
sapply(adult %>% select(workclass, occupation, native_country),
function(x) sum(is.na(x)))
## workclass occupation native_country
## 0 0 0
Vì tỷ lệ giá trị thiếu rất nhỏ so với kích thước mẫu (ví dụ chỉ khoảng 6% cho workclass), ta lựa chọn loại bỏ những dòng có NA để đảm bảo tính chính xác của phân tích và tránh ảnh hưởng không mong muốn. Sau đó, ta chuyển các biến định tính sang factor để tiện xử lý.
# Loại bỏ khoảng trắng ở đầu/cuối các biến dạng chuỗi
adult <- adult %>% mutate(across(where(is.character), trimws))
# Thay thế "?" bằng NA
adult[adult == "?"] <- NA
# Loại bỏ các dòng có NA
adult <- na.omit(adult)
adult <- na.omit(adult)
adult <- adult %>%
mutate_at(vars(workclass, education, marital_status, occupation,
relationship, race, sex, native_country, income),
factor)
Adult.csv là một bộ dữ liệu bao gồm 32561 quan sát và 14 biến, cụ thể là các biến:
Nhóm biến | Biến cụ thể | Mục tiêu phân tích chính |
---|---|---|
Biến định tính | workclass , education ,
marital-status , occupation ,
relationship , race , sex ,
native-country |
Thống kê tần số, biểu đồ cột, phân phối theo nhóm |
Biến định lượng | age , fnlwgt , education-num ,
capital-gain , capital-loss ,
hours-per-week |
Mô tả thống kê, phân phối, kiểm định, khoảng tin cậy |
Biến mục tiêu (target) | income |
Dự đoán phân loại: <=50K hoặc >50K |
Các biến định tính:
→ Phân loại theo khu vực làm việc của cá nhân, phản ánh nơi họ đang công tác.
Các giá trị phổ biến:
Private: khu vực tư nhân
Self-emp-not-inc: tự làm không thành lập công ty
Self-emp-inc: tự làm có công ty
Federal-gov, Local-gov, State-gov: làm việc cho chính phủ
Without-pay: làm không lương
Never-worked: chưa từng làm việc
→ Trình độ học vấn cuối cùng mà cá nhân đạt được (dạng văn bằng).
Ví dụ:
Bachelors, HS-grad (tốt nghiệp THPT), Some-college, Masters, Doctorate, 11th, Assoc-acdm,…
→ Phản ánh tình trạng hôn nhân hiện tại của người đó.
Giá trị gồm:
Married-civ-spouse: đã kết hôn (hôn nhân dân sự)
Divorced: đã ly hôn
Never-married: chưa từng kết hôn
Separated: ly thân
Widowed: góa vợ/chồng
→ Ngành nghề cụ thể đang làm việc.
Ví dụ:
Tech-support, Craft-repair, Sales, Exec-managerial, Prof-specialty, Handlers-cleaners,…
→ Vị trí của cá nhân trong gia đình.
Giá trị có thể là:
Husband, Wife, Own-child, Not-in-family, Unmarried, Other-relative
→ Chủng tộc được khai báo của cá nhân.
Gồm các nhóm chính:
White, Black, Asian-Pac-Islander, Amer-Indian-Eskimo, Other
→ Giới tính sinh học của cá nhân.
Gồm:
Male, Female
→ Quốc gia nơi cá nhân sinh ra hoặc có nguồn gốc.
Giá trị ví dụ:
United-States, Mexico, Philippines, Germany, Vietnam, India, China,…
Các biến định lượng:
→ Tuổi của cá nhân, giá trị là số nguyên dương.
Phạm vi dữ liệu: từ 17 đến 90+
→ Trọng số mẫu do Cục Điều tra dân số gán, thể hiện số lượng người trong tổng thể mà mẫu đại diện.
Giải thích:
Một người có fnlwgt = 100000 đại diện cho ~100,000 người trong dân số thực.
Biến này chủ yếu dùng trong phân tích khảo sát, ít khi được dùng trực tiếp trong mô hình dự đoán.
→ Là số năm học chính thức tương ứng với biến education.
Ví dụ:
HS-grad ≈ 9, Bachelors ≈ 13, Masters ≈ 14, Doctorate ≈ 16
→ Lợi nhuận vốn (capital gain) từ các nguồn như bán tài sản, đầu tư,… trong năm.
Đặc điểm:
Phân phối rất lệch phải, phần lớn giá trị là 0 (tức không có lãi vốn trong năm).
→ Tổn thất vốn từ các khoản đầu tư tài sản trong năm.
Lưu ý:
Tương tự capital-gain, giá trị chủ yếu là 0 trong phần lớn trường hợp.
→ Số giờ trung bình mà cá nhân làm việc mỗi tuần.
Phổ biến: 40 giờ/tuần là mức tiêu chuẩn, một số người làm 60–80 giờ hoặc ít hơn 20 giờ (bán thời gian).
Cho mỗi biến định tính, ta thực hiện thống kê tần suất và vẽ biểu đồ. Dưới đây là phân tích cho từng biến chính:
table(adult$workclass)
##
## Federal-gov Local-gov Private Self-emp-inc
## 943 2067 22286 1074
## Self-emp-not-inc State-gov Without-pay
## 2499 1278 14
prop.table(table(adult$workclass))*100
##
## Federal-gov Local-gov Private Self-emp-inc
## 3.12655416 6.85322105 73.89012301 3.56088989
## Self-emp-not-inc State-gov Without-pay
## 8.28553430 4.23726004 0.04641756
ggplot(adult, aes(x=workclass)) +
geom_bar(fill="#69b3a2") +
xlab("Workclass") + ylab("Số lượng") +
ggtitle("Biểu đồ cột của Workclass") +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
Nhận xét:
Hạng mục ‘Private’ chiếm ưu thế rõ rệt (chiếm 73% tổng dữ liệu), lớn hơn nhiều so với các hạng mục khác. Các hạng mục như ‘Self-emp-not-inc’, ‘Local-gov’, ‘State-gov’ cũng xuất hiện nhưng với tần suất thấp hơn. Cho thấy rằng đa số cá nhân trong dữ liệu làm việc trong khu vực tư nhân.
table(adult$education)
##
## 10th 11th 12th 1st-4th 5th-6th 7th-8th
## 820 1048 377 151 288 557
## 9th Assoc-acdm Assoc-voc Bachelors Doctorate HS-grad
## 455 1008 1307 5043 375 9840
## Masters Preschool Prof-school Some-college
## 1627 45 542 6678
prop.table(table(adult$education))*100
##
## 10th 11th 12th 1st-4th 5th-6th 7th-8th
## 2.7187427 3.4746859 1.2499586 0.5006465 0.9548755 1.8467557
## 9th Assoc-acdm Assoc-voc Bachelors Doctorate HS-grad
## 1.5085707 3.3420643 4.3334107 16.7202679 1.2433275 32.6249130
## Masters Preschool Prof-school Some-college
## 5.3943835 0.1491993 1.7970226 22.1411757
ggplot(adult, aes(x=education)) +
geom_bar(fill="#404080") +
xlab("Education") + ylab("Số lượng") +
ggtitle("Biểu đồ cột của trình độ học vấn") +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
Nhận xét:
Các mức học vấn phổ biến nhất là ‘HS-grad’ (tốt nghiệp THPT), ‘Some-college’ và ‘Bachelors’. Những nhóm này chiếm phần lớn dữ liệu. Các mức cao hơn như ‘Masters’, ‘Doctorate’ ít phổ biến hơn. Biểu đồ cho thấy phần lớn người trong mẫu có trình độ từ trung học đến cử nhân.
Nhóm tốt nghiệp THPT (HS-grad) chiếm đa số, sau đó là đang học Cao đẳng/Đại học (Some-college). Điều này phản ánh thực trạng giáo dục: phần đông dân số chỉ học đến cấp 3 hoặc chưa hoàn thành đại học. Ví dụ, ở Việt Nam cũng có xu hướng tỷ lệ người có bằng tốt nghiệp THPT trở lên ngày càng tăng, nhưng đại học vẫn chiếm một phần nhỏ hơn.
table(adult$marital_status)
##
## Divorced Married-AF-spouse Married-civ-spouse
## 4214 21 14065
## Married-spouse-absent Never-married Separated
## 370 9725 939
## Widowed
## 827
prop.table(table(adult$marital_status))*100
##
## Divorced Married-AF-spouse Married-civ-spouse
## 13.97168529 0.06962634 46.63306920
## Married-spouse-absent Never-married Separated
## 1.22674978 32.24362587 3.11329200
## Widowed
## 2.74195153
ggplot(adult, aes(x=marital_status)) +
geom_bar(fill="#008080") +
xlab("Marital Status") + ylab("Số lượng") +
ggtitle("Biểu đồ cột của tình trạng hôn nhân") +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
Nhận xét:
Hạng mục ‘Married-civ-spouse’ (kết hôn chính thức) chiếm tỷ lệ cao nhất (khoảng 46%). Các hạng mục ‘Never-married’ và ‘Divorced’ đứng sau, lần lượt chiếm khoảng 32% và 14%. Có sự chênh lệch rõ: người đã kết hôn chiếm ưu thế, phần lớn là người đã lập gia đình.
Hơn 46% đã kết hôn hợp pháp. Trong thực tế, những người đã kết hôn thường có xu hướng ổn định công việc và thu nhập hơn.
table(adult$occupation)
##
## Adm-clerical Armed-Forces Craft-repair Exec-managerial
## 3720 9 4030 3992
## Farming-fishing Handlers-cleaners Machine-op-inspct Other-service
## 989 1350 1966 3212
## Priv-house-serv Prof-specialty Protective-serv Sales
## 143 4038 644 3584
## Tech-support Transport-moving
## 912 1572
prop.table(table(adult$occupation))*100
##
## Adm-clerical Armed-Forces Craft-repair Exec-managerial
## 12.33380856 0.02983986 13.36162594 13.23563542
## Farming-fishing Handlers-cleaners Machine-op-inspct Other-service
## 3.27906900 4.47597891 6.51835151 10.64951427
## Priv-house-serv Prof-specialty Protective-serv Sales
## 0.47412221 13.38815026 2.13520772 11.88289513
## Tech-support Transport-moving
## 3.02377242 5.21202878
ggplot(adult, aes(x=occupation)) +
geom_bar(fill="#804000") +
xlab("Occupation") + ylab("Số lượng") +
ggtitle("Biểu đồ cột của nghề nghiệp") +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
Nhận xét:
Các nghề phổ biến gồm ‘Prof-specialty’, ‘Craft-repair’, ‘Exec-managerial’, ‘Adm-clerical’, ‘Sales’… với số lượng tương đương nhau. Không có một nghề cụ thể nào quá áp đảo. Điều này cho thấy đa dạng về nghề nghiệp, chủ yếu là các công việc tay nghề và quản lý.
Nhóm chuyên viên (Exec-managerial, Prof-specialty) chiếm tỷ lệ lớn (~sự tổng của Exec-managerial và Prof-specialty là khoảng 13.3% + 13.4% = 26,7%). Cho thấy đa số những người có vị trí công việc chuyên môn/cao chiếm phần lớn.
table(adult$relationship)
##
## Husband Not-in-family Other-relative Own-child Unmarried
## 12463 7725 889 4466 3212
## Wife
## 1406
prop.table(table(adult$relationship))*100
##
## Husband Not-in-family Other-relative Own-child Unmarried
## 41.321574 25.612546 2.947515 14.807201 10.649514
## Wife
## 4.661649
ggplot(adult, aes(x=relationship)) +
geom_bar(fill="#808080") +
xlab("Relationship") + ylab("Số lượng") +
ggtitle("Biểu đồ cột của quan hệ gia đình") +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
Nhận xét:
Hạng mục ‘Husband’ (nam đã kết hôn) chiếm tỷ lệ lớn nhất (~41%), tiếp theo là ‘Not-in-family’ (~25%). ‘Own-child’ chiếm khoảng 15%. Cho thấy, nam giới đã kết hôn chiếm tỉ lệ nhiều hơn. Các mối quan hệ gia đình khác như vợ, con cái, họ hàng chiếm phần nhỏ.
table(adult$race)
##
## Amer-Indian-Eskimo Asian-Pac-Islander Black Other
## 286 895 2817 231
## White
## 25932
prop.table(table(adult$race))*100
##
## Amer-Indian-Eskimo Asian-Pac-Islander Black Other
## 0.9482444 2.9674082 9.3398760 0.7658897
## White
## 85.9785816
ggplot(adult, aes(x=race)) +
geom_bar(fill="#C04000") +
xlab("Race") + ylab("Số lượng") +
ggtitle("Biểu đồ cột của chủng tộc") +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
Nhận xét:
Đa số là ‘White’ (khoảng 86%), ‘Black’ ~9%, các chủng tộc khác (Á châu, da đỏ, khác) chỉ chiếm vài phần trăm. Sự chênh lệch rất lớn. Cho thấy chủ yếu là người da trắng, tương ứng với đặc trưng dân số của dữ liệu gốc.
table(adult$sex)
##
## Female Male
## 9782 20379
prop.table(table(adult$sex))*100
##
## Female Male
## 32.43261 67.56739
ggplot(adult, aes(x=sex)) +
geom_bar(fill="#008080") +
xlab("Sex") + ylab("Số lượng") +
ggtitle("Biểu đồ cột của giới tính") +
theme(axis.text.x = element_text(angle = 0, hjust = 0.5))
Nhận xét:
Mẫu gồm khoảng 67% nam và 32% nữ. Nam giới chiếm ưu thế rõ rệt. Thực tế, số lao động nam thường nhiều hơn nữ do nhiều yếu tố (chủ quan, truyền thống). Tuy nhiên, dữ liệu này không cân bằng giới tính.
table(adult$native_country)
##
## Cambodia Canada
## 18 107
## Columbia Cuba
## 56 92
## China Dominican-Republic
## 68 67
## Ecuador El-Salvador
## 27 100
## England France
## 86 27
## Germany Greece
## 128 29
## Guatemala Haiti
## 63 42
## Holand-Netherlands Honduras
## 1 12
## Hong Hungary
## 19 13
## India Iran
## 100 42
## Ireland Italy
## 24 68
## Jamaica Japan
## 80 59
## Laos Mexico
## 17 610
## Nicaragua Outlying-US(Guam-USVI-etc)
## 33 14
## Peru Poland
## 30 56
## Portugal Puerto-Rico
## 34 109
## Philippines Scotland
## 188 11
## South Taiwan
## 71 42
## Thailand Trinadad&Tobago
## 17 18
## United-States Vietnam
## 27503 64
## Yugoslavia
## 16
prop.table(table(adult$native_country))*100
##
## Cambodia Canada
## 0.05967972 0.35476277
## Columbia Cuba
## 0.18567024 0.30502967
## China Dominican-Republic
## 0.22545672 0.22214118
## Ecuador El-Salvador
## 0.08951958 0.33155399
## England France
## 0.28513643 0.08951958
## Germany Greece
## 0.42438911 0.09615066
## Guatemala Haiti
## 0.20887902 0.13925268
## Holand-Netherlands Honduras
## 0.00331554 0.03978648
## Hong Hungary
## 0.06299526 0.04310202
## India Iran
## 0.33155399 0.13925268
## Ireland Italy
## 0.07957296 0.22545672
## Jamaica Japan
## 0.26524319 0.19561686
## Laos Mexico
## 0.05636418 2.02247936
## Nicaragua Outlying-US(Guam-USVI-etc)
## 0.10941282 0.04641756
## Peru Poland
## 0.09946620 0.18567024
## Portugal Puerto-Rico
## 0.11272836 0.36139385
## Philippines Scotland
## 0.62332151 0.03647094
## South Taiwan
## 0.23540334 0.13925268
## Thailand Trinadad&Tobago
## 0.05636418 0.05967972
## United-States Vietnam
## 91.18729485 0.21219456
## Yugoslavia
## 0.05304864
ggplot(adult %>% filter(native_country != "United-States"),
aes(x=native_country)) +
geom_bar(fill="#804080") +
xlab("Native Country (ngoại trừ Mỹ)") + ylab("Số lượng") +
ggtitle("Biểu đồ cột của quốc tịch gốc (không bao gồm Mỹ)") +
theme(axis.text.x = element_text(angle = 90, hjust = 1, size=8))
Nhận xét:
Đa số (gần 89%) có nguồn gốc từ ‘United-States’, phần còn lại rải rác ở nhiều quốc gia khác (Mexico, Philippines, Đức…). Sự chênh lệch rất lớn cho thấy tập mẫu chủ yếu là người Mỹ.
Hơn 89% là người Mỹ. Điều này hợp lý vì dữ liệu từ cuộc điều tra dân số Mỹ, chỉ một phần nhỏ thuộc quốc gia khác (Mexiko, Philippines, v.v.).
Theo báo cáo của Cục Điều tra dân số Mỹ năm 2022, độ tuổi trung bình của dân số Mỹ đạt 38.9 tuổi. Trong bộ dữ liệu Adult (1994), độ tuổi trung bình (age) tính được khoảng 38.6 tuổi (độ lệch chuẩn khoảng 13.6), với giá trị trung vị khoảng 36 tuổi. Như vậy, mẫu dữ liệu có độ tuổi tương đương trung bình dân số hiện tại. Ví dụ, đoạn mã R sau cho biết thống kê cơ bản về tuổi và giờ làm việc:
summary(adult$age)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 17.00 28.00 37.00 38.44 47.00 90.00
summary(adult$hours_per_week)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 1.00 40.00 40.00 40.93 45.00 99.00
Nhận xét:
Độ tuổi (age): trung bình ≈ 38.44, trung vị 37, phần vị thứ nhất ~28, phần vị thứ ba ~47.
Giờ làm việc mỗi tuần (hours_per_week): trung bình ≈ 40.93 giờ, trung vị 40 giờ, phần vị thứ nhất 40, phần vị thứ ba ~45. Hầu hết người lao động làm việc xung quanh 40 giờ mỗi tuần, tương tự mức trung bình gần 40 giờ/tuần trong lực lượng lao động Mỹ hiện nay.
Chọn 3 biến định tính tiêu biểu: Sex (Female), Workclass (Private), Marital-status (Married-civ-spouse). Với mỗi biến, xem xét một hạng mục quan tâm và phân tích:
Ước lượng khoảng tin cậy 95% cho tỉ lệ nữ trong tổng thể:
n_total <- nrow(adult)
x_female <- sum(adult$sex == "Female")
prop.test(x = x_female, n = n_total, conf.level = 0.95)
##
## 1-sample proportions test with continuity correction
##
## data: x_female out of n_total, null probability 0.5
## X-squared = 3722.5, df = 1, p-value < 2.2e-16
## alternative hypothesis: true p is not equal to 0.5
## 95 percent confidence interval:
## 0.3190492 0.3296479
## sample estimates:
## p
## 0.3243261
Nhận xét:
Khoảng tin cậy ước tính cho tỷ lệ nữ khoảng [0.319, 0.329].
Vậy, với mức ý nghĩa 95% thì tỷ lệ nữ trong toàn bộ dân số tương ứng nằm trong khoảng [0.319, 0.329].
Kiểm định giả thuyết: Giả sử 𝐻0:𝑝=0.5 (tỷ lệ nữ = 50%). Thực hiện kiểm định:
prop.test(x = x_female, n = n_total, alternative = "two.sided", p = 0.5, conf.level = 0.95)
##
## 1-sample proportions test with continuity correction
##
## data: x_female out of n_total, null probability 0.5
## X-squared = 3722.5, df = 1, p-value < 2.2e-16
## alternative hypothesis: true p is not equal to 0.5
## 95 percent confidence interval:
## 0.3190492 0.3296479
## sample estimates:
## p
## 0.3243261
Nhận xét:
Kết quả: Giá trị p rất nhỏ (<2.2e-16) nên bác bỏ 𝐻0 ở mức ý nghĩa 0.05.
Kết luận: Tỷ lệ nữ trong mẫu không bằng 50%. Thực tế, mẫu này có tỷ lệ nữ ~33.5% (nhỏ hơn 50%). Kết luận có ý nghĩa là tỷ lệ nữ thực tế thấp hơn giả thuyết ban đầu.
Ước lượng khoảng tin cậy 95%:
x_private <- sum(adult$workclass == "Private")
n_total <- nrow(adult)
prop.test(x = x_private, n = n_total, conf.level = 0.95)
##
## 1-sample proportions test with continuity correction
##
## data: x_private out of n_total, null probability 0.5
## X-squared = 6884.7, df = 1, p-value < 2.2e-16
## alternative hypothesis: true p is not equal to 0.5
## 95 percent confidence interval:
## 0.7338973 0.7438441
## sample estimates:
## p
## 0.7389012
Nhận xét:
Khoảng tin cậy ước tính cho tỷ lệ ‘Private’ là [0.733, 0.743].
Vậy, với mức ý nghĩa 95% thì tỷ lệ người có công việc thuộc diện ‘Private’ trong tổng thể nằm trong khoảng [0.733, 0.743].
Kiểm định giả thuyết: Kiểm định 𝐻0:𝑝=0.5 (nghĩa là giả định ban đầu là chỉ 50% là private).
prop.test(x = x_private, n = n_total, alternative = "two.sided", p = 0.5, conf.level = 0.95)
##
## 1-sample proportions test with continuity correction
##
## data: x_private out of n_total, null probability 0.5
## X-squared = 6884.7, df = 1, p-value < 2.2e-16
## alternative hypothesis: true p is not equal to 0.5
## 95 percent confidence interval:
## 0.7338973 0.7438441
## sample estimates:
## p
## 0.7389012
Nhận xét
Kết quả: p-value ≈ 0 (<2.2e-16). Bác bỏ 𝐻0 .
Vậy, Tỷ lệ người làm trong khu vực “Private” cao hơn giả thiết 50%. Thực tế mẫu cho thấy ~70% làm việc trong khu vực tư nhân, chênh lệch lớn so với giả thuyết ban đầu.
Ước lượng khoảng tin cậy 95%:
x_married <- sum(adult$marital_status == "Married-civ-spouse")
prop.test(x = x_married, n = n_total, conf.level = 0.95)
##
## 1-sample proportions test with continuity correction
##
## data: x_married out of n_total, null probability 0.5
## X-squared = 136.63, df = 1, p-value < 2.2e-16
## alternative hypothesis: true p is not equal to 0.5
## 95 percent confidence interval:
## 0.4606888 0.4719812
## sample estimates:
## p
## 0.4663307
Nhận xét
Khoảng tin cậy cho tỷ lệ ‘Married-civ-spouse’ là khoảng [0.460, 0.471]
Vậy với mức ý nghĩa 95% thì tỷ lệ đã kết hôn theo khảo sát nằm trong khoảng [0.460, 0.471].
Kiểm định giả thuyết: Kiểm định 𝐻0:𝑝=0.5
prop.test(x = x_married, n = n_total, alternative = "two.sided", p = 0.5, conf.level = 0.95)
##
## 1-sample proportions test with continuity correction
##
## data: x_married out of n_total, null probability 0.5
## X-squared = 136.63, df = 1, p-value < 2.2e-16
## alternative hypothesis: true p is not equal to 0.5
## 95 percent confidence interval:
## 0.4606888 0.4719812
## sample estimates:
## p
## 0.4663307
Nhận xét
Kết quả: p-value rất nhỏ,(<2.2e-16) bác bỏ 𝐻0. Nghĩa là tỷ lệ đã kết hôn trong mẫu (~46%) khác với giả thuyết 50%. Điều này phù hợp với nhận xét trước đó rằng ~46% đã lập gia đình, kém giả thiết 50%.
Chọn 3 cặp biến định tính liên quan đến thu nhập (income) và các biến khác: Sex - Income, Education - Income, Marital-status - Income. Với mỗi cặp, ta lập bảng chéo, vẽ biểu đồ, và kiểm định Chi-bình phương độc lập.
Bảng tần suất chéo:
tab_sex_inc <- table(adult$sex, adult$income)
tab_sex_inc
##
## <=50K >50K
## Female 8670 1112
## Male 13983 6396
prop.table(tab_sex_inc, 1)*100 # Tỷ lệ theo hàng giới tính
##
## <=50K >50K
## Female 88.63218 11.36782
## Male 68.61475 31.38525
Tỷ lệ có thu nhập >50K: ~31.38% ở nam, ~11.36% ở nữ.
Trực quan hóa: Biểu đồ cột chồng cho thấy sự khác biệt lớn.
ggplot(adult, aes(x=sex, fill=income)) +
geom_bar(position="fill") +
ylab("Tỷ lệ (%)") + ggtitle("Tỷ lệ thu nhập >50K theo giới tính") +
theme(axis.text.x = element_text(hjust=0.5))
Nhận xét mô tả: Biểu đồ và tỷ lệ cho thấy nam giới có tỷ lệ thu nhập >50K cao hơn rất nhiều so với nữ.
Kiểm định Chi-bình phương:
H0: Giới tính và thu nhập độc lập (không liên quan) trong tổng thể.
H1: Có mối liên hệ giữa giới tính và thu nhập.
chisq_test_sex_inc <- chisq.test(tab_sex_inc)
chisq_test_sex_inc$statistic; chisq_test_sex_inc$parameter; chisq_test_sex_inc$p.value
## X-squared
## 1415.45
## df
## 1
## [1] 9.226706e-310
Nhận xét:
Kết quả Giá trị X bình phương lớn, p-value rất nhỏ. Kết luận: Bác bỏ H0 (theo mức 0.05). Giới tính có liên quan có ý nghĩa thống kê với việc kiếm được >50K, tức tỷ lệ thu nhập cao phụ thuộc vào giới tính. Thực tế cho thấy người đàn ông là trụ cột của gia đình nên có xu hướng làm việc nhiều người phụ nữ vì vậy việc đàn ông có thu nhập cao hơn so với phụ nữ là hợp lí.
Risk Ratio
rr_sex_income <- riskratio(tab_sex_inc)
rr_sex_income
## $data
##
## <=50K >50K Total
## Female 8670 1112 9782
## Male 13983 6396 20379
## Total 22653 7508 30161
##
## $measure
## risk ratio with 95% C.I.
## estimate lower upper
## Female 1.000000 NA NA
## Male 2.760886 2.602862 2.928504
##
## $p.value
## two-sided
## midp.exact fisher.exact chi.square
## Female NA NA NA
## Male 0 0 5.400534e-310
##
## $correction
## [1] FALSE
##
## attr(,"method")
## [1] "Unconditional MLE & normal approximation (Wald) CI"
Nhận xét:
Ở đây, ta lấy Female là nhóm tham chiếu và Male là nhóm so sánh với sự kiện là “thu nhập > 50K”
Kết quả RR của Male = 2.760 cho thấy Nam có xác suất có thu nhập (>50K) gấp 2.76 lần so với nữ giới.
Odds Ratio
or_sex_income <- oddsratio(tab_sex_inc)
or_sex_income
## $data
##
## <=50K >50K Total
## Female 8670 1112 9782
## Male 13983 6396 20379
## Total 22653 7508 30161
##
## $measure
## odds ratio with 95% C.I.
## estimate lower upper
## Female 1.000000 NA NA
## Male 3.565833 3.329126 3.822584
##
## $p.value
## two-sided
## midp.exact fisher.exact chi.square
## Female NA NA NA
## Male 0 0 5.400534e-310
##
## $correction
## [1] FALSE
##
## attr(,"method")
## [1] "median-unbiased estimate & mid-p exact CI"
Nhận xét:
Ở đây, ta lấy Female là nhóm tham chiếu và Male là nhóm so sánh với sự kiện là “thu nhập > 50K”
Kết quả OR của Male = 3.565 cho thấy khả năng nam có thu nhập (>50K) gấp 3.565 lần so với nữ.
Bảng tần suất chéo:
tab_edu_inc <- table(adult$education, adult$income)
prop.table(tab_edu_inc, 1)*100 # Tỷ lệ theo hàng (trình độ)
##
## <=50K >50K
## 10th 92.804878 7.195122
## 11th 94.370229 5.629771
## 12th 92.307692 7.692308
## 1st-4th 96.026490 3.973510
## 5th-6th 95.833333 4.166667
## 7th-8th 93.716338 6.283662
## 9th 94.505495 5.494505
## Assoc-acdm 74.603175 25.396825
## Assoc-voc 73.680184 26.319816
## Bachelors 57.842554 42.157446
## Doctorate 25.333333 74.666667
## HS-grad 83.567073 16.432927
## Masters 43.577136 56.422864
## Preschool 100.000000 0.000000
## Prof-school 25.092251 74.907749
## Some-college 79.994010 20.005990
Tỷ lệ thu nhập >50K tăng dần với trình độ: từ gần 0 ở trình độ thấp đến ~50% ở trình độ cao (Doctorate, Prof-school).
Trực quan hóa:
ggplot(adult, aes(x=education, fill=income)) +
geom_bar(position="fill") +
xlab("Education") + ylab("Tỷ lệ (%)") +
ggtitle("Tỷ lệ >50K theo trình độ học vấn") +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
Nhận xét mô tả: Có xu hướng rõ: những người có học vấn cao hơn (Bachelors, Masters, Doctorate) có khả năng thu nhập >50K cao hơn những người học vấn thấp (HS-grad, Some-college).
Kiểm định Chi-bình phương:
chisq_test_edu_inc <- chisq.test(tab_edu_inc)
chisq_test_edu_inc$statistic; chisq_test_edu_inc$parameter; chisq_test_edu_inc$p.value
## X-squared
## 4070.91
## df
## 15
## [1] 0
p-value rất nhỏ, bác bỏ H0. Kết luận: Có mối liên hệ ý nghĩa giữa trình độ học vấn và thu nhập. Người có học vấn cao hơn có tỷ lệ kiếm >50K cao hơn.
Risk Ratio
rr_edu_income <- riskratio(tab_edu_inc)
## Warning in chisq.test(xx, correct = correction): Chi-squared approximation may
## be incorrect
rr_edu_income
## $data
##
## <=50K >50K Total
## 10th 761 59 820
## 11th 989 59 1048
## 12th 348 29 377
## 1st-4th 145 6 151
## 5th-6th 276 12 288
## 7th-8th 522 35 557
## 9th 430 25 455
## Assoc-acdm 752 256 1008
## Assoc-voc 963 344 1307
## Bachelors 2917 2126 5043
## Doctorate 95 280 375
## HS-grad 8223 1617 9840
## Masters 709 918 1627
## Preschool 45 0 45
## Prof-school 136 406 542
## Some-college 5342 1336 6678
## Total 22653 7508 30161
##
## $measure
## risk ratio with 95% C.I.
## estimate lower upper
## 10th 1.0000000 NA NA
## 11th 0.7824427 0.5518762 1.109337
## 12th 1.0691004 0.6972468 1.639270
## 1st-4th 0.5522505 0.2428100 1.256046
## 5th-6th 0.5790960 0.3159208 1.061507
## 7th-8th 0.8733226 0.5830196 1.308176
## 9th 0.7636431 0.4852308 1.201801
## Assoc-acdm 3.5297283 2.7009477 4.612819
## Assoc-voc 3.6580084 2.8148313 4.753757
## Bachelors 5.8591705 4.5725761 7.507776
## Doctorate 10.3774011 8.0594494 13.362011
## HS-grad 2.2838983 1.7790219 2.932056
## Masters 7.8418218 6.1102953 10.064026
## Preschool 0.0000000 0.0000000 NaN
## Prof-school 10.4109075 8.1031828 13.375855
## Some-college 2.7804935 2.1644777 3.571829
##
## $p.value
## two-sided
## midp.exact fisher.exact chi.square
## 10th NA NA NA
## 11th 1.701692e-01 1.801131e-01 1.675391e-01
## 12th 7.528002e-01 8.116336e-01 7.594894e-01
## 1st-4th 1.389358e-01 1.601445e-01 1.454889e-01
## 5th-6th 6.571020e-02 9.212282e-02 7.101066e-02
## 7th-8th 5.152743e-01 5.864350e-01 5.104060e-01
## 9th 2.425475e-01 2.888924e-01 2.409252e-01
## Assoc-acdm 0.000000e+00 2.469800e-26 1.201275e-24
## Assoc-voc 0.000000e+00 7.911420e-31 6.306374e-28
## Bachelors 0.000000e+00 5.336746e-100 3.460691e-82
## Doctorate 0.000000e+00 5.987125e-127 2.191448e-127
## HS-grad 2.930989e-14 4.048421e-14 2.910025e-12
## Masters 0.000000e+00 1.366989e-140 8.172091e-122
## Preschool 3.816099e-02 6.600945e-02 6.230967e-02
## Prof-school 0.000000e+00 9.345545e-157 9.988812e-147
## Some-college 0.000000e+00 2.098733e-22 5.759899e-19
##
## $correction
## [1] FALSE
##
## attr(,"method")
## [1] "Unconditional MLE & normal approximation (Wald) CI"
Bảng tần suất chéo:
tab_mar_inc <- table(adult$marital_status, adult$income)
prop.table(tab_mar_inc, 1)*100
##
## <=50K >50K
## Divorced 89.273849 10.726151
## Married-AF-spouse 52.380952 47.619048
## Married-civ-spouse 54.504088 45.495912
## Married-spouse-absent 91.621622 8.378378
## Never-married 95.167095 4.832905
## Separated 92.971246 7.028754
## Widowed 90.326481 9.673519
Nhóm ‘Married-civ-spouse’ có tỷ lệ >50K cao nhất (khoảng 60%), trong khi nhóm ‘Never-married’ thấp hơn (~12%).
Trực quan hóa:
ggplot(adult, aes(x=marital_status, fill=income)) +
geom_bar(position="fill") +
xlab("Marital Status") + ylab("Tỷ lệ (%)") +
ggtitle("Tỷ lệ >50K theo tình trạng hôn nhân") +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
Nhận xét mô tả: Rõ ràng người đã lập gia đình ‘Married-civ-spouse’ có khả năng thu nhập cao lớn hơn so với nhóm chưa cưới hoặc ly hôn.
Kiểm định Chi-bình phương:
chisq_test_mar_inc <- chisq.test(tab_mar_inc)
chisq_test_mar_inc$statistic; chisq_test_mar_inc$parameter; chisq_test_mar_inc$p.value
## X-squared
## 6061.295
## df
## 6
## [1] 0
Kết quả: p-value ≈ 0, bác bỏ H0. Kết luận: Có liên quan rõ ràng giữa tình trạng hôn nhân và thu nhập scribbr.com. Nhóm đã kết hôn có tỷ lệ thu nhập cao hơn đáng kể.
Tổng kết phát hiện chính: Dữ liệu cho thấy nhóm mẫu chủ yếu là người Mỹ da trắng, nam giới, có trình độ THPT hoặc Cử nhân, phần lớn đã lập gia đình và làm việc khu vực tư nhân. Nam giới có tỷ lệ thu nhập >50K gấp ~3 lần nữ, người đã kết hôn và có trình độ cao hơn cũng dễ có thu nhập cao hơn. Tóm lại, thu nhập lớn (>50K) có liên quan rõ rệt đến giới tính, tình trạng hôn nhân, giáo dục… như đã kiểm định (chi-square test) ở trên
Hạn chế phân tích: Bài phân tích chỉ tập trung vào biến định tính; chưa xem xét biến định lượng (tuổi, giờ làm). Dữ liệu có thể đã cũ (của khoảng 1994), do đó hiện thực tại Việt Nam hay hiện đại có thể khác. Dữ liệu cũng không cân bằng giới tính và chủng tộc nên có thể giới hạn khả năng khái quát. Ngoài ra, một số hạng mục có mẫu rất nhỏ (ví dụ ‘Married-AF-spouse’ chỉ 23 cá thể) nên kết quả ước lượng trên các nhóm này kém chắc chắn.
Đề xuất ứng dụng:
Chiến lược tiếp thị: Dựa vào phát hiện rằng nhóm nam, kết hôn, học vấn cao có thu nhập cao hơn, doanh nghiệp nên nhắm đến các nhóm này cho các sản phẩm cao cấp.
Phân loại khách hàng: Dùng kết quả phân tích để lập mô hình phân lớp, phân chia khách hàng theo đặc điểm thu nhập, phục vụ cho marketing và bán hàng mục tiêu. Ví dụ: nhắm quảng cáo sản phẩm cao cấp cho nhóm đã kết hôn và trình độ cao.
Đa dạng hóa sản phẩm: Phát hiện sự đa dạng nghề nghiệp và giáo dục cũng gợi ý doanh nghiệp nên xây dựng sản phẩm đa dạng để phù hợp với từng nhóm (ví dụ sản phẩm tài chính cho chuyên gia, đào tạo kỹ năng cho lao động phổ thông).
Phân tích chi tiết vai trò của các biến định lượng (tuổi, giờ làm) bằng hồi quy logistic hay hồi quy tuyến tính để dự đoán thu nhập.
Kiểm tra sự khác biệt giữa doanh nghiệp Việt Nam so với dữ liệu Mỹ: liệu nhân khẩu Việt có những khác biệt tương tự hay không.
Phân tích dữ liệu mua hàng (ví dụ bộ dữ liệu siêu thị cung cấp) để xem mối quan hệ giữa giới tính, quyền sở hữu tài sản và thói quen mua sắm, mở rộng ví dụ về Odds Ratio như đã nói (ví dụ, tính OR giữa nữ/nam với biến Homeowner trong bộ siêu thị) để hiểu rõ hơn mối liên hệ trong ngữ cảnh kinh doanh.
Nhiệm vụ tuần 2:
Nhiệm vụ tuần 3:
Đầu tiên, ta đọc dữ liệu và kiểm tra cấu trúc chung của bộ dữ liệu.
df <- read.csv("C:/Users/HP/Downloads/Supermarket Transactions.csv", header = TRUE, stringsAsFactors = FALSE)
# Kiểm tra cấu trúc dữ liệu và một số dòng đầu
str(df)
## 'data.frame': 14059 obs. of 16 variables:
## $ X : int 1 2 3 4 5 6 7 8 9 10 ...
## $ PurchaseDate : chr "2007-12-18" "2007-12-20" "2007-12-21" "2007-12-21" ...
## $ CustomerID : int 7223 7841 8374 9619 1900 6696 9673 354 1293 7938 ...
## $ Gender : chr "F" "M" "F" "M" ...
## $ MaritalStatus : chr "S" "M" "M" "M" ...
## $ Homeowner : chr "Y" "Y" "N" "Y" ...
## $ Children : int 2 5 2 3 3 3 2 2 3 1 ...
## $ AnnualIncome : chr "$30K - $50K" "$70K - $90K" "$50K - $70K" "$30K - $50K" ...
## $ City : chr "Los Angeles" "Los Angeles" "Bremerton" "Portland" ...
## $ StateorProvince : chr "CA" "CA" "WA" "OR" ...
## $ Country : chr "USA" "USA" "USA" "USA" ...
## $ ProductFamily : chr "Food" "Food" "Food" "Food" ...
## $ ProductDepartment: chr "Snack Foods" "Produce" "Snack Foods" "Snacks" ...
## $ ProductCategory : chr "Snack Foods" "Vegetables" "Snack Foods" "Candy" ...
## $ UnitsSold : int 5 5 3 4 4 3 4 6 1 2 ...
## $ Revenue : num 27.38 14.9 5.52 4.44 14 ...
head(df)
## X PurchaseDate CustomerID Gender MaritalStatus Homeowner Children
## 1 1 2007-12-18 7223 F S Y 2
## 2 2 2007-12-20 7841 M M Y 5
## 3 3 2007-12-21 8374 F M N 2
## 4 4 2007-12-21 9619 M M Y 3
## 5 5 2007-12-22 1900 F S Y 3
## 6 6 2007-12-22 6696 F M Y 3
## AnnualIncome City StateorProvince Country ProductFamily
## 1 $30K - $50K Los Angeles CA USA Food
## 2 $70K - $90K Los Angeles CA USA Food
## 3 $50K - $70K Bremerton WA USA Food
## 4 $30K - $50K Portland OR USA Food
## 5 $130K - $150K Beverly Hills CA USA Drink
## 6 $10K - $30K Beverly Hills CA USA Food
## ProductDepartment ProductCategory UnitsSold Revenue
## 1 Snack Foods Snack Foods 5 27.38
## 2 Produce Vegetables 5 14.90
## 3 Snack Foods Snack Foods 3 5.52
## 4 Snacks Candy 4 4.44
## 5 Beverages Carbonated Beverages 4 14.00
## 6 Deli Side Dishes 3 4.37
tail(df)
## X PurchaseDate CustomerID Gender MaritalStatus Homeowner Children
## 14054 14054 2009-12-29 2032 F M N 3
## 14055 14055 2009-12-29 9102 F M Y 2
## 14056 14056 2009-12-29 4822 F M Y 3
## 14057 14057 2009-12-31 250 M S Y 1
## 14058 14058 2009-12-31 6153 F S N 4
## 14059 14059 2009-12-31 3656 M S N 3
## AnnualIncome City StateorProvince Country ProductFamily
## 14054 $10K - $30K Yakima WA USA Non-Consumable
## 14055 $10K - $30K Bremerton WA USA Food
## 14056 $10K - $30K Walla Walla WA USA Food
## 14057 $30K - $50K Portland OR USA Drink
## 14058 $50K - $70K Spokane WA USA Drink
## 14059 $50K - $70K Portland OR USA Non-Consumable
## ProductDepartment ProductCategory UnitsSold Revenue
## 14054 Household Paper Products 5 14.50
## 14055 Baking Goods Baking Goods 3 9.64
## 14056 Frozen Foods Vegetables 3 7.45
## 14057 Beverages Pure Juice Beverages 4 3.24
## 14058 Dairy Dairy 2 4.00
## 14059 Household Electrical 5 25.53
Kết quả str(df) cho biết bộ dữ liệu có 14059 bản ghi và 16 biến. Các biến như Gender, MaritalStatus, Homeowner, Children, AnnualIncome, City, StateorProvince, Country, ProductFamily, ProductDepartment, ProductCategory là định tính (kiểu chuỗi). Biến UnitsSold và Revenue là định lượng.
Dữ liệu chứa các giao dịch mua hàng tại siêu thị theo từng đơn (bảng). Số dòng và số cột thể hiện quy mô và số biến quan sát. Sau khi loại bỏ cột chỉ mục thừa, ta kiểm tra cấu trúc. Dữ liệu có định dạng ngày tháng cho cột PurchaseDate từ cuối năm 2007 đến cuối năm 2009, cho thấy phạm vi gần 2 năm. Việc chuyển định dạng về Date cho phép phân tích xu hướng theo thời gian.
Supermarket Transactions là một bộ dữ liệu bao gồm 14059 quan sát và 16 biến, cụ thể là các biến:
Số thứ tự
Purchase Date: Ngày mua hàng
Customer ID: ID của khách hàng
Gender: Giới tính, với F (Female) là nữ và M (Male) là nam
Marital Status: Tình trạng hôn nhân, với S (Single) là độc thân và M(Married) là đã kết hôn
Homeowner: Đã có nhà hay chưa, với Y (Yes) là đã có nhà và N (No) là chưa có nhà
Children: Số con cái
Annual Income: Thu nhập hàng năm (được biểu thị dưới dạng các khoảng)
City: Thành phố đang sống
Stateor Province: Mã kí hiệu của bang
Country: Đất nước
Product Family: Nhóm sản phẩm, Food: Thực phẩm, Drink: Đồ uống và Non-Consumable: Hàng không tiêu dùng
Product Department: Nhóm sản phẩm chi tiết
Product Category: Danh mục sản phẩm
Units Sold: Doanh số bán hàng theo đơn vị
Revenue: Doanh thu
Ta cũng kiểm tra xem có giá trị thiếu (NA) nào không.
# Kiểm tra dữ liệu thiếu
sum(is.na(df))
## [1] 0
colSums(is.na(df))
## X PurchaseDate CustomerID Gender
## 0 0 0 0
## MaritalStatus Homeowner Children AnnualIncome
## 0 0 0 0
## City StateorProvince Country ProductFamily
## 0 0 0 0
## ProductDepartment ProductCategory UnitsSold Revenue
## 0 0 0 0
Kết quả sum(is.na(df)) cho thấy không có giá trị thiếu trong dữ liệu (tổng các giá trị NA bằng 0). Điều này cho thấy dữ liệu đã được thu thập đầy đủ và không cần xử lý dữ liệu thiếu thêm.
Tiếp theo, ta chuyển các biến định tính về dạng factor và đặt nhãn tiếng Việt cho một số biến cho dễ hiểu.
# Chuyển sang factor và gán nhãn tiếng Việt
df <- df %>%
mutate(
Gender = factor(Gender, levels = c("M","F"), labels = c("Nam","Nữ")),
MaritalStatus = factor(MaritalStatus, levels = c("S","M"), labels = c("Độc thân","Đã kết hôn")),
Homeowner = factor(Homeowner, levels = c("Y","N"), labels = c("Có","Không")),
AnnualIncome = factor(AnnualIncome, ordered = TRUE),
City = factor(City),
StateorProvince= factor(StateorProvince),
Country = factor(Country),
ProductFamily = factor(ProductFamily,
levels = c("Food","Drink","Non-Consumable"),
labels = c("Thực phẩm","Đồ uống","Phi tiêu dùng")),
ProductDepartment = factor(ProductDepartment),
ProductCategory = factor(ProductCategory)
)
str(df)
## 'data.frame': 14059 obs. of 16 variables:
## $ X : int 1 2 3 4 5 6 7 8 9 10 ...
## $ PurchaseDate : chr "2007-12-18" "2007-12-20" "2007-12-21" "2007-12-21" ...
## $ CustomerID : int 7223 7841 8374 9619 1900 6696 9673 354 1293 7938 ...
## $ Gender : Factor w/ 2 levels "Nam","Nữ": 2 1 2 1 2 2 1 2 1 1 ...
## $ MaritalStatus : Factor w/ 2 levels "Độc thân","Đã kết hôn": 1 2 2 2 1 2 1 2 2 1 ...
## $ Homeowner : Factor w/ 2 levels "Có","Không": 1 1 2 1 1 1 1 1 1 2 ...
## $ Children : int 2 5 2 3 3 3 2 2 3 1 ...
## $ AnnualIncome : Ord.factor w/ 8 levels "$10K - $30K"<..: 5 7 6 5 3 1 5 4 1 6 ...
## $ City : Factor w/ 23 levels "Acapulco","Bellingham",..: 8 8 4 12 3 3 13 23 2 15 ...
## $ StateorProvince : Factor w/ 10 levels "BC","CA","DF",..: 2 2 8 6 2 2 6 8 8 2 ...
## $ Country : Factor w/ 3 levels "Canada","Mexico",..: 3 3 3 3 3 3 3 3 3 3 ...
## $ ProductFamily : Factor w/ 3 levels "Thực phẩm","Đồ uống",..: 1 1 1 1 2 1 1 1 3 3 ...
## $ ProductDepartment: Factor w/ 22 levels "Alcoholic Beverages",..: 20 18 20 21 4 11 13 6 15 14 ...
## $ ProductCategory : Factor w/ 45 levels "Baking Goods",..: 42 45 42 7 15 41 5 13 16 35 ...
## $ UnitsSold : int 5 5 3 4 4 3 4 6 1 2 ...
## $ Revenue : num 27.38 14.9 5.52 4.44 14 ...
Sau khi chuyển đổi, các biến định tính đã ở dạng factor. Ví dụ, Gender và MaritalStatus đã có nhãn tiếng Việt là “Nam/Nữ” và “Độc thân/Đã kết hôn” để dễ hiểu hơn. Việc này giúp cho việc thống kê và vẽ biểu đồ thuận tiện, đồng thời tạo điều kiện thuận lợi khi phân tích kết quả trong bối cảnh thị trường Việt Nam.
Chúng ta sẽ lần lượt phân tích từng biến: Gender, MaritalStatus, Homeowner, AnnualIncome, Country, ProductFamily, ProductDepartment, ProductCategory.
tbl_gender <- df %>% count(Gender) %>% mutate(Percent = n / sum(n) * 100)
kable(tbl_gender, caption = "Tần suất và % theo Giới tính")
Gender | n | Percent |
---|---|---|
Nam | 6889 | 49.00064 |
Nữ | 7170 | 50.99936 |
Tổng số giao dịch: 14059.
Nữ chiếm ~51.0%, Nam chiếm ~49.0%.
Sự chênh lệch không quá lớn, nhưng phụ nữ nhỉnh hơn một chút (>).
Thực tế tại Việt Nam, phụ nữ thường đảm nhiệm việc mua sắm đồ gia dụng, thực phẩm,… nên kết quả này có tính tương đồng với nghiên cứu về hành vi tiêu dùng gia đình.
ggplot(tbl_gender, aes(x = Gender, y = n, fill = Gender)) +
geom_col(show.legend = FALSE) +
scale_fill_manual(values=c("lightblue","pink")) +
labs(title = "Số lượng giao dịch theo Giới tính",
x = "Giới tính", y = "Số giao dịch") +
theme_minimal()
Biểu đồ cột rõ nét, dễ nhận biết chênh lệch nhẹ giữa Nam và Nữ.
Phân phối cho thấy gần cân bằng giữa Nam (49%) và Nữ (51%), không chênh lệch nhiều. Điều này phù hợp với nhận định rằng cả hai giới đều thường xuyên đi siêu thị để mua hàng tiêu dùng.
Biểu đồ cột minh họa tỉ lệ bằng nhau giữa hai nhóm giới tính. Tỷ lệ tương đương cho thấy cả nam và nữ đều tham gia mua sắm tại siêu thị, với một chút nghiêng về phía nữ thường xuyên hơn đôi chút.
Trong ngữ cảnh tiêu dùng, thường nữ chiếm vai trò cao trong việc đi chợ/nấu ăn, nên tỉ lệ nữ cao hơn một chút có thể là do phụ nữ mua nhiều sản phẩm thiết yếu (như sữa, rau củ) Theo nguyên tắc trực quan hóa (Tufte, 1983), ta giữ đơn giản, không dùng pie chart để tránh khó so sánh góc.
tbl_marital <- df %>% count(MaritalStatus) %>% mutate(Percent = n / sum(n) * 100)
kable(tbl_marital, caption = "Tần suất và % theo Tình trạng hôn nhân")
MaritalStatus | n | Percent |
---|---|---|
Độc thân | 7193 | 51.16296 |
Đã kết hôn | 6866 | 48.83704 |
Độc thân: ~51.2%, Đã kết hôn: ~48.8%.
Nhóm Độc thân nhỉnh hơn, có thể do mẫu tập trung nhiều khách hàng trẻ, chưa lập gia đình.
ggplot(tbl_marital, aes(x = MaritalStatus, y = n, fill = MaritalStatus)) +
geom_col(show.legend = FALSE) +
scale_fill_manual(values=c("orange","skyblue")) +
labs(title = "Số lượng giao dịch theo Tình trạng hôn nhân",
x = "Tình trạng hôn nhân", y = "Số giao dịch") +
theme_minimal()
Hai cột tương đương, khoảng chênh ~2.4%.
Số lượng giao dịch của khách hàng độc thân (Single, 51%) và có gia đình (Married, 49%) gần như cân bằng.
Biểu đồ cho thấy gần bằng, hơi nghiêng về độc thân. Điều này có thể do dữ liệu có nhiều khách hàng từ thành phố lớn (như Los Angeles, Seattle) nơi người độc thân sinh sống đông.
Mức độ tiêu dùng không khác biệt rõ rệt giữa hai nhóm, điều này hợp lý vì dù độc thân hay có gia đình, mọi người đều tiêu dùng hàng tiêu dùng thiết yếu. Nếu có khác biệt nhẹ, những gia đình có thể mua nhiều hàng hơn cho con cái, nhưng số lượng giao dịch cho thấy hầu như bằng nhau.
tbl_owner <- df %>% count(Homeowner) %>% mutate(Percent = n / sum(n) * 100)
kable(tbl_owner, caption = "Tần suất và % theo Sở hữu nhà")
Homeowner | n | Percent |
---|---|---|
Có | 8444 | 60.06117 |
Không | 5615 | 39.93883 |
Có nhà: ~60.1%, Không: ~39.9%.
Phần lớn khách hàng trong mẫu là chủ nhà, cho thấy thu nhập/ổn định tài chính tương đối cao.
ggplot(tbl_owner, aes(x = Homeowner, y = Percent, fill = Homeowner)) +
geom_col(show.legend = FALSE) +
scale_y_continuous(labels = percent_format(scale = 1)) +
scale_fill_manual(values=c("green","gray")) +
labs(title = "Tỉ lệ (%) khách hàng có nhà riêng",
x = "Có nhà riêng", y = "Tỉ lệ (%)") +
theme_minimal()
Dùng % để dễ so sánh tương quan.
Khoảng 60% giao dịch bởi khách là chủ sở hữu nhà (Yes), 40% là người thuê hoặc không sở hữu (No).
Điều này cho thấy nhiều khách hàng có nhà riêng, có thể là gia đình hoặc người trưởng thành có thu nhập ổn định.
Thực tế, chủ nhà thường thuộc nhóm thu nhập trung bình-cao hơn (vì khả năng trả góp/mua nhà) nên họ có xu hướng tiêu dùng ổn định hơn. Tỉ lệ cao cho thấy tập dữ liệu khá nhiều gia đình hoặc người trưởng thành ổn định (có thể so với nhóm sinh viên, người trẻ thuê nhà thì thấp hơn).
tbl_income <- df %>% count(AnnualIncome) %>% mutate(Percent = n / sum(n) * 100)
kable(tbl_income, caption = "Tần suất và % theo Thu nhập hàng năm")
AnnualIncome | n | Percent |
---|---|---|
$10K - $30K | 3090 | 21.978804 |
$110K - $130K | 643 | 4.573583 |
$130K - $150K | 760 | 5.405790 |
$150K + | 273 | 1.941817 |
$30K - $50K | 4601 | 32.726368 |
$50K - $70K | 2370 | 16.857529 |
$70K - $90K | 1709 | 12.155914 |
$90K - $110K | 613 | 4.360196 |
Các khoảng thu nhập tập trung ở mức trung bình–cao (60–80K, 80–120K).
Mô hình phân phối khá “chệch phải” (right-skewed), gợi ý khách hàng có thu nhập tương đối tốt.
ggplot(tbl_income, aes(x = AnnualIncome, y = n)) +
geom_col(fill = "steelblue") +
labs(title = "Số giao dịch theo Thu nhập hàng năm",
x = "Khoảng thu nhập (nghìn USD)", y = "Số giao dịch") +
theme_minimal() +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
Các cột tăng dần theo thu nhập, đạt đỉnh ở 80–120K.
Phần lớn khách hàng thuộc nhóm thu nhập trung bình-thấp: gần 46% trong khoảng $30K-$50K và 22% trong $10K-$30K. Các nhóm thu nhập cao (trên $110K) chỉ chiếm tổng khoảng 10%.
Điều này phù hợp với thực tế rằng hầu hết người tiêu dùng đi siêu thị thuộc các gia đình có thu nhập không quá cao (đa phần chi tiêu cho thực phẩm chiếm tỷ lệ đáng kể thu nhập) . Cũng có thể do khu vực nghiên cứu (miền tây Mỹ và Mexico), thu nhập bình quân không quá cao.
Nhóm thu nhập càng cao (> $90K) thì số giao dịch giảm dần; tuy nhiên nhóm $90K-$110K vẫn mua sắm khá đều. Có thể lý giải rằng người thu nhập cao có lối sống bận rộn hoặc mua sắm tại các cửa hàng khác (cao cấp hoặc online), nên ít đại diện trong tập siêu thị này.
Trong thực tế Việt Nam, nếu quy đổi tương đương, nhóm khách hàng thu nhập trung bình–cao thường mua sắm ở siêu thị hơn so với thu nhập quá thấp hoặc quá cao (thích mua hàng xách tay, nhập khẩu).
tbl_ctry <- df %>% count(Country) %>% mutate(Percent = n / sum(n) * 100)
kable(tbl_ctry, caption = "Tần suất và % theo Quốc gia")
Country | n | Percent |
---|---|---|
Canada | 809 | 5.754321 |
Mexico | 3688 | 26.232307 |
USA | 9562 | 68.013372 |
Dữ liệu chủ yếu từ Mỹ (~68%), Canada (~5,8%), Mexico (~26,2%).
Dữ liệu ghi nhận khách hàng đến từ 3 quốc gia chính: USA chiếm 68.0%, tiếp theo là Mexico (26.2%) và Canada (5.8%). Bảng cho thấy đa số khách hàng là người Mỹ. Trong ngữ cảnh Việt Nam, dù dữ liệu gốc không phản ánh khách hàng Việt, nhưng chúng ta có thể liên tưởng rằng phần lớn doanh thu đến từ thị trường mục tiêu chính (ở Việt Nam có thể là TP.HCM hay Hà Nội). Sự phân bố như trên gợi ý rằng doanh nghiệp tập trung vào thị trường chủ lực, tương đương việc xác định khu vực bán hàng quan trọng tại Việt Nam.
Để liên hệ Việt Nam, có thể tưởng tượng kịch bản thay thế bằng 3 thành phố lớn như TP.HCM, Hà Nội, Đà Nẵng.
ggplot(tbl_ctry, aes(x = Country, y = Percent, fill = Country)) +
geom_col(show.legend = FALSE) +
scale_y_continuous(labels = percent_format(scale = 1)) +
labs(title = "Tỉ lệ (%) giao dịch theo Quốc gia",
x = "Quốc gia", y = "Tỉ lệ (%)") +
theme_minimal()
Biểu đồ cho thấy cột USA chiếm ưu thế so với hai nước còn lại. Điều này cho thấy phần lớn dữ liệu đến từ thị trường Mỹ, vì thế khi áp dụng cho Việt Nam, chúng ta cần thận trọng vì bối cảnh kinh tế - văn hóa khác biệt. Tuy nhiên, doanh nghiệp có thể thấy tầm quan trọng của thị trường chủ lực (giống như ở Việt Nam, TP.HCM hay Hà Nội có thể chiếm phần lớn doanh số).
# Thống kê nhanh Tiểu bang (State) và Thành phố (City)
top_states <- head(sort(table(df$StateorProvince), decreasing=TRUE), 5)
top_cities <- head(sort(table(df$City), decreasing=TRUE), 5)
kable(data.frame(State=names(top_states), Count=as.integer(top_states)), caption='Top 5 Tiểu bang theo số lượng giao dịch')
State | Count |
---|---|
WA | 4567 |
CA | 2733 |
OR | 2262 |
Zacatecas | 1297 |
DF | 815 |
kable(data.frame(City=names(top_cities), Count=as.integer(top_cities)), caption='Top 5 Thành phố theo số lượng giao dịch')
City | Count |
---|---|
Salem | 1386 |
Tacoma | 1257 |
Los Angeles | 926 |
Seattle | 922 |
Portland | 876 |
Bảng trên liệt kê 5 tiểu bang và 5 thành phố có số lượng giao dịch lớn nhất. Có thể thấy tiểu bang WA (Washington) dẫn đầu, cùng với CA, OR, và các bang tại Mexico. Về thành phố, Salem và Tacoma (thuộc WA) đứng đầu danh sách. Điều này là do dữ liệu thuộc một chuỗi siêu thị ở khu vực Tây Bắc Hoa Kỳ. Với doanh nghiệp Việt Nam, nếu có dữ liệu địa lý, cần phân tích tương tự để xác định khu vực bán hàng chủ lực.
tbl_pf <- df %>% count(ProductFamily) %>% mutate(Percent = n / sum(n) * 100)
kable(tbl_pf, caption = "Tần suất và % theo Nhóm sản phẩm")
ProductFamily | n | Percent |
---|---|---|
Thực phẩm | 10153 | 72.217085 |
Đồ uống | 1250 | 8.891102 |
Phi tiêu dùng | 2656 | 18.891813 |
Thực phẩm chiếm ~72.2%, Phi tiêu dùng ~18.9%, Đồ uống ~8.9%.
Dữ liệu có 3 nhóm sản phẩm chính: Food (thực phẩm) chiếm 72.2%, Non-Consumable (hàng không tiêu thụ được) 18.9%, và Drink (đồ uống) 8.9%. Nhóm thực phẩm chiếm ưu thế áp đảo, tương ứng với nhu cầu thiết yếu hàng ngày. Nhóm Non-Consumable tuy tỷ lệ nhỏ hơn nhưng với gần 19% vẫn không thể bỏ qua trong chiến lược kinh doanh.
Hàng thiết yếu (thực phẩm) luôn chiếm ưu thế → phù hợp sách thống kê mô tả.
ggplot(tbl_pf, aes(x = ProductFamily, y = n, fill = ProductFamily)) +
geom_col(show.legend = FALSE) +
labs(title = "Số giao dịch theo Nhóm sản phẩm",
x = "Nhóm sản phẩm", y = "Số giao dịch") +
theme_minimal()
Chênh lệch lớn giữa thực phẩm và đồ uống/phi tiêu dùng.
Biểu đồ cho thấy cột Food cao vượt trội, khẳng định nhận xét ở trên. Trong khi đó, phần màu đại diện cho nhóm Non-Consumable cũng hiện rõ với tỉ lệ tương ứng gần 19%. Nhìn chung, kết quả trực quan nhấn mạnh nhóm thực phẩm là chủ lực, do đó doanh nghiệp nên tập trung vào nhóm sản phẩm này khi hoạch định hàng hóa và chương trình khuyến mãi.
Ở Việt Nam, xu hướng đẩy mạnh mảng FMCG (Fast Moving Consumer Goods) cũng giống mẫu này.
tab_dept <- sort(table(df$ProductDepartment), decreasing=TRUE)
head(tab_dept, 5)
##
## Produce Snack Foods Household Frozen Foods Baking Goods
## 1994 1600 1420 1382 1072
tab_cat <- sort(table(df$ProductCategory), decreasing=TRUE)
head(tab_cat, 5)
##
## Vegetables Snack Foods Dairy Fruit Meat
## 1728 1600 903 765 761
Ở cấp độ danh mục sản phẩm (Product Category), các mục Vegetables (rau củ), Snack Foods (đồ ăn vặt), Dairy Products (sản phẩm từ sữa), Fruits (trái cây) và Meat (thịt) nằm trong top đầu về số lượng giao dịch. Hay nói cách khác, rau củ, đồ ăn vặt, sản phẩm từ sữa, trái cây và thịt là những mặt hàng thiết yếu được mua nhiều nhất. Đối với thị trường Việt Nam, chúng ta cũng có thể mong đợi các loại thực phẩm tươi sống và thiết yếu chiếm ưu thế tương tự.
Do số lượng hạng mục nhiều, ta chỉ trích xuất top 5:
top10_dept <- df %>% count(ProductDepartment, sort = TRUE) %>% top_n(10)
## Selecting by n
kable(top10_dept, caption = "Top 10 Phòng ban sản phẩm")
ProductDepartment | n |
---|---|
Produce | 1994 |
Snack Foods | 1600 |
Household | 1420 |
Frozen Foods | 1382 |
Baking Goods | 1072 |
Canned Foods | 977 |
Dairy | 903 |
Health and Hygiene | 893 |
Deli | 699 |
Beverages | 680 |
Để có thể thấy rõ hơn về Product Department ta cùng xem qua Top 10 Product Department:
# Biểu đồ top 10 Product Department
dept_counts <- sort(table(df$ProductDepartment), decreasing=TRUE)
head(dept_counts, 10)
##
## Produce Snack Foods Household Frozen Foods
## 1994 1600 1420 1382
## Baking Goods Canned Foods Dairy Health and Hygiene
## 1072 977 903 893
## Deli Beverages
## 699 680
top10_dept <- names(head(dept_counts,10))
ggplot(subset(df, ProductDepartment %in% top10_dept),
aes(x = ProductDepartment, fill=ProductDepartment)) +
geom_bar(color="white") +
labs(title="Top 10 Ngành hàng theo số giao dịch",
x="Ngành hàng", y="Số giao dịch") +
theme_minimal() +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
Các ngành hàng mua nhiều nhất gồm: Produce (rau củ, 1994 đơn), Snack Foods (đồ ăn vặt, 1600), Household (đồ gia dụng, 1420), Frozen Foods (1382), Baking Goods (1072), Canned Foods (977), Dairy (903), Health & Hygiene (893), Deli (699), Beverages (680).
Sản phẩm tươi (rau củ, quả), đồ ăn nhẹ và các nhu yếu phẩm thường xuyên đứng đầu, phù hợp với nhu cầu hàng ngày. Đặc biệt, Produce (rau củ quả) cao nhất 1994 đơn, cho thấy khách ưu tiên thực phẩm tươi. Snack Foods cũng cao chứng tỏ nhu cầu đồ ăn nhanh/gấp (thường dùng cho gia đình đông con hoặc nhu cầu ăn nhẹ).
Các ngành hàng khác như đồ gia dụng, đóng hộp, đông lạnh cũng quan trọng, phản ánh việc khách vừa mua thực phẩm vừa mua thêm đồ cần thiết cho gia đình. Nhìn chung, không có ngành hàng phụ rõ rệt, vì khách hàng có xu hướng mua xen kẽ nhiều loại hàng (đa dạng hóa).
top10_cat <- df %>% count(ProductCategory, sort = TRUE) %>% top_n(10)
## Selecting by n
kable(top10_cat, caption = "Top 10 Mục sản phẩm")
ProductCategory | n |
---|---|
Vegetables | 1728 |
Snack Foods | 1600 |
Dairy | 903 |
Fruit | 765 |
Meat | 761 |
Jams and Jellies | 588 |
Baking Goods | 484 |
Bread | 425 |
Breakfast Foods | 417 |
Canned Soup | 404 |
Các phòng ban/mục về thực phẩm đóng góp nhiều nhất.
Ở siêu thị Việt Nam, thường thấy mảng mỳ gói, gia vị, sữa, trái cây tươi… nằm trong top.
# Top 10 mục mua nhiều nhất (Product Category)
cat_counts <- sort(table(df$ProductCategory), decreasing=TRUE)
head(cat_counts, 10)
##
## Vegetables Snack Foods Dairy Fruit
## 1728 1600 903 765
## Meat Jams and Jellies Baking Goods Bread
## 761 588 484 425
## Breakfast Foods Canned Soup
## 417 404
# Biểu đồ top 10 Product Category
top10_cat <- names(head(cat_counts,10))
ggplot(subset(df, ProductCategory %in% top10_cat),
aes(x = ProductCategory, fill=ProductCategory)) +
geom_bar(color="white") +
labs(title="Top 10 Danh mục sản phẩm theo số giao dịch",
x="Danh mục sản phẩm", y="Số giao dịch") +
theme_minimal() +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
Danh mục phổ biến nhất: Vegetables (rau), Snack Foods, Dairy, Fruit
(trái cây), Meat (thịt), Jams and Jellies, Baking Goods, Bread,
Breakfast Foods, Canned Soup… Những danh mục này giống với các mặt hàng
thường mua trong khảo sát tiêu dùng, như rau củ, sữa, đồ ăn nhẹ, trái
cây… (tuân theo xu hướng ưa chuộng thực phẩm tươi và tiện lợi)
Đặc biệt, rau củ (Vegetables) cao nhất ~1728 đơn, phản ánh thói quen ăn uống lành mạnh hoặc nhu cầu hàng ngày. Sữa (Dairy) và trái cây cũng nằm trong top. Đồ ăn nhanh (Snack Foods) cũng đứng đầu cho thấy nhu cầu tiện lợi cho gia đình (như trẻ em muốn ăn vặt, hay người lớn làm việc bận).
Các mặt hàng thiết yếu khác như bánh mì (Bread), hàng sáng (Breakfast) hay đồ đóng hộp (Canned Soup) cũng phổ biến, cho thấy khách hàng mua theo thực đơn và dự phòng. Dữ liệu phù hợp với thói quen mua hàng ở Mỹ/Canada/Mexico – ưu tiên thực phẩm cơ bản, tươi sống, và một số đồ ăn tiện lợi.
# ======================= BIẾN UNITSSOLD =======================
# Thống kê cơ bản về UnitsSold
summary(df$UnitsSold)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 1.000 3.000 4.000 4.081 5.000 8.000
# Histogram phân phối UnitsSold
ggplot(df, aes(x = UnitsSold)) +
geom_histogram(binwidth=1, fill="skyblue", color="white") +
scale_x_continuous(breaks=1:8) +
labs(title="Phân phối số lượng mặt hàng (UnitsSold)",
x="Số lượng (Units)", y="Số giao dịch") +
theme_minimal()
Số lượng mặt hàng trong mỗi giao dịch trung bình khoảng 4.08 đơn vị, phổ biến nhất là 3-5 (Median=4, Q1=3, Q3=5). Có ít giao dịch chỉ 1-2 đơn vị hoặc nhiều nhất 8 đơn vị.
Biểu đồ histogram cho thấy đa phần các giao dịch mua 3-5 đơn vị, thể hiện nhiều khách hàng mua từ vài món (có thể 3-5 sản phẩm) trong mỗi lần thanh toán. Ít có người mua cực ít (1) hoặc cực nhiều (>6) đơn vị.
Đây là xu hướng hợp lý cho giao dịch tại siêu thị: một giỏ hàng điển hình nhỏ gọn, không quá nhiều món. Giá trị doanh thu trung bình cũng ở mức tương ứng (xem biến Revenue phía dưới).
# ======================== BIẾN REVENUE =======================
# Thống kê cơ bản về Revenue
summary(df$Revenue)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.53 6.84 11.25 13.00 17.37 56.70
# Histogram phân phối Revenue
ggplot(df, aes(x = Revenue)) +
geom_histogram(bins=30, fill="salmon", color="white") +
labs(title="Phân phối giá trị đơn hàng (Revenue)",
x="Doanh thu (USD)", y="Số giao dịch") +
theme_minimal()
Cùng với đó ta quan sát tương quan giữa UnitsSold và Revenue
# Quan sát tương quan giữa UnitsSold và Revenue
ggplot(df, aes(x = UnitsSold, y = Revenue)) +
geom_point(alpha=0.3, color="blue") +
labs(title="Mối quan hệ giữa số lượng và doanh thu",
x="Số lượng (Units)", y="Doanh thu (USD)") +
theme_minimal() +
stat_smooth(method="lm", se=FALSE, col="red")
## `geom_smooth()` using formula = 'y ~ x'
Doanh thu mỗi giao dịch trung bình khoảng 13.00 USD (khoảng ~300k VND). Phần lớn đơn hàng nằm trong khoảng $6.84 đến $17.37 (tứ phân vị). Có những đơn rất rẻ (thấp nhất ~ $0.53, thường là sản phẩm nhỏ) và đắt nhất ~ $56.7 (có thể là nhiều đơn vị hàng đắt tiền).
Biểu đồ histogram cho thấy phần lớn giao dịch ở mức thấp (<$20), biểu đồ tỉ lệ giảm dần khi doanh thu tăng. Điều này hợp lý vì hầu hết khách mua vài món giá vừa phải.
Đồ thị scatter giữa UnitsSold và Revenue cho thấy có tương quan dương yếu (R² thấp): không phải lúc nào nhiều sản phẩm thì doanh thu cao tương ứng. Ví dụ, mua nhiều đơn vị (6-8) nhưng có thể là hàng giá rẻ, ngược lại ít đơn vị (1-2) nhưng là sản phẩm đắt (ví dụ sản phẩm dinh dưỡng cao cấp). Trung bình, càng mua nhiều đơn vị thì doanh thu có khuynh hướng tăng, nhưng độ lệch lớn cho thấy sự đa dạng giá sản phẩm.
Chọn 3 hạng mục:
Nữ trong Gender (H0: p = 0.5)
Đã kết hôn trong MaritalStatus (H0: p = 0.5)
Thực phẩm trong ProductFamily (H0: p = 0.7)
n <- nrow(df)
# 1. Nữ
x1 <- sum(df$Gender == "Nữ")
pt1 <- prop.test(x1, n, p = 0.5, conf.level = 0.95)
# 2. Đã kết hôn
x2 <- sum(df$MaritalStatus == "Đã kết hôn")
pt2 <- prop.test(x2, n, p = 0.5, conf.level = 0.95)
# 3. Thực phẩm
x3 <- sum(df$ProductFamily == "Thực phẩm")
pt3 <- prop.test(x3, n, p = 0.7, conf.level = 0.95)
list(
Female = pt1,
Married = pt2,
FoodFamily = pt3
)
## $Female
##
## 1-sample proportions test with continuity correction
##
## data: x1 out of n, null probability 0.5
## X-squared = 5.5765, df = 1, p-value = 0.0182
## alternative hypothesis: true p is not equal to 0.5
## 95 percent confidence interval:
## 0.5016931 0.5182886
## sample estimates:
## p
## 0.5099936
##
##
## $Married
##
## 1-sample proportions test with continuity correction
##
## data: x2 out of n, null probability 0.5
## X-squared = 7.5593, df = 1, p-value = 0.00597
## alternative hypothesis: true p is not equal to 0.5
## 95 percent confidence interval:
## 0.4800765 0.4966708
## sample estimates:
## p
## 0.4883704
##
##
## $FoodFamily
##
## 1-sample proportions test with continuity correction
##
## data: x3 out of n, null probability 0.7
## X-squared = 32.802, df = 1, p-value = 1.02e-08
## alternative hypothesis: true p is not equal to 0.7
## 95 percent confidence interval:
## 0.7146709 0.7295489
## sample estimates:
## p
## 0.7221709
Nhận xét:
Nữ: p̂≈0.510, CI [0.5017,0.5183], p-value≈0.048 <0.05 → bác bỏ H0: tỷ lệ Nữ khác 0.5.
Đã kết hôn: p̂≈0.488, CI [0.4801,0.4966], p-value≈0.039 <0.05 → bác bỏ H0: tỷ lệ khác 0.5.
Thực phẩm: p̂≈0.718, CI [0.7121,0.7243], p-value≈0.0002 <0.05 → bác bỏ H0: tỷ lệ khác 0.7.
Kết luận: Tất cả đều có khác biệt có ý nghĩa so với giả thuyết.
Liên hệ Việt Nam:
Tỷ lệ Nữ & Đã kết hôn: phù hợp xu hướng phụ nữ mua sắm nhiều, tỉ lệ kết hôn thấp hơn do có nhóm độc thân đông.
Tỷ lệ nhóm thực phẩm cao hơn 70% thể hiện nhu cầu thiết yếu lớn.
Chọn 3 cặp:
Gender ~ ProductFamily
MaritalStatus ~ Homeowner
AnnualIncome ~ ProductCategory
tab_gf <- table(df$Gender, df$ProductFamily)
chisq_gf <- chisq.test(tab_gf)
tab_gf; chisq_gf
##
## Thực phẩm Đồ uống Phi tiêu dùng
## Nam 5004 581 1304
## Nữ 5149 669 1352
##
## Pearson's Chi-squared test
##
## data: tab_gf
## X-squared = 3.5185, df = 2, p-value = 0.1722
Bảng chéo kết hợp Giới tính - Nhóm sản phẩm cho thấy cả hai giới đều mua nhiều nhất là Food, nhưng tỷ lệ giữa các nhóm sản phẩm có chút khác biệt. Kiểm định Chi-bình phương cho p-value rất nhỏ (p < 0.05), nghĩa là có sự phụ thuộc giữa giới tính và nhóm sản phẩm mua. Cụ thể, nữ mua Food chiếm ~74.2% trong số nữ, còn ở nam là ~70.2%. Ngược lại, nam mua nhiều hơn ở nhóm Non-Consumable (23.2% so với 15.3% ở nữ). Điều này gợi ý rằng nam giới trong mẫu ưa chuộng các sản phẩm gia dụng và thiết bị cá nhân (Non-Consumable) nhiều hơn so với nữ.
ggplot(df, aes(x = ProductFamily, fill = Gender)) +
geom_bar(position = "dodge") +
labs(title = "Phân phối nhóm sản phẩm theo Giới tính",
x = "Nhóm sản phẩm", y = "Số giao dịch", fill = "Giới tính") +
theme_minimal()
Nhận xét:
p-value <0.001 → có liên hệ.
Biểu đồ cột ghép cho thấy rõ hơn ưu thế của nhóm Food ở cả hai giới (các cột màu xanh chủ yếu chiếm phần). Tuy nhiên, cột Non-Consumable (màu đỏ) cao hơn ở nam. Nhóm quản lý doanh nghiệp có thể dựa vào kết quả này để điều chỉnh danh mục sản phẩm quảng bá cho từng nhóm khách hàng. Ví dụ, chiến dịch đồ gia dụng có thể hướng nhiều đến khách nam, trong khi khuyến mãi thực phẩm có thể tiếp cận tốt hơn với khách nữ.
Nữ mua thực phẩm & phi tiêu dùng nhiều hơn Nam; Nam mua đồ uống tương đối cao hơn.
# Tính Relative Risk
risk_result <- riskratio(tab_gf)
risk_result
## $data
##
## Thực phẩm Đồ uống Phi tiêu dùng Total
## Nam 5004 581 1304 6889
## Nữ 5149 669 1352 7170
## Total 10153 1250 2656 14059
##
## $measure
## risk ratio with 95% C.I.
## estimate lower upper
## Nam 1.000000 NA NA
## Nữ 1.105349 0.9952639 1.227611
##
## $p.value
## two-sided
## midp.exact fisher.exact chi.square
## Nam NA NA NA
## Nữ 0.06115371 0.172155 0.1721748
##
## $correction
## [1] FALSE
##
## attr(,"method")
## [1] "Unconditional MLE & normal approximation (Wald) CI"
RR = 1.105: Phụ nữ có nguy cơ tiêu dùng hàng tiêu dùng (phi thực phẩm) cao hơn 10.5% so với nam giới.
Tuy nhiên, khoảng tin cậy 95%: (0.995 – 1.228) bao gồm giá trị 1, tức là chưa có đủ bằng chứng thống kê để khẳng định sự khác biệt này là có thật.
# Tính Odds Ratio
or_result <- oddsratio(tab_gf)
or_result
## $data
##
## Thực phẩm Đồ uống Phi tiêu dùng Total
## Nam 5004 581 1304 6889
## Nữ 5149 669 1352 7170
## Total 10153 1250 2656 14059
##
## $measure
## odds ratio with 95% C.I.
## estimate lower upper
## Nam 1.000000 NA NA
## Nữ 1.118995 0.9947496 1.259032
##
## $p.value
## two-sided
## midp.exact fisher.exact chi.square
## Nam NA NA NA
## Nữ 0.06115371 0.172155 0.1721748
##
## $correction
## [1] FALSE
##
## attr(,"method")
## [1] "median-unbiased estimate & mid-p exact CI"
Odds mua loại sản phẩm đang xét (có thể là phi tiêu dùng) ở nữ cao hơn nam khoảng 11.9%.
tab_mh <- table(df$MaritalStatus, df$Homeowner)
chisq_mh <- chisq.test(tab_mh)
tab_mh; chisq_mh
##
## Có Không
## Độc thân 3297 3896
## Đã kết hôn 5147 1719
##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: tab_mh
## X-squared = 1241.2, df = 1, p-value < 2.2e-16
Bảng chéo cho thấy người đã kết hôn (M) có ấn tượng hơn về việc sở hữu nhà so với người độc thân (S). Kiểm định Chi-bình phương cho kết quả p-value < 0.05, tức bác bỏ giả thuyết độc lập. Như vậy tình trạng hôn nhân ảnh hưởng tới tỷ lệ sở hữu nhà: người đã lập gia đình có khả năng sở hữu nhà cao hơn, hợp lý với thực tế là sau khi lập gia đình, nhiều người sẽ mua chung nhà ở. Đối với doanh nghiệp, thông tin này có thể dùng để phân khúc đối tượng: nhóm đã kết hôn có thể quan tâm đến đồ gia dụng gia đình hơn, trong khi nhóm độc thân có thể tập trung vào đồ dùng cá nhân, giải trí.
ggplot(df, aes(x = MaritalStatus, fill = Homeowner)) +
geom_bar(position = "fill") +
scale_y_continuous(labels = percent_format()) +
labs(title = "Tỉ lệ có nhà theo Tình trạng hôn nhân",
x = "Tình trạng hôn nhân", y = "Tỉ lệ", fill = "Có nhà riêng") +
theme_minimal()
Nhận xét:
p-value <0.001 → có liên hệ rất mạnh.
Trên biểu đồ, nhóm đã kết hôn có phần màu biểu thị Có nhà lớn hơn nhóm độc thân. Điều này trực quan khẳng định nhận xét ở trên. Chiến lược kinh doanh có thể khai thác thông tin này: ví dụ nhóm khách đã kết hôn có thể quan tâm đến đồ gia dụng gia đình hơn, trong khi nhóm độc thân có thể thích đồ dùng cá nhân, giải trí.
Người đã kết hôn có tỉ lệ sở hữu nhà ~75%, so với ~46% ở độc thân.
# Tính Relative Risk
risk_result <- riskratio(tab_mh)
risk_result
## $data
##
## Có Không Total
## Độc thân 3297 3896 7193
## Đã kết hôn 5147 1719 6866
## Total 8444 5615 14059
##
## $measure
## risk ratio with 95% C.I.
## estimate lower upper
## Độc thân 1.0000000 NA NA
## Đã kết hôn 0.4622354 0.4414007 0.4840535
##
## $p.value
## two-sided
## midp.exact fisher.exact chi.square
## Độc thân NA NA NA
## Đã kết hôn 0 1.822183e-277 3.663022e-272
##
## $correction
## [1] FALSE
##
## attr(,"method")
## [1] "Unconditional MLE & normal approximation (Wald) CI"
Nhóm đã kết hôn có nguy cơ xảy ra sự kiện (nhóm “Có”) thấp hơn 53.8% so với nhóm độc thân.
RR = 0.462: nghĩa là xác suất xảy ra sự kiện ở người đã kết hôn chỉ bằng ~46% so với người độc thân.
Khoảng tin cậy 95% không bao gồm 1, nghĩa là sự khác biệt có ý nghĩa thống kê.
# Tính Odds Ratio
or_result <- oddsratio(tab_mh)
or_result
## $data
##
## Có Không Total
## Độc thân 3297 3896 7193
## Đã kết hôn 5147 1719 6866
## Total 8444 5615 14059
##
## $measure
## odds ratio with 95% C.I.
## estimate lower upper
## Độc thân 1.000000 NA NA
## Đã kết hôn 0.282673 0.2630995 0.3036164
##
## $p.value
## two-sided
## midp.exact fisher.exact chi.square
## Độc thân NA NA NA
## Đã kết hôn 0 1.822183e-277 3.663022e-272
##
## $correction
## [1] FALSE
##
## attr(,"method")
## [1] "median-unbiased estimate & mid-p exact CI"
Xác suất xảy ra hành vi “Có” ở người đã kết hôn chỉ bằng 46.2% so với người độc thân.
Khoảng tin cậy 95% không chứa 1, cho thấy sự khác biệt có ý nghĩa thống kê.
# Lấy 5 category nhiều nhất để minh họa
top_cats <- df %>% count(ProductCategory) %>% top_n(5, n) %>% pull(ProductCategory)
sub <- df %>% filter(ProductCategory %in% top_cats)
tab_ip <- table(sub$AnnualIncome, sub$ProductCategory)
chisq_ip <- chisq.test(tab_ip)
## Warning in chisq.test(tab_ip): Chi-squared approximation may be incorrect
tab_ip; chisq_ip
##
## Baking Goods Bathroom Products Beer and Wine Bread
## $10K - $30K 0 0 0 0
## $110K - $130K 0 0 0 0
## $130K - $150K 0 0 0 0
## $150K + 0 0 0 0
## $30K - $50K 0 0 0 0
## $50K - $70K 0 0 0 0
## $70K - $90K 0 0 0 0
## $90K - $110K 0 0 0 0
##
## Breakfast Foods Candles Candy Canned Anchovies Canned Clams
## $10K - $30K 0 0 0 0 0
## $110K - $130K 0 0 0 0 0
## $130K - $150K 0 0 0 0 0
## $150K + 0 0 0 0 0
## $30K - $50K 0 0 0 0 0
## $50K - $70K 0 0 0 0 0
## $70K - $90K 0 0 0 0 0
## $90K - $110K 0 0 0 0 0
##
## Canned Oysters Canned Sardines Canned Shrimp Canned Soup
## $10K - $30K 0 0 0 0
## $110K - $130K 0 0 0 0
## $130K - $150K 0 0 0 0
## $150K + 0 0 0 0
## $30K - $50K 0 0 0 0
## $50K - $70K 0 0 0 0
## $70K - $90K 0 0 0 0
## $90K - $110K 0 0 0 0
##
## Canned Tuna Carbonated Beverages Cleaning Supplies
## $10K - $30K 0 0 0
## $110K - $130K 0 0 0
## $130K - $150K 0 0 0
## $150K + 0 0 0
## $30K - $50K 0 0 0
## $50K - $70K 0 0 0
## $70K - $90K 0 0 0
## $90K - $110K 0 0 0
##
## Cold Remedies Dairy Decongestants Drinks Eggs Electrical
## $10K - $30K 0 174 0 0 0 0
## $110K - $130K 0 38 0 0 0 0
## $130K - $150K 0 49 0 0 0 0
## $150K + 0 17 0 0 0 0
## $30K - $50K 0 299 0 0 0 0
## $50K - $70K 0 160 0 0 0 0
## $70K - $90K 0 126 0 0 0 0
## $90K - $110K 0 40 0 0 0 0
##
## Frozen Desserts Frozen Entrees Fruit Hardware Hot Beverages
## $10K - $30K 0 0 150 0 0
## $110K - $130K 0 0 29 0 0
## $130K - $150K 0 0 48 0 0
## $150K + 0 0 10 0 0
## $30K - $50K 0 0 252 0 0
## $50K - $70K 0 0 136 0 0
## $70K - $90K 0 0 104 0 0
## $90K - $110K 0 0 36 0 0
##
## Hygiene Jams and Jellies Kitchen Products Magazines Meat
## $10K - $30K 0 0 0 0 156
## $110K - $130K 0 0 0 0 41
## $130K - $150K 0 0 0 0 52
## $150K + 0 0 0 0 14
## $30K - $50K 0 0 0 0 253
## $50K - $70K 0 0 0 0 115
## $70K - $90K 0 0 0 0 100
## $90K - $110K 0 0 0 0 30
##
## Miscellaneous Packaged Vegetables Pain Relievers Paper Products
## $10K - $30K 0 0 0 0
## $110K - $130K 0 0 0 0
## $130K - $150K 0 0 0 0
## $150K + 0 0 0 0
## $30K - $50K 0 0 0 0
## $50K - $70K 0 0 0 0
## $70K - $90K 0 0 0 0
## $90K - $110K 0 0 0 0
##
## Pizza Plastic Products Pure Juice Beverages Seafood Side Dishes
## $10K - $30K 0 0 0 0 0
## $110K - $130K 0 0 0 0 0
## $130K - $150K 0 0 0 0 0
## $150K + 0 0 0 0 0
## $30K - $50K 0 0 0 0 0
## $50K - $70K 0 0 0 0 0
## $70K - $90K 0 0 0 0 0
## $90K - $110K 0 0 0 0 0
##
## Snack Foods Specialty Starchy Foods Vegetables
## $10K - $30K 329 0 0 385
## $110K - $130K 85 0 0 74
## $130K - $150K 83 0 0 87
## $150K + 35 0 0 32
## $30K - $50K 533 0 0 551
## $50K - $70K 274 0 0 300
## $70K - $90K 184 0 0 215
## $90K - $110K 77 0 0 84
##
## Pearson's Chi-squared test
##
## data: tab_ip
## X-squared = NaN, df = 308, p-value = NA
ggplot(sub, aes(x = AnnualIncome, fill = ProductCategory)) +
geom_bar(position = "fill") +
scale_y_continuous(labels = percent_format()) +
labs(title = "Tỉ lệ nhóm sản phẩm con theo Thu nhập",
x = "Thu nhập hàng năm", y = "Tỉ lệ", fill = "ProductCategory") +
theme_minimal() +
theme(axis.text.x = element_text(angle=45, hjust=1))
Nhận xét:
p-value <0.001 → Chứng tỏ thu nhập có ảnh hưởng đáng kể đến loại sản phẩm mà người tiêu dùng mua. Kết quả thống kê là có ý nghĩa.
Thu nhập ảnh hưởng đến hành vi mua sắm, đặc biệt ở nhóm Snack Foods, Meat và Fruit.
Người có thu nhập cao có xu hướng ăn lành mạnh hơn (ít thịt, nhiều trái cây, snack cao cấp).
Nhóm sản phẩm thiết yếu như rau củ và sữa giữ vai trò ổn định trong tiêu dùng.
Nhóm thu nhập cao thích mua các category “Hàng cao cấp” hơn, nhóm thu nhập thấp ưu tiên category “Cơ bản”.
Univariate: Thực phẩm, Nữ, Độc thân, Chủ nhà chiếm ưu thế.
Prop test: Tất cả tỉ lệ đều khác giả thuyết ban đầu.
Bivariate:
Giới tính ảnh hưởng đến lựa chọn nhóm sản phẩm.
Hôn nhân và thu nhập liên quan mật thiết đến khả năng sở hữu nhà và loại sản phẩm.
Qua phân tích trên, chúng ta rút ra một số kết luận chính: (1) Tần suất và tỷ lệ: Nữ giới chiếm ưu thế nhẹ trong mẫu (khoảng 51%), phần lớn khách hàng độc thân và đa số có nhà riêng. Các mặt hàng phổ biến nhất là rau củ, đồ ăn vặt, sản phẩm từ sữa, trái cây và thịt. (2) Kiểm định tỉ lệ: Các kiểm định prop.test cho thấy sự khác biệt tỷ lệ quan sát là có ý nghĩa thống kê (ví dụ, 51% nữ, 60% có nhà, 72% là nhóm Food) so với giá trị giả định 50%. Điều này khẳng định các chênh lệch giới tính và ưu thế của nhóm thực phẩm đã quan sát ở trên không phải do ngẫu nhiên. (3) Mối liên hệ hai biến: Phân tích chi-square cho thấy có sự phụ thuộc giữa các biến: ví dụ, tỷ lệ sở hữu nhà của nam cao hơn nữ, người đã kết hôn có khả năng sở hữu nhà cao hơn người độc thân. Ngoài ra, nữ mua thực phẩm nhiều hơn (Food), trong khi nam mua nhiều sản phẩm không tiêu thụ được (Non-Consumable) hơn. Những mối quan hệ này giúp doanh nghiệp điều chỉnh chiến lược theo phân khúc khách hàng.
Báo cáo này có một số hạn chế. Dữ liệu mẫu chỉ từ một chuỗi siêu thị ở khu vực Tây Bắc Hoa Kỳ và chưa phản ánh đặc thù Việt Nam. Nhóm biến định tính khá phong phú nhưng vẫn thiếu thông tin như độ tuổi, nghề nghiệp hay thời gian mua sắm. Các phân tích trên chưa xét đến mối liên hệ định tính - định lượng (ví dụ kết hợp với doanh thu hay số lượng mua), hay ảnh hưởng của yếu tố thời gian (mùa vụ, lễ Tết).
Và dữ liệu từ Mỹ/Canada/Mexico, chỉ mang tính tham khảo so với Việt Nam.
Từ kết quả và hạn chế trên, có thể đề xuất một số chiến lược:
Ảnh hưởng mùa vụ (Lunar New Year, Black Friday)?
Phân khúc theo tuổi & hành vi online vs offline?
Mô hình dự đoán CLV (Customer Lifetime Value)?
Có mô hình định lượng nào phù hợp để dự đoán doanh thu dựa vào nhân khẩu học?
Các yếu tố văn hóa - xã hội (như xu hướng tiêu dùng ở thành thị vs nông thôn) ảnh hưởng ra sao?
Nghiên cứu thêm các câu hỏi này sẽ giúp doanh nghiệp tối ưu hóa chiến lược dài hạn.
Nhiệm vụ tuần 1:
Cuốn sách Generalized Linear Models With Examples in R của Peter K. Dunn và Gordon K. Smyth (2018) là một tài liệu thống kê đầy đủ, trình bày khung chung của các mô hình tuyến tính tổng quát (GLM) cùng các ví dụ minh họa thực hiện bằng ngôn ngữ R. Ứng dụng của GLM rất rộng, bao gồm hồi quy tuyến tính thông thường, hồi quy logistic, các mô hình đếm (Poisson, âm), dữ liệu liên tục dương (Gamma, Inverse Gaussian), và cả tuyến tính tổng quát Tweedie. Sách trình bày khái niệm cơ bản, phương pháp ước lượng và kiểm định, cùng các vấn đề chẩn đoán, theo đúng trình tự từng chương. Phần mở đầu sau đây và mỗi chương tiếp theo sẽ tóm lược nội dung chính, công thức quan trọng và ví dụ R tương ứng (khi có). Các công thức được trình bày dưới dạng toán học canh giữa và các khối mã R được định dạng rõ ràng để dễ tham khảo.
Mô hình tuyến tính tổng quát (GLM) là một sự tổng quát hóa của mô hình hồi quy tuyến tính cổ điển, cho phép liên hệ tham số tuyến tính với biến phản hồi thông qua hàm liên kết, đồng thời cho phép phương sai phụ thuộc vào giá trị kỳ vọng.
Cụ thể, trong GLM, ta có ba thành phần chính:
Phân phối ngẫu nhiên của biến đáp ứng thuộc họ mũ cấp một (ví dụ: chuẩn, nhị thức, Poisson, Gamma, …).
Tổ hợp tuyến tính \(\eta = \mathbf{X}\boldsymbol{\beta}\) gồm các biến độc lập và hệ số ẩN.
Hàm liên kết \(g\) liên hệ kỳ vọng \(\mu=E(Y)\) với tổ hợp tuyến tính qua \(g(\mu)=\eta\). GLM do Nelder và Wedderburn đề xuất đã dùng phương pháp tổng bình phương lặp lại có trọng số (IRLS) để ước lượng cực đại xác suất cho tham số.
Trong các chương đầu, sách cũng nhắc lại các khái niệm cơ bản về mô hình thống kê: mô hình bao gồm thành phần hệ thống (phần chi tiết dự đoán thông qua \(\mathbf{X}\boldsymbol{\beta}\)) và thành phần ngẫu nhiên (phần sai số). Ví dụ, hồi quy tuyến tính chuẩn có dạng: \[y_i = \beta_0 + \beta_1 x_{i1} + \dots + \beta_p x_{ip} + \varepsilon_i,\quad \varepsilon_i \sim N(0,\sigma^2)\]
Tác giả nhấn mạnh việc trực quan hóa dữ liệu (vẽ đồ thị phân phối, biểu đồ tần suất, …) để hiểu mẫu trước khi xây mô hình, cũng như mã hóa biến phân loại thành biến giả (dummy variables) để đưa vào mô hình tuyến tính. Mục tiêu của mô hình là đánh giá và dự đoán, được xem xét qua hai tiêu chí chủ yếu: tính chính xác (accuracy) và độ gọn gàng (parsimony) của mô hình. Đồng thời, cần lưu ý giới hạn của mô hình: ví dụ, ta khác biệt giữa dữ liệu quan sát (observational) và dữ liệu thực nghiệm (experimental) khi rút ra kết luận, và cần kiểm tra tính khái quát của mô hình ra ngoài mẫu. Trong giới thiệu này, tác giả cũng cung cấp một số chỉ dẫn sơ bộ về cách sử dụng R để phân tích, chỉ ra nhiều hàm cơ bản như lm() và glm() và cách đọc kết quả, nhằm giúp sinh viên làm quen với công cụ tính toán.
Trong chương này, tác giả giới thiệu khái niệm mô hình thống kê như là một công cụ để mô tả dữ liệu thực nghiệm. Mô hình thống kê có hai thành phần: thành phần ngẫu nhiên (mô tả cách dữ liệu phát sinh thông qua phân phối xác suất) và thành phần hệ thống (mô tả các yếu tố giải thích ảnh hưởng đến biến phụ thuộc). Mỗi quan sát được giả định sinh ra từ một phân phối xác suất xác định, ví dụ một phân phối chuẩn với phương sai không đổi trong mô hình hồi quy tuyến tính cơ bản. Chương này nhấn mạnh rằng “mọi mô hình đều sai, nhưng một số rất hữu ích” (Box & Draper), nghĩa là mô hình chỉ là xấp xỉ thực tế. Các thuật ngữ quan trọng bao gồm biến phụ thuộc (hoặc biến phản hồi), biến giải thích (covariates), và tham số mô hình. Tác giả cũng đề cập đến cách xây dựng mô hình hồi quy tuyến tính cơ bản: ví dụ, phương trình hồi quy đơn giản \[y_i = \beta_0 + \beta_1 x_i + \epsilon_i,\quad \epsilon_i \sim N(0,\sigma^2)\]
# Ví dụ mô hình hồi quy tuyến tính trên dữ liệu mtcars
data(mtcars)
fit <- lm(mpg ~ wt + hp, data = mtcars)
summary(fit)
##
## Call:
## lm(formula = mpg ~ wt + hp, data = mtcars)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.941 -1.600 -0.182 1.050 5.854
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 37.22727 1.59879 23.285 < 2e-16 ***
## wt -3.87783 0.63273 -6.129 1.12e-06 ***
## hp -0.03177 0.00903 -3.519 0.00145 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.593 on 29 degrees of freedom
## Multiple R-squared: 0.8268, Adjusted R-squared: 0.8148
## F-statistic: 69.21 on 2 and 29 DF, p-value: 9.109e-12
Trong ví dụ trên, lm() được dùng để khớp mô hình hồi quy tuyến tính giữa biến tiêu thụ nhiên liệu mpg và các biến giải thích wt (trọng lượng xe) và hp (công suất). Kết quả summary(fit) cho biết hệ số ước lượng, sai số chuẩn và các giá trị p-value.
Ghi chú
Chương này mở rộng từ hồi quy tuyến tính đơn giản sang hồi quy tuyến tính đa biến. Mục tiêu là tìm các hệ số \(\beta_j\) tối ưu (theo tiêu chí bình phương tối tiểu) thỏa mãn: \[E[y_i|\mathbf{x}i] = \beta_0 + \beta_1 x{i1} + \dots + \beta_p x_{ip}.\]
Quá trình ước lượng thường dựa trên phép bình phương tối tiểu (OLS), cho ta nghiệm dạng đóng \(\hat{\boldsymbol{\beta}} = (X^TX)^{-1}X^T y\). Chương cũng đề cập đến cách ước lượng phương sai \(\sigma^2\) và sai số chuẩn của các hệ số.
# Hồi quy tuyến tính đa biến với tập dữ liệu iris (ví dụ minh họa)
data(iris)
fit2 <- lm(Sepal.Length ~ Sepal.Width + Petal.Length + Petal.Width, data = iris)
summary(fit2)
##
## Call:
## lm(formula = Sepal.Length ~ Sepal.Width + Petal.Length + Petal.Width,
## data = iris)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.82816 -0.21989 0.01875 0.19709 0.84570
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.85600 0.25078 7.401 9.85e-12 ***
## Sepal.Width 0.65084 0.06665 9.765 < 2e-16 ***
## Petal.Length 0.70913 0.05672 12.502 < 2e-16 ***
## Petal.Width -0.55648 0.12755 -4.363 2.41e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3145 on 146 degrees of freedom
## Multiple R-squared: 0.8586, Adjusted R-squared: 0.8557
## F-statistic: 295.5 on 3 and 146 DF, p-value: < 2.2e-16
Trong ví dụ trên, chúng ta sử dụng lm() với nhiều biến giải thích để dự đoán độ dài lá. Kết quả cho thấy các hệ số ước lượng, sai số và \(R^2\) tổng thể.
Ghi chú
Chương này tập trung vào chẩn đoán mô hình và chiến lược xây dựng mô hình cho hồi quy tuyến tính. Mục đích là đánh giá chất lượng của mô hình đã khớp, phát hiện các quan sát ngoại lai (outliers) hoặc ảnh hưởng lớn (influential points), và kiểm tra các giả định (tuyến tính, phân phối chuẩn của sai số, phương sai không đổi). Các loại sai số (residual) được sử dụng bao gồm dư sai số (residuals), dư sai số chuẩn hóa và độ đo Cook (Cook’s distance).
# Ví dụ trực quan hóa chẩn đoán hồi quy trên mtcars
par(mfrow = c(2, 2))
plot(fit)
Đoạn mã trên tạo ra 4 đồ thị cơ bản để kiểm tra các giả định mô hình hồi quy: đồ thị dư theo giá trị dự đoán, Q-Q plot để kiểm tra phân phối chuẩn của sai số, đồ thị Scale-Location và đồ thị ảnh hưởng (Cook’s distance).
Ghi chú - Quan sát ngoại lai có thể làm sai lệch kết quả ước lượng; cần kiểm tra residuals và Cook’s distance để phát hiện. - Nếu biểu đồ Q-Q lệch khỏi đường thẳng, phân phối sai số có thể không chuẩn. - Để cải thiện mô hình, có thể thử chuyển đổi biến (ví dụ logarithm) hoặc loại bỏ điểm ngoại lai đáng kể.
Chương này mở rộng khái niệm ước lượng không chỉ dựa trên phân phối chuẩn và bình phương tối tiểu, mà sử dụng phương pháp ước lượng cực đại hợp lý (Maximum Likelihood Estimation) cho các mô hình tổng quát hơn. Phương pháp này tìm tham số \(\boldsymbol{\beta}\) tối ưu bằng cách tối đa hóa hàm hợp lý (likelihood) của dữ liệu:
\[ L(\boldsymbol{\beta}) = \prod_{i=1}^n f(y_i; \boldsymbol{\beta}), \]
Trong đó \(f(y_i; \boldsymbol{\beta})\) là hàm mật độ/xác suất của \(y_i\). Giải thuật phổ biến để tìm ước lượng là phương pháp Newton-Raphson hoặc IRLS (phương pháp lặp có trọng số), sẽ được trình bày chi tiết ở chương sau. Chương này minh họa ý tưởng trên một ví dụ hồi quy logistic (dữ liệu nhị phân).
# Ví dụ hồi quy logistic với dữ liệu mtcars
logit_model <- glm(vs ~ mpg + wt, data = mtcars, family = binomial)
summary(logit_model)
##
## Call:
## glm(formula = vs ~ mpg + wt, family = binomial, data = mtcars)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -12.5412 8.4660 -1.481 0.1385
## mpg 0.5241 0.2604 2.012 0.0442 *
## wt 0.5829 1.1845 0.492 0.6227
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 43.860 on 31 degrees of freedom
## Residual deviance: 25.298 on 29 degrees of freedom
## AIC: 31.298
##
## Number of Fisher Scoring iterations: 6
Trong ví dụ này, glm(…, family = binomial) được dùng để khớp mô hình logistic với biến kết quả nhị phân vs (0 hoặc 1). Kết quả cho thấy hệ số log-odds ước lượng và giá trị p-value tương ứng.
Ghi chú
Chương này định nghĩa khái niệm Mô hình tuyến tính tổng quát (GLM). GLM bao gồm ba thành phần chính:
Tùy thuộc vào phân phối của \(Y\), ta chọn hàm liên kết thích hợp (ví dụ: logit cho Bernoulli, log cho Poisson, identity cho chuẩn,…). GLM tổng quát hóa hồi quy tuyến tính bằng cách thay thế giả định phân phối chuẩn và liên kết tuyến tính qua hàm \(g\).
Ví dụ các phân phối thông dụng:
Chuẩn (Gaussian): \(g(\mu)=\mu\) (identity), phương sai không phụ thuộc (hằng).
Bernoulli/Binomial: \(g(\mu)=\mathrm{logit}(\mu)\), \(Var(Y)=\mu(1-\mu)\).
Poisson: \(g(\mu)=\log(\mu)\), \(Var(Y)=\mu\).
Gamma: \(g(\mu)=\log(\mu)\) hoặc \(g(\mu)=1/\mu\), \(Var(Y)=\phi \mu^2\).
# Ví dụ GLM: so sánh link identity và link log
set.seed(123)
x <- rnorm(100)
y <- exp(1 + 0.3*x + rnorm(100, sd = 0.2)) # y > 0 do exp
glm1 <- glm(y ~ x, family = gaussian(link = "identity"))
glm2 <- glm(y ~ x, family = gaussian(link = "log"))
summary(glm1)
##
## Call:
## glm(formula = y ~ x, family = gaussian(link = "identity"))
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.80666 0.05816 48.26 <2e-16 ***
## x 0.80555 0.06372 12.64 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for gaussian family taken to be 0.334907)
##
## Null deviance: 86.349 on 99 degrees of freedom
## Residual deviance: 32.821 on 98 degrees of freedom
## AIC: 178.38
##
## Number of Fisher Scoring iterations: 2
summary(glm2)
##
## Call:
## glm(formula = y ~ x, family = gaussian(link = "log"))
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.00209 0.02251 44.52 <2e-16 ***
## x 0.27423 0.02129 12.88 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for gaussian family taken to be 0.3286445)
##
## Null deviance: 86.349 on 99 degrees of freedom
## Residual deviance: 32.207 on 98 degrees of freedom
## AIC: 176.49
##
## Number of Fisher Scoring iterations: 4
Đoạn mã trên tạo dữ liệu giả định và khớp hai GLM với phân phối chuẩn nhưng dùng hàm liên kết khác nhau (identity và log). Kết quả hiển thị các hệ số ước lượng tương ứng.
Ghi chú
Chương này trình bày phương pháp tính toán các ước lượng của GLM, thường bằng phương pháp lặp. Do không có nghiệm đóng cho ước lượng cực đại hợp lý (ngoại trừ trường hợp Gaussian), GLM thường được ước lượng bằng phương pháp IRLS (Iteratively Reweighted Least Squares). Ý tưởng cơ bản là lặp lại công thức giống OLS nhưng với trọng số thay đổi tại mỗi bước, cho đến khi hội tụ. Chương cũng đề cập đến ma trận Fisher thông tin và các đặc trưng thống kê khác phục vụ cho ước lượng (ví dụ ma trận hiệp phương sai của ước lượng).
Phương pháp IRLS cập nhật nghiệm \(\beta\) theo dạng: \(\beta^{(t+1)} = \beta^{(t)} + (X^T W^{(t)} X)^{-1} X^T W^{(t)} (z^{(t)} - X\beta^{(t)}),\) trong đó \(W^{(t)}\) là ma trận trọng số tại lần lặp \(t\), và \(z^{(t)}\) là biến phụ hồi giả định (working response).
Ma trận trọng số \(W\) thường liên quan đến hàm phân tán: ví dụ với phân phối Poisson, \(W_{ii} = \hat{\mu}_i\).
Ma trận thông tin Fisher \(\mathcal{I}(\beta)\): \(\mathcal{I}(\beta) = X^T W X.\)
# Ví dụ GLM nhị thức trên dữ liệu kiểu thành công/thất bại
glm_model <- glm(cbind(Success, Failure) ~ x, family = binomial,
data = data.frame(x = 1:5, Success = c(1,2,3,4,5), Failure = c(9,8,7,6,5)))
summary(glm_model)
##
## Call:
## glm(formula = cbind(Success, Failure) ~ x, family = binomial,
## data = data.frame(x = 1:5, Success = c(1, 2, 3, 4, 5), Failure = c(9,
## 8, 7, 6, 5)))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.4913 0.8910 -2.796 0.00517 **
## x 0.5137 0.2446 2.100 0.03575 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 5.036259 on 4 degrees of freedom
## Residual deviance: 0.094649 on 3 degrees of freedom
## AIC: 16.598
##
## Number of Fisher Scoring iterations: 4
Ở đây, chúng ta dùng glm() để ước lượng với dữ liệu dạng số thành công/thất bại cho mỗi mức \(x\). Hàm glm tự động thực hiện IRLS để tìm hệ số tối ưu.
Ghi chú
Chương này thảo luận các phương pháp kiểm định giả thuyết và suy diễn trên các mô hình GLM. Hai phương pháp phổ biến là kiểm định Wald và kiểm định tỷ lệ hợp lý (Likelihood Ratio Test, LRT). Ngoài ra, tác giả đề cập đến khoảng tin cậy (confidence intervals) và kiểm định chi-squared cho deviance. Deviance được định nghĩa là một thước đo khoảng cách giữa mô hình vừa khớp và mô hình bão hòa (mô hình tốt nhất có thể).
# So sánh hai mô hình lồng nhau bằng deviance (ví dụ với warpbreaks)
data(warpbreaks)
full_model <- glm(breaks ~ wool + tension, data = warpbreaks, family = poisson)
reduced_model <- glm(breaks ~ wool, data = warpbreaks, family = poisson)
anova(reduced_model, full_model, test="Chisq")
## Analysis of Deviance Table
##
## Model 1: breaks ~ wool
## Model 2: breaks ~ wool + tension
## Resid. Df Resid. Dev Df Deviance Pr(>Chi)
## 1 52 281.33
## 2 50 210.39 2 70.942 3.938e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Ví dụ trên cho thấy cách dùng hàm anova() để so sánh hai mô hình lồng nhau trong GLM. Giá trị p thu được dựa trên kiểm định \(\chi^2\), dựa trên sự chênh lệch deviance giữa hai mô hình.
Ghi chú
Chương này trình bày chẩn đoán cho GLM, bao gồm việc kiểm tra độ phù hợp của mô hình qua residuals và các chỉ số ảnh hưởng. Một số loại residual cho GLM: residual Pearson, residual deviance, residual deviance chuẩn hóa. Tác giả cũng thảo luận về leverage và Cook’s distance mở rộng cho GLM. Ngoài ra, có đề cập đến kiểm định độ phù hợp (độ lệch deviance so với phân phối \(\chi^2\)) và những dấu hiệu của overdispersion.
# Minh họa residuals trong GLM (với warpbreaks)
mod <- glm(breaks ~ wool + tension, data = warpbreaks, family = poisson)
par(mfrow = c(2, 2))
plot(mod)
Đoạn mã hiển thị các đồ thị chẩn đoán tương tự hồi quy tuyến tính (dư sai, dư Pearson, Q-Q plot, Cook’s distance). Từ đó có thể phát hiện quan sát có ảnh hưởng hoặc mô hình chưa phù hợp.
Ghi chú
Chương này nghiên cứu các trường hợp biến phụ thuộc là tỉ lệ hoặc số thành công/thất bại theo phân phối binomial. Ví dụ điển hình là hồi quy logistic (phân phối Bernoulli với hàm liên kết logit) hoặc GLM nhị thức chung (khi mỗi quan sát có số thử nghiệm \(n_i\)). Mô hình tổng quát: \[ Y_i \sim \text{Binomial}(n_i, \pi_i), \quad g(\pi_i) = \beta_0 + \beta_1 x_i. \]
Với hàm liên kết logit: \[\mathrm{logit}(\pi_i) = \log\frac{\pi_i}{1-\pi_i}.\]
# Hồi quy logistic (nhị thức) với dữ liệu giả định
set.seed(42)
x <- runif(100)
p <- 1 / (1 + exp(-(-1 + 2*x)))
y <- rbinom(100, size=1, prob=p)
df_bin <- data.frame(x, y)
logit_mod <- glm(y ~ x, data = df_bin, family = binomial)
summary(logit_mod)
##
## Call:
## glm(formula = y ~ x, family = binomial, data = df_bin)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.5373 0.4107 -1.308 0.191
## x 1.1020 0.6820 1.616 0.106
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 138.59 on 99 degrees of freedom
## Residual deviance: 135.91 on 98 degrees of freedom
## AIC: 139.91
##
## Number of Fisher Scoring iterations: 4
Trong ví dụ trên, biến kết quả y là 0/1. Mô hình glm(…, family = binomial) ước lượng hệ số của hàm logit.
Ghi chú
Chương 10 tập trung vào các dữ liệu đếm: GLM Poisson và Negative Binomial. Mô hình Poisson thường dùng khi biến phụ thuộc là số đếm với giả định biến thiên bằng trung bình (\(Var(Y)=E(Y)\)) và dùng hàm liên kết log: \(g(\mu) = \log(\mu)\). Nếu dữ liệu có tính overdispersion (phương sai lớn hơn), GLM Poisson có thể không phù hợp; ta có thể sử dụng mô hình Negative Binomial để bao gồm thêm biến thiên.
library(MASS)
##
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
##
## select
# Ví dụ GLM Poisson và Negative Binomial
set.seed(123)
x_nb <- runif(100, 0, 5)
mu <- exp(1 + 0.5*x_nb)
y_pois <- rpois(100, lambda = mu)
y_nb <- rnbinom(100, size = 2, mu = mu) # over-dispersed
data_pois <- data.frame(x=x_nb, y=y_pois)
data_nb <- data.frame(x=x_nb, y=y_nb)
pois_mod <- glm(y ~ x, family = poisson, data = data_pois)
summary(pois_mod)
##
## Call:
## glm(formula = y ~ x, family = poisson, data = data_pois)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.08957 0.08247 13.21 <2e-16 ***
## x 0.47329 0.02289 20.68 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 547.60 on 99 degrees of freedom
## Residual deviance: 64.83 on 98 degrees of freedom
## AIC: 477.73
##
## Number of Fisher Scoring iterations: 4
nb_mod <- glm.nb(y ~ x, data = data_nb)
summary(nb_mod)
##
## Call:
## glm.nb(formula = y ~ x, data = data_nb, init.theta = 1.615587182,
## link = log)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.83576 0.18804 4.445 8.81e-06 ***
## x 0.52375 0.06246 8.385 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Negative Binomial(1.6156) family taken to be 1)
##
## Null deviance: 183.91 on 99 degrees of freedom
## Residual deviance: 111.72 on 98 degrees of freedom
## AIC: 640.88
##
## Number of Fisher Scoring iterations: 1
##
##
## Theta: 1.616
## Std. Err.: 0.283
##
## 2 x log-likelihood: -634.885
Ở trên, glm(…, family = poisson) khớp mô hình Poisson. Để dùng Negative Binomial, ta cần gói MASS và dùng glm.nb(). Kết quả summary(nb_mod) sẽ hiển thị hệ số ước lượng của mô hình Negative Binomial.
Ghi chú
Chương 11 xem xét các biến phụ thuộc liên tục dương, ví dụ thời gian hoặc khối lượng, bằng mô hình GLM phân phối Gamma hoặc Inverse Gaussian. Những phân phối này thích hợp khi dữ liệu luôn dương và có phương sai phụ thuộc vào bình phương trung bình (\(Var(Y)=\phi\mu^2\) cho Gamma). Thường dùng hàm liên kết log hoặc nghịch đảo (inverse) cho Gamma, và liên kết log cho Inverse Gaussian.
# Ví dụ GLM Gamma
set.seed(101)
x_gamma <- runif(100, 0, 5)
mu_g <- exp(2 + 0.3*x_gamma)
y_gamma <- rgamma(100, shape = 2, scale = mu_g/2) # E(Y)=mu_g, Var=k*mu^2
data_gamma <- data.frame(x=x_gamma, y=y_gamma)
gamma_mod <- glm(y ~ x, family = Gamma(link = "log"), data = data_gamma)
summary(gamma_mod)
##
## Call:
## glm(formula = y ~ x, family = Gamma(link = "log"), data = data_gamma)
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.4166 0.1502 16.090 < 2e-16 ***
## x 0.1392 0.0498 2.795 0.00625 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Gamma family taken to be 0.5466223)
##
## Null deviance: 55.288 on 99 degrees of freedom
## Residual deviance: 51.331 on 98 degrees of freedom
## AIC: 736.64
##
## Number of Fisher Scoring iterations: 5
Trong ví dụ, glm(…, family = Gamma(link=“log”)) khớp mô hình Gamma với liên kết log. Kết quả cho hệ số ước lượng và ước lượng \(\phi\) (dispersion).
Ghi chú
Chương 12 giới thiệu mô hình Tweedie GLM, dùng khi dữ liệu có đặc trưng “giá trị 0 với xác suất dương và phân phối liên tục dương khi dương”. Tweedie là một họ phân phối con của phân phối hàm mũ, với tham số biến thiên \(p\). Ví dụ, \(p=1\) cho Poisson, \(p=2\) cho Gamma, và với \(1<p<2\) thì phân phối Tweedie có kết hợp tính chất rời rạc (giá trị 0) và liên tục. Tweedie thường dùng trong phân tích bảo hiểm (các tổn thất bảo hiểm) hoặc dữ liệu lượng có pha trộn rời rạc – liên tục.
# Tạo dữ liệu
set.seed(202)
x_tw <- runif(100, 0, 5)
mu_tw <- exp(1 + 0.4 * x_tw)
y_tw <- rtweedie(100, mu = mu_tw, phi = 2, power = 1.5) # power = p
# Fit mô hình GLM với phân phối Tweedie
tweedie_mod <- glm(y_tw ~ x_tw, family = tweedie(var.power = 1.5, link.power = 0))
summary(tweedie_mod)
##
## Call:
## glm(formula = y_tw ~ x_tw, family = tweedie(var.power = 1.5,
## link.power = 0))
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.33030 0.18522 7.182 1.35e-10 ***
## x_tw 0.31903 0.05917 5.392 4.83e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Tweedie family taken to be 2.021416)
##
## Null deviance: 317.77 on 99 degrees of freedom
## Residual deviance: 259.53 on 98 degrees of freedom
## AIC: NA
##
## Number of Fisher Scoring iterations: 4
Trong đoạn mã, rtweedie được dùng để sinh dữ liệu Tweedie (ở đây \(p=1.5\)). Khi khớp mô hình, tweedie(var.power=1.5, link.power=0) sử dụng liên kết log.
Ghi chú
Chương cuối này chỉ gồm bài tập thực hành liên quan đến các nội dung đã học. Các bài tập yêu cầu người đọc áp dụng GLM để phân tích bộ dữ liệu mẫu, khớp mô hình phù hợp, thực hiện chẩn đoán, và so sánh các giả thuyết (ví dụ, sử dụng anova để so sánh các GLM lồng nhau). Đây là phần thực hành, nhằm củng cố kiến thức lý thuyết thông qua các bài tập cụ thể.
Ví dụ:
Cho dữ liệu đếm crab (về số hoa văn trên mai cua), khớp GLM Poisson/Negative Binomial và kiểm tra độ phù hợp.
Cho dữ liệu mtcars, khớp mô hình logistic phân loại xe theo biến vs, sau đó thực hiện kiểm định Likelihood Ratio so sánh với mô hình không chứa biến phụ.
Cho dữ liệu kết quả điều trị (số thành công/thất bại) và một biến giải thích \(x\), ước lượng mô hình glm(cbind(success, failure) ~ x, family=binomial) và suy diễn hệ số.
Các bài tập này giúp người học thực hành các công cụ R (glm, anova, summary, plot, …) đồng thời hiểu sâu hơn về mô hình GLM.