Mục tiêu học tập

• Nắm rõ từng giai đoạn trong quy trình xử lý của DESeq2
• Xác định nguồn biến thiên và kiểm tra hệ số chuẩn hóa (size factors)
• Đánh giá độ phân tán gen và vai trò của nó trong mô hình
• Hiểu tại sao độ phân tán là yếu tố then chốt trong phân tích biểu hiện khác biệt

Quy trình xử lý của DESeq2

Khi tạo đối tượng DESeqDataSet và gọi hàm DESeq(), toàn bộ quy trình thống kê đã được thực thi:

dds_obj <- DESeqDataSetFromTximport(tx_data, colData = sample_meta, design = ~ condition)
dds_result <- DESeq(dds_obj)

Dưới đây là các bước chính diễn ra bên trong:

1. Hệ số chuẩn hóa (Size Factors)

DESeq2 tự động tính toán hệ số chuẩn hóa để điều chỉnh sự chênh lệch về chiều sâu giải trình tự giữa các mẫu. Phương pháp mặc định là sử dụng trung vị tỷ lệ (median-of-ratios).

Ví dụ với dữ liệu MOV10:

# Xem hệ số chuẩn hóa
sizeFactors(dds_obj)

Irrel_kd_1 Irrel_kd_2 Irrel_kd_3 Mov10_kd_2 Mov10_kd_3 Mov10_oe_1 Mov10_oe_2 
 1.1149694  0.9606733  0.7492240  1.5633640  0.9359695  1.2262649  1.1405026 
Mov10_oe_3 
 0.6542030 

Hệ số lớn hơn 1 tương ứng với mẫu có tổng số reads cao hơn — khi chuẩn hóa, ta chia số đếm thô cho hệ số này để cân bằng.

# Tổng reads thô mỗi mẫu
colSums(counts(dds_obj))

# Tổng reads sau chuẩn hóa
colSums(counts(dds_obj, normalized = TRUE))

Lưu ý: Giá trị sau chuẩn hóa không hoàn toàn bằng nhau do DESeq2 còn điều chỉnh theo thành phần RNA — tránh bị ảnh hưởng bởi một số ít gen biến thiên mạnh.

2. Độ phân tán theo gen (Gene-wise Dispersion)

Trong dữ liệu RNA-seq, phương sai không bằng trung bình — đặc biệt, gen có mức biểu hiện thấp thường có phương sai dao động quanh trung bình nhưng rất thất thường, trong khi gen biểu hiện cao có phương sai ổn định hơn nhưng lớn hơn trung bình.

Độ phân tán (dispersion) đo lường mức độ phương sai lệch khỏi kỳ vọng Poisson (nơi phương sai = trung bình). Trong thực tế sinh học, dispersion luôn nhỏ hơn 1.

Công thức mô hình:

phương_sai = trung_bình + trung_bình² × độ_phân_tán

Suy ra:

độ_phân_tán = (phương_sai - trung_bình) / trung_bình²

DESeq2 ước lượng độ phân tán riêng cho từng gen dựa trên trung bình biểu hiện và phương sai quan sát được từ các mẫu lặp lại.

3. Khớp đường cong chuẩn (Fitting the Dispersion Curve)

Vì số mẫu lặp lại thường ít (3–6), ước lượng độ phân tán cho từng gen riêng lẻ không đáng tin cậy. DESeq2 giải quyết bằng cách khớp một đường cong xu hướng chung — giả định rằng các gen có mức biểu hiện tương đương nên có độ phân tán tương tự.

Đường cong này (thường là màu đỏ) đại diện cho kỳ vọng về độ phân tán tại mỗi mức biểu hiện trung bình.

4. Co rút độ phân tán (Shrinkage)

DESeq2 "co" các ước lượng độ phân tán riêng lẻ về phía đường cong chuẩn để tăng độ chính xác. Mức độ co phụ thuộc vào:

• Khoảng cách giữa ước lượng gen và đường cong
• Số lượng mẫu lặp lại (nhiều mẫu → ít co hơn)

Việc co giúp giảm tỷ lệ dương tính giả — đặc biệt quan trọng khi kích thước mẫu nhỏ. Tuy nhiên, các gen có độ phân tán cực cao (do nhiễu kỹ thuật hoặc sinh học bất thường) sẽ không bị co, nhằm tránh làm sai lệch kết quả.

Cảnh báo: Nếu đồ thị phân tán trông như đám mây rời rạc, hoặc độ phân tán tăng khi biểu hiện tăng — đó là dấu hiệu dữ liệu không phù hợp mô hình. Cần kiểm tra lại chất lượng mẫu, ô nhiễm (vd: mtRNA), hoặc mẫu ngoại lai.

5. Kiểm tra với dữ liệu MOV10

# Vẽ đồ thị phân tán
plotDispEsts(dds_obj)

Ở dữ liệu này, ta thấy:

• Độ phân tán giảm rõ rệt khi biểu hiện trung bình tăng → phù hợp kỳ vọng
• Các điểm đen (ước lượng thô) nằm rải quanh đường cong → mô hình hợp lý
• Có sự co mạnh do chỉ có 2–3 mẫu lặp lại mỗi nhóm → nên tăng số mẫu lặp (≥4) để cải thiện độ tin cậy