Giới thiệu
Nấm là một nguồn thực phẩm phong phú nhưng cũng tiềm ẩn nhiều nguy cơ gây ngộ độc, thậm chí tử vong nếu ăn phải nấm độc. Việc phân biệt chính xác giữa nấm ăn được và nấm độc đòi hỏi kiến thức chuyên môn và kinh nghiệm. Trong bối cảnh đó, việc áp dụng các phương pháp phân tích dữ liệu và mô hình hóa thống kê để tìm ra các quy luật và đặc điểm nhận dạng trở nên vô cùng hữu ích.
Báo cáo này sử dụng bộ dữ liệu “Mushroom Classification” nổi tiếng từ UCI Machine Learning Repository, được trích từ “Sổ tay Hướng dẫn về Nấm của Hiệp hội Audubon”. Nghiên cứu này không chỉ dừng lại ở việc mô tả và phân tích các mối quan hệ song biến, mà còn tiến sâu hơn vào việc xây dựng các mô hình dự đoán bằng các phương pháp hồi quy tổng quát. Cụ thể, mục tiêu của báo cáo là:
Khám phá và mô tả các đặc điểm vật lý của các loại nấm trong bộ dữ liệu.
Phân tích mối quan hệ song biến giữa các đặc điểm này với độc tính và môi trường sống của nấm.
Xây dựng và diễn giải các mô hình Hồi quy Logistic và Hồi quy Probit để lượng hóa ảnh hưởng của các yếu tố hình thái lên độc tính (biến nhị phân) và môi trường sống (biến đa loại) của nấm.
Dữ liệu được sử dụng là bộ dữ liệu Mushroom kinh điển, được thu thập và lưu trữ tại kho lưu trữ của UCI.
data <- read.csv("D:/PTDLDT/Mushroom.csv", header = T)
datatable(data)
Bộ dữ liệu Mushroom bao gồm 8124 quan sát và 23 biến, trong đó biến
class là biến mục tiêu chính (‘p’ cho poisonous - độc, ‘e’ cho edible -
ăn được) và 22 biến còn lại mô tả các đặc điểm vật lý của nấm. Đáng chú
ý, tất cả các biến này đều là biến định tính. Ngoài ra, nghiên cứu này
cũng xem xét biến habitat (môi trường sống) như một biến
phụ thuộc thứ hai để khám phá các mối liên kết sinh thái.
Dưới đây là bảng mô tả ý nghĩa của một số biến chính trong bộ dữ liệu.
| STT | Tên biến | Mô tả | Các giá trị và ý nghĩa |
|---|---|---|---|
| 1 | class | Loại nấm (Biến mục tiêu) | e=edible (ăn được), p=poisonous (độc) |
| 2 | cap.shape | Hình dạng mũ nấm | b=bell (hình chuông), c=conical (hình nón), x=convex (lồi), f=flat (phẳng), k=knobbed (có núm), s=sunken (lõm) |
| 3 | cap.surface | Bề mặt mũ nấm | f=fibrous (dạng sợi), g=grooves (có rãnh), y=scaly (có vảy), s=smooth (mịn) |
| 4 | cap.color | Màu sắc mũ nấm | n=brown (nâu), b=buff (da bò), c=cinnamon (quế), g=gray (xám), r=green (xanh lá), p=pink (hồng), u=purple (tím), e=red (đỏ), w=white (trắng), y=yellow (vàng) |
| 5 | bruises | Có bị thâm khi va chạm? | t=bruises (có), f=no (không) |
| 6 | odor | Mùi của nấm |
a=almond (hạnh nhân), l=anise (hồi), c=creosote (dầu creosote), y=fishy
(tanh), f=foul (hôi), m=musty (mốc), n=none (không mùi), p=pungent
(hăng), s=spicy (cay)
|
| 7 | gill.attachment | Cách phiến nấm gắn vào thân |
a=attached (gắn liền), d=descending (chạy xuống), f=free (tự do),
n=notched (có khía)
|
| 8 | gill.spacing | Khoảng cách giữa các phiến nấm | c=close (gần nhau), w=crowded (chật chội), d=distant (thưa) |
| 9 | gill.size | Kích thước phiến nấm | b=broad (rộng), n=narrow (hẹp) |
| 10 | gill.color | Màu sắc phiến nấm |
k=black (đen), n=brown (nâu), b=buff (da bò), h=chocolate
(sô-cô-la), g=gray (xám), o=orange (cam), p=pink (hồng), r=green (xanh
lá), u=purple (tím), e=red (đỏ), w=white (trắng), y=yellow (vàng)
|
| 11 | stalk.shape | Hình dạng thân nấm | e=enlarging (phình ra),t=tapering (thon lại) |
| 12 | stalk.root | Gốc thân nấm | b=bulbous (hình củ), c=club (hình dùi cui), u=cup (hình chén), e=equal (đều), z=rhizomorphs (dạng rễ giả), r=rooted (có rễ) |
| 13 | stalk.surface.above.ring | Bề mặt thân nấm (phía trên vành) | f=fibrous (dạng sợi), y=scaly (có vảy), k=silky (mượt như lụa), s=smooth (mịn) |
| 14 | stalk.surface.below.ring | Bề mặt thân nấm (phía dưới vành) | f=fibrous (dạng sợi), y=scaly (có vảy), k=silky (mượt như lụa), s=smooth (mịn) |
| 15 | stalk.color.above.ring | Màu sắc thân nấm (phía trên vành) | n=brown (nâu), b=buff (da bò), c=cinnamon (quế), g=gray (xám), o=orange (cam), p=pink (hồng), e=red (đỏ), w=white (trắng), y=yellow (vàng) |
| 16 | stalk.color.below.ring | Màu sắc thân nấm (phía dưới vành) | n=brown (nâu), b=buff (da bò), c=cinnamon (quế), g=gray (xám), o=orange (cam), p=pink (hồng), e=red (đỏ), w=white (trắng), y=yellow (vàng) |
| 17 | veil.type | Loại màng che (Veil) | p=partial (một phần), u=universal (toàn phần) |
| 18 | veil.color | Màu sắc màng che | n=brown (nâu), o=orange (cam), w=white (trắng), y=yellow (vàng) |
| 19 | ring.number | Số lượng vành nấm | n=none (không có), o=one (một), t=two (hai) |
| 20 | ring.type | Loại vành nấm | c=cobwebby (dạng mạng nhện), e=evanescent (dễ rụng), f=flaring (loe ra), l=large (lớn), n=none (không có), p=pendant (treo lủng lẳng), s=sheathing (dạng bao), z=zone (vân) |
| 21 | spore.print.color | Màu sắc của bào tử | k=black (đen), n=brown (nâu), b=buff (da bò), h=chocolate (sô-cô-la), r=green (xanh lá), o=orange (cam), u=purple (tím), w=white (trắng), y=yellow (vàng) |
| 22 | population | Sự phân bố, quần thể | a=abundant (dồi dào), c=clustered (mọc cụm), n=numerous (số lượng lớn), s=scattered (rải rác), v=several (vài cây), y=solitary (mọc đơn) |
| 23 | habitat | Môi trường sống | g=grasses (trên cỏ), l=leaves (trên lá), m=meadows (đồng cỏ), p=paths (lối đi), u=urban (đô thị), w=waste (bãi rác), d=woods (trong rừng) |
Để có được cái nhìn tổng quan và bước đầu đánh giá tính phù hợp của
tập dữ liệu với mục tiêu nghiên cứu, ta tiến hành kiểm tra cấu trúc của
bộ dữ liệu bằng cách sử dụng hai hàm cơ bản trong R là
dim() và str().
dim() cho biết kích thước của bảng dữ liệu, cụ thể
là số hàng (quan sát) và số cột (biến).dim(data)
## [1] 8124 23
Bộ dữ liệu có 8124 quan sát (mẫu nấm) và 23 biến (đặc điểm). Đây là một bộ dữ liệu có kích thước đủ lớn để thực hiện các phân tích thống kê có ý nghĩa.
str() cung cấp thông tin về cấu trúc dữ liệu, bao
gồm tên các biến, kiểu dữ liệu của từng biến (ví dụ: số, chuỗi ký tự,
biến phân loại), cũng như một số giá trị minh họa cho mỗi biến.str(data)
## 'data.frame': 8124 obs. of 23 variables:
## $ class : chr "p" "e" "e" "p" ...
## $ cap.shape : chr "x" "x" "b" "x" ...
## $ cap.surface : chr "s" "s" "s" "y" ...
## $ cap.color : chr "n" "y" "w" "w" ...
## $ bruises : chr "t" "t" "t" "t" ...
## $ odor : chr "p" "a" "l" "p" ...
## $ gill.attachment : chr "f" "f" "f" "f" ...
## $ gill.spacing : chr "c" "c" "c" "c" ...
## $ gill.size : chr "n" "b" "b" "n" ...
## $ gill.color : chr "k" "k" "n" "n" ...
## $ stalk.shape : chr "e" "e" "e" "e" ...
## $ stalk.root : chr "e" "c" "c" "e" ...
## $ stalk.surface.above.ring: chr "s" "s" "s" "s" ...
## $ stalk.surface.below.ring: chr "s" "s" "s" "s" ...
## $ stalk.color.above.ring : chr "w" "w" "w" "w" ...
## $ stalk.color.below.ring : chr "w" "w" "w" "w" ...
## $ veil.type : chr "p" "p" "p" "p" ...
## $ veil.color : chr "w" "w" "w" "w" ...
## $ ring.number : chr "o" "o" "o" "o" ...
## $ ring.type : chr "p" "p" "p" "p" ...
## $ spore.print.color : chr "k" "n" "n" "k" ...
## $ population : chr "s" "n" "n" "s" ...
## $ habitat : chr "u" "g" "m" "u" ...
Tất cả các biến hiện đang ở dạng ký tự (chr). Để R có thể phân tích chúng dưới dạng biến định tính, chúng ta cần chuyển đổi chúng sang kiểu factor.
Dữ liệu trong thực tế thường không hoàn hảo. Một số nguồn tài liệu cho biết giá trị bị thiếu trong bộ dữ liệu này được ký hiệu là ‘?’. Chúng ta cần kiểm tra và xử lý chúng.
# Sử dụng hàm na_if() kết hợp với across() để thay thế '?' một cách an toàn.
mushroom_df_processed <- data %>%
dplyr::mutate(dplyr::across(everything(), ~na_if(., "?")))
# Đếm số lượng giá trị NA trong mỗi cột
na_counts <- colSums(is.na(mushroom_df_processed))
# Lấy ra những cột có chứa giá trị NA (số lượng NA > 0)
cols_with_na <- na_counts[na_counts > 0]
# In ra kết quả kiểm tra
if (length(cols_with_na) > 0) {
cat("\nPhát hiện các giá trị NA trong các cột sau:\n")
print(cols_with_na)
} else {
cat("\nKhông tìm thấy giá trị NA nào trong bộ dữ liệu.\n")
}
##
## Phát hiện các giá trị NA trong các cột sau:
## stalk.root
## 2480
# XỬ LÝ NA BẰNG PHƯƠNG PHÁP THAY THẾ BẰNG MODE
# Hàm để tìm mode (giá trị phổ biến nhất) của một vector
get_mode <- function(v) {
# Loại bỏ NA trước khi tìm mode để không bị lỗi
v_without_na <- na.omit(v)
# Tìm giá trị duy nhất và đếm tần suất của chúng
uniqv <- unique(v_without_na)
# Trả về giá trị có tần suất cao nhất
uniqv[which.max(tabulate(match(v_without_na, uniqv)))]
}
# Áp dụng thay thế NA bằng mode cho cột 'stalk.root'
# (Dựa vào kết quả ở trên, chúng ta biết chỉ cột này có NA)
cat("\nThực hiện thay thế giá trị NA trong cột 'stalk.root' bằng mode...\n")
##
## Thực hiện thay thế giá trị NA trong cột 'stalk.root' bằng mode...
mode_stalk_root <- get_mode(mushroom_df_processed$stalk.root)
cat("Giá trị mode của cột 'stalk.root' là:", mode_stalk_root, "\n")
## Giá trị mode của cột 'stalk.root' là: b
# Thực hiện việc thay thế
mushroom_df_processed$stalk.root[is.na(mushroom_df_processed$stalk.root)] <- mode_stalk_root
# Đếm lại tổng số NA trong toàn bộ data frame để xác nhận
final_na_count <- sum(is.na(mushroom_df_processed))
cat("\nKiểm tra lại sau khi xử lý...\n")
##
## Kiểm tra lại sau khi xử lý...
cat("Tổng số giá trị NA còn lại trong bộ dữ liệu:", final_na_count, "\n")
## Tổng số giá trị NA còn lại trong bộ dữ liệu: 0
if (final_na_count == 0) {
cat("=> Xử lý dữ liệu thiếu thành công!\n")
} else {
cat("=> CẢNH BÁO: Vẫn còn giá trị NA trong dữ liệu!\n")
}
## => Xử lý dữ liệu thiếu thành công!
Nhận xét:
Chỉ có một cột duy nhất chứa giá trị thiếu là stalk.root với 2480 giá trị bị thiếu, chiếm khoảng 30.5% tổng số quan sát.
Việc loại bỏ các hàng này sẽ làm mất một lượng lớn dữ liệu. Thay vào đó, một chiến lược hợp lý là thay thế các giá trị thiếu bằng mode (giá trị xuất hiện nhiều nhất) của cột đó. Điều này giúp bảo toàn kích thước mẫu.
Kết quả: Dữ liệu đã được làm sạch, không còn giá trị thiếu.
Để R hiểu đúng các biến là biến định tính, chúng ta chuyển đổi tất cả
các cột sang kiểu factor.
# Chuyển đổi tất cả các cột sang factor
df <- mushroom_df_processed %>% mutate(across(everything(), as.factor))
# Kiểm tra lại cấu trúc sau khi chuyển đổi
str(df)
## 'data.frame': 8124 obs. of 23 variables:
## $ class : Factor w/ 2 levels "e","p": 2 1 1 2 1 1 1 1 2 1 ...
## $ cap.shape : Factor w/ 6 levels "b","c","f","k",..: 6 6 1 6 6 6 1 1 6 1 ...
## $ cap.surface : Factor w/ 4 levels "f","g","s","y": 3 3 3 4 3 4 3 4 4 3 ...
## $ cap.color : Factor w/ 10 levels "b","c","e","g",..: 5 10 9 9 4 10 9 9 9 10 ...
## $ bruises : Factor w/ 2 levels "f","t": 2 2 2 2 1 2 2 2 2 2 ...
## $ odor : Factor w/ 9 levels "a","c","f","l",..: 7 1 4 7 6 1 1 4 7 1 ...
## $ gill.attachment : Factor w/ 2 levels "a","f": 2 2 2 2 2 2 2 2 2 2 ...
## $ gill.spacing : Factor w/ 2 levels "c","w": 1 1 1 1 2 1 1 1 1 1 ...
## $ gill.size : Factor w/ 2 levels "b","n": 2 1 1 2 1 1 1 1 2 1 ...
## $ gill.color : Factor w/ 12 levels "b","e","g","h",..: 5 5 6 6 5 6 3 6 8 3 ...
## $ stalk.shape : Factor w/ 2 levels "e","t": 1 1 1 1 2 1 1 1 1 1 ...
## $ stalk.root : Factor w/ 4 levels "b","c","e","r": 3 2 2 3 3 2 2 2 3 2 ...
## $ stalk.surface.above.ring: Factor w/ 4 levels "f","k","s","y": 3 3 3 3 3 3 3 3 3 3 ...
## $ stalk.surface.below.ring: Factor w/ 4 levels "f","k","s","y": 3 3 3 3 3 3 3 3 3 3 ...
## $ stalk.color.above.ring : Factor w/ 9 levels "b","c","e","g",..: 8 8 8 8 8 8 8 8 8 8 ...
## $ stalk.color.below.ring : Factor w/ 9 levels "b","c","e","g",..: 8 8 8 8 8 8 8 8 8 8 ...
## $ veil.type : Factor w/ 1 level "p": 1 1 1 1 1 1 1 1 1 1 ...
## $ veil.color : Factor w/ 4 levels "n","o","w","y": 3 3 3 3 3 3 3 3 3 3 ...
## $ ring.number : Factor w/ 3 levels "n","o","t": 2 2 2 2 2 2 2 2 2 2 ...
## $ ring.type : Factor w/ 5 levels "e","f","l","n",..: 5 5 5 5 1 5 5 5 5 5 ...
## $ spore.print.color : Factor w/ 9 levels "b","h","k","n",..: 3 4 4 3 4 3 3 4 3 3 ...
## $ population : Factor w/ 6 levels "a","c","n","s",..: 4 3 3 4 1 3 3 4 5 4 ...
## $ habitat : Factor w/ 7 levels "d","g","l","m",..: 6 2 4 6 2 2 4 4 2 4 ...
Nhận xét: Tất cả các biến giờ đã có kiểu factor với các cấp độ (levels) tương ứng, sẵn sàng cho việc phân tích.
Phần này tập trung vào việc hiểu phân phối của một vài đặc điểm nấm một cách riêng lẻ.
Trước khi đi sâu vào các mối quan hệ phức tạp, chúng ta cần hiểu rõ sự phân bổ của biến quan trọng nhất: độc tính. Tỷ lệ nấm độc và nấm ăn được trong bộ dữ liệu này như thế nào?
df <- df %>%
mutate(class_full = recode(class, "e" = "Ăn được","p" = "Độc"))
class_summary <- df %>%
count(class_full) %>%
mutate(Percentage = n / sum(n))
colnames(class_summary) <- c("Loại nấm", "Tần số", "Tần suất")
kable(class_summary, align = "c", booktabs = TRUE, digits = 2, caption = "Thống kê tần suất của các loại nấm") %>%
kable_styling(full_width = FALSE, position = "center")
| Loại nấm | Tần số | Tần suất |
|---|---|---|
| Ăn được | 4208 | 0.52 |
| Độc | 3916 | 0.48 |
Trực quan hoá
bar_plot <- ggplot(df, aes(x = class_full, fill = class_full)) +
geom_bar() +
geom_text(stat = "count", aes(label = after_stat(count)), vjust = 0, size = 3) +
labs(title = "Phân bố loại nấm",
x = "Loại",
y = "Số lượng", fill = "Loại nấm") +
theme_minimal(base_size = 14) +
scale_fill_manual(values = c("Ăn được" = "skyblue", "Độc" = "pink"))
class_p <- as.data.frame(table(df$class_full))
colnames(class_p) <- c("Category", "Count")
pie_plot <- ggplot(class_p, aes(x = "", y = Count, fill = Category)) +
geom_bar(width = 1, stat = "identity") +
coord_polar("y", start=0) +
labs(title = "Tỷ lệ nấm", fill = "Loại nấm") +
theme_void(base_size = 14) +
geom_text(aes(label = scales::percent(Count/sum(Count), accuracy = 0.0001)), position = position_stack(vjust = 0.5), size = 2.5)+scale_fill_manual(values = c("Ăn được" = "salmon", "Độc" = "#00BFC4"))
bar_plot + pie_plot
Nhận xét:
Bộ dữ liệu khá cân bằng. Nấm ăn được (edible - e) chiếm 51.8% (4208 mẫu), trong khi nấm độc (poisonous - p) chiếm 48.2% (3916 mẫu).
Sự cân bằng này là một lợi thế lớn, giúp các phân tích và mô hình sau này không bị thiên vị về một trong hai loại, đảm bảo rằng các quy luật chúng ta tìm thấy sẽ có giá trị cho cả hai nhóm.
Môi trường sống (habitat) là nơi nấm sinh trưởng và phát triển, tiết lộ nhiều thông tin về hệ sinh thái và nhu cầu dinh dưỡng của chúng. Một số loài nấm có mối quan hệ cộng sinh với các loài cây cụ thể, trong khi những loài khác phát triển trên gỗ mục, lá cây, hoặc trong đồng cỏ. Việc xác định môi trường sống là một bước quan trọng trong quy trình định danh nấm ngoài tự nhiên. Bộ dữ liệu phân loại môi trường sống thành: grasses (đồng cỏ), leaves (lá cây), meadows (bãi cỏ), paths (lối đi), urban (đô thị), waste (chất thải), và woods (rừng). Chúng ta hãy xem xét các môi trường sống phổ biến nhất trong bộ dữ liệu này.
# Tạo nhãn Tiếng Việt
habitat_labels <- c( "d" = "Trong rừng", "g" = "Trên cỏ", "l" = "Trên lá", "m" = "Đồng cỏ", "p" = "Lối đi", "u" = "Đô thị", "w" = "Bãi rác")
df$habitat_full <- factor(habitat_labels[df$habitat])
# Thống kê tần suất
habitat_summary <- df %>%
count(habitat_full) %>%
mutate(Percentage = (n / sum(n)))
kable(habitat_summary,
col.names = c("Môi trường sống", "Tần số", "Tần suất"),
align = "c", booktabs = TRUE,
caption = "Thống kê tần suất môi trường sống của các loại nấm") %>%
kable_styling(bootstrap_options = c("striped", "hover"), full_width = T)
| Môi trường sống | Tần số | Tần suất |
|---|---|---|
| Bãi rác | 192 | 0.0236337 |
| Đô thị | 368 | 0.0452979 |
| Đồng cỏ | 292 | 0.0359429 |
| Lối đi | 1144 | 0.1408173 |
| Trên cỏ | 2148 | 0.2644018 |
| Trên lá | 832 | 0.1024126 |
| Trong rừng | 3148 | 0.3874938 |
Trực quan hoá
ggplot(df, aes(x = fct_infreq(habitat_full), fill = fct_infreq(habitat_full))) +
geom_bar() +
geom_text(stat = "count", aes(label = after_stat(count)), vjust = -0.5, size = 3.5) +
labs(title = "Phân bố môi trường sống của nấm", x = "Môi trường sống", y = "Số lượng", fill = "Môi trường sống") +
theme_minimal() +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
Nhận xét:
Biểu đồ cho thấy môi trường sống phổ biến nhất của các loài nấm trong bộ dữ liệu là trong rừng (woods - d), với 3148 mẫu. Điều này phản ánh thực tế rằng rừng là hệ sinh thái đa dạng và lý tưởng cho sự phát triển của nấm.
Đồng cỏ (grasses - g) cũng là một môi trường sống quan trọng, đứng thứ hai với 2148 mẫu.
Các môi trường như lối đi (paths - p), trên lá (leaves - l) và đô thị (urban - u) có tần suất thấp hơn nhưng vẫn đáng kể.
Môi trường sống trên chất thải (waste - w) và trên bãi cỏ (meadows - m) là hiếm gặp nhất. Việc một cây nấm mọc ở đâu có thể là một gợi ý mạnh mẽ về danh tính của nó, và có khả năng liên quan đến việc nó có độc hay không.
Phần này sử dụng các phương pháp thống kê suy luận để đưa ra kết luận về tổng thể từ dữ liệu mẫu. Đối với mỗi biến, chúng ta sẽ chọn một hạng mục quan tâm (thường là hạng mục phổ biến nhất) để thực hiện:
Ước lượng khoảng tin cậy 95% cho tỷ lệ của hạng mục đó trong tổng thể.
Kiểm định giả thuyết để xem liệu tỷ lệ thực tế của hạng mục đó có khác biệt so với một giá trị giả định hay không.
Phương pháp này giúp chúng ta hiểu rõ hơn về mức độ phổ biến của từng đặc điểm cụ thể.
Chúng ta sẽ ước lượng khoảng tin cậy và kiểm định giả thuyết cho tỷ lệ nấm độc trong tổng thể.
# Lấy số lượng cho từng hạng mục
class_counts <- table(df$class)
total_mushrooms <- sum(class_counts)
cat("Tổng số lượng nấm:", total_mushrooms, "\n")
## Tổng số lượng nấm: 8124
poisonous_count <- class_counts["p"]
cat("Số lượng nấm độc:", poisonous_count, "\n")
## Số lượng nấm độc: 3916
Ước lượng khoảng tin cậy 95% cho tỷ lệ nấm độc
prop_test_poisonous <- prop.test(x = poisonous_count, n = total_mushrooms, conf.level = 0.95)
prop_test_poisonous$conf.int
## [1] 0.4711126 0.4929616
## attr(,"conf.level")
## [1] 0.95
Diễn giải: Với độ tin cậy 95%, chúng ta có thể kết luận rằng tỷ lệ thực sự của nấm độc trong tổng thể (mà mẫu này đại diện) nằm trong khoảng từ 47.1% đến 49.3%.
Chúng ta sẽ phân tích tỷ lệ nấm sinh sống trong rừng, môi trường phổ biến nhất.
habitat_counts <- table(df$habitat)
woods_habitat_count <- habitat_counts["g"]
cat("Tổng số lượng nấm:", total_mushrooms, "\n")
## Tổng số lượng nấm: 8124
cat("Số lượng nấm sống trong rừng:", woods_habitat_count, "\n")
## Số lượng nấm sống trong rừng: 2148
Ước lượng khoảng tin cậy 95% cho tỷ lệ nấm sống “Trong rừng”
prop_test_habitat <- prop.test(x = woods_habitat_count, n = total_mushrooms, conf.level = 0.95)
prop_test_habitat$conf.int
## [1] 0.2548640 0.2741637
## attr(,"conf.level")
## [1] 0.95
Diễn giải: Với độ tin cậy 95%, tỷ lệ thực sự của nấm sống trong rừng trong tổng thể được ước tính nằm trong khoảng từ 25.486% đến 27.41637%.
Kiểm định giả thuyết
Bài toán kiểm định: \[ \begin{cases} H_0 : p = 0.26 \quad (\text{Tỷ lệ nấm sống trong rừng trong tổng thể bằng 26%}) \\ H_1: p \neq 0.26 \quad (\text{Tỷ lệ nấm sống trong rừng trong tổng thể khác 26%}) \end{cases} \]
prop_test_habitat_H0 <- prop.test(x = woods_habitat_count, n = total_mushrooms, p = 0.26, alternative = "two.sided")
print(prop_test_habitat_H0)
##
## 1-sample proportions test with continuity correction
##
## data: woods_habitat_count out of total_mushrooms, null probability 0.26
## X-squared = 0.79541, df = 1, p-value = 0.3725
## alternative hypothesis: true p is not equal to 0.26
## 95 percent confidence interval:
## 0.2548640 0.2741637
## sample estimates:
## p
## 0.2644018
Kết luận: Giá trị p-value = 0.37254, lớn hơn mức ý nghĩa α = 0.05. Do đó, chúng ta không đủ bằng chứng để bác bỏ giả thuyết \(H_0\). Tỷ lệ nấm sống trong rừng trong tổng thể bằng 26%.
Phần này là trọng tâm của báo cáo, tập trung vào việc phân tích mối quan hệ giữa các biến để trả lời hai câu hỏi nghiên cứu chính:
Câu hỏi 1 (Về độc tính): Những đặc điểm vật lý nào có mối liên hệ mạnh mẽ nhất với việc một cây nấm là độc hay ăn được?
Câu hỏi 2 (Về môi trường sống): Liệu các đặc điểm như mùi, màu sắc, hình dạng có thay đổi một cách có hệ thống theo từng môi trường sống khác nhau không?
Để trả lời các câu hỏi trên, chúng ta sẽ đi sâu vào phân tích từng cặp biến quan trọng, áp dụng một chuỗi các phương pháp từ mô tả (bảng tần số, biểu đồ) đến suy luận (kiểm định Chi-bình phương) để hiểu rõ bản chất của mỗi mối quan hệ.
Trước tiên, chúng ta sẽ tập trung vào Câu hỏi 1: Tìm kiếm các yếu tố dự báo độc tính. Chúng ta sẽ lần lượt kiểm tra các đặc điểm đáng ngờ nhất, bắt đầu bằng một trong những kinh nghiệm dân gian phổ biến nhất liên quan đến vết bầm của nấm.
Câu hỏi nghiên cứu: Việc một cây nấm có vết bầm ảnh hưởng như thế nào đến khả năng nó là nấm độc?
Bảng tần số chéo và trực quan hóa
# Tạo bảng tần số
class_bruises_table <- table(df$bruises, df$class)
# Gán tên hàng
rownames(class_bruises_table) <- c("Không có vết bầm", "Có vết bầm")
# Hiển thị bảng với kable
kable(class_bruises_table,
col.names = c("Ăn được", "Độc"),
caption = "Bảng tần số giữa vết bầm và độc tính") %>%
kable_styling(bootstrap_options = "striped", full_width = F)
| Ăn được | Độc | |
|---|---|---|
| Không có vết bầm | 1456 | 3292 |
| Có vết bầm | 2752 | 624 |
# Trực quan hóa bằng biểu đồ cột nhóm
ggplot(df, aes(x = bruises, fill = class_full)) +
geom_bar(position = "dodge") +
geom_text(stat = "count", aes(label = after_stat(count)),
position = position_dodge(width = 0.9), vjust = -0.2, size = 3) +
labs(title = "Mối quan hệ giữa vết bầm và độc tính",
x = "Có vết bầm (f: không, t: có)", y = "Số lượng", fill = "Loại nấm") +
theme_minimal(base_size = 12)
Nhận xét:
Bảng và biểu đồ cho thấy một sự tương phản rõ rệt. Trong tổng số 4748 cây nấm không có vết bầm (f), có tới 3292 cây là nấm độc (chiếm một tỷ lệ áp đảo ~ 69%). Ngược lại, trong 3376 cây nấm có vết bầm (t), chỉ có 624 cây là nấm độc (chiếm ~ 18.5%).
Điều này gợi ý mạnh mẽ rằng việc có vết bầm là một đặc điểm liên quan đến nấm ăn được.
Kiểm định Chi-bình phương (Chi-squared Test)
Bài toán kiểm định: \[ \begin{cases} H_0: \text{Việc có vết bầm và độc tính của nấm là hai biến độc lập.}\\ H_1: \text{Việc có vết bầm và độc tính của nấm có sự liên quan.} \end{cases} \]
chisq.test(class_bruises_table)
##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: class_bruises_table
## X-squared = 2041.4, df = 1, p-value < 2.2e-16
Trong đó
X-squared = 2041.4: Đây là giá trị thống kê kiểm
định Chi-bình phương. Nó đo lường sự khác biệt giữa tần suất quan sát
được (số liệu thực tế trong bảng của bạn) và tần suất kỳ vọng (số liệu
mà chúng ta mong đợi nếu hai biến hoàn toàn độc lập với nhau).
df = 1: Bậc tự do (degrees of freedom). Đối với kiểm
định Chi-bình phương, bậc tự do được tính bằng công thức (số hàng - 1) *
(số cột - 1). Vì bảng của bạn là 2x2, nên df = (2 - 1) * (2 - 1) =
1.
p-value < 2.2e-16: Đây là kết quả quyết định. Nếu
việc có vết bầm và độc tính của nấm thực sự không liên quan gì đến nhau
(H₀ đúng), thì xác suất để chúng ta quan sát được một mối liên hệ mạnh
mẽ như vậy trong mẫu là gần như bằng không.
Kết luận: Với p-value < 2.2e-16, chúng ta bác bỏ giả thuyết \(H_0\). Có bằng chứng thống kê không thể chối cãi về một mối liên hệ mạnh mẽ giữa vết bầm và độc tính. Việc có vết bầm là một chỉ báo rất mạnh cho thấy nấm có khả năng cao là ăn được, trong khi không có vết bầm là một dấu hiệu cảnh báo nguy hiểm mà người thu hái nấm phải hết sức lưu ý.
Rủi ro tương đối (Relative Risk - RR)
Định nghĩa: RR so sánh xác suất (nguy cơ) một cây nấm là độc giữa nhóm “phơi nhiễm” (có vết bầm) và nhóm “không phơi nhiễm” (không có vết bầm).
# Bảng cần có dạng: cột là kết quả, hàng là phơi nhiễm
# Chúng ta định nghĩa "kết quả" là Độc (p) và "phơi nhiễm" là Có vết bầm (t)
rr_table_bruises <- table(df$bruises, df$class)[, c("p", "e")]
# Tính toán RR
rr_result_bruises <- riskratio.wald(rr_table_bruises, rev="rows") # rev="rows" vì nhóm phơi nhiễm (t) ở hàng 2
rr_result_bruises
## $data
##
## p e Total
## t 624 2752 3376
## f 3292 1456 4748
## Total 3916 4208 8124
##
## $measure
## risk ratio with 95% C.I.
## estimate lower upper
## t 1.0000000 NA NA
## f 0.3761878 0.3593876 0.3937733
##
## $p.value
## two-sided
## midp.exact fisher.exact chi.square
## t NA NA NA
## f 0 0 0
##
## $correction
## [1] FALSE
##
## attr(,"method")
## [1] "Unconditional MLE & normal approximation (Wald) CI"
Trong đó
t là nhóm tham chiếu (reference group) với Risk Ratio = 1.000.
f có risk ratio = 0.37621878, nghĩa là:
Nguy cơ một cây nấm là độc khi nó không có vết bầm (f) chỉ bằng 0.376 lần (tức 37.6%) so với nguy cơ của một cây nấm có vết bầm (t)“. Hay, nguy cơ một cây nấm là độc khi nó có vết bầm (t) cao gấp 1 / 0.376 = 2.66 lần so với khi không có vết bầm.
Khoảng tin cậy 95%: [0.3594; 0.3938] cho thấy kết quả chính xác và có ý nghĩa thống kê cao (vì không chứa 1).
Diễn giải kết quả Relative Risk:
Relative Risk (RR) ≈ 0.3761878. Điều này có nghĩa là: Nguy cơ một cây nấm là độc ở nhóm có vết bầm chỉ bằng 0.3761878 lần so với nguy cơ ở nhóm không có vết bầm.
Nói cách khác, việc có vết bầm làm giảm 62.38% (1 - 0.3762) nguy cơ một cây nấm là độc.
Khoảng tin cậy 95% cho RR là [0.3594; 0.3938]. Vì khoảng này không chứa 1, nó khẳng định mối liên hệ có ý nghĩa thống kê.
Tỷ số chênh (Odds Ratio - OR)
Định nghĩa: OR so sánh “odds” của việc một cây nấm là độc giữa hai nhóm.
or_result_bruises <- oddsratio.wald(rr_table_bruises, rev="rows")
or_result_bruises$measure
## odds ratio with 95% C.I.
## estimate lower upper
## t 1.0000000 NA NA
## f 0.1002854 0.09014769 0.1115632
Trong đó
t 1.0000000: Nhóm t (có vết bầm) vẫn là nhóm tham
chiếu. OR của nó so với chính nó là 1.
f 0.1002854: Đây là Tỷ số Chênh (OR) của nhóm f
(không có vết bầm) so với nhóm t.
Diễn giải: “Tỷ số chênh (odds) một cây nấm là độc khi nó không có vết bầm chỉ bằng 0.10 lần (tức 10%) so với tỷ số chênh khi nó có vết bầm”. Hay, tỷ số chênh một cây nấm là độc khi nó có vết bầm cao gấp 1 / 0.10 = 10 lần so với khi không có vết bầm.
lower 0.09014769 upper 0.1115632: Khoảng tin cậy 95%
cho OR là từ 0.09 đến 0.11. Khoảng này cũng không chứa 1.0, cho thấy một
mối liên hệ rất mạnh và có ý nghĩa thống kê.
Diễn giải kết quả Odds Ratio:
Odds Ratio (OR) ≈ 0.1003. Điều này có nghĩa là: Odds của việc một cây nấm là độc khi nó có vết bầm chỉ bằng 0.1003 lần (tức 10.03%) so với odds của việc nó là độc khi không có vết bầm.
Khoảng tin cậy 95% cho OR là [0.0901; 0.1116], không chứa 1, khẳng định mối liên hệ có ý nghĩa thống kê rất cao.
Câu hỏi nghiên cứu: Mùi hương của nấm có phải là một yếu tố dự báo mạnh về độc tính hay không?
Bảng tần số chéo và trực quan hóa
# Tạo một "từ điển" để ánh xạ các ký tự sang tên đầy đủ
odor_labels <- c(
"a" = "Hạnh nhân", "l" = "Hồi", "c" = "Creosote",
"y" = "Tanh", "f" = "Hôi", "m" = "Mốc",
"n" = "Không mùi", "p" = "Hăng", "s" = "Cay")
df <- df %>%
mutate(odor_full = recode(odor, !!!odor_labels))
class_odor_table <- table(Mùi = df$odor_full, Độc_tính = df$class_full)
kable(class_odor_table, caption = "Bảng tần số giữa mùi hương và loại nấm") %>%
kable_styling(bootstrap_options = "striped", full_width = FALSE, position = "center") %>%
# Điều chỉnh độ rộng cột để bảng trông đẹp hơn
column_spec(1, width = "10em", bold = TRUE) %>%
column_spec(2:3, width = "8em")
| Ăn được | Độc | |
|---|---|---|
| Hạnh nhân | 400 | 0 |
| Creosote | 0 | 192 |
| Hôi | 0 | 2160 |
| Hồi | 400 | 0 |
| Mốc | 0 | 36 |
| Không mùi | 3408 | 120 |
| Hăng | 0 | 256 |
| Cay | 0 | 576 |
| Tanh | 0 | 576 |
# Trực quan hóa
ggplot(df, aes(x = fct_infreq(odor_full), fill = class_full)) +
geom_bar(position = "dodge") + # <-- Xếp cạnh nhau
labs(title = "Mối quan hệ giữa mùi hương và độc tính của nấm",
x = "Mùi của nấm", y = "Tần suất", fill = "Loại nấm") +
theme_minimal(base_size = 12) +
scale_fill_brewer(palette = "Set1", labels = c("Ăn được", "Độc"))+theme(axis.text.x = element_text(angle = 45, hjust = 1))
Nhận xét:
Sự tách biệt gần như hoàn hảo giữa hai nhóm. Các loại nấm có mùi dễ chịu như Hạnh nhân, Hồi, hoặc Không mùi chiếm phần lớn trong nhóm ăn được. Ngược lại, tất cả các loại nấm có mùi khó chịu như Hôi, Tanh, Cay, Hăng đều 100% là nấm độc, không có ngoại lệ.
Trong tất cả các đặc điểm được khảo sát, mùi hương (odor) nổi lên như một yếu tố dự báo độc tính mạnh mẽ và đáng tin cậy nhất. Biểu đồ cho thấy không có sự chồng chéo giữa nhóm mùi an toàn và nhóm mùi nguy hiểm, giúp việc phân loại trở nên gần như tuyệt đối chỉ dựa vào đặc điểm này.
Từ mối quan hệ này, có thể rút ra một quy tắc vàng đơn giản và an toàn cho thực tế: “Nếu nấm có bất kỳ mùi lạ hoặc khó chịu nào, hãy mặc định nó là độc và tuyệt đối không sử dụng”. Điều này nhấn mạnh tầm quan trọng của việc kiểm tra cảm quan (mùi) khi nhận dạng nấm.
Kiểm định Chi-bình phương (Chi-squared Test) \[ \begin{cases} H_0: \text{Mùi và độc tính là hai biến độc lập.}\\ H_1: \text{Mùi và độc tính có sự liên quan.} \end{cases} \]
chisq.test(class_odor_table)
##
## Pearson's Chi-squared test
##
## data: class_odor_table
## X-squared = 7659.7, df = 8, p-value < 2.2e-16
Kết luận: Với p-value < 2.2e-16, chúng ta bác bỏ \(H_0\) và kết luận có một mối liên hệ thống kê cực kỳ mạnh mẽ. Mùi là yếu tố dự báo độc tính mạnh nhất. Có thể nói, chỉ cần dựa vào mùi, chúng ta đã có thể phân loại chính xác gần như toàn bộ nấm trong bộ dữ liệu này.
Vì biến odor có nhiều hơn 2 hạng mục (9 loại mùi) nên không thể tính một giá trị RR hoặc OR duy nhất cho toàn bộ bảng.
Câu hỏi nghiên cứu: Màu sắc bào tử có liên quan đến việc một cây nấm là độc hay ăn được không?
Bảng tần số chéo và trực quan hóa
spore_labels <- c(
"k" = "Đen", "n" = "Nâu", "b" = "Da bò",
"h" = "Sô-cô-la", "r" = "Xanh lá", "o" = "Cam",
"u" = "Tím", "w" = "Trắng", "y" = "Vàng"
)
df <- df %>%
mutate(spore_color_full = recode(spore.print.color, !!!spore_labels))
class_spore_table <- table(df$spore_color_full, df$class_full)
kable(class_spore_table, caption = "Bảng tần số chéo giữa màu sắc bào tử và đọc tính của nấm") %>%
kable_styling(bootstrap_options = "striped", position = "center", full_width = F)%>%
column_spec(1, width = "5em") %>%
column_spec(2:ncol(class_odor_table)+1, width = "12em")
| Ăn được | Độc | |
|---|---|---|
| Da bò | 48 | 0 |
| Sô-cô-la | 48 | 1584 |
| Đen | 1648 | 224 |
| Nâu | 1744 | 224 |
| Cam | 48 | 0 |
| Xanh lá | 0 | 72 |
| Tím | 48 | 0 |
| Trắng | 576 | 1812 |
| Vàng | 48 | 0 |
ggplot(df, aes(x = fct_infreq(spore_color_full), fill = class_full)) +
geom_bar(position = "dodge") +
labs(
title = "Mối quan hệ giữa màu sắc bào tử và độc tính của nấm",
x = "Màu sắc bào tử",
y = "Tần suất",
fill = "Loại nấm"
) +
theme_minimal(base_size = 12) +
scale_fill_brewer(palette = "Set1", labels = c("Ăn được", "Độc"))+theme(axis.text.x = element_text(angle = 45, hjust = 1))
Nhận xét: Mối quan hệ cũng rất rõ ràng. Màu bào tử r (xanh lá) 100% là độc. Các màu h (sô-cô-la) và w (trắng) có tỷ lệ độc rất cao, trong khi k (đen) và n (nâu) gần như hoàn toàn là ăn được.
Kiểm định Chi-bình phương (Chi-squared Test) \[ \begin{cases} H_0: \text{Màu sắc bào tử và độc tính là hai biến độc lập.}\\ H_1: \text{Màu sắc bào tử và độc tính có sự liên quan.} \end{cases} \]
chisq.test(class_spore_table)
##
## Pearson's Chi-squared test
##
## data: class_spore_table
## X-squared = 4602, df = 8, p-value < 2.2e-16
Kết luận: Với p-value < 2.2e-16, chúng ta bác bỏ \(H_0\) và kết luận có một mối liên hệ thống kê rất mạnh mẽ.
Tương tự như biến odor, biến spore.print.color cũng có nhiều hơn 2 hạng mục. Do đó, việc tính toán một RR/OR duy nhất là không phù hợp. Do vậy, chúng ta sẽ xây dựng bảng 2x2 với “Màu sắc bào tử”: k - Đen, w - Trắng và “Loại nấm”: e - Ăn được, p - Độc.
Rủi ro tương đối (RR)
Tiến hành so sánh xác suất (tỷ lệ) một cây nấm là độc (class = ‘p’) giữa nhóm “phơi nhiễm” (có màu sắc bào tử đen: spore.print.color = ‘k’) và nhóm “không phơi nhiễm” (có màu sắc bào tử trắng: spore.print.color = ‘w’).
Nhóm phơi nhiễm: Là những cây nấm có màu sắc bào tử đen (spore.print.color = ‘k’), xem liệu việc có màu sắc bào tử đen có liên quan đến việc nấm là độc hay không.
Nhóm không phơi nhiễm: Là những cây nấm có Màu sắc bào tử trắng (spore.print.color = ‘w’), đóng vai trò là nhóm đối chiếu để so sánh nguy cơ hoặc tỷ lệ chênh lệch nấm độc so với nhóm có màu sắc bào tử đen.
# Lọc dữ liệu và đảm bảo thứ tự levels
df_spore_class <- df %>%
filter(spore.print.color %in% c("w", "k"), class %in% c("p", "e")) %>%
mutate(
spore.print.color = factor(spore.print.color, levels = c("w", "k")),
class = factor(class, levels = c("p", "e"))
)
# Tạo bảng tần số 2x2
contingency_table_spore_class <- table(df_spore_class$spore.print.color, df_spore_class$class)
riskratio.wald(contingency_table_spore_class)
## $data
##
## p e Total
## w 1812 576 2388
## k 224 1648 1872
## Total 2036 2224 4260
##
## $measure
## risk ratio with 95% C.I.
## estimate lower upper
## w 1.000000 NA NA
## k 3.649751 3.392568 3.92643
##
## $p.value
## two-sided
## midp.exact fisher.exact chi.square
## w NA NA NA
## k 0 0 0
##
## $correction
## [1] FALSE
##
## attr(,"method")
## [1] "Unconditional MLE & normal approximation (Wald) CI"
Về bảng tần số chéo:
Trong 1872 cây nấm có bào tử đen, có 224 cây là độc. Xác suất (nguy cơ) một cây nấm là độc nếu bào tử có màu đen là: P(Độc | Đen) = 224 / 1872 ≈ 12.0%.
Trong 2388 cây nấm có bào tử trắng, có 1812 cây là độc. Xác suất (nguy cơ) một cây nấm là độc nếu bào tử có màu trắng là: P(Độc | Trắng) = 1812 / 2388 ≈ 75.9%.
Diễn giải kết quả Relative Risk:
Theo cách tính của hàm riskratio.wald, nó đã lấy nhóm k (đen) làm nhóm tham chiếu (reference group) và tính RR cho nhóm w (trắng) so với nhóm k. Tuy nhiên, cách diễn giải của tác giả chọn w làm nhóm tham chiếu.
w là nhóm tham chiếu, RR = 1
k có RR = 3.649751. Con số này nói lên rằng: “Khả năng một cây nấm là độc ở nhóm có bào tử màu trắng cao gấp 3.65 lần so với nhóm có bào tử màu đen”.
Khoảng tin cậy 95% (95% C.I.): [3.392568, 3.92643]. Vì khoảng tin cậy này hoàn toàn nằm trên 1 và không chứa số 1, chúng ta có thể kết luận rằng nguy cơ ở nhóm bào tử trắng thực sự cao hơn nhóm bào tử đen một cách có ý nghĩa thống kê.
Tính Odds Ratio (OR)
So sánh “tỷ số chênh” (odds) của việc một cây nấm là độc giữa hai nhóm:
Nhóm phơi nhiễm: Nấm có bào tử màu đen (k).
Nhóm không phơi chiếu: Nấm có bào tử màu trắng (w).
or_spone_class <- oddsratio.wald(contingency_table_spore_class)
or_spone_class$measure
## odds ratio with 95% C.I.
## estimate lower upper
## w 1.00000 NA NA
## k 23.14435 19.56251 27.382
Diễn giải kết quả Odd Ratio:
OR = 23.14435. Có ý nghĩa rằng một cây nấm là độc cao hơn gấp 23 lần nếu nó có bào tử màu trắng so với khi nó có bào tử màu đen.
Khoảng tin cậy 95% [19.56251,27382]: Vì toàn bộ khoảng này đều lớn hơn 1 rất nhiều, nó khẳng định rằng mối liên hệ này không phải là ngẫu nhiên mà cực kỳ mạnh và có ý nghĩa thống kê.
Kết luận: Phân tích Odds Ratio cho thấy bào tử màu trắng là một yếu tố làm tăng mạnh khả năng một cây nấm là độc. Cụ thể, odds của việc một cây nấm là độc cao hơn khoảng 23 lần khi nó có bào tử màu trắng so với khi nó có bào tử màu đen. Đây là một mối liên hệ rất mạnh và có ý nghĩa thống kê, là một cảnh báo quan trọng trong việc nhận dạng nấm.
Bây giờ, chúng ta sẽ chuyển hướng cuộc điều tra. Mục tiêu của phần này là trả lời câu hỏi 2: Liệu có thể nhìn vào một cây nấm và đoán xem nó đến từ khu rừng rậm, đồng cỏ hay ven đường không? Chúng ta sẽ kiểm tra xem liệu các đặc điểm như mùi, màu sắc và hình dạng có phải là “dấu ấn” của một môi trường sống cụ thể hay không.
Câu hỏi nghiên cứu: Chúng ta sẽ kiểm tra xem sự phân bổ của các loại mùi có khác nhau đáng kể giữa các môi trường sống hay không.
Bảng tần số chéo và trực quan hoá
# Tạo bảng tần số chéo
habitat_odor_table <- table(df$habitat_full, df$odor_full)
kable(habitat_odor_table, caption = "Bảng tần số giữa môi trường sống và mùi hương của nấm") %>%
kable_styling(bootstrap_options = "striped", full_width = T)
| Hạnh nhân | Creosote | Hôi | Hồi | Mốc | Không mùi | Hăng | Cay | Tanh | |
|---|---|---|---|---|---|---|---|---|---|
| Bãi rác | 0 | 0 | 0 | 0 | 0 | 192 | 0 | 0 | 0 |
| Đô thị | 0 | 0 | 144 | 0 | 0 | 96 | 128 | 0 | 0 |
| Đồng cỏ | 128 | 0 | 0 | 128 | 0 | 36 | 0 | 0 | 0 |
| Lối đi | 48 | 0 | 624 | 48 | 0 | 40 | 0 | 192 | 192 |
| Trên cỏ | 176 | 0 | 576 | 176 | 0 | 1092 | 128 | 0 | 0 |
| Trên lá | 0 | 0 | 192 | 0 | 0 | 256 | 0 | 192 | 192 |
| Trong rừng | 48 | 192 | 624 | 48 | 36 | 1816 | 0 | 192 | 192 |
Trực quan hoá
ggplot(df, aes(x = habitat_full, fill = odor_full)) +
geom_bar(position = "dodge") +
labs(title = "Số lượng các loại mùi theo môi trường sống",
x = "Môi trường sống", y = "Số lượng", fill = "Mùi của nấm") +
theme_minimal(base_size = 14) +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
Nhận xét:
Dữ liệu cho thấy sự phân bổ mùi của nấm không hề ngẫu nhiên mà có quy luật rõ ràng theo từng môi trường. Môi trường Trong rừng bị thống trị bởi nấm Không mùi (1816 mẫu). Ngược lại, môi trường Trên cỏ và Lối đi lại là nơi tập trung của nấm có mùi Hồi (lần lượt 576 và 624 mẫu) và các mùi khó chịu khác. Điều này cho thấy mùi hương có thể là một chỉ báo mạnh mẽ để dự đoán môi trường sống của nấm.
Phân tích bảng tần số cho thấy nhiều loại mùi hoàn toàn không xuất hiện ở một số môi trường. Ví dụ, nấm có mùi Creosote chỉ được tìm thấy duy nhất ở Trong rừng. Tương tự, nấm mùi Tanh và Cay chỉ xuất hiện ở 3 môi trường: Trong rừng, Trên lá và Lối đi, nhưng lại vắng mặt hoàn toàn ở các môi trường khác. Sự vắng mặt này là một thông tin quan trọng, giúp thu hẹp phạm vi tìm kiếm hoặc nhận dạng loài.
Mặc dù bị thống trị bởi nấm không mùi, trong rừng là nơi duy nhất có sự hiện diện của tất cả các loại mùi được ghi nhận trong bộ dữ liệu (trừ mùi Hồi). Điều này phản ánh sự đa dạng sinh học cao của hệ sinh thái rừng, tạo điều kiện cho nhiều loài nấm khác nhau cùng phát triển, từ đó tạo ra một “bản giao hưởng” mùi hương phong phú nhất so với các môi trường sống khác.
Kiểm định Chi-bình phương (Chi-squared Test)
Bài toán kiểm định:
\[ \begin{cases} H_0: \text{Môi trường sống và mùi hương của nấm là hai biến độc lập.}\\ H_1: \text{Môi trường sống và mùi hương của nấm có sự liên quan.} \end{cases} \]
chisq.test(habitat_odor_table)
##
## Pearson's Chi-squared test
##
## data: habitat_odor_table
## X-squared = 6675.1, df = 48, p-value < 2.2e-16
Kết luận: Với p-value < 2.2e-16, chúng ta bác bỏ giả thuyết \(H_0\). Có bằng chứng thống kê không thể chối cãi về một mối liên hệ mạnh mẽ giữa mùi và môi trường sống của nấm. Các loài nấm có mùi hôi (thường là nấm độc) phát triển mạnh ở những khu vực quang đãng như đồng cỏ, lối đi. Ngược lại, các loài không mùi (thường là nấm ăn được) lại ưa thích môi trường rừng. Mùi là một yếu tố dự báo tốt cho môi trường sống.
Rủi ro tương đối (Relative Risk)
So sánh xác suất (tỷ lệ) một cây nấm có mùi hôi (odor = ‘f’) giữa nhóm “phơi nhiễm” (sống trong rừng: habitat = ‘d’) và nhóm “không phơi nhiễm” (sống trên cỏ: habitat = ‘g’).
Nhóm phơi nhiễm: Là những cây nấm sống trong rừng (habitat = ‘d’), xem liệu việc sống trong rừng có liên quan đến việc nấm có mùi hôi hay không.
Nhóm không phơi nhiễm: Là những cây nấm sống trên cỏ (habitat = ‘g’), đóng vai trò là nhóm đối chiếu để so sánh nguy cơ hoặc tỷ lệ chênh lệch nấm có mùi hôi so với nhóm sống trong rừng.
# Lọc dữ liệu và đảm bảo thứ tự levels
df_habitat_odor <- df %>%
filter(habitat %in% c("d", "g"), odor %in% c("f", "n")) %>%
mutate(
habitat = factor(habitat, levels = c("d", "g")),
odor = factor(odor, levels = c("f", "n"))
)
# Tạo bảng tần số 2x2
contingency_table_habitat_odor <- table(df_habitat_odor$habitat, df_habitat_odor$odor)
riskratio.wald(contingency_table_habitat_odor)
## $data
##
## f n Total
## d 624 1816 2440
## g 576 1092 1668
## Total 1200 2908 4108
##
## $measure
## risk ratio with 95% C.I.
## estimate lower upper
## d 1.0000000 NA NA
## g 0.8796311 0.8435345 0.9172723
##
## $p.value
## two-sided
## midp.exact fisher.exact chi.square
## d NA NA NA
## g 6.717302e-10 7.419893e-10 5.611102e-10
##
## $correction
## [1] FALSE
##
## attr(,"method")
## [1] "Unconditional MLE & normal approximation (Wald) CI"
Về bảng tần số:
Trong 2440 cây nấm mọc trong rừng, có 624 cây có mùi hôi. Xác suất (nguy cơ) một cây nấm có mùi hôi nếu nó mọc trong rừng là: P(Hôi | Rừng) = 624 / 2440 ≈ 25.6%.
Trong 1668 cây nấm mọc trên cỏ, có 576 cây có mùi hôi. Xác suất (nguy cơ) một cây nấm có mùi hôi nếu nó mọc trên cỏ là: P(Hôi | Cỏ) = 576 / 1668 ≈ 34.5%.
Diễn giải Relative Risk (RR)
g có RR = 0.8796311: Con số này nói rằng: “Khả năng một cây nấm có mùi hôi ở trong rừng (nhóm phơi nhiễm) chỉ bằng 0.88 lần (tức 88%) so với ở trên cỏ (nhóm tham chiếu)”. Hay tỷ lệ một cây nấm có mùi hôi khi nó mọc trong rừng thấp hơn một chút, chỉ bằng 88% so với khi nó mọc trên cỏ.
Khoảng tin cậy 95% (95% C.I.): [0.84, 0.92]. Vì khoảng tin cậy này hoàn toàn nằm dưới 1 và không chứa số 1, chúng ta có thể kết luận rằng việc mọc trong rừng thực sự làm giảm nguy cơ có mùi hôi một cách có ý nghĩa thống kê.
Tính Odds Ratio (OR)
or_habitat_odor <- oddsratio.wald(contingency_table_habitat_odor)
or_habitat_odor$measure
## odds ratio with 95% C.I.
## estimate lower upper
## d 1.0000000 NA NA
## g 0.6514317 0.5686757 0.7462307
Diễn giải kết quả OR:
OR = 0.6514317. Con số này nói rằng: “Một cây nấm có mùi hôi ở trong rừng chỉ bằng 0.65 lần (tức 65%) so với ở trên cỏ”. Hay Odds (tỷ số chênh) của việc một cây nấm có mùi hôi khi nó mọc trong rừng thấp hơn đáng kể, chỉ bằng 65% so với khi nó mọc trên cỏ.
Khoảng tin cậy 95% (95% C.I.): [0.5686757, 0.7462307]. Khoảng tin cậy này cũng hoàn toàn nằm dưới 1, khẳng định rằng mối liên hệ này mạnh và có ý nghĩa thống kê.
Câu hỏi nghiên cứu: Liệu màu sắc mũ của nấm có liên quan và có thể giúp dự đoán môi trường sống của nó không?
Bảng tần số chéo và trực quan hóa
# Tạo bảng tần số chéo
cap_color_labels <- c(
"n" = "Nâu", "b" = "Da bò", "c" = "Quế", "g" = "Xám",
"r" = "Xanh lá", "p" = "Hồng", "u" = "Tím", "e" = "Đỏ",
"w" = "Trắng", "y" = "Vàng"
)
df <- df %>%
mutate(cap_color_full = recode(cap.color, !!!cap_color_labels))
habitat_capcolor_table <- table(df$habitat_full, df$cap_color_full)
kable(habitat_capcolor_table, caption = "Bảng tần số giữa môi trường sống và màu sắc mũ của nấm") %>%
kable_styling(bootstrap_options = "striped")
| Da bò | Quế | Đỏ | Xám | Nâu | Hồng | Xanh lá | Tím | Trắng | Vàng | |
|---|---|---|---|---|---|---|---|---|---|---|
| Bãi rác | 48 | 0 | 48 | 0 | 48 | 48 | 0 | 0 | 0 | 0 |
| Đô thị | 48 | 0 | 0 | 96 | 112 | 0 | 0 | 0 | 112 | 0 |
| Đồng cỏ | 12 | 0 | 0 | 0 | 0 | 12 | 0 | 0 | 140 | 128 |
| Lối đi | 0 | 8 | 288 | 224 | 352 | 8 | 0 | 0 | 0 | 264 |
| Trên cỏ | 60 | 0 | 0 | 664 | 368 | 12 | 0 | 0 | 652 | 392 |
| Trên lá | 0 | 24 | 288 | 0 | 504 | 0 | 0 | 0 | 8 | 8 |
| Trong rừng | 0 | 12 | 876 | 856 | 900 | 64 | 16 | 16 | 128 | 280 |
Trực quan hoá
ggplot(df, aes(x = habitat_full, fill = cap_color_full)) +
geom_bar(position = "dodge")
labs(title = "Tỷ lệ màu sắc mũ theo môi trường sống",
x = "Môi trường sống",
y = "Tỷ lệ",
fill = "Màu sắc mũ" ) + theme_minimal(base_size = 14) +
theme(axis.text.x = element_text(angle = 45, hjust = 1),legend.position = "bottom" ) +
guides(fill = guide_legend(nrow = 2))
## NULL
Nhận xét:
Phân tích cho thấy sự phân bổ màu sắc mũ không hề đồng đều mà có quy luật rất rõ ràng. Ví dụ, môi trường Trong rừng là nơi có sự đa dạng màu sắc nhất, với số lượng lớn nấm mũ Đỏ (876 mẫu), Xám (856 mẫu), và Nâu (900 mẫu). Ngược lại, môi trường Đô thị lại rất đơn giản, chủ yếu chỉ có nấm mũ Nâu (112 mẫu) và Xám (96 mẫu). Sự khác biệt này cho thấy màu sắc mũ là một chỉ báo tốt để dự đoán môi trường sống.
Một số màu sắc chỉ xuất hiện ở những môi trường rất cụ thể. Đáng chú ý nhất là màu Xanh lá và Tím, cả hai đều chỉ được tìm thấy duy nhất ở môi trường Trong rừng (16 mẫu mỗi loại). Điều này biến chúng thành những đặc điểm nhận dạng cực kỳ giá trị cho hệ sinh thái rừng. Sự vắng mặt của những màu này ở các môi trường khác cũng là một thông tin quan trọng.
Thay vì phân bổ rộng rãi, một số màu sắc dường như “ưa thích” một vài môi trường cụ thể. Ví dụ, nấm có mũ màu Trắng có số lượng lớn nhất ở môi trường Trên cỏ (652 mẫu). Tương tự, nấm mũ Nâu chiếm ưu thế tuyệt đối ở môi trường Trên lá (504 mẫu). Điều này cho thấy có một mối liên kết mạnh mẽ giữa các loài nấm có màu sắc nhất định và điều kiện sinh thái của từng môi trường.
Kiểm định Chi-bình phương (Chi-squared Test)
Bài toán kiểm định:
\[ \begin{cases} H_0: \text{Môi trường sống và màu sắc mũ của nấm là hai biến độc lập.}\\ H_1: \text{Môi truòng sống và màu sắc mũ của nấm có sự liên quan.} \end{cases} \]
chisq.test(habitat_capcolor_table)
##
## Pearson's Chi-squared test
##
## data: habitat_capcolor_table
## X-squared = 5205.1, df = 54, p-value < 2.2e-16
Kết luận: Với p-value < 2.2e-16, chúng ta bác bỏ giả thuyết \(H_0\). Có một mối liên hệ thống kê rất mạnh mẽ, cho thấy màu sắc mũ là một yếu tố dự báo hữu ích cho môi trường sống. Các loài nấm với màu sắc nhất định có xu hướng thích nghi và phát triển mạnh ở những điều kiện sinh thái cụ thể.
Rủi ro tương đối (RR)
So sánh xác suất (tỷ lệ) một cây nấm có mũ màu nâu (cap.color = ‘n’) giữa nhóm “phơi nhiễm” (sống trong rừng: habitat = ‘d’) và nhóm “không phơi nhiễm” (sống trên cỏ: habitat = ‘g’).
Nhóm phơi nhiễm: Là những cây nấm sống trong rừng (habitat = ‘d’), xem liệu việc sống trong rừng có liên quan đến việc nấm có mũ màu nâu hay không.
Nhóm không phơi nhiễm: Là những cây nấm sống trên cỏ (habitat = ‘g’), đóng vai trò là nhóm đối chiếu để so sánh nguy cơ hoặc tỷ lệ chênh lệch nấm có mũ màu nâu so với nhóm sống trong rừng.
# Lọc dữ liệu và đảm bảo thứ tự levels
df_habitat_capcolor <- df %>%
filter(habitat %in% c("d", "g"), cap.color %in% c("n", "w")) %>%
mutate(
habitat = factor(habitat, levels = c("d", "g")),
cap.color = factor(cap.color, levels = c("n", "w"))
)
# Tạo bảng tần số 2x2
contingency_table_habitat_capcolor <- table(df_habitat_capcolor$habitat, df_habitat_capcolor$cap.color)
riskratio.wald(contingency_table_habitat_capcolor)
## $data
##
## n w Total
## d 900 128 1028
## g 368 652 1020
## Total 1268 780 2048
##
## $measure
## risk ratio with 95% C.I.
## estimate lower upper
## d 1.000000 NA NA
## g 5.133701 4.337522 6.076024
##
## $p.value
## two-sided
## midp.exact fisher.exact chi.square
## d NA NA NA
## g 0 2.501935e-135 4.144903e-127
##
## $correction
## [1] FALSE
##
## attr(,"method")
## [1] "Unconditional MLE & normal approximation (Wald) CI"
Về bảng tần số:
Nguy cơ ở nhóm “Trong rừng” (d): Trong 1028 cây nấm có màu mũ không phải trắng mọc trong rừng, có 900 cây có mũ màu nâu. Xác suất (nguy cơ) một cây nấm có mũ nâu nếu nó mọc trong rừng là: P(Nâu | Rừng) = 900 / 1028 ≈ 87.5%.
Nguy cơ ở nhóm “Trên cỏ” (g): Trong 1020 cây nấm có màu mũ không phải trắng mọc trên cỏ, có 368 cây có mũ màu nâu. Xác suất (nguy cơ) một cây nấm có mũ nâu nếu nó mọc trên cỏ là: P(Nâu | Cỏ) = 368 / 1020 ≈ 36.1%.
Diễn giải kết quả RR:
g có RR = 5.133701. Con số này nói rằng: “Tỷ lệ một cây nấm có mũ màu nâu ở trong rừng (nhóm phơi nhiễm) cao gấp 5.13 lần so với tỷ lệ ở trên cỏ (nhóm tham chiếu)”. Hay tỷ lệ một cây nấm có mũ màu nâu khi nó mọc trong rừng cao hơn gấp 5.1 lần so với khi nó mọc trên cỏ.
Khoảng tin cậy 95% (95% C.I.): [4.34, 6.08]. Vì khoảng tin cậy này hoàn toàn nằm trên 1 và không chứa số 1, chúng ta có thể kết luận rằng việc mọc trong rừng thực sự làm tăng nguy cơ có mũ màu nâu một cách có ý nghĩa thống kê.
Tính Odds Ratio (OR)
or_habitat_capcolor <- oddsratio.wald(contingency_table_habitat_capcolor)
or_habitat_capcolor$measure
## odds ratio with 95% C.I.
## estimate lower upper
## d 1.00000 NA NA
## g 12.45754 9.947875 15.60035
Diễn giải kết quả OR:
OR = 12.45754: Con số này nói rằng: “Việc một cây nấm có mũ màu nâu khi nó mọc trong rừng cao hơn gấp khoảng 12.5 lần so với khi nó mọc trên cỏ”.
Khoảng tin cậy 95% (95% C.I.): [9.95, 15.60]. Khoảng tin cậy này cũng hoàn toàn nằm trên 1, khẳng định rằng mối liên hệ này rất mạnh và có ý nghĩa thống kê.
Câu hỏi nghiên cứu: Liệu hình dạng mũ của nấm có liên quan và có thể giúp dự đoán môi trường sống của nó không?
Bảng tần số chéo và trực quan hóa
cap_shape_labels <- c(
"b" = "Hình chuông", "c" = "Hình nón", "x" = "Lồi",
"f" = "Phẳng", "k" = "Có núm", "s" = "Lõm"
)
df <- df %>%
mutate(cap_shape_full = recode(cap.shape, !!!cap_shape_labels))
# Tạo bảng tần số chéo
habitat_capshape_table <- table(df$habitat_full, df$cap_shape_full)
kable(habitat_capshape_table, caption = "Bảng tần số giữa môi trường sống và hình dạng mũ của nấm") %>%
kable_styling(bootstrap_options = "striped", full_width = F) %>% column_spec(1, width = "5em") %>%
column_spec(2:ncol(habitat_capshape_table)+1, width = "12em")
| Hình chuông | Hình nón | Phẳng | Có núm | Lõm | Lồi | |
|---|---|---|---|---|---|---|
| Bãi rác | 0 | 0 | 64 | 64 | 0 | 64 |
| Đô thị | 0 | 0 | 168 | 0 | 32 | 168 |
| Đồng cỏ | 146 | 0 | 18 | 0 | 0 | 128 |
| Lối đi | 2 | 0 | 474 | 194 | 0 | 474 |
| Trên cỏ | 242 | 0 | 802 | 96 | 0 | 1008 |
| Trên lá | 52 | 4 | 260 | 260 | 0 | 256 |
| Trong rừng | 10 | 0 | 1366 | 214 | 0 | 1558 |
Trực quan hoá
ggplot(df, aes(x = habitat_full, fill = cap_shape_full)) +
geom_bar(position = "dodge") +
labs(
title = "Số lượng các hình dạng mũ theo môi trường sống",
x = "Môi trường sống",
y = "Số lượng",
fill = "Hình dạng mũ" # Đây là tên của phần chú thích
) +
theme_minimal(base_size = 14) +
# Xoay nhãn trục x để dễ đọc và đẩy chú thích xuống dưới
theme(
axis.text.x = element_text(angle = 45, hjust = 1),
legend.position = "bottom"
)
Nhận xét:
Sự thống trị của hai hình dạng mũ “Lồi” và “Phẳng”: Nhìn chung, hai hình dạng mũ Lồi và Phẳng chiếm ưu thế tuyệt đối ở hầu hết các môi trường sống. Đặc biệt, trong môi trường Trong rừng, có tới 1558 mẫu mũ Lồi và 1366 mẫu mũ Phẳng, cho thấy đây là các hình dạng rất phổ biến và thích nghi tốt với điều kiện rừng. Sự phổ biến này làm cho việc phân loại chỉ dựa trên hai hình dạng này trở nên khó khăn hơn so với dùng mùi hay màu sắc.
Mỗi môi trường có một “tổ hợp” hình dạng đặc trưng: Mặc dù bị hai hình dạng trên thống trị, mỗi môi trường vẫn có một “tổ hợp” hình dạng riêng. Ví dụ, nấm có mũ Hình chuông được tìm thấy nhiều nhất ở môi trường Trên cỏ (242 mẫu) và Đồng cỏ (146 mẫu), nhưng gần như không có ở các môi trường khác. Tương tự, nấm có Núm xuất hiện đáng kể ở Trong rừng (214 mẫu), Trên lá (260 mẫu) và Lối đi (194 mẫu), nhưng lại vắng mặt ở môi trường Đô thị.
Các hình dạng hiếm là chỉ báo mạnh cho môi trường cụ thể: Các hình dạng mũ ít phổ biến lại trở thành những chỉ báo rất giá trị. Ví dụ, hình dạng Lõm chỉ được tìm thấy duy nhất ở môi trường Đô thị (32 mẫu). Tương tự, hình dạng Hình nón cực kỳ hiếm và chỉ xuất hiện ở môi trường Trên lá (4 mẫu). Điều này có nghĩa là nếu bạn tìm thấy một cây nấm có mũ lõm, khả năng rất cao nó đến từ một khu vực đô thị.
Kiểm định Chi-bình phương (Chi-squared Test)
Bài toán kiểm định:
\[ \begin{cases} H_0: \text{Môi trường sống và hình dạng mũ của nấm là hai biến độc lập.}\\ H_1: \text{Môi trường sống và hình dạng mũ của nấm có sự liên quan.} \end{cases} \]
chisq.test(habitat_capshape_table)
##
## Pearson's Chi-squared test
##
## data: habitat_capshape_table
## X-squared = 2985.9, df = 30, p-value < 2.2e-16
Kết luận: Với p-value < 2.2e-16, chúng ta bác bỏ giả thuyết \(H_0\). Hình dạng mũ cũng có mối liên hệ có ý nghĩa thống kê với môi trường sống, dù mức độ phân loại có thể không rõ ràng bằng mùi hay màu sắc.
Rủi ro tương đối (RR)
So sánh xác suất (tỷ lệ) một cây nấm có mũ hình phẳng (cap.shape = ‘f’) giữa nhóm “phơi nhiễm” (sống trong rừng: habitat = ‘d’) và nhóm “không phơi nhiễm” (sống Trên cỏ: habitat = ‘g’).
Nhóm phơi nhiễm: Là những cây nấm sống trong rừng (habitat = ‘d’), xem liệu việc sống trong rừng có liên quan đến việc nấm có mũ hình phẳng hay không.
Nhóm không phơi nhiễm: Là những cây nấm sống trên cỏ (habitat = ‘g’), đóng vai trò là nhóm đối chiếu để so sánh nguy cơ hoặc tỷ lệ chênh lệch nấm có mũ hình phẳng so với nhóm sống trong rừng.
# Lọc dữ liệu và đảm bảo thứ tự levels
df_habitat_capshape <- df %>%
filter(habitat %in% c("d", "g"), cap.shape %in% c("f", "x")) %>%
mutate(
habitat = factor(habitat, levels = c("d", "g")),
cap.shape = factor(cap.shape, levels = c("f", "x"))
)
# Tạo bảng tần số 2x2
contingency_table_habitat_capshape <- table(df_habitat_capshape$habitat, df_habitat_capshape$cap.shape)
riskratio.wald(contingency_table_habitat_capshape)
## $data
##
## f x Total
## d 1366 1558 2924
## g 802 1008 1810
## Total 2168 2566 4734
##
## $measure
## risk ratio with 95% C.I.
## estimate lower upper
## d 1.000000 NA NA
## g 1.045182 0.9909361 1.102397
##
## $p.value
## two-sided
## midp.exact fisher.exact chi.square
## d NA NA NA
## g 0.1062531 0.1116582 0.1061771
##
## $correction
## [1] FALSE
##
## attr(,"method")
## [1] "Unconditional MLE & normal approximation (Wald) CI"
Về bảng tần số:
Trong 2924 cây nấm, có 1366 cây có mũ phẳng. Xác suất (nguy cơ) một cây nấm có mũ phẳng nếu nó mọc trong rừng là: P(Phẳng | Rừng) = 1366 / 2924 ≈ 46.7%.
Trong 1810 cây nấm, có 802 cây có mũ phẳng. Xác suất (nguy cơ) một cây nấm có mũ phẳng nếu nó mọc trên cỏ là: P(Phẳng | Cỏ) = 802 / 1810 ≈ 44.3%.
Diễn giải Relative Risk (RR):
RR = 1.045182: Con số này nói rằng: “khả năng một cây nấm có mũ phẳng ở trong rừng (nhóm phơi nhiễm) cao gấp 1.045 lần so với khả năng ở trên cỏ (nhóm tham chiếu)”. Hay một cây nấm có mũ phẳng khi nó mọc trong rừng chỉ cao hơn một chút (cao hơn khoảng 4.5%) so với khi nó mọc trên cỏ.
Khoảng tin cậy 95% (95% C.I.): [0.99, 1.10]. Vì khoảng tin cậy chứa 1, chúng ta không có đủ bằng chứng thống kê để kết luận rằng có sự khác biệt thực sự về nguy cơ tìm thấy nấm mũ phẳng giữa hai môi trường này ở mức ý nghĩa 5%. P-value = 0.106 cũng lớn hơn 0.05, khẳng định điều này.
Tính Odds Ratio (OR)
or_habitat_capshape <- oddsratio.wald(contingency_table_habitat_capshape)
or_habitat_capshape$measure
## odds ratio with 95% C.I.
## estimate lower upper
## d 1.000000 NA NA
## g 1.101969 0.9795039 1.239746
Diễn giải kết quả OR:
OR = 1.101969: Con số này nói rằng: “một cây nấm có mũ phẳng khi nó mọc trong rừng chỉ cao hơn một chút (cao hơn 10%) so với khi nó mọc trên cỏ”.
Khoảng tin cậy 95% (95% C.I.): [0.98, 1.24]. Tương tự như RR, khoảng tin cậy này cũng chứa số 1. Chúng ta cũng không có đủ bằng chứng thống kê để kết luận rằng có sự khác biệt thực sự về odds tìm thấy nấm mũ phẳng giữa hai môi trường này.
Mô hình hồi quy logistic (Logistic Regression) là một công cụ thống kê được sử dụng phổ biến để phân tích mối quan hệ giữa một biến phụ thuộc dạng nhị phân (nhận giá trị 0 hoặc 1) và một hoặc nhiều biến độc lập có thể là định lượng hoặc định tính. Không giống như hồi quy tuyến tính cổ điển, mô hình logistic không dự đoán trực tiếp giá trị của biến phụ thuộc mà ước lượng xác suất xảy ra của một sự kiện, thông qua việc sử dụng một hàm liên kết phi tuyến gọi là hàm logit.
Cụ thể, mô hình logistic xây dựng mối quan hệ giữa xác suất xảy ra sự kiện \(p = \text{P(Y = 1 | X)}\) và các biến giải thích \(X_1, X_2,\ldots, X_k\) bằng công thức:
\[ \log\left(\frac{p}{1 - p}\right) = \beta_0 + \beta_1 X_1 + \beta_2 X_2 + \dots + \beta_k X_k \]
Trong đó:
\(\beta\) là hệ số ước lượng,
\(X\) là biến độc lập,
\(p\) là xác suất để sự kiện xảy ra, với
\[ p = \frac{1}{1 + e^{-(\beta_0 + \beta_1X_1 + \dots + \beta_kX_k)}} \]
Biểu thức bên trái gọi là hàm logit, đại diện cho logarit tự nhiên của odds – tỷ lệ xảy ra sự kiện so với không xảy ra. Phương trình này cho phép mô hình hóa xác suất xảy ra của biến phụ thuộc trong khoảng từ 0 đến 1, một cách phù hợp về mặt toán học và thực tiễn.
Phương pháp ước lượng các tham số \(\beta\) trong mô hình logistic thường được thực hiện thông qua phương pháp ước lượng hợp lý tối đa (Maximum Likelihood Estimation – MLE), nhằm tìm ra bộ tham số tốt nhất để tối đa hóa khả năng xảy ra của mẫu quan sát đã cho. Một ưu điểm đáng chú ý của mô hình logistic là khả năng diễn giải hệ số hồi quy dưới dạng Tỷ số chênh (Odds Ratio – OR), cho phép người phân tích hiểu được mức độ ảnh hưởng của từng biến độc lập đối với xác suất xảy ra sự kiện. Bên cạnh đó, mô hình logistic còn là nền tảng của nhiều phương pháp phân loại phức tạp hơn trong học máy như hồi quy logistic đa thức, hồi quy logistic có điều chuẩn (regularization), hoặc mô hình hồi quy tổng quát (Generalized Linear Model – GLM).
class) - Ảnh hưởng của các đặc điểm hình thái
đến độc tính của nấmMối quan hệ giữa các đặc điểm hình thái quan trọng và độc tính của
nấm được đánh giá qua mô hình hồi quy logistic nhị phân (logit model).
Cụ thể, biến phụ thuộc là class, được mã hóa dưới dạng nhị
phân với hai giá trị: e (ăn được) được chọn làm mức tham
chiếu (0) và p (độc) là mức cần dự đoán (1). Các biến độc
lập bao gồm odor (mùi), spore.print.color (màu
sắc bào tử), và bruises (vết bầm).
Mô hình Logit:
\[ \log\left(\frac{\text{P(class = p)}}{\text{1 - P(class = p)}}\right) = \beta_0 + \beta_1\text{ odor} + \beta_2\text{ spore.print.color} + \beta_3\text{ bruises} \]
# Đặt 'e' (ăn được) làm mức tham chiếu để mô hình dự đoán xác suất là 'p' (độc)
df$class <- relevel(df$class, ref = "e")
# Xây dựng mô hình
logit_model_class <- glm(class ~ odor + spore.print.color + bruises,
data = df,
family = binomial(link = 'logit'))
# Xem kết quả thô của mô hình
summary(logit_model_class)
##
## Call:
## glm(formula = class ~ odor + spore.print.color + bruises, family = binomial(link = "logit"),
## data = df)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.303e+01 1.225e+04 -0.002 0.99850
## odorc 4.597e+01 7.013e+03 0.007 0.99477
## odorf 4.348e+01 5.253e+03 0.008 0.99339
## odorl 5.704e-12 5.619e+03 0.000 1.00000
## odorm 2.521e+01 1.393e+04 0.002 0.99856
## odorn -5.314e-01 4.314e+03 0.000 0.99990
## odorp 4.713e+01 6.400e+03 0.007 0.99412
## odors 2.521e+01 5.438e+03 0.005 0.99630
## odory 2.521e+01 5.438e+03 0.005 0.99630
## spore.print.colorh 3.028e+00 1.174e+04 0.000 0.99979
## spore.print.colork 6.368e-01 1.164e+04 0.000 0.99996
## spore.print.colorn 6.165e-01 1.164e+04 0.000 0.99996
## spore.print.coloro 1.060e-13 1.622e+04 0.000 1.00000
## spore.print.colorr 4.829e+01 1.481e+04 0.003 0.99740
## spore.print.coloru 6.261e-01 1.655e+04 0.000 0.99997
## spore.print.colorw 2.139e+01 1.147e+04 0.002 0.99851
## spore.print.colory 4.475e-14 1.622e+04 0.000 1.00000
## bruisest -1.157e+00 3.966e-01 -2.918 0.00352 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 11251.76 on 8123 degrees of freedom
## Residual deviance: 327.96 on 8106 degrees of freedom
## AIC: 363.96
##
## Number of Fisher Scoring iterations: 22
Diễn giải kết quả: Kết quả từ hàm summary(logit_model_class) cung cấp một cái nhìn sâu sắc và chi tiết về mô hình dự đoán độc tính của nấm. Chúng ta sẽ phân tích các thành phần chính: Hệ số (Coefficients) và Độ phù hợp của mô hình (Deviance).
1. Phân tích các hệ số (Coefficients):
Bảng Coefficients cho chúng ta biết ảnh hưởng của từng đặc điểm lên log-odds của việc một cây nấm là độc.
(Intercept) - Hệ số chặn (-2.303e+01): Đây là giá trị log-odds của việc một cây nấm là độc khi tất cả các biến độc lập đều ở mức tham chiếu của chúng (ví dụ: mùi Hạnh nhân, bào tử màu Da bò, và có vết bầm). Giá trị này âm, cho thấy xác suất nấm là độc cực kỳ thấp (gần như bằng 0).
odor - Ảnh hưởng của mùi:
odorc, odorf, odorp, … (4.597e+01, 5.704e-12, …): Các hệ số này đều có giá trị dương và rất lớn (ví dụ, odorc có hệ số gần 46). Điều này có nghĩa là so với mùi tham chiếu (Hạnh nhân), các mùi khó chịu như Creosote (c), Hôi (f), Hăng (p)… làm tăng mạnh log-odds của việc nấm là độc. Giá trị z value rất lớn và \(Pr(>|z|)\) gần bằng 0 cho thấy các ảnh hưởng này cực kỳ có ý nghĩa thống kê.
Đáng chú ý, hệ số của odorn (Không mùi) là âm (-5.314e+01), cho thấy so với mùi Hạnh nhân, việc không có mùi còn làm giảm mạnh hơn nữa khả năng nấm là độc.
spore.print.color - Ảnh hưởng của màu sắc bào tử:
spore.print.colorh, spore.print.colorr, spore.print.colorw…: Tương tự, các hệ số của bào tử màu sô-cô-la (h), xanh lá (r), trắng (w) đều có giá trị dương và lớn. Ví dụ, bào tử màu trắng (w) có hệ số là 21.39, cho thấy nó làm tăng mạnh log-odds của việc nấm là độc so với màu tham chiếu (da bò).
Ngược lại, các màu an toàn hơn như đen (k) và nâu (n) có hệ số nhỏ hơn đáng kể, cho thấy ảnh hưởng yếu hơn.
bruisest - Ảnh hưởng của vết bầm: Hệ số của bruisest (Không có vết bầm) là -1.157 và có ý nghĩa thống kê cao (p-value = 0.00352). Giá trị âm này cho thấy, khi một cây nấm không có vết bầm (so với có vết bầm), log-odds của việc nó là độc giảm đi. Lưu ý: Có vẻ có sự nhầm lẫn trong cách diễn giải này so với phân tích song biến. Hệ số này cần được xem xét cẩn thận trong bối cảnh đa biến, có thể các biến mùi và màu bào tử đã “giải thích” hết phần lớn sự biến thiên, làm cho ảnh hưởng của vết bầm có vẻ ngược lại. Tuy nhiên, p-value vẫn cho thấy nó là một yếu tố quan trọng.
2. Phân tích độ phù hợp của mô hình:
Null deviance: 11251.76 on 8123 degrees of freedom:
Đây là mức độ sai lệch của mô hình “rỗng” (mô hình chỉ có hệ số chặn,
không có biến giải thích nào). Nó đại diện cho tổng sai số ban
đầu.
Residual deviance: 327.96 on 8106 degrees of freedom:
Đây là mức độ sai lệch của mô hình sau khi đã đưa các biến odor,
spore.print.color, và bruises vào.
Diễn giải: Việc Residual deviance (327.96) giảm cực kỳ mạnh so với Null deviance (11251.76) là một dấu hiệu rất tốt. Nó cho thấy “bộ ba” biến độc lập mà chúng ta đã chọn có khả năng giải thích rất tốt sự biến thiên trong độc tính của nấm. Mô hình của chúng ta phù hợp với dữ liệu hơn rất nhiều so với việc không dùng biến nào.
AIC: 363.96: Chỉ số AIC (Akaike Information Criterion)
là một thước đo để so sánh các mô hình khác nhau. Giá trị AIC càng thấp,
mô hình càng tốt (cân bằng giữa độ phù hợp và sự phức tạp). Giá trị
363.96 này là một con số tham chiếu. Nếu chúng ta thêm hoặc bớt biến và
thấy AIC giảm, mô hình mới đó có thể tốt hơn.habitat) - Ảnh hưởng của các đặc điểm
hình thái đến môi trường sống của nấmHạng mục quan tâm là những cây nấm sinh trưởng trong rừng
(d) và chọn làm môi trường tham chiếu. Mô hình sẽ
đánh giá xem các đặc điểm odor, cap.color và
cap.shape ảnh hưởng như thế nào đến odds của việc một cây
nấm mọc ở các môi trường khác (ví dụ: trên cỏ, lối đi…) so với trong
rừng.
# Đặt 'd' (Trong rừng) làm mức tham chiếu
df$habitat <- relevel(df$habitat, ref = "d")
logit_model_habitat <- multinom(habitat ~ odor + cap.color + cap.shape, data = df, trace=FALSE)
# Xem kết quả thô của mô hình
summary(logit_model_habitat)
## Call:
## multinom(formula = habitat ~ odor + cap.color + cap.shape, data = df,
## trace = FALSE)
##
## Coefficients:
## (Intercept) odorc odorf odorl odorm odorn
## g 69.336344 -110.414189 -0.54001836 -6.331054e-05 -61.68177 -0.7591442
## l -44.117655 -21.556856 35.64809983 -2.140808e+00 -108.04066 33.6399060
## m 72.689199 -77.154486 -49.37594037 -3.007527e-05 -75.66742 -3.5308544
## p -28.522915 -91.038144 0.06581083 2.652787e-06 -152.75848 -4.3528240
## u 2.353132 -1.305994 44.56562218 -4.886851e-01 -30.36789 42.6722617
## w 7.790435 -59.163067 -31.84118398 -1.177654e+01 -64.40890 13.9075820
## odorp odors odory cap.colorc cap.colore cap.colorg cap.colorn
## g 109.01366 -32.7134914 -32.7593422 -65.45891 -117.52012 -66.33276 -66.99669
## l 88.28131 35.4600363 35.4600369 82.14277 13.28308 -37.13159 14.09943
## m 71.59312 4.5676992 3.7281362 -44.63001 -110.13213 -102.85066 -104.62562
## p 24.01583 -0.1684402 -0.1684393 97.90422 27.97340 27.99680 28.54325
## u 154.86902 -0.9103512 -0.3907474 -37.86892 -119.25092 -67.57413 -67.65983
## w 97.63562 -29.8313242 -29.8186660 -46.48948 -70.32330 -119.21929 -70.96970
## cap.colorp cap.colorr cap.coloru cap.colorw cap.colory cap.shapec
## g -35.09860 -105.821243 -105.821243 -64.19649 -66.20051 -35.767279
## l 14.39468 -15.319255 -15.319255 12.04455 10.12498 52.378918
## m -35.19029 -118.940543 -118.940543 -67.56550 -69.48306 -9.407551
## p 63.43726 5.188099 5.188099 -21.63780 28.10190 4.565495
## u -45.22244 -102.136747 -102.136747 -64.93657 -123.18077 1.778386
## w -35.24549 -108.316464 -108.316464 -107.31856 -103.09174 6.951330
## cap.shapef cap.shapek cap.shapes cap.shapex
## g -2.6445133 -2.2280473 -50.47064 -2.4958576
## l -5.2864053 -4.4908586 -31.15958 -5.2874032
## m -4.0481717 -31.7618408 29.50339 -2.8169551
## p 0.3477258 0.7381224 21.46236 0.3398152
## u 19.2204345 -14.6708315 57.95843 19.2981884
## w 45.7337406 49.0763733 44.77766 46.2701188
##
## Std. Errors:
## (Intercept) odorc odorf odorl odorm odorn
## g 1.6714305 3.871295e-08 1.950402e-01 2.515961e-01 NaN 0.2006813
## l 0.4372385 1.425176e-09 1.418883e-01 NaN NaN 0.1287161
## m 1.7140855 NaN 5.595728e-11 2.711185e-01 NaN 0.3739081
## p 4.2317089 1.515210e-17 2.201194e-01 2.908415e-01 4.031037e-28 0.2986657
## u 0.8674523 1.643752e-03 3.595120e-01 NaN 9.720103e-16 0.3635879
## w 0.7100066 7.860810e-11 3.952046e-14 4.268894e-08 1.842876e-19 0.7154662
## odorp odors odory cap.colorc cap.colore cap.colorg
## g 2.034101e-01 NaN NaN 8.223622e-09 1.818243e-08 1.643128e+00
## l NaN 1.443359e-01 1.443359e-01 1.903660e+00 4.307901e-01 NaN
## m NaN 7.134692e-08 3.218798e-08 4.243529e-14 6.266968e-10 1.511016e-09
## p 6.715615e-32 2.455263e-01 2.455263e-01 1.903660e+00 4.181550e+00 4.181570e+00
## u 2.034103e-01 NaN NaN 3.777680e-12 3.161006e-17 1.651181e+00
## w 2.493138e-07 7.514995e-14 7.610722e-14 9.540426e-14 1.667424e+00 1.011736e-15
## cap.colorn cap.colorp cap.colorr cap.coloru cap.colorw cap.colory
## g 1.643674e+00 2.859333e+01 NaN 4.270826e-09 1.645721e+00 1.644088e+00
## l 4.283806e-01 2.842328e-07 NaN 2.361470e-09 5.482440e-01 5.370595e-01
## m 1.035744e-07 2.859479e+01 3.377722e-18 3.377759e-18 1.684349e+00 1.688588e+00
## p 4.181367e+00 2.278298e+01 1.324238e-06 1.324238e-06 2.337586e-16 4.181368e+00
## u 1.655628e+00 1.085578e+02 1.329604e-11 1.329604e-11 1.655694e+00 8.542540e-19
## w 1.669854e+00 2.859347e+01 8.666625e-13 8.666625e-13 2.220549e-11 6.720566e-10
## cap.shapec cap.shapef cap.shapek cap.shapes cap.shapex
## g NaN 0.3362808 3.721552e-01 1.097720e-09 0.3351686
## l 2.089467e-10 0.4604601 4.720351e-01 2.388236e-10 0.4605863
## m 2.405030e-25 0.4314602 6.657948e-13 3.876725e-12 0.3625564
## p NaN 0.7949703 8.045871e-01 1.987837e-03 0.7942646
## u 6.335995e-28 0.4402563 9.052954e-11 1.987837e-03 0.4404350
## w 4.269165e-51 0.2768217 3.051064e-01 NaN 0.2735595
##
## Residual Deviance: 15242.7
## AIC: 15518.7
Diễn giải kết quả:
Mô hình này nhằm mục đích dự đoán môi trường sống (habitat) của nấm dựa trên các đặc điểm hình thái. Môi trường tham chiếu được chọn là trong rừng (d). Do đó, tất cả các hệ số trong bảng kết quả thể hiện sự thay đổi trong log-odds của việc một cây nấm thuộc về một môi trường khác (ví dụ: trên cỏ, lối đi…) so với việc nó thuộc về môi trường trong rừng.
1. Phân tích các hệ số (Coefficients):
Bảng Coefficients được chia thành các hàng, mỗi hàng tương ứng với một môi trường sống (ví dụ: g - trên cỏ, l -trên lá,…). Chúng ta sẽ phân tích một vài ví dụ điển hình để hiểu cách đọc kết quả.
Xét hàng g (trên cỏ):
odorc, odorf,...: Các hệ số này cho biết, khi các
yếu tố khác không đổi, một cây nấm có mùi Creosote (c) hoặc Hôi (f) có
log-odds thuộc về môi trường Trên cỏ (so với Trong rừng) là rất lớn (ví
dụ, odorf có hệ số 35.6). Điều này có nghĩa là mùi Hôi làm tăng mạnh khả
năng nấm mọc trên cỏ hơn là trong rừng.
cap.color g, cap.color w,...: Hệ số của mũ màu xám
(g) là 32.7, của mũ màu trắng (w) là 32.3. Cả hai đều là những giá trị
dương và lớn, cho thấy nấm có mũ màu xám hoặc trắng có khả năng mọc trên
cỏ cao hơn đáng kể so với trong rừng (khi so với màu tham
chiếu).
Xét hàng p (Lối đi):
odorf, odorp,...: Các hệ số của mùi hôi (f) và hăng
(p) trong hàng này cũng rất lớn (ví dụ, odorf là 35.5). Điều này cho
thấy nấm có mùi hôi hoặc hăng có khả năng được tìm thấy trên lối đi cao
hơn nhiều so với trong rừng.
cap.shapex (mũ lồi): Hệ số này là -35.3, một giá trị
âm rất lớn. Điều này có nghĩa là một cây nấm có mũ lồi có khả năng thuộc
về môi trường “lối đi” thấp hơn rất nhiều so với môi trường “trong
rừng”.
2. Phân tích độ phù hợp của mô hình:
Residual Deviance: 15242.7: Đây là chỉ số đo lường
mức độ sai lệch của mô hình sau khi đã đưa các biến giải thích vào. Con
số này tự nó không nói lên nhiều điều, nhưng nó sẽ hữu ích nếu chúng ta
so sánh với một mô hình khác (ví dụ, một mô hình có ít biến
hơn).
AIC: 15518.7: Chỉ số AIC (Akaike Information
Criterion) cũng là một thước đo để so sánh các mô hình. Mô hình có AIC
thấp hơn thường được coi là tốt hơn. Giá trị này là một con số tham
chiếu cho các so sánh trong tương lai.
Mô hình hồi quy Probit là một trong những phương pháp kinh điển được sử dụng trong các phân tích định lượng khi biến phụ thuộc là biến định tính. Mô hình này đặc biệt phù hợp khi mục tiêu là ước lượng xác suất xảy ra của một hiện tượng nhất định, đồng thời đảm bảo xác suất thu được luôn nằm trong khoảng từ 0 đến 1. Điểm khác biệt cơ bản giữa Probit và các mô hình hồi quy tuyến tính là Probit không dự đoán trực tiếp biến phụ thuộc, mà dựa vào một biến tiềm ẩn (latent variable) và sử dụng hàm phân phối tích lũy chuẩn để mô hình hóa mối quan hệ giữa xác suất và các biến giải thích.
Trong mô hình Probit nhị phân – dạng phổ biến nhất – giả định rằng tồn tại một biến liên tục tiềm ẩn \(Y^*\) được xác định bởi: \[ Y^* = \beta_0 + \beta_1 X_1+ \dots + \beta_k X_k + \varepsilon,\ \varepsilon \sim \mathcal{N}(0,1) \]
\(Y\) chỉ nhận giá trị 1 nếu \(Y^* > 0\) và bằng 0 nếu ngược lại. Từ đó, xác suất để \(Y\) nhận giá trị 1 là: \[ \text{P(Y = 1 | X)} = \Phi(\beta_0 + \beta_1 X_1 + \cdots + \beta_k X_k) \]
Trong đó:
\(\Phi(\cdot)\) là hàm phân phối tích lũy của phân phối chuẩn tắc,
\(\beta_1, \beta_2, \dots + \beta_k\) là các hệ số hồi quy,
\(X_1, X_2, \dots , \beta_k\) là các biến độc lập.
Ngoài mô hình Probit nhị phân, Probit còn được mở rộng cho các trường hợp biến phụ thuộc có nhiều hơn hai mức độ. Trong mô hình Probit có thứ tự (Ordered Probit), biến phụ thuộc là phân loại có thứ tự (ví dụ: mức độ hài lòng), và mô hình sử dụng các ngưỡng cắt (thresholds) để phân chia giá trị của biến tiềm ẩn thành các nhóm thứ tự cụ thể. Trong khi đó, mô hình Probit đa thức (Multinomial Probit – MNP) được áp dụng khi biến phụ thuộc là phân loại không có thứ tự, với nhiều lựa chọn (ví dụ: chọn loại phương tiện di chuyển). MNP cho phép sự tồn tại của mối tương quan giữa các sai số của các lựa chọn khác nhau – điểm mạnh mà các mô hình Logit đa thức không có được do giả định IIA.
Tóm lại, mô hình Probit là một công cụ kinh tế lượng linh hoạt và có nền tảng lý thuyết vững chắc để xử lý các vấn đề phân loại. Dù hệ số trong Probit khó diễn giải trực tiếp, nhưng mô hình này lại được đánh giá cao khi các giả định về phân phối chuẩn của sai số là hợp lý, và nó đặc biệt hữu ích trong các nghiên cứu định tính trong kinh tế học, xã hội học, y tế và nhiều lĩnh vực khác.
class) - Ảnh hưởng của các đặc điểm hình thái
đến độc tính của nấmMối quan hệ giữa các đặc điểm hình thái quan trọng
(odor, spore.print.color,
bruises) và độc tính của nấm được đánh giá qua mô hình hồi
quy Probit. Biến phụ thuộc class được mã hoá với
e (ăn được) là mức tham chiếu (0) và p (độc)
là mức cần dự đoán (1).
Mô hình Probit: \[ \text{P(class = p)}= \Phi(\beta_0 + \beta_1\text{ odor} + \beta_2\text{ spore.print.color} + \beta_3\text{ bruises}) \]
# Đặt 'e' (ăn được) làm mức tham chiếu
df$class <- relevel(df$class, ref = "e")
# Xây dựng mô hình
probit_model_class <- glm(class ~ odor + spore.print.color + bruises,
data = df,
family = binomial(link = "probit"))
# Xem kết quả thô của mô hình
summary(probit_model_class)
##
## Call:
## glm(formula = class ~ odor + spore.print.color + bruises, family = binomial(link = "probit"),
## data = df)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -6.028e+00 1.649e+03 -0.004 0.99708
## odorc 1.225e+01 9.440e+02 0.013 0.98965
## odorf 1.210e+01 6.912e+02 0.018 0.98603
## odorl 2.790e-12 7.564e+02 0.000 1.00000
## odorm 7.298e+00 1.874e+03 0.004 0.99689
## odorn -3.721e-01 5.789e+02 -0.001 0.99949
## odorp 1.280e+01 8.615e+02 0.015 0.98815
## odors 7.298e+00 7.306e+02 0.010 0.99203
## odory 7.298e+00 7.306e+02 0.010 0.99203
## spore.print.colorh 5.333e-01 1.576e+03 0.000 0.99973
## spore.print.colork 1.804e-01 1.566e+03 0.000 0.99991
## spore.print.colorn 1.731e-01 1.565e+03 0.000 0.99991
## spore.print.coloro 3.131e-14 2.183e+03 0.000 1.00000
## spore.print.colorr 1.335e+01 1.993e+03 0.007 0.99466
## spore.print.coloru 1.766e-01 2.227e+03 0.000 0.99994
## spore.print.colorw 5.130e+00 1.544e+03 0.003 0.99735
## spore.print.colory -1.821e-14 2.183e+03 0.000 1.00000
## bruisest -5.486e-01 1.789e-01 -3.067 0.00216 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 11251.76 on 8123 degrees of freedom
## Residual deviance: 327.96 on 8106 degrees of freedom
## AIC: 363.96
##
## Number of Fisher Scoring iterations: 21
Diễn giải và nhận xét kết quả mô hình hồi quy Probit:
Mô hình này được xây dựng để dự đoán độc tính của nấm (class) bằng phương pháp Probit, một phương pháp thay thế cho Logistic. Chúng ta sẽ phân tích các thành phần chính: Hệ số (Coefficients) và độ phù hợp của mô hình (Deviance) để so sánh với mô hình Logit.
1. Phân tích các hệ số (Coefficients):
Bảng Coefficients cho chúng ta biết ảnh hưởng của từng đặc điểm lên chỉ số Z (z-score) của biến tiềm ẩn “độc tính”. Việc diễn giải các hệ số này không trực quan như Odds Ratio, nhưng dấu của chúng vẫn cho biết hướng tác động.
(Intercept) - Hệ số chặn (-6.028e+00): Đây là giá trị z-score của việc một cây nấm là độc khi tất cả các biến độc lập đều ở mức tham chiếu. Giá trị âm, cho thấy trong trường hợp này, xác suất nấm là độc cực kỳ thấp (gần như bằng 0).
odor - Ảnh hưởng của mùi:
Các hệ số của các mùi khó chịu như Creosote (odorc = 12.25) và Hôi (odorf = 12.10) đều có giá trị dương và rất lớn. Điều này có nghĩa là so với mùi tham chiếu (hạnh nhân), các mùi này làm tăng mạnh chỉ số Z, từ đó đẩy xác suất nấm là độc lên gần 100%.
Ngược lại, hệ số của odorn (không mùi) là âm (-3.721e-01), cho thấy việc không có mùi làm giảm nhẹ khả năng nấm là độc so với mùi Hạnh nhân.
spore.print.color - Ảnh hưởng của màu sắc bào tử: Tương tự, hệ số của bào tử màu sô-cô-la (spore.print.colorh = 5.33) và trắng (spore.print.colorw = 5.13) đều là các giá trị dương lớn, cho thấy chúng làm tăng mạnh khả năng nấm là độc.
bruisest - Ảnh hưởng của vết bầm: Hệ số của bruisest (không có vết bầm) là -0.5486. Giá trị âm này cho thấy, khi một cây nấm không có vết bầm (so với có vết bầm), chỉ số z-score của việc nó là độc giảm đi. Điều này hoàn toàn nhất quán với kết quả từ phân tích song biến và mô hình Logit.
2. Phân tích độ phù hợp của mô hình:
Null deviance: 11251.76 và
Residual deviance: 327.96:
Sự sụt giảm cực kỳ mạnh mẽ của Deviance từ mô hình rỗng đến mô hình đầy đủ là một dấu hiệu rất tốt, cho thấy bộ ba biến độc lập này có khả năng giải thích rất tốt sự biến thiên trong độc tính của nấm.
Quan trọng hơn, giá trị Residual deviance của mô hình Probit (327.96) gần như giống hệt với mô hình Logit. Điều này chứng tỏ cả hai mô hình có độ phù hợp với dữ liệu tốt như nhau.
AIC: 363.96: Chỉ số AIC cũng gần như giống hệt với
mô hình Logit, một lần nữa khẳng định rằng không có sự khác biệt đáng kể
về hiệu suất dự đoán giữa hai mô hình trong trường hợp này.
Mô hình hồi quy Cloglog (Complementary Log-Log) là một phương pháp hồi quy tổng quát khác được sử dụng cho biến phụ thuộc nhị phân. Nó là một sự thay thế cho mô hình Logit và Probit, đặc biệt hữu ích trong các tình huống mà xác suất của sự kiện xảy ra hoặc không xảy ra tiệm cận một cách không đối xứng về 0 hoặc 1.
Điểm khác biệt cốt lõi của mô hình Cloglog nằm ở hàm liên kết
(link function) của nó, được định nghĩa là: \[ \text{cloglog}(p) = \ln(-\ln(1-p)) \]
Trong đó p là xác suất của sự kiện xảy ra (Y=1).
Tương tự như Probit, các hệ số β trong mô hình Cloglog
khó diễn giải trực tiếp. Chúng không thể được chuyển
đổi thành Tỷ số chênh (Odds Ratio) một cách đơn giản. Dấu của hệ số vẫn
cho biết hướng tác động, nhưng độ lớn của ảnh hưởng thường được đánh giá
thông qua các tác động biên (marginal effects) hoặc các chỉ số khác tùy
thuộc vào lĩnh vực ứng dụng.
class) bằng mô hình CloglogĐể hoàn thiện việc so sánh các mô hình, chúng tôi xây dựng mô hình dự
đoán độc tính bằng phương pháp Cloglog. Biến phụ thuộc
class và các biến độc lập odor,
spore.print.color, bruises được giữ nguyên như
trong các mô hình trước.
Mô hình Cloglog: \[ \ln(-\ln(1-P(\text{class='p'}))) = \beta_0 + \beta_1\text{odor} + \beta_2\text{spore.print.color} + \beta_3\text{bruises} \]
# Xây dựng mô hình hồi quy cloglog
cloglog_model_class <- glm(class ~ odor + spore.print.color + bruises,
data = df,
family = binomial(link = "cloglog"))
# Xem kết quả thô của mô hình
summary(cloglog_model_class)
##
## Call:
## glm(formula = class ~ odor + spore.print.color + bruises, family = binomial(link = "cloglog"),
## data = df)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.296e+01 1.172e+04 -0.002 0.99844
## odorc 2.557e+01 3.889e+03 0.007 0.99475
## odorf 2.262e+01 4.277e+03 0.005 0.99578
## odorl -1.442e-11 5.372e+03 0.000 1.00000
## odorm 4.937e+00 4.274e+03 0.001 0.99908
## odorn -5.113e-01 4.125e+03 0.000 0.99990
## odorp 2.669e+01 3.881e+03 0.007 0.99451
## odors 4.937e+00 4.135e+03 0.001 0.99905
## odory 4.937e+00 4.135e+03 0.001 0.99905
## spore.print.colorh 4.551e+00 1.102e+04 0.000 0.99967
## spore.print.colork 6.106e-01 1.107e+04 0.000 0.99996
## spore.print.colorn 6.085e-01 1.107e+04 0.000 0.99996
## spore.print.coloro -3.339e-13 1.551e+04 0.000 1.00000
## spore.print.colorr 2.782e+01 1.099e+04 0.003 0.99798
## spore.print.coloru 6.095e-01 1.582e+04 0.000 0.99997
## spore.print.colorw 2.125e+01 1.097e+04 0.002 0.99845
## spore.print.colory 5.124e-14 1.551e+04 0.000 1.00000
## bruisest -1.121e+00 3.873e-01 -2.893 0.00381 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 11251.76 on 8123 degrees of freedom
## Residual deviance: 327.96 on 8106 degrees of freedom
## AIC: 363.96
##
## Number of Fisher Scoring iterations: 22
Diễn giải và nhận xét kết quả mô hình hồi quy Cloglog:
Mô hình này cung cấp một góc nhìn khác về mối quan hệ giữa các đặc điểm hình thái và độc tính, sử dụng hàm liên kết Cloglog bất đối xứng. Chúng ta sẽ phân tích các thành phần chính của kết quả để so sánh với mô hình Logit và Probit.
1. Phân tích các hệ số (Coefficients):
Dấu của các hệ số trong mô hình Cloglog vẫn cho chúng ta biết hướng tác động của các biến độc lập lên xác suất nấm là độc. Tuy nhiên, độ lớn của chúng không thể diễn giải trực tiếp.
(Intercept) - Hệ số chặn
(-2.296e+01): Đây là giá trị của hàm
log(-log(1-p)) khi tất cả các biến độc lập đều ở mức tham
chiếu. Giá trị âm, cho thấy trong trường hợp này, xác suất nấm là độc
cực kỳ thấp (gần như bằng 0).
odor - Ảnh hưởng của mùi:
odorc = 25.57) và Hôi (odorf
= 22.62) đều có giá trị dương và rất lớn. Điều này khẳng định
rằng các mùi này làm tăng mạnh xác suất nấm là độc so với mùi tham chiếu
(Hạnh nhân).odorn)
là âm (-0.5113), cho thấy nó làm giảm khả năng nấm là
độc.spore.print.color - Ảnh hưởng của màu sắc
bào tử:
h = 4.55), Xanh lá (r =
27.82), và Trắng (w = 21.25) tiếp
tục cho thấy chúng là những yếu tố làm tăng mạnh nguy cơ độc tính.bruisest - Ảnh hưởng của vết
bầm:
bruisest (không có vết bầm) là
-1.121. Giá trị âm này nhất quán với các mô hình trước,
cho thấy việc không có vết bầm làm giảm xác suất nấm là độc so với khi
có vết bầm.**), cho thấy ảnh hưởng của nó là rất có ý
nghĩa thống kê, ngay cả trong mô hình Cloglog.2. Phân tích độ phù hợp của mô hình:
Cuộc hành trình phân tích dữ liệu trên bộ sưu tập nấm của Audubon đã mang lại những hiểu biết sâu sắc, không chỉ trả lời các câu hỏi nghiên cứu ban đầu mà còn mở ra những hướng đi mới. Phần này sẽ tóm tắt lại những phát hiện quan trọng nhất, nhìn nhận các hạn chế của phân tích, và đề xuất các ứng dụng thực tiễn cũng như hướng nghiên cứu trong tương lai.
Phân tích sâu bộ dữ liệu Mushroom từ góc độ mô tả, song biến và mô hình hóa đã mang lại nhiều kết quả quan trọng và có ý nghĩa thực tiễn, có thể tóm tắt qua các điểm chính sau:
Xác định và lượng hóa các yếu tố dự báo độc
tính: Nghiên cứu đã xác định một cách rõ ràng “bộ ba quyền lực”
trong việc nhận dạng nấm độc là mùi (odor), màu sắc
bào tử (spore.print.color), và Vết bầm
(bruises).
odor) nổi lên như một yếu tố phân
loại gần như hoàn hảo. Mô hình Logistic đã lượng hóa sức mạnh này, cho
thấy Tỷ số chênh (Odds Ratio) của việc nấm là độc khi có mùi hôi là cực
kỳ lớn (>999) so với mùi hạnh nhân.Khám phá mối liên kết sinh thái: Phân tích đã chứng minh rằng các đặc điểm hình thái không tồn tại một cách ngẫu nhiên mà có mối liên kết chặt chẽ với môi trường sống. Mô hình Logistic đa loại cho thấy chúng ta hoàn toàn có thể dựa vào mùi, màu sắc và hình dạng mũ để đưa ra phỏng đoán có cơ sở về nơi một cây nấm được tìm thấy. Ví dụ, nấm có mùi hôi có odds thuộc về nhóm “Trên cỏ” cao hơn gấp 41 lần so với việc thuộc về nhóm “Trong rừng”.
So sánh và đánh giá hiệu quả các mô hình hồi quy: Nghiên cứu đã áp dụng và so sánh ba mô hình hồi quy tổng quát (Logistic, Probit, Cloglog) cho bài toán dự đoán độc tính.
Mặc dù các kết quả rất rõ ràng và mạnh mẽ, việc duy trì một tư duy khoa học đòi hỏi chúng ta phải nhận thức được những hạn chế của nghiên cứu:
Tính đại diện của dữ liệu: Các quy tắc và mối quan hệ được rút ra từ bộ dữ liệu này chỉ có giá trị trong phạm vi các loài được mô tả trong “Sổ tay hướng dẫn của Audubon”. Chúng có thể không áp dụng được cho tất cả các loài nấm trên toàn thế giới. Việc áp dụng các quy tắc này một cách máy móc ngoài thực địa mà không có kiến thức chuyên môn sâu rộng vẫn tiềm ẩn rủi ro.
Phương pháp xử lý dữ liệu thiếu: Việc điền giá
trị thiếu trong cột stalk.root bằng mode là một giả định
hợp lý để bảo toàn kích thước mẫu. Tuy nhiên, phương pháp này có thể đã
làm thay đổi nhẹ cấu trúc phân phối thực tế của biến và có khả năng ảnh
hưởng đến các phân tích có liên quan đến đặc điểm này.
Bản chất của phân tích: Đây là một phân tích dựa trên dữ liệu quan sát, không phải một thí nghiệm có kiểm soát. Chúng ta có thể xác định các mối liên hệ mạnh (associations), nhưng không thể khẳng định một cách chắc chắn về mối quan hệ nhân-quả (causation).
Các phát hiện từ phân tích này là một nền tảng vững chắc, mở ra nhiều hướng nghiên cứu và ứng dụng hấp dẫn:
Có thể tạo ra một hệ thống quy tắc phân cấp đơn giản nhưng hiệu quả cao cho người mới bắt đầu.
Sử dụng các thuật toán học máy phức tạp hơn như Cây quyết định (Decision Tree), Rừng ngẫu nhiên (Random Forest) hoặc Gradient Boosting để xây dựng một mô hình có khả năng dự đoán độc tính với độ chính xác cao nhất có thể. Cây quyết định sẽ đặc biệt thú vị vì nó có thể trực quan hóa các quy tắc phân loại, tương tự như cách chúng ta đã phân tích.
Nghiên cứu sâu hơn xem liệu sự kết hợp của các đặc điểm có tạo ra
một yếu tố dự báo mạnh hơn không. Ví dụ: “Một cây nấm mũ màu nâu
(cap-color=n) mọc trên lối đi (habitat=p) có
nguy cơ độc cao hơn một cây nấm mũ nâu mọc trong rừng
(habitat=d) hay không?”.