Kiểm định xu thế với Mann-Kendall & Sen's Slope
Tiêu chuẩn phân tích xu hướng môi trường, thủy văn và khí hậu từ USGS & WMO
Vì sao Hồi quy OLS thất bại với dữ liệu thực địa?
Khi theo dõi lượng mưa, nhiệt độ hay nồng độ ô nhiễm theo thời gian, phản xạ đầu tiên thường là vẽ đường thẳng xu hướng bằng Hồi quy bình phương tối thiểu (OLS). Tuy nhiên, chuỗi số liệu môi trường thực tế hiếm khi lý tưởng như sách giáo khoa.
3 nguyên nhân khiến OLS đưa ra kết luận thiếu chính xác:
- Ngoại lai cực trị (Outliers): Chỉ cần 1 năm bão lịch sử hoặc 1 trận hạn hán bất thường, đường thẳng OLS sẽ bị bẻ cong hoàn toàn.
- Không tuân theo phân phối chuẩn: Số liệu dòng chảy, lượng mưa thường bị lệch phải nặng hoặc chứa nhiều giá trị 0.
- Ép buộc tuyến tính: Thực tế thiên nhiên biến đổi theo dạng đơn điệu (Monotonic - tăng liên tục hoặc giảm liên tục) chứ không bắt buộc phải thẳng hàng tăm tắp.
| Tiêu chí so sánh | Hồi quy Tuyến tính (OLS) | Mann-Kendall & Sen's Slope |
|---|---|---|
| Phương pháp luận | Tham số (Parametric) | Phi tham số (Non-parametric) |
| Yêu cầu phân phối | Phần dư chuẩn, phương sai đồng nhất | Tự do phân phối (Phân phối bất kỳ) |
| Độ nhạy với Outliers | Rất cao (Điểm phá vỡ $\approx 0\%$) | Bền bỉ (Điểm phá vỡ đạt tới $\approx 29.3\%$) |
| Dạng xu hướng phát hiện | Đường thẳng cố định ($y = ax + b$) | Xu hướng đơn điệu (Tăng/giảm đều, cong hay thẳng đều nhận) |
Nguyên lý chấm điểm từng cặp dữ liệu
Được thiết lập bởi Mann (1945) và Kendall (1975), kiểm định này chỉ đặt một câu hỏi trực giác: "Càng tiến về các năm sau, số liệu có xu hướng cao hơn các năm trước đó hay không?"
1. Hệ giả thuyết thống kê
- $H_0$ (Giả thuyết không): Không có xu hướng đơn điệu (dữ liệu độc lập và phân phối ngẫu nhiên theo thời gian).
- $H_1$ (Giả thuyết đối): Dữ liệu tồn tại một xu hướng đơn điệu tăng hoặc giảm thực sự.
2. Cách tính thống kê điểm số $S$
Cho chuỗi thời gian gồm $n$ mốc quan trắc: $x_1, x_2, \dots, x_n$. Ta so sánh tất cả các cặp $(x_k, x_j)$ thỏa mãn mốc thời gian sau đứng sau mốc trước ($j > k$):
Hàm dấu $\text{sgn}(\theta)$ hoạt động như một trọng tài ghi điểm:
- $\mathbf{+1}$ nếu $x_j - x_k > 0$ (thời điểm sau cao hơn thời điểm trước $\to$ tín hiệu tăng).
- $\mathbf{0}$ nếu $x_j - x_k = 0$ (hai thời điểm bằng nhau).
- $\mathbf{-1}$ nếu $x_j - x_k < 0$ (thời điểm sau tụt xuống $\to$ tín hiệu giảm).
Tổng số cặp so sánh có thể ghép được là $N = \frac{n(n-1)}{2}$.
3. Phương sai $\text{Var}(S)$ & Xử lý mẫu hòa (Tied values)
Khi $n \ge 10$, thống kê $S$ tiệm cận chuẩn với kỳ vọng $E(S) = 0$. Phương sai của $S$ được tính theo công thức chuẩn hóa của Gilbert (1987) có tính đến các nhóm dữ liệu bị bằng nhau:
(Trong đó $g$ là số nhóm có giá trị trùng lặp, $t_p$ là số lượng phần tử trùng nhau trong nhóm thứ $p$. Nếu toàn bộ số liệu không có số nào trùng, cụm trừ phía sau bằng 0).
4. Chuẩn hóa $Z_{MK}$ & Đưa ra quyết định
Để kiểm định ý nghĩa thống kê, $S$ được chuẩn hóa thành giá trị $Z_{MK}$ kèm hiệu chỉnh liên tục ($\pm 1$):
Quy tắc quyết định ở mức ý nghĩa $\alpha = 0.05$ (độ tin cậy 95%):
• Nếu $|Z_{MK}| > 1.96$ (tương đương $p\text{-value} < 0.05$): Bác bỏ $H_0$, kết luận **xu hướng có ý nghĩa thống kê**.
• Dấu của $Z_{MK}$ cho biết chiều hướng: $Z_{MK} > 0$ là tăng, $Z_{MK} < 0$ là giảm.
Ước lượng độ dốc Sen (Theil-Sen Estimator)
Kiểm định MK chỉ khẳng định xu hướng có tồn tại hay không. Để trả lời câu hỏi: "Mỗi năm đại lượng này tăng hoặc giảm bình quân bao nhiêu đơn vị?", ta dùng phương pháp phi tham số của Sen (1968).
Quy trình 3 bước:
- Tính tốc độ thay đổi giữa toàn bộ các cặp điểm ($j > k$):
$$Q_i = \frac{x_j - x_k}{j - k}$$
- Với $n$ điểm dữ liệu, ta thu được danh sách gồm $N = \frac{n(n-1)}{2}$ giá trị độ dốc $Q_i$.
- Độ dốc Sen ($\beta$) chính là Trung vị (Median) của toàn bộ danh sách $Q_i$:
$$\beta = \text{Median}(Q_1, Q_2, \dots, Q_N)$$
Điểm phá vỡ (Breakdown Point) đạt $29.3\%$: Nghĩa là gần $30\%$ quan trắc trong chuỗi có thể là ngoại lai cực đoan (như máy hỏng, lũ quét đột biến) mà độ dốc ước lượng vẫn phản ánh đúng thực tế của chuỗi nền.
Thực hành từng bước với chuỗi 5 năm
Giả sử theo dõi nhiệt độ mùa hè qua 5 năm liên tiếp: $x = [20, 21, 21, 24, 26]^\circ\text{C}$ với $n=5$.
Số cặp so sánh: $N = \frac{5 \times 4}{2} = 10$ cặp.
| Cặp $(k, j)$ | Mốc năm | Hiệu số $(x_j - x_k)$ | Điểm $\text{sgn}$ | Độ dốc cặp $Q = \frac{x_j - x_k}{j - k}$ |
|---|---|---|---|---|
| (1, 2) | Năm 1 $\to$ Năm 2 | $21 - 20 = 1$ | $+1$ | $1 / 1 = 1.00$ |
| (1, 3) | Năm 1 $\to$ Năm 3 | $21 - 20 = 1$ | $+1$ | $1 / 2 = 0.50$ |
| (1, 4) | Năm 1 $\to$ Năm 4 | $24 - 20 = 4$ | $+1$ | $4 / 3 \approx 1.33$ |
| (1, 5) | Năm 1 $\to$ Năm 5 | $26 - 20 = 6$ | $+1$ | $6 / 4 = 1.50$ |
| (2, 3) | Năm 2 $\to$ Năm 3 | $21 - 21 = 0$ | $0$ | $0 / 1 = 0.00$ |
| (2, 4) | Năm 2 $\to$ Năm 4 | $24 - 21 = 3$ | $+1$ | $3 / 2 = 1.50$ |
| (2, 5) | Năm 2 $\to$ Năm 5 | $26 - 21 = 5$ | $+1$ | $5 / 3 \approx 1.67$ |
| (3, 4) | Năm 3 $\to$ Năm 4 | $24 - 21 = 3$ | $+1$ | $3 / 1 = 3.00$ |
| (3, 5) | Năm 3 $\to$ Năm 5 | $26 - 21 = 5$ | $+1$ | $5 / 2 = 2.50$ |
| (4, 5) | Năm 4 $\to$ Năm 5 | $26 - 24 = 2$ | $+1$ | $2 / 1 = 2.00$ |
1. Tính thống kê $S$:
$$S = (+1) + (+1) + (+1) + (+1) + 0 + (+1) + (+1) + (+1) + (+1) + (+1) = \mathbf{+9}$$
2. Tính độ dốc Sen:
Sắp xếp 10 giá trị $Q$ theo thứ tự tăng dần:
Trung vị là giá trị trung bình giữa phần tử thứ 5 và thứ 6: $$\beta = \frac{1.50 + 1.50}{2} = \mathbf{1.5^\circ\text{C}/\text{năm}}$$
Xem xét cẩn thận 2 vấn đề sau đây:
1. Hiện tượng Tự tương quan chuỗi (Serial Autocorrelation)
Nếu số liệu hôm nay phụ thuộc vào hôm qua (ngày nắng thì xác suất ngày mai vẫn nắng), giả định độc lập của kiểm định MK bị phá vỡ. Hiện tượng này làm thu hẹp phương sai $\text{Var}(S)$, dẫn đến dương tính giả (kết luận có xu hướng trong khi thực tế không có).
Giải pháp: Áp dụng thuật toán Trend-Free Pre-Whitening (TFPW) của Yue et al. (2002) hoặc kiểm định Modified Mann-Kendall (Hamed & Rao, 1998) để triệt tiêu hệ số tự tương quan trễ bậc 1 trước khi chạy MK.
2. Tính chu kỳ và mùa vụ (Seasonality)
Với dữ liệu thu thập theo tháng hoặc quý, nhiệt độ mùa hè luôn cao hơn mùa đông. Nếu lấy tháng 7 so sánh với tháng 1 trước đó sẽ làm sai lệch bản chất kiểm định.
Giải pháp: Áp dụng Seasonal Mann-Kendall Test (Hirsch et al., 1982 - chuẩn USGS): Thuật toán chỉ thực hiện so sánh các tháng cùng tên giữa các năm với nhau (tháng 1 năm nay chỉ so với tháng 1 các năm trước).
VÍ DỤ MINH HOẠ
Chọn các kịch bản thực tế bên dưới hoặc chỉnh sửa trực tiếp chuỗi số liệu:
Hãy chú ý kịch bản "Ngoại lai cực đoan": Đường OLS (màu đỏ) bị điểm dị thường kéo dốc ngược lên trên, trong khi Sen's Slope (màu xanh lá) giữ nguyên xu thế nền của dữ liệu.
Câu hỏi: Tại sao kiểm định Mann-Kendall lại được gọi là kiểm định "phi tham số"?
Bạn đã nắm được phần lý thuyết? Chuyển đến công cụ thực hành, chạy dữ liệu của bạn ngay trên trình duyệt