Trực quan hóa Dữ liệu Khuyết với gói VIM

Dữ liệu khuyết là một vấn đề trong công việc hằng ngày của mọi nhà khoa học dữ liệu. Chúng ta có nên nội suy (impute) không? Nếu có, phương pháp nào phù hợp? Hay có thể đơn giản loại bỏ các quan sát có điểm dữ liệu bị khuyết? Để trả lời những câu hỏi này, bạn cần biết cơ chế đằng sau dữ liệu khuyết là gì. Việc phát hiện bằng các kiểm định thống kê khá phức tạp và đôi khi chỉ đưa ra các nhận định mơ hồ. Ngược lại, các công cụ trực quan hóa rất dễ sử dụng và không chỉ giúp phát hiện cơ chế dữ liệu khuyết mà còn mang lại hiểu biết về các khía cạnh khác của chất lượng dữ liệu. Trong hướng dẫn này, chúng ta sẽ giới thiệu một tập các phương pháp vẽ đồ thị có sẵn trong gói VIM để cho thấy chúng có thể giúp bạn nắm chắc các mẫu hình về cách dữ liệu bị khuyết như thế nào.

Cơ chế Dữ liệu Khuyết

Có thể có nhiều lý do khiến một tập dữ liệu không đầy đủ. Việc tìm hiểu các nguyên nhân có thể dẫn đến dữ liệu khuyết là điều quan trọng, vì điều này ảnh hưởng đến cách chúng ta xử lý vấn đề. Chẳng hạn, nếu lo ngại về không hồi đáp trong các khảo sát liên quan đến việc làm, ta có thể kỳ vọng những người giàu nhất và nghèo nhất không tiết lộ thu nhập trong bảng hỏi, nghĩa là các điểm dữ liệu khuyết không phân bố đều trong toàn bộ tập dữ liệu. Trong trường hợp này, chỉ đơn giản loại bỏ các quan sát không đầy đủ trước khi phân tích sẽ dẫn đến kết quả chệch. Mặt khác, ví dụ nếu dữ liệu bị thiếu do thiết bị thu thập dữ liệu gặp trục trặc, có thể vị trí của các điểm bị khuyết trong tập dữ liệu là hoàn toàn ngẫu nhiên.

Có ba kiểu mẫu phân biệt theo đó dữ liệu có thể bị khuyết. Chúng thường được gọi là cơ chế dữ liệu khuyết.

• Khuyết hoàn toàn ngẫu nhiên (MCAR)
  Với MCAR, không có mẫu hình có hệ thống nào trong vị trí các điểm dữ liệu bị khuyết trong tập dữ liệu: chúng xảy ra hoàn toàn ngẫu nhiên. Về mặt hình thức, điều này có nghĩa xác suất một quan sát bị thiếu không phụ thuộc vào giá trị của các biến khác cũng như giá trị của chính nó. Trong trường hợp này, loại bỏ các quan sát không đầy đủ sẽ không gây sai lệch cho kết quả phân tích sau đó.
  Ví dụ: Nhiệt độ được cảm biến đo liên tục, thu thập dữ liệu và gửi qua Internet tới cơ sở dữ liệu. Vì lý do không rõ, kết nối Internet đôi khi bị gián đoạn.

• Khuyết ngẫu nhiên có điều kiện (MAR)
  Với MAR, xác suất một quan sát cụ thể bị thiếu vẫn độc lập với giá trị của chính nó, nhưng phụ thuộc vào giá trị của các biến khác. Trong trường hợp này, loại bỏ các quan sát không đầy đủ khiến mẫu kém đại diện hơn.
  Ví dụ: Một số dữ liệu vào ban đêm bị thiếu do công tác bảo trì cảm biến, vốn luôn diễn ra qua đêm.

• Khuyết không ngẫu nhiên (MNAR)
  Với MNAR, xác suất một quan sát bị thiếu phụ thuộc vào chính các giá trị không quan sát được của nó. Trong trường hợp này, một lần nữa, loại bỏ dữ liệu không đầy đủ dẫn đến phân tích sai lệch.
  Ví dụ: Cảm biến bị đóng băng ở -20 độ C và không đo nhiệt độ thấp hơn giá trị này.

Trong thực tế, rất khó xác định cơ chế dữ liệu khuyết nào trong ba cơ chế trên áp dụng cho một trường hợp cụ thể. Để thấy điều này, hãy vẽ một số dữ liệu MNAR. Chúng ta sẽ dùng các gói ggplot2 và gridExtra cần được nạp trước. Ta cũng sẽ nạp hai gói khác để dùng sau: VIM cho phân tích trực quan dữ liệu khuyết và dplyr cho một số tiền xử lý dữ liệu. Hãy nhớ install.packages() chúng trước nếu chúng vẫn chưa có trong không gian làm việc của bạn.

library("VIM")
library("dplyr")
library("ggplot2")
library("gridExtra")
Được hỗ trợ bởi

Hãy tạo một tập dữ liệu với hai biến phân phối chuẩn, độc lập nhau, x và y, trong đó y có một số giá trị bị khuyết ở phía đuôi phải của phân phối, tức là ở các giá trị cao nhất. Lưu ý rằng ta cũng giữ lại các giá trị thực và đặt chúng thành thiếu ở một biến riêng.

set.seed(2)
mnar_data <- data.frame(x = rnorm(100), y = rnorm(100)) %>%
    mutate(y_miss = ifelse(y > 1, y, NA),
           y = ifelse(is.na(y_miss), y, NA),
           x_miss = ifelse(is.na(y), x, NA))
Được hỗ trợ bởi

Bây giờ ta sẽ vẽ hai biểu đồ phân tán của x theo y: một biểu đồ gồm cả các điểm mà ta giả định y bị thiếu, và một biểu đồ chỉ với dữ liệu quan sát được.

grid.arrange(
  # Plot highlighting missing data points
  ggplot(mnar_data, aes(x, y)) +
    geom_point(size = 4, alpha = 0.6) +
    geom_point(aes(x, y_miss), col = "red", size = 4, alpha = 0.6) +
    geom_hline(aes(yintercept = y_miss), col = "red", alpha = 0.6,
               linetype = "dashed") +
    ylim(NA, max(mnar_data$y_miss, na.rm = TRUE)) +
    ggtitle("Where the data points are missing"),
  # What data scientist can see
  ggplot(mnar_data, aes(x, y)) +
    geom_point(size = 4, alpha = 0.6) +
    geom_vline(aes(xintercept = x_miss), col = "red", alpha = 0.6,
               linetype = "dashed") +
    ylim(NA, max(mnar_data$y_miss, na.rm = TRUE)) +
    ggtitle("What the data scientist can see"),
  # Arrange the two plots side-by-side
  ncol = 2
)
Được hỗ trợ bởi

Trong khi từ biểu đồ đầu tiên có thể thấy rõ dữ liệu thực sự là MNAR, những gì nhà khoa học dữ liệu nhìn thấy chỉ là tọa độ x của các điểm bị khuyết. Có thể rút ra kết luận gì từ biểu đồ thứ hai? Dữ liệu có phải MCAR không? Có thể có, vị trí các giá trị bị thiếu dường như phân bố khá đồng đều theo trục x. Hay là MAR? Quan sát kỹ hơn, có vẻ có nhiều giá trị y bị thiếu khi x âm hơn khi x dương, nên khả năng này cũng có thể xảy ra. Cuối cùng, MNAR không bao giờ có thể bị loại trừ, vì để biết dữ liệu có MNAR hay không, theo định nghĩa, chúng ta cần quan sát các giá trị không quan sát được.

May mắn thay, ta không hoàn toàn bế tắc. Các công cụ trực quan hóa dữ liệu có thể cung cấp một số định hướng về các mẫu hình trong dữ liệu khuyết. Giờ ta sẽ phân tích một số đồ thị hữu ích để phát hiện các mẫu hình đó dựa trên tập dữ liệu biopics từ gói fivethirtyeight, chứa một số thông tin về các phim tiểu sử.

Làm sạch Dữ liệu

Trước khi chuyển sang phần vẽ đồ thị, dữ liệu cần được tiền xử lý. Khối mã dưới đây sẽ thực hiện việc đó. Lời gọi select ở đầu pipeline dplyr trích xuất các biến quan tâm và mutate sau đó gom các chủng tộc liên quan lại với nhau để đủ số quan sát cho mỗi nhóm. Sau đó, ta điều chỉnh kiểu của một số biến: cần factor cho biến chuỗi và integer cho biến logic. Cuối cùng, các biến được đổi tên ngắn gọn hơn để biểu đồ dễ đọc.

data(biopics, package = "fivethirtyeight")

biopics <- biopics %>%
  select(country, year_release, box_office, number_of_subjects,
         type_of_subject, subject_race, person_of_color, subject_sex) %>%
  mutate(subject_race = ifelse(grepl("^Hispanic", subject_race), "Hispanic",
                               subject_race),
         subject_race = ifelse(grepl("^African", subject_race), "African",
                               subject_race),
         subject_race = ifelse(grepl("^Middle", subject_race), "Mid Eastern",
                               subject_race),
         subject_race = ifelse(subject_race %in% c("White", "Asian", "African",
                                                   "Hispanic", "Mid Eastern",
                                                   "Multi racial", NA),
                               subject_race, "other")) %>%
  mutate(country = as.factor(country),
         type_of_subject = as.factor(type_of_subject),
         subject_race = as.factor(subject_race),
         subject_sex = as.factor(subject_sex),
         person_of_color = as.integer(person_of_color)) %>%
      as.data.frame()

colnames(biopics) <- c("country", "year", "earnings", "sub_num",
                       "sub_type", "sub_race", "non_white", "sub_sex")
Được hỗ trợ bởi

Kết quả là một tập dữ liệu gồm tám biến:

• country - quốc gia hoặc các quốc gia sản xuất phim,
• year - năm phát hành,
• earnings - doanh thu gộp tại phòng vé Hoa Kỳ,
• sub_num - số nhân vật được thể hiện trong phim,
• sub_type - nghề nghiệp của nhân vật hoặc lý do được ghi nhận,
• sub_race - chủng tộc của nhân vật,
• non_white - biến giả cho biết là người da màu,
• sub_sex - giới tính của nhân vật.

Biểu đồ Tổng hợp

Câu hỏi cơ bản đầu tiên cần đặt ra là: ở những biến nào có quan sát bị thiếu, và bao nhiêu? Biểu đồ tổng hợp là công cụ hữu ích để trả lời các câu hỏi này. Dòng lệnh một dòng dưới đây là tất cả những gì bạn cần.

aggr(biopics, numbers = TRUE, prop = c(TRUE, FALSE))
Được hỗ trợ bởi

Ta chỉ định đối số numbers = TRUE để hiển thị con số trên đỉnh các cột. Đối số prop là giá trị logic cho biết có dùng tỷ lệ giá trị khuyết và các tổ hợp thay vì tổng số lượng hay không. Ta đặt TRUE cho biểu đồ đầu và FALSE cho biểu đồ thứ hai.

Có thể thấy giá trị khuyết chỉ xuất hiện ở hai biến: chúng chiếm hơn 40% ở earnings và khoảng 25% ở sub_race. Trong biểu đồ các tổ hợp ở bên phải, lưới thể hiện tất cả các tổ hợp giá trị khuyết (đỏ) và quan sát được (xanh) có trong dữ liệu. Có 317 quan sát đầy đủ, và trong 77 hàng cả hai biến đều bị khuyết.

Việc nhìn đồng thời hai biểu đồ có thể hơi rối, nên hãy dồn cùng thông tin đó vào một trực quan hóa đơn bằng dòng lệnh sau.

aggr(biopics, combined = TRUE, numbers = TRUE)
Được hỗ trợ bởi

Đối số combined cho phép gộp hai biểu đồ thành một. Ở đây, các cột ngang bên phải lưới cho thấy tần suất của các tổ hợp tương ứng, trong khi các cột dọc phía trên thể hiện tỷ lệ giá trị khuyết ở mỗi biến. Cộng các giá trị của mọi tổ hợp có thể cho earnings, tức 0.32 và 0.10, ta thấy tỷ lệ tổng thể giá trị khuyết ở biến này, vốn ta ước lượng là hơn 40% từ biểu đồ cột trước, thực tế là 42%.

Spinogram và Spineplot

Cho đến giờ, ta đã có cái nhìn tổng quan về dữ liệu khuyết. Đã đến lúc xem kỹ hơn các tương tác giữa từng biến cụ thể. Spinogram và spineplot cho phép nghiên cứu tỷ lệ giá trị khuyết ở một biến theo các giá trị khác nhau của biến còn lại. Nếu biến còn lại là số, ta có spinogram; khi nó là phân loại, trực quan hóa này gọi là spineplot.

Cả hai có thể được tạo bằng hàm spineMiss(), nhận vào một dataframe hai cột. Biến được chỉ định đầu tiên là biến dùng để chia dữ liệu. Nó được ánh xạ lên trục ngang của đồ thị. Biến thứ hai là biến mà ta quan tâm đến mẫu hình dữ liệu khuyết. Để thấy trong thực tế, hãy tạo một spineplot. Ta chỉ định biến phân loại sub_race trước, theo sau là earnings, có thể hiểu là: tỷ lệ giá trị khuyết trong earnings cho mỗi nhóm của sub_race là bao nhiêu?

spineMiss(biopics[, c("sub_race", "earnings")])
Được hỗ trợ bởi

Độ rộng tương đối của các cột cho từng nhóm trong sub_race phản ánh tần suất của nhóm đó trong tập dữ liệu: chẳng hạn, trong phần lớn phim, nhân vật chính là người da trắng. Bên trong mỗi cột, tỷ lệ khuyết của earnings được hiển thị, trong khi cột bóng mờ ở bên phải thể hiện tỷ lệ này cho toàn bộ tập dữ liệu. Có vẻ khi nhân vật chính là người châu Phi, ta có khả năng cao nhất có đầy đủ thông tin doanh thu.

Bằng cách đảo thứ tự các biến như trong dòng mã dưới đây, ta có thể tạo một spinogram để trả lời câu hỏi theo hướng ngược lại: tỷ lệ dữ liệu khuyết trong sub_race là bao nhiêu với các giá trị khác nhau của earnings?

spineMiss(biopics[, c("earnings", "sub_race")])
Được hỗ trợ bởi

Vì earnings là biến số, các giá trị của nó được chia vào các bin, trong đó độ rộng các bin tương ứng với phân phối của chính biến đó. Spinogram cho thấy earnings lệch trái rất mạnh: chỉ vài bộ phim đạt lợi nhuận cao nhất. Thú vị là, với những phim bom tấn đó, ta có khả năng cao nhất bị thiếu chủng tộc của nhân vật, được phản ánh qua cột đỏ cao nhất ở các giá trị lớn nhất của earnings.

Biểu đồ Mosaic

Spinogram và spineplot cho phép nghiên cứu tương tác giữa hai biến. Ý tưởng này được tổng quát hóa cho nhiều biến dưới dạng một biểu đồ mosaic. Biểu đồ này là tập hợp các ô, trong đó mỗi ô tương ứng với một tổ hợp cụ thể các nhóm (đối với biến phân loại) hoặc các bin (đối với biến số) của hai biến trở lên. Bên trong mỗi ô, tỷ lệ điểm dữ liệu bị khuyết ở một biến khác được hiển thị. Về nguyên tắc, có thể tạo biểu đồ mosaic cho bất kỳ số lượng biến nào, nhưng nếu quá nhiều sẽ trở nên rối. Ngoài ra, chúng dễ đọc hơn nếu các biến phân chia là factor và không có quá nhiều mức.

Hãy xem tỷ lệ khuyết trong earnings khi chia theo sub_sex và US_movie. Biến sau, được tạo trong khối mã dưới đây, là biến logic cho biết liệu Hoa Kỳ có tham gia sản xuất bộ phim hay không. Ta cung cấp dataframe ba biến làm đối số đầu tiên cho mosaicMiss(). Đối số plotvars đặt thành vector gồm phần tử 1 và 2 nghĩa là ta muốn chia dữ liệu theo hai cột đầu của dataframe này. Đặt highlight là 3 nghĩa là ta muốn hiển thị tỷ lệ khuyết của biến ở cột thứ ba trong các ô.

biopics <- biopics %>%
    mutate(US_movie = ifelse(grepl("US", country), TRUE, FALSE))
mosaicMiss(biopics[, c("sub_sex", "US_movie", "earnings")], highlight = 3,
           plotvars = 1:2, miss.labels = FALSE)
Được hỗ trợ bởi

Một lần nữa, kích thước các ô tương ứng với tần suất xuất hiện của các tổ hợp trong tập dữ liệu. Chẳng hạn, ô lớn nhất ở dưới cùng bên phải cho biết hầu hết các phim có nhân vật chính là nam và được sản xuất, ít nhất một phần, ở Mỹ.

Đối với các phim có nhân vật chính là nam, thông tin doanh thu có khả năng bị thiếu nhiều hơn ở các phim không thuộc Mỹ. Mặt khác, nếu nhân vật chính là nữ, các phim của Mỹ có tỷ lệ thiếu doanh thu nhỉnh hơn một chút. Tuy nhiên, các khác biệt này dường như quá nhỏ để coi là có ý nghĩa.

Boxplot Song song

Một kiểu trực quan hóa khác là boxplot song song. Ý tưởng là chia tập dữ liệu thành hai tập con: một tập chỉ gồm các giá trị quan sát được của một biến không đầy đủ và tập kia chỉ gồm các giá trị bị khuyết của biến đó. Với cả hai tập con, vẽ boxplot của một biến số được chọn. Điều này cho phép bạn kiểm tra xem phân phối của biến được chọn có bị ảnh hưởng bởi tình trạng khuyết của biến dùng để chia hay không.

Để tạo boxplot song song, ta có thể dùng hàm pbox(). Ta chỉ truyền một đối số: một dataframe hai cột. Boxplot sẽ được vẽ cho biến ở cột đầu, còn biến ở cột thứ hai được dùng để chia. Trong ví dụ dưới đây, ta lấy log của earnings vì chúng lệch rất mạnh và nếu không biến đổi, boxplot sẽ giống như một đường thẳng.

biopics <- biopics %>%
  mutate(log_earnings = log(earnings))
pbox(biopics[, c("log_earnings", "sub_race")])
Được hỗ trợ bởi

Hộp trắng ở bên trái cho thấy phân phối tổng thể của log_earnings, trong khi hộp xanh và đỏ lần lượt cho thấy phân phối của nó cho các tập con có giá trị quan sát được và bị khuyết trong sub_race. Độ rộng tương đối của các hộp phản ánh kích thước các tập con: hộp xanh rộng hơn nghĩa là có nhiều giá trị quan sát được hơn giá trị khuyết trong sub_race. Ngoài ra, hai hộp trông khá giống nhau và cũng giống với hộp trắng tổng thể. Điều này gợi ý rằng việc thiếu thông tin chủng tộc không ảnh hưởng đến phân phối doanh thu.

Biểu đồ Tọa độ Song song

Ta đã xem các biến đơn lẻ và tương tác của chúng. Giờ hãy chuyển sang phân tích tất cả các biến trong tập dữ liệu cùng lúc. Trực quan hóa cho phép điều này là biểu đồ tọa độ song song. Trong biểu đồ này, mỗi biến được biến đổi về cùng một thang đo và được biểu diễn bằng một trục song song. Với biến phân loại, thang đo trên trục được chia thành một tập các điểm cách đều, mỗi điểm cho một nhóm. Các giá trị khuyết được đặt phía trên các trục dọc, ngoài vùng vẽ. Mỗi đường trên biểu đồ tương ứng với một quan sát trong dữ liệu, màu của đường cho biết dữ liệu khuyết ở biến được chọn.

Ta có thể tạo biểu đồ tọa độ song song bằng hàm parcoordMiss(). Mặc định, hàm dùng tất cả các biến trong dataframe được truyền làm đối số đầu tiên. Đặt highlight là earnings sẽ tô màu khác cho các đường đại diện các quan sát có giá trị khuyết ở earnings. Đặt alpha-blending là 0.6 làm các đường hơi trong suốt, tăng độ dễ đọc của biểu đồ.

# Remove variables created for the previous plots
biopics <- biopics %>%
  select(- US_movie, - log_earnings)

parcoordMiss(biopics, highlight = 'earnings', alpha = 0.6)
Được hỗ trợ bởi

Trong biểu đồ trên, các đường đỏ đậm biểu thị các quan sát có giá trị khuyết ở earnings. Những quan sát được tô sáng này có vẻ hành xử khác với phần còn lại. Đặc biệt, chỉ một vài trong số đó tương ứng với mức cao thứ hai của country, vốn bao gồm rất nhiều quan sát ở phần còn lại! Hóa ra đó là các phim US/UK (để thấy điều này, chỉ cần chạy levels(biopics$country) và tìm giá trị áp chót!).

Hơn nữa, quan sát nhanh trục year cho thấy một đoạn nào đó khoảng hai phần ba chiều cao trục có ít đường đỏ đậm cắt qua. Có vẻ đã có một giai đoạn không lâu trước đây dữ liệu doanh thu đầy đủ hơn. Điều này gợi ý rằng cả country và year có thể giúp giải thích phân bố giá trị khuyết trong earnings!

Biểu đồ Ma trận

Công cụ trực quan cuối cùng ta xem xét là biểu đồ ma trận. Nó trực quan hóa mọi ô của ma trận dữ liệu bằng các hình chữ nhật. Dữ liệu quan sát được hiển thị theo thang màu xám-đen liên tục (màu càng đậm, giá trị càng cao), trong khi các giá trị khuyết được tô đỏ. Thực hành tốt là sắp xếp dữ liệu theo một trong các biến không đầy đủ - điều này giúp biểu đồ dễ diễn giải hơn. Khối mã dưới đây làm điều đó, sắp xếp theo earnings.

matrixplot(biopics, sortby = c('earnings'))
Được hỗ trợ bởi

Biểu đồ xác nhận một số phát hiện trước đó: các quan sát có giá trị khuyết ở earnings có xu hướng thuộc về lâu trước đây (giá trị year thấp) và cũng có điểm thấp ở country - điều này được thể hiện qua màu nhạt hơn của hai biến này trong các hàng mà earnings có màu đỏ. Ngoài ra, có vẻ có ít phim với nhân vật không phải người da trắng khi earnings bị khuyết.

Kết luận

Tóm lại, có vẻ dữ liệu khuyết trong tập biopics không phải MCAR. Vị trí các điểm dữ liệu khuyết trong earnings và sub_race giúp giải thích lẫn nhau, như spineplot và spinogram cho thấy. Hơn nữa, biểu đồ ma trận và biểu đồ tọa độ song song gợi ý rằng country, year và non_white cũng có thể giúp giải thích phân bố giá trị khuyết trong earnings. Do đó, không phải là ý tưởng hay khi loại bỏ các quan sát không đầy đủ, vì điều này rất có thể sẽ đưa sai lệch vào suy luận của chúng ta.

Suy ngẫm cuối

Bạn đã đi được một chặng đường dài! Giờ bạn đã biết các cơ chế dữ liệu khuyết là gì và chúng khác nhau thế nào, và quan trọng hơn, cách suy luận về chúng bằng nhiều công cụ trực quan hóa. Làm tốt lắm! Giờ bạn đã sẵn sàng để phân tích các tập dữ liệu không đầy đủ của riêng mình!

Nếu bạn quan tâm học thêm về R, hãy tham gia khóa Data Visualization with ggplot2 (Part 1) của DataCamp và xem Hướng dẫn Công thức R của chúng tôi.