Giới thiệu dự án

Đại dịch COVID-19 bùng phát trên quy mô toàn cầu đã đặt ra những thách thức chưa từng có đối với hệ thống y tế công cộng, an sinh xã hội và sự phát triển kinh tế toàn cầu. Theo dữ liệu từ Tổ chức Y tế Thế giới (WHO) và kho dữ liệu mở Our World in Data (OWID), hàng trăm triệu ca nhiễm và hàng triệu ca tử vong đã được ghi nhận với diễn biến biến động phức tạp theo từng giai đoạn và biến thể (Alpha, Delta, Omicron). Việc theo dõi, tổng hợp và mô hình hóa toán học đối với chuỗi dữ liệu thời gian (Time-series Data) trở thành tiền đề cốt lõi để các nhà hoạch định chính sách đưa ra quyết định phòng chống dịch kịp thời.

Tuy nhiên, việc khai thác dữ liệu dịch tễ học trong thực tế đối mặt với nhiều rào cản kỹ thuật:

  • Dữ liệu nhiễu và sai lệch định dạng: Tồn tại các giá trị âm bất thường (do điều chỉnh số liệu hồi tố), dữ liệu khuyết thiếu (NA) do hệ thống báo cáo không đồng bộ theo ngày cuối tuần.
  • Sự pha trộn giữa thực thể địa lý và dữ liệu tổng hợp: Tập dữ liệu thô chứa cả các mã định danh khu vực (OWID_HIC, OWID_ASI, OWID_WRL) lẫn quốc gia độc lập, gây sai lệch nghiêm trọng cho các phép toán gộp (aggregation).
  • Độ lệch phân phối cực lớn (High Skewness): Phân bố số ca nhiễm mới và tử vong có phương sai rất cao, khiến các chỉ số đo lường xu hướng trung tâm truyền thống như kỳ vọng toán học (Mean) mất đi tính đại diện.

Đồ án "Phân Tích Dữ Liệu Covid-19: Thống Kê và Kết Quả Nghiên Cứu" (Mã học phần: CO1007 - Cấu trúc Rời rạc cho Khoa học Máy tính, Khoa Khoa học và Kỹ thuật Máy tính, Trường Đại học Bách Khoa - ĐHQG TP.HCM) được phát triển nhằm giải quyết triệt để các bài toán xử lý và phân tích số liệu trên.

+-------------------------------------------------------------------------+
|                              MỤC TIÊU DỰ ÁN                             |
+-------------------------------------------------------------------------+
| 1. Xây dựng Data Pipeline tự động chuẩn hóa 163,090+ bản ghi dữ liệu.    |
| 2. Triển khai bộ công cụ thống kê mô tả: Tứ phân vị, IQR, Boxplot, SD.  |
| 3. Phát triển thuật toán lọc nhiễu 7-day Simple Moving Average (SMA).   |
| 4. Mô hình hóa tương quan Pearson và Hồi quy tuyến tính bằng OLS Matrix.|
| 5. Phát hiện chuỗi bùng phát bằng giải thuật Run-Length Encoding (RLE). |
+-------------------------------------------------------------------------+

Phạm vi và giới hạn nghiên cứu

  • Nguồn dữ liệu gốc: Kho dữ liệu dịch tễ học chính thức OWID với 163,090 quan sát và 6 thuộc tính nền tảng (iso_code, continent, location, date, new_cases, new_deaths).
  • Không gian phân tích chuyên sâu: Bộ ba quốc gia/vùng lãnh thổ thuộc Châu Đại Dương gồm Australia, New Caledonia và New Zealand trong khung thời gian định kỳ (tháng 02, 05, 06, 09) và giai đoạn cao điểm 2020 – 2022.
  • Giới hạn kỹ thuật: Phân tích tập trung vào các mối quan hệ tuyến tính và thống kê tần suất rời rạc, chưa tích hợp các mạng nơ-ron hồi quy phi tuyến tính (RNN/LSTM) hoặc mô hình dịch tễ động học vi phân (SIR/SEIR).

Phân tích và thiết kế giải pháp

Phân tích hiện trạng

Trước khi xây dựng hệ thống xử lý chuyên dụng bằng R, việc phân tích dữ liệu dịch tễ thường phụ thuộc vào bảng tính thủ công hoặc các nền tảng BI đóng gói sẵn.

Tiêu chí đánh giá Bảng tính truyền thống (MS Excel) Nền tảng BI đóng gói (Tableau / PowerBI) Pipeline Lập trình R Chuyên dụng (Đồ án)
Khả năng xử lý dữ liệu lớn Hạn chế (< 1 triệu dòng, dễ tràn RAM) Tốt, nhưng khó can thiệp thuật toán sâu Tối ưu hóa bộ nhớ thông qua Vectorization
Tính linh hoạt của thuật toán Thấp, phụ thuộc công thức dựng sẵn Trung bình, giới hạn bởi hàm DAX/Custom Tuyệt đối: Tự cài đặt đại số ma trận, RLE
Xử lý Outlier & Phân vị Rườm rà, khó tự động hóa chuỗi Có sẵn biểu đồ nhưng khó tùy biến biên IQR Tự động hóa tính toán $Q_1, Q_2, Q_3$, lọc ngoại lai
Làm mịn dữ liệu chuỗi thời gian Phức tạp khi dịch chuyển trượt dòng Hỗ trợ hàm cơ bản Triển khai cửa sổ trượt linh hoạt ($k=7$)
Khả năng tái lập (Reproducibility) Rất thấp, dễ lỗi thao tác thủ công Trung bình Rất cao: 100% tái lập qua source script

Yêu cầu hệ thống theo mô hình MoSCoW

  • Must have: Module làm sạch dữ liệu (loại bỏ 13 mã thực thể phi quốc gia, khử giá trị âm bằng trị tuyệt đối, định dạng lại chuẩn ISO Date), tính toán chính xác $Q_1, Q_2, Q_3$, IQR, phương sai, độ lệch chuẩn.
  • Should have: Module khử nhiễu chu kỳ tuần bằng trung bình động 7 ngày (7-day SMA), biểu đồ phân phối tích lũy và Boxplot trực quan hóa đa chiều.
  • Could have: Tự xây dựng hàm giải hệ phương trình chuẩn tắc đại số tuyến tính cho hồi quy OLS mà không dựa hoàn toàn vào hàm đóng gói sẵn.
  • Won't have: Giao diện Web Client tương tác thời gian thực (ưu tiên báo cáo phân tích tự động hóa).

Thiết kế hệ thống

Kiến trúc hệ thống được thiết kế theo luồng dữ liệu một chiều (Unidirectional Pipeline Architecture), tách biệt giữa tầng tiền xử lý, tầng tính toán thống kê rời rạc và tầng trực quan hóa nâng cao.

flowchart TD
    A["Dữ liệu thô OWID (owid-covid-data.csv)"] --> B["Data Ingestion & Type Casting (Date format, pacman, dplyr)"]
    B --> C["Data Sanitization (Xử lý âm abs(), Lọc NOT.COUNTRY, Filter Missing NA)"]
    C --> D["Tầng Thống kê Mô tả (Descriptive Engine)"]
    C --> E["Tầng Xử lý Chuỗi Thời gian (Time-series Engine)"]
    C --> F["Tầng Tương quan & Dự báo (Inference Engine)"]
    
    D --> D1["Mean, Median (Q2), Q1, Q3, IQR"]
    D --> D2["Outlier Detection (1.5 x IQR Rule)"]
    D --> D3["Box-and-Whisker Visualization"]
    
    E --> E1["7-day Simple Moving Average (SMA)"]
    E --> E2["Run-Length Encoding (RLE Streak Analysis)"]
    E --> E3["Cumulative Frequency & Relative Curves"]
    
    F --> F1["Pearson Correlation Coefficient (r)"]
    F --> F2["OLS Normal Equation Matrix Inversion: w = (X^T X)^-1 X^T y"]
    F --> F3["Giao thoa đợt bùng phát (Temporal Intersection)"]

Bảng thông số công nghệ (Technology Stack)

Thành phần / Thư viện Phiên bản khuyến nghị Vai trò kiến trúc trong hệ thống
R Core Runtime >= 4.1.2 Môi trường thông dịch và tính toán số học chính xác cao
pacman 0.5.1 Quản lý vòng đời và tự động cài đặt / nạp các package phụ thuộc
dplyr / tidyverse 1.0.8 / 1.3.1 Thao tác ngữ nghĩa dữ liệu (group_by, mutate, summarise)
ggplot2 & cowplot 3.3.5 / 1.1.1 Dựng đồ họa thống kê đa lớp, trực quan hóa Boxplot, Scatter
matlib 0.9.5 Thực hiện tính ma trận nghịch đảo phục vụ phương pháp OLS
scales 1.1.1 Chuẩn hóa tỉ lệ trục tọa độ và định dạng phần trăm (percent)

Phương pháp nghiên cứu và cơ sở toán học

  1. Thước đo xu hướng trung tâm và độ biến thiên: $$\text{Mean}(D) = \bar{x} = \frac{1}{k}\sum_{i=1}^k x_i, \quad V(D) = \sigma^2 = \frac{1}{k}\sum_{i=1}^k (x_i - \bar{x})^2, \quad \text{Std}(D) = \sigma = \sqrt{V(D)}$$
  2. Xác định Tứ phân vị và Ngoại lai (Outlier Detection): $$Q_1(D) = \text{Median}({x \mid x < \text{Median}(D)}), \quad Q_3(D) = \text{Median}({x \mid x > \text{Median}(D)})$$ $$\text{IQR} = Q_3 - Q_1 \implies \text{Outlier} \notin [Q_1 - 1.5 \times \text{IQR}, ; Q_3 + 1.5 \times \text{IQR}]$$
  3. Mô hình Hồi quy Tuyến tính Bình phương Cực tiểu (OLS Normal Equation): Tìm vector trọng số $w = [w_0, w_1]^T$ tối thiểu hóa hàm mất mát $L(F) = \sum_{i=1}^n (w_0 + w_1 x_i - y_i)^2$: $$X w - y \perp \text{Col}(X) \iff X^T(Xw - y) = 0 \iff w = (X^T X)^{-1} X^T y$$
  4. Hệ số tương quan tuyến tính Pearson ($r$): $$r_{xy} = \frac{\sum_{i=1}^n (x_i - \bar{x})(y_i - \bar{y})}{\sqrt{\sum_{i=1}^n (x_i - \bar{x})^2 \sum_{i=1}^n (y_i - \bar{y})^2}}$$

Implementation và kết quả

Quy trình phát triển và Code Snippets cốt lõi

Quy trình tiền xử lý được chuẩn hóa nhằm chuyển đổi 163,090 dòng dữ liệu thô sang 153,316 quan sát quốc gia hợp lệ, loại bỏ các mã tổng hợp OWID:

# Tiền xử lý: Chuẩn hóa kiểu ngày tháng và khử nhiễu số âm
COVID.DATA$date <- as.Date(COVID.DATA$date, format="%m/%d/%Y")
COVID.DATA$new_cases <- abs(COVID.DATA$new_cases)
COVID.DATA$new_deaths <- abs(COVID.DATA$new_deaths)

# Tách lọc các thực thể địa lý không phải quốc gia
NOT.COUNTRY <- unique(filter(COVID.DATA, continent == "")$iso_code)
COVID.DATA <- COVID.DATA[!COVID.DATA$iso_code %in% NOT.COUNTRY, ]

1. Thuật toán xác định tứ phân vị và lọc Outlier

Hàm tính toán được tối ưu hóa cho cấu trúc mảng rời rạc, tuân thủ đúng định nghĩa phân vị thực nghiệm:

Q1 <- function(attr) {
  sorted_data <- sort(attr, decreasing = FALSE)
  return(median(sorted_data[1 : floor(length(sorted_data) / 2)]))
}

Q3 <- function(attr) {
  sorted_data <- sort(attr, decreasing = FALSE)
  return(median(sorted_data[ceiling(length(sorted_data) / 2 + 1) : length(sorted_data)]))
}

Outliers <- function(attr) {
  q1_val <- Q1(attr)
  q3_val <- Q3(attr)
  iqr_val <- q3_val - q1_val
  return(sum(attr < (q1_val - 1.5 * iqr_val) | attr > (q3_val + 1.5 * iqr_val)))
}

2. Thuật toán Hồi quy Tuyến tính bằng Đại số Tuyến tính Ma trận

Hiện thực hóa trực tiếp công thức phương trình chuẩn tắc OLS $w = (X^T X)^{-1} X^T y$:

linreg <- function(x, y) {
  # Khởi tạo ma trận thiết kế X với cột bias toàn 1
  X_mat <- cbind(rep(1, length(x)), x)
  # Tính toán vector hệ số w thông qua ma trận nghịch đảo
  w <- matlib::inv(t(X_mat) %*% X_mat) %*% (t(X_mat) %*% y)
  return(w)
}

3. Thuật toán xác định chuỗi bùng phát dịch bằng Run-Length Encoding (RLE)

get_duration <- function(data, k, attr) {
  # Mã hóa độ dài vệt thỏa mãn ngưỡng bùng phát k
  streak <- rle(tidyr::replace_na(data[[attr]], 0) >= k)
  duration <- data.frame(Begin = as.Date(character()), End = as.Date(character()), Streak = integer())
  current_date <- as.Date(data$date[1])
  
  for (id in 1:length(streak$lengths)) {
    if (streak$values[id]) {
      end_date <- current_date + streak$lengths[id]
      duration <- rbind(duration, data.frame(
        Begin = current_date, 
        End = end_date, 
        Streak = streak$lengths[id]
      ))
    }
    current_date <- current_date + streak$lengths[id]
  }
  return(duration)
}

Kiểm thử và Đánh giá Thực nghiệm

Tập dữ liệu đã trải qua các bài kiểm thử tính toàn vẹn (Data Integrity Tests) nhằm đảm bảo không còn giá trị rỗng hoặc giá trị âm gây suy biến ma trận hồi quy.

Total Processed Records: 153,316 observations
Total Evaluated Countries: 225 nations across 6 continents
Africa: 38,647 records | Asia: 35,528 records | Europe: 36,375 records
North America: 24,438 records | South America: 9,335 records | Oceania: 8,993 records

Kết quả Đạt được

Bảng tổng hợp thống kê mô tả ca nhiễm và tử vong mới

Quốc gia / Vùng lãnh thổ Chỉ số Min $Q_1$ Median ($Q_2$) $Q_3$ Max Trung bình ($\bar{x}$) Độ lệch chuẩn ($\sigma$) Số lượng Outliers
Australia Ca nhiễm 0 10.0 26.0 646.0 175,271 6,445.94 22,051.40 127
Ca tử vong 0 0.0 0.0 7.0 155 6.82 15.79 60
New Caledonia Ca nhiễm 0 0.0 0.0 1.0 3,343 62.01 279.89 151
Ca tử vong 0 0.0 0.0 1.0 16 1.08 3.35 33
New Zealand Ca nhiễm 0 1.0 4.0 23.0 2,336 26.89 132.41 132
Ca tử vong 0 0.0 0.0 0.0 4 0.08 0.35 41

Đánh giá Tương quan và Chu kỳ Bùng phát Dịch

  1. Hiện tượng lệch phân phối (Skewness): Tại Australia, $Q_2 = 26$ nhưng Mean lên tới $6,445.94$ với độ lệch chuẩn khổng lồ $\sigma = 22,051.40$. Điều này chứng minh giá trị kỳ vọng bị kéo lệch hoàn toàn bởi làn sóng lây nhiễm biến thể Omicron cuối năm 2021 – đầu năm 2022. Trung vị ($Q_2$) là đại lượng phản ánh trung tâm chính xác hơn kỳ vọng toán học trong bài toán này.
  2. Hiệu quả của bộ lọc 7-day Simple Moving Average: Khi chưa làm mịn, hệ số tương quan Pearson giữa ca nhiễm và ca tử vong biến động không đồng nhất (nhiều tháng nhận giá trị NA hoặc $r \approx 0$). Sau khi áp dụng 7-day SMA, hệ số tương quan tuyến tính tăng vọt lên mức có ý nghĩa thống kê cao:
    • Australia (Tháng 09/2021): $r = 0.852$ (Tương quan thuận rất chặt chẽ).
    • New Caledonia (Tháng 09/2021): $r = 0.760$.
    • New Zealand (Tháng 05/2020): $r = 0.869$.
  3. Phân tích đợt bùng phát đỉnh điểm ($k \ge 1,000$ ca nhiễm/ngày):
    • Australia: Chuỗi bùng phát kéo dài kỷ lục 97 ngày liên tiếp (từ ngày 15/11/2021 đến 20/02/2022).
    • New Caledonia: Chuỗi đạt cực đại 5 ngày (14/02/2022 đến 19/02/2022).
    • New Zealand: Chuỗi đạt cực đại 5 ngày (15/02/2022 đến 20/02/2022).
    • Khoảng thời gian giao thoa bùng phát lớn nhất giữa cả 3 quốc gia: Chính xác từ 15/02/2022 đến 20/02/2022 (5 ngày trùng lặp đỉnh dịch toàn khu vực).

Đổi mới và đóng góp

  • Tự động hóa toàn diện quy trình đại số ma trận: Thay vì phụ thuộc vào các hàm cấp cao có sẵn như lm(), đồ án đã hiện thực trực tiếp công thức nghịch đảo ma trận OLS chuẩn tắc $w = (X^T X)^{-1} X^T y$ tích hợp cùng thư viện matlib, giúp sinh viên và người nghiên cứu hiểu sâu bản chất tối ưu hóa vi tích phân trong không gian đa chiều.
  • Giải thuật phát hiện vệt bùng phát (Streak Detection) dựa trên RLE: Tận dụng triệt để cấu trúc Run-Length Encoding để quét chuỗi thời gian tuyến tính $O(N)$, xác định chính xác các khoảng thời gian liên tục thỏa mãn điều kiện ngưỡng dịch tễ mà không cần lặp lồng phức tạp.
  • Kỹ thuật xử lý bất đối xứng dữ liệu y tế: Đề xuất phương pháp kết hợp linh hoạt giữa kiểm định đồ họa Boxplot không có Outlier (outline = FALSE) cùng đường tần số tích lũy dạng bậc thang để so sánh mức độ đồng đều trong việc thu thập dữ liệu giữa các lục địa.

Ứng dụng thực tế và triển khai

+-------------------------------------------------------------------------+
|                  KỊCH BẢN ỨNG DỤNG THỰC TẾ (USE CASES)                  |
+-------------------------------------------------------------------------+
| 1. Giám sát dịch tễ học cấp quốc gia: Đưa ra cảnh báo sớm khi chỉ số   |
|    7-day SMA vượt qua ngưỡng phương sai lịch sử 2 sigma.               |
| 2. Hoạch định nguồn lực y tế: Ước lượng số giường ICU và máy thở dự phòng|
|    dựa trên mô hình hồi quy OLS tương quan giữa ca nhiễm và tử vong.   |
| 3. Đánh giá hiệu quả Vaccine: So sánh hệ số góc w1 của đường hồi quy    |
|    giữa các làn sóng dịch trước và sau khi tiêm chủng diện rộng.       |
+-------------------------------------------------------------------------+

Hướng dẫn triển khai mã nguồn (Deployment Instructions)

# 1. Clone repository chứa mã nguồn và dữ liệu
git clone https://github.com/owid/covid-19-data.git
cd /L03-NHOM2-9562/Data

# 2. Cài đặt R Runtime (>= 4.1.2) và chạy script phân tích chính
Rscript -e 'if(!require(pacman)) install.packages("pacman"); pacman::p_load(dplyr, tidyverse, cowplot, scales, matlib, here)'
Rscript -e 'source(here::here("function.R"))'

Hạn chế và hướng phát triển

  • Giới hạn mô hình tuyến tính đơn biến: Giả định hồi quy OLS bậc nhất giữa ca nhiễm và ca tử vong chưa phản ánh độ trễ sinh học (Time-lag Delay thông thường từ 14 đến 21 ngày kể từ khi nhiễm bệnh đến khi chuyển nặng).
  • Ảnh hưởng của biến số ẩn (Confounding Variables): Chưa tích hợp các yếu tố về tỷ lệ bao phủ vaccine, mật độ dân số, cơ cấu độ tuổi và năng lực xét nghiệm PCR vào mô hình tương quan.
  • Hướng nâng cấp đề xuất:
    • Tích hợp mô hình chuỗi thời gian ARIMA/SARIMA và kiến trúc Deep Learning LSTM để dự báo dài hạn.
    • Xây dựng mô hình phân tích hồi quy đa biến (Multivariate Regression) kết hợp phân tích thành phần chính (PCA) để giảm chiều không gian dữ liệu dịch tễ.

Đối tượng hưởng lợi

  • Sinh viên ngành Khoa học Máy tính & Dữ liệu: Nắm vững phương pháp ứng dụng cấu trúc rời rạc, xác suất thống kê và đại số tuyến tính vào việc giải quyết bài toán dữ liệu thực tế bằng ngôn ngữ R.
  • Kỹ sư dữ liệu (Data Engineers): Tham khảo kiến trúc Pipeline tiền xử lý, làm sạch và chuẩn hóa dữ liệu chuỗi thời gian quy mô lớn.
  • Nhà nghiên cứu Y tế Công cộng: Sử dụng bộ chỉ số thống kê mô tả ($Q_1, Q_2, Q_3$, SMA, RLE Outbreak Detection) để định lượng chính xác mức độ khốc liệt của từng làn sóng dịch bệnh.

Câu hỏi thường gặp

1. Tại sao phải dùng Trung vị ($Q_2$) thay cho Trung bình cộng (Mean) khi đánh giá ca nhiễm COVID-19?

Do số liệu ca nhiễm có độ phân tán cực kỳ lớn và xuất hiện nhiều giá trị ngoại lai (Outliers) trong các đợt bùng phát. Giá trị kỳ vọng (Mean) rất nhạy cảm với các Outlier này và sẽ bị kéo lệch lên cao, không phản ánh đúng mức độ lây nhiễm trung tâm của đa số các ngày trong năm.

2. Ý nghĩa cốt lõi của việc tính toán đường trung bình động 7 ngày (7-day SMA) là gì?

7-day SMA giúp triệt tiêu yếu tố "nhiễu chu kỳ tuần" (các ngày cuối tuần xét nghiệm và báo cáo ít hơn ngày trong tuần), làm phẳng các đột biến ngắn hạn và giúp các mô hình hồi quy tuyến tính đạt hệ số tương quan $r$ chính xác hơn.

3. Làm thế nào để giải quyết trường hợp ma trận $X^T X$ không khả nghịch trong bài toán OLS?

Trong bài toán thực nghiệm trên, ma trận thiết kế $X$ gồm cột vector hằng số 1 và vector đặc trưng độc lập tuyến tính $x$, đảm bảo $\det(X^T X) \neq 0$. Nếu xảy ra suy biến đa cộng tuyến, phương pháp chính quy hóa Ridge Regression ($X^T X + \lambda I$) sẽ được áp dụng.

4. Thuật toán RLE phát hiện đợt bùng phát hoạt động như thế nào?

RLE chuyển đổi một chuỗi boolean (TRUE/FALSE dựa trên điều kiện $x_i \ge k$) thành các cặp (value, length). Từ đó, hệ thống lọc ra tất cả các đoạn liên tục có value == TRUE và trích xuất trực tiếp ngày bắt đầu, ngày kết thúc và độ dài đợt dịch.

5. Yêu cầu phần cứng tối thiểu để thực thi trọn vẹn source code của đồ án là gì?

Môi trường yêu cầu tối thiểu CPU 2 Cores, 4GB RAM khả dụng và không gian ổ đĩa trống 500MB để lưu trữ file dữ liệu CSV và các biểu đồ xuất bản định dạng phân giải cao từ ggplot2.

Kết luận

Đồ án "Phân Tích Dữ Liệu Covid-19: Thống Kê và Kết Quả Nghiên Cứu" đã xây dựng thành công một quy trình toàn diện từ tiếp nhận, tiền xử lý, chuẩn hóa cấu trúc đến phân tích thống kê chuyên sâu trên tập dữ liệu dịch tễ học thực tế với hơn 150,000 bản ghi. Bằng việc kết hợp chặt chẽ giữa nền tảng lý thuyết Cấu trúc rời rạc, Xác suất thống kê, Đại số tuyến tính ma trận và năng lực lập trình R chuyên nghiệp, nhóm nghiên cứu đã làm sáng tỏ các quy luật biến động của đại dịch tại khu vực Châu Đại Dương, chứng minh vai trò không thể thay thế của Khoa học dữ liệu trong việc giải quyết các thách thức y tế toàn cầu.