MỤC LỤC
MỞ ĐẦU .1
Chƣơng 1. TỔNG QUAN .3
1.1. Cơ sở khoa học của dự báo khí hậu mùa .3
1.2. Các nghiên cứu trên thế giới.8
1.3. Các nghiên cứu ở trong nước.13
Chƣơng 2. MÔ HÌNH ARIMA VÀ SỐ LIỆU SỬ DỤNG .21
2.1. Giới thiệu cấu trúc của mô hình ARIMA .21
2.1.1. Mô hình tự hồi quy trung bình trượt ARIMA .22
2.1.2. Mô hình động thái ARIMAX.23
2.2. Phương pháp áp dụng mô hình ARIMA và ARIMAX đối với bài toán dự
báo mưa mùa.24
2.2.1. Xác định tính ổn định ngẫu nhiên của chuỗi thời gian .25
2.2.2. Nhận dạng cấu trúc của mô hình .28
2.2.3. Xác định các tham số của mô hình.32
2.2.4. Kiểm định mô hình .35
2.2.5. Phần mềm thống kê SAS đối với mô hình ARIMA và ARIMAX .36
2.3. Các nguồn số liệu được sử dụng.36
2.3.1. Số liệu quan trắc mưa từ các trạm khí tượng.37
2.3.2. Số liệu về các chỉ số khí hậu .37
2.3.3. Số liệu về số vết đen mặt trời (Sunspot Number) .38
2.3.4. Xử lý số liệu.38
Chƣơng 3. KẾT QUẢ VÀ NHẬN XÉT.41
3.1. Xây dựng mô hình dự báo mưa vụ đông xuân bằng mô hình ARIMA .41
3.1.1. Xác định tính ổn định của chuỗi lượng mưa vụ Đông xuân.41
3.1.2. Nhận dạng mô hình ARIMA .43
3.1.3. Xác định các thông và kiểm định mô hình ARIMA .44
3.2. Xây dựng mô hình dự báo lượng mưa vụ đông xuân bằng mô hình động
thái ARIMAX .46
3.2.1. Xác định tính ổn định của chuỗi nhân tố dự báo .46
81 trang |
Chia sẻ: mimhthuy20 | Lượt xem: 989 | Lượt tải: 1
Bạn đang xem trước 20 trang tài liệu Luận văn Nghiên cứu ứng dụng mô hình arima để dự báo lượng mưa vụ đông xuân ở một số tỉnh vùng đồng bằng bắc bộ, để xem tài liệu hoàn chỉnh bạn click vào nút DOWNLOAD ở trên
h trượt bậc q (q = 1, 2,) và được ký hiệu là
ARIMA(p,d,q). Dạng tổng quát của mô hình ARIMA(p,d,q) có thể được viết như
sau [20]:
Wt = µ + p1Wt-l + p2Wt-2 ++ ppWt-p - q1at-1 - q2at-2 -- qqat-q + at (2.1)
Trong đó:
Wt = Δ
d
yt
d là bậc sai phân, μ là hằng số
Với d = 0 Wt = yt ; với d = 1 Δyt = yt - yt-1;
yt, yt-l, yt-2, , yt-p là giá trị quan trắc ở các bước thời gian t, t-1, t-2,, t-p
at, at-1, at-2, , at-q là sai số ngẫu nghiên (giữa giá trị thực và giá trị tính toán)
ở các bước thời gian t, t-1, t-2,, t-q;
p1, p2, , pp ; q1, q3, ..., qq là các tham số hồi quy.
Phương trình 2.1 cũng có thể viết gọn lại thông qua phép toán dịch chuyển lùi
23
tt a
Bp
Bq
W
)(
)(
hoặc p(B)(1-B)d yt = μ + q(B)at (2.2)
Trong đó:
yt, at như đã trình bày ở trên
B là phép tính dịch chuyển lùi: BWt = Wt-1 hay B
k
Wt = Wt-k
p(B) = (1 – p1B – p2B
2
– – ppB
p) là phép toán tự hồi quy
q(B) = (1 – q1B – q2B
2
– – qqB
q) là phép toán trung bình trượt
2.1.2. Mô hình động thái ARIMAX
Mô hình động thái ARIMAX có sự khác biệt cơ bản so với mô hình tự hồi
quy trung bình trượt ARIMA là ngoài việc xem xét quá trình tự hồi quy trung bình
trượt của chuỗi yếu tố dự báo, nó còn cho phép xem xét ảnh hưởng của các chuỗi
thời gian khác tác động đến yếu tố dự báo, chuỗi tác động (biến độc lập) được gọi là
chuỗi nhập, chuỗi bị tác động (biến phụ thuộc) được gọi là chuỗi xuất.
Giả sử ta có các chuỗi độc lập Xit (i = 1,2m; t=1,2n) và chuỗi phụ thuộc
Yt (t=1,2n), khi đó mô hình động thái ARIMAX được viết dưới dạng tổng quát
như sau:
tti
ki
m
i
r
i
s
i
t a
Bp
Bq
XB
BS
BU
Y
)(
)(
)(
)(
,
1
(2.3)
Trong đó:
Yt là giá trị quan trắc ở các bước thời gian t; μ là hằng số;
B là phép toán dịch chuyển lùi theo quy tắc : BXt = Xt-1 , B
k
Xt = Xt-k ;
s
isi1i0 B U B U U)( BU
s
i ;
s
isi1i0 BS BS S)( BS
s
i là
những trọng số động thái của chuỗi độc lập thứ i;
k là thời điểm tác động của chuỗi độc lập thứ i tại thời điểm t = k;
p(B) = (1 – p1B – p2B
2
– – ppB
p
); q(B) = (1 – q1B – q2B
2
– – qqB
q
) là phép
toán tự hồi quy và trung bình trượt của chuỗi phụ thuộc;
24
at, là sai số ngẫu nghiên (giữa giá trị thực và giá trị tính toán).
Lưu ý : các chuỗi Xit và Yt trong công thức 2.3 phải là các chuỗi có tính ổn
định ngẫu nhiên, nếu chuỗi không ổn định, sẽ cần phải thông qua bước sai phân để
đưa chuỗi về dạng ổn định ngẫu nhiên.
Tóm lại: Bản chất của các mô hình ARIMA và ARIMAX là mô hình ngẫu
nhiên. Việc phân tích chuỗi thời gian trong các mô hình này bắt buộc phải chấp
nhận một giả thiết hết sức cơ bản là tính ổn định của các quá trình ngẫu nhiên, tính
ổn định ở đây có nghĩa là các đặc trưng thống kê (hay phân phối xác suất) không
thay đổi theo thời gian. Trong thực tế nhiều quá trình ngẫu nhiên có tính ổn định
trong một khoảng thời gian gián đoạn hữu hạn nào đó có thể coi là ổn định. Ví dụ
chuỗi tổng lượng mưa tháng là chuỗi không dừng, còn chuỗi tổng lượng mưa năm
có thể coi là dừng vì khi đó qui luật bên trong năm bị loại trừ. Các chuỗi không
dừng có thể trở thành dừng nhờ một số phép biến đổi sai phân. Lợi thế cơ bản của
các mô hình này là cho phép dự báo với độ chính xác nhất định, mặc dù chưa hiểu
rõ bản chất của các quá trình tác động từ các nhân tố dự báo đến yếu tố dự báo.
2.2. Phƣơng pháp áp dụng mô hình ARIMA và ARIMAX đối với bài
toán dự báo mƣa mùa
Trong mục 2.1 đã trình bày các dạng tổng quát của mô hình ARIMA và
ARIMAX, nó có thể bao gồm nhiều thành phần tham gia vào mô hình như: thành
phần tự hồi quy, thành phần trung bình trượt, thành phần sai phân, thành phần ảnh
hưởng của các chuỗi nhập khác (các chuỗi nhân tố dự báo), trong mỗi thành phần
lại có các thành phần con khác nhau. Bài toán cần giải quyết ở đây là đưa ra được
phương pháp xác định các thành phần có ý nghĩa về mặt thống kê để tham gia vào
mô hình dự báo mưa hạn mùa. Đây là bài toán khá phức tạp, độ chính xác của mô
hình dự báo không chỉ phụ thuộc vào các chuỗi nhập, chuỗi xuất mà còn phụ thuộc
việc lựa chọn chính xác các thành phần tham gia vào mô hình dự báo. Để giải quyết
bài toán này luận văn đã thực hiện theo các bước sau:
25
1) Áp dụng phương pháp thống kê sai phân để xác định tính ổn định ngẫu nhiên
của các chuỗi dữ liệu tham gia vào mô hình ARIMA và ARIMAX;
2) Kế thừa phương pháp Box Jenkin đối với mô hình ARIMA và phương pháp
Box Tao đối với mô hình ARIMAX trong việc nhận dạng các thành phần tự
hồi quy, thành phần trung bình trượt và thành phần ảnh hưởng của các chuỗi
nhập đến chuỗi lượng mưa thông qua việc xem xét sự biến đổi các hàm tự
tương quan, tự tương quan riêng phần và tương quan chéo;
3) Sử dụng phương pháp bình phương tối thiểu trong việc xác định các tham số
trong mô hình ARIMA và ARIMAX;
4) Áp dụng các phương pháp kiểm nghiệm giả thiết thống kê trong khí hậu để
chọn lựa các tham số có đủ độ tin cậy thống kê tham gia trong mô hình
ARIMA và ARIMAX;
5) Sử dụng công cụ phần mềm thống kê SAS để tính toán các đặc trưng của
chuỗi thời gian và các tham số trong mô hình ARIMA và ARIMAX.
Sau đây sẽ trình bày cụ thể từng nội dung này:
2.2.1. Xác định tính ổn định ngẫu nhiên của chuỗi thời gian
Chuỗi thời gian là chuỗi số liệu được sắp xếp theo trình tự thời gian. Nếu
một chuỗi thời gian có giá trị trung bình và phương sai không đổi theo thời gian thì
chuỗi đó được xem là ổn định ngẫu nhiên (chuỗi có tính dừng) hay nói một cách
khác cụ thể hơn đó là một chuỗi thời gian không có xu thế, không có chu kỳ, mà chỉ
dao động xung quanh kỳ vọng của nó.
Một chuỗi quan trắc khí hậu trung bình tháng thường bao gồm 3 thành phần:
1) thành phần ngẫu nhiên là sự tăng lên hay giảm đi thường xen kẽ nhau, góp phần
làm cho các trị số khí hậu dao động xung quanh một giá trị nào đó. Giá trị đó có thể
là trung bình số học, nếu chuỗi không có thành phần chu kỳ và xu thế. 2) Thành
phần chu kỳ là những biến đổi của chuỗi lặp lại nhiều lần sau những khoảng thời
gian nhất định nào đó. Mối tương quan giữa các thành phần trong một chu kỳ
thường đạt trị số lớn nhất. 3) Thành phần xu thế là biểu hiện xu hướng tăng hoặc
26
giảm theo thời gian của các thành phần trong chuỗi, trị số đầu của xu thế là cực tiểu
hoặc cực đại và trị số cuối của xu thế là cực đại hoặc cực tiểu. Biểu đồ minh họa 3
thành phần này được trình bày trong hình 2.1.
Hình 2.1. Các thành phần trong chuỗi quan trắc khí hậu [10]
Để loại bỏ thành phần xu thế và chu kỳ nhằm đưa các chuỗi quan trắc về
dạng ổn định ngẫu nhiên, thường sử dụng phép lọc sai phân, phép lọc Loga, phép
lọc căn thức...[10]. Trong luận văn này chúng tôi chọn phép lọc sai phân, cụ thể như
sau:
- Đối với việc loại bỏ thành phần xu thế: sử dụng phép biến đổi sai phân bậc 1
hoặc bậc 2. Sai phân bậc 1 là chênh lệch giữa 2 giá trị kề nhau trong chuỗi.
ΔYt = Yt - Yt-1 (2.5)
Trong đó: ΔYt là giá trị của sai phân bậc 1
Yt và Yt-1
là các thời đoạn trước và thời đoạn sau đó.
a) b)
c) d)
27
Nếu sai phân bậc 1 vẫn còn thể hiện xu thế thì thực hiện tiếp sai phân bậc 2.
Sai phân bậc 2 chính là sai phân của sai phân bậc 1:
Δ2(Yt) = ΔYt - ΔYt-1 = (Yt - Yt-1) - (Yt-1 - Yt-2) (2.6)
Nếu sai phân bậc 2 chưa đạt được tính dừng ta có thể tiếp tục lấy sai phân
bậc 3 hoặc cao hơn.
- Đối với việc loại bỏ thành phần mùa và chu kỳ: Sai phân mùa là chênh lệch giá
trị của hai quan trắc cách nhau khoảng thời gian L, L có thể là một năm, hai năm
hay một mùa Ví dụ : nếu là số liệu tổng lượng mưa tháng, ta có L =12. Do đó sai
phân mùa bậc 1 có tính mùa là:
ΔYt = Yt - Yt-L = Yt - Yt-12 (2.7)
Cũng có thể lấy sai phân bậc 2 của sai phân mùa bậc 1 khi chuỗi chưa đạt
được độ ổn định:
- Kiểm tra tính ổn định ngẫu nhiên của chuỗi
Trong thực hành, để kiểm tra các chuỗi thời gian tham ra trong mô hình
ARIMA hoặc ARIMAX đã đạt tiêu chuẩn ổn định ngẫu nhiên hay chưa, thường dựa
vào hàm tự tương quan. Theo Quenouille đã chứng minh chuỗi thời gian được xem
là ổn định ngẫu nhiên khi hầu hết hệ số tự tương quan của chuỗi (rk ) thỏa mản biểu
thức giới hạn tin cậy (2.8) và tiến dần về 0, ngoại trừ một số bước trễ như bước
mùa, vụ, chu kỳ , nằm ngoài khoảng này [6, 37].
Biểu thức giới hạn tin cậy có thể viết dưới dạng sau:
rnkrn StrSt 1,2/1,2/ (2.8)
n
Sr
1
(2.9)
Trong đó: Sr là sai số chuẩn của các hệ số tự tương quan rk ; 1,2/ nt là điểm
phần trăm =0.05 của phân bố Student với n-1 bậc tự do.
Hình 2.2 minh họa chuỗi dữ liệu tổng lượng mưa tháng trước khi sai phân
và sau khi sai phân. Trên hình này, phần hình A, A‟ ở phía trên ứng với trường hợp
28
chuỗi chưa ổn định, phần hình B, B‟ ứng với chuỗi sau khi sai phân và được xem là
chuỗi ổn định ngẫu nhiên.
Hình 2.2. Minh họa diễn biến của chuỗi lượng mưa tháng và hàm tự tương
quan đối với trạm Hà Nội trước khí sai phân (A,A‟) và sau khi sai phân (B,B‟).
2.2.2. Nhận dạng cấu trúc của mô hình
Sau khi đã loại bỏ được các thành phần chu kỳ, xu thế của chuỗi thời gian, sẽ
tiến hành nhận dạng cấu trúc của mô hình. Box - Jenkin đã đưa ra phương pháp
nhận dạng cấu trúc của mô hình ARIMA thông qua việc xem xét sự biến đổi của
hàm tự tương quan (Autocorrelation function - ACF) và tự tương quan riêng phần
(Part autocorrelation function - PAFC) để xác định các thành phần tự hồi quy (AR)
và thành phần trung bình trượt (MA). Đối với mô hình ARIMAX, Box Tao đã đưa
ra một số dáng điệu chính của hàm tương quan chéo (Cross correlation function -
CCF) để xác định mức độ ảnh hưởng (hàm truyền) của các chuỗi nhập đến chuỗi
yếu tố dự báo. Định nghĩa và thuật toán để tính các hàm ACF, PACF và CCF được
trình bày chi tiết trong tài liệu ARIMA [37].
29
Xuất phát từ bản chất của hàm ACF và PACF, Box - Jenkin đã đưa ra một số
dạng biểu đồ thường gặp đối với hàm ACF và PACF, tương ứng với nó là các dạng
của mô hình ARIMA nhằm hỗ trợ cho việc nhận dạng cấu trúc của mô hình, các
dạng biểu đồ này được trình bày trong hình 2.3. Chi tiết về cách áp dụng biểu đồ
này được trình bày dưới đây:
Nếu biểu đồ hàm ACF có dạng nhỏ dần theo các bước trễ thời gian và biểu
đồ hàm PACF chỉ có giá trị khác 0 tại bước thời gian t-1, sau đó giảm đột
ngột về 0, (cụm từ „khác 0‟ hay „bằng 0‟ ở đây được hiểu theo thuật ngữ
thống kê, nếu các giá trị này nằm trong khoảng từ - 1,2/ nt *Sr đến + 1,2/ nt * Sr
được xem là bằng 0, ngoài khoảng này được xem là khác 0 ) thì có một thông
số tự hồi qui (p=1) được chọn, mô hình có dạng ARIMA(1,0,0). Ngược lại,
khi biểu đồ hàm PACF tắt dần, hàm ACF có giá trị khác 0 bước t-1, sau đó
giảm đột ngột về 0, trong trường hợp này mô hình có dạng ARIMA(0,0,1).
Dáng điệu của hàm ACF và PACF trong các trường hợp này được minh họa
trên hình (Hình 2.3a).
Tương tự như trên, nếu biểu đồ hàm ACF tắt dần, hàm PACF có giá trị khác
0 ở các bước thời gian t-1, t-2, sau đó giảm đột ngột về 0 thì mô hình có hai
thông số tự hồi qui (p=2), mô hình có dạng ARIMA(2,0,0). Ngược lại, nếu
hàm PACF có dạng tắt dần, hàm ACF có giá trị khác 0 ở các bước thời gian
t-1, t-2, sau đó giảm đột ngột về 0 thì mô hình có hai thông số trung bình
trượt (q=2), mô hình có dạng ARIMA(0,0,2). Đồ thị của hạn ACF và PACF
trong các trường hợp này được minh họa trên hình (Hình 2.3b).
Khi biểu đồ hàm ACF có dạng tắt dần và có giá trị khác 0 ở các bước thời
gian t-1, t-2t-p, tương tự hàm hàm PACF có dạng tắt dần và có giá trị khác
0 ở các bước thời gian t-1, t-2t-q, trong trường hợp này cả 2 thành phần
AR và MA đều có trong mô hình, dạng của mô hình trong trường hợp này sẽ
là ARIMA(p,0,q). Đồ thị của hạn ACF và PACF trong các trường hợp này
được minh họa trên hình (Hình 2.3c).
30
(2.3a)
(2.3b)
(2.3c)
Hình 2.3 Một số dạng chính của hàm ACF và PACF tưng ứng với các dạng
mô hình ARIMA khác nhau [20]
Đối với việc xác định ảnh hưởng của từng chuỗi nhập đến yếu tố dự báo
trong mô hình ARIMAX, Box Tao cũng sử dụng phương pháp trực quan để xem xét
sự biến đổi của hàm tương quan chéo (Cross correlation function - CCF), từ đó đưa
31
ra hàm truyền tương ứng của chuỗi nhập X tham gia trong mô hình ARIMAX. Nội
dung của phương pháp có thể bao gồm 4 dạng chính sau:
Hình 2.4. Một số dạng chính của hàm tương quan chéo giữa biến nhập (X)
và biến phụ thuộc (Y) tưng ứng với các dạng mô hình ARIMA khác nhau [20]
1) Nếu hàm tương quan chéo (CCF) giữa biến độc lập X và biến phụ thuộc Y
có giá trị „khác 0‟ tại bước thời gian (t-b), sau đó giảm đột ngột, các bước
(2.4A)
(2.4B)
(2.4C)
(2.4D)
32
thời gian khác đều có giá trị „bằng 0‟. Khi đó hàm truyền của biến X tham
gia vào mô hình ARIMAX sẽ có dạng U0Xt-b, (hình 2.4A).
2) Nếu hàm CCF có giá trị „khác 0‟ tại bước thời gian (t-b) và (t-b-1), sau đó
giảm đột ngột, các bước thời gian khác đều có giá trị „bằng 0‟. Khi đó hàm
truyền của biến X tham gia vào mô hình ARIMAX sẽ có dạng (U0 +
U1B)Xt-b, (hình 2.4B).
3) Nếu hàm CCF có giá trị „khác 0‟ tại bước thời gian (t-b), sau đó có xu thế
giảm dần nhưng vẫn „khác 0‟ ở bước (t-b-1), còn các bước thời gian khác
có giá trị „bằng 0‟. Khi đó hàm truyền của biến X tham gia vào mô hình
ARIMAX sẽ có dạng sau : (hình 2.4C).
4) Nếu hàm CCF có giá trị „khác 0‟ và tại bước thời gian (t-b-1) và (t-b), đạt
cao nhất tại (t-b-1), sau đó có xu thế giảm dần nhưng vẫn „khác 0‟ ở bước
(t-b-2), còn các bước thời gian khác có giá trị „bằng 0‟. Khi đó hàm truyền
của biến X tham gia vào mô hình ARIMAX sẽ có dạng sau : (hình 2.4D).
2.2.3. Xác định các tham số của mô hình
Như đã trình bày ở trên (mục 2.2.2), tùy thuộc vào đặc tính và mối quan hệ của
các chuỗi nhập và chuỗi xuất, mô hình ARIMA và mô hình ARIMAX có thể có một
trong các dạng chính sau:
Khi chỉ có thành phần AR(p), mô hình sẽ có dạng:
yt = p1yt-l +p2yt-2 + +ppyt-p + at (2.10)
Khi chỉ có thành phần MA(q), mô hình sẽ có dạng:
yt = µ + q1at-1 + q2at-2 ++ qqat-q + at (2.11)
Khi có đủ cả 2 thành phần AR và MA, mô hình sẽ có dạng:
yt = µ - p1yt-l - p2yt-2 -- ppyt-p + q1at-1 + q2at-2 ++ qqat-q + at (2.12)
btX
BS
BUU
1
10
1
btX
BS
U
1
0
1
33
Khi có đủ tất cả các thành phần AR, MA và thành phần động thái X1, X2,,
Xm mô hình ARIMAX có thể có dạng:
yt = µ - p1yt-l - p2yt-2 -- ppyt-p + q1at-1 + q2at-2 ++ qqat-q +
(Ui1B+ Ui2B
2
++ UihB
h
)Xit-b + at;
(2.13)
trong đó: (i = 1,2,,m; h = 1,2,3,).
Việc xác định các tham số hồi quy trong các mô hình này được dựa theo
nguyên tắc bình phương tối thiểu, trong đó các hệ số hồi quy được xác định sao cho
tổng bình phương độ lệch giữa giá trị thực và giá trị mô phỏng là nhỏ nhất. Chi tiết
về phương pháp này được trình bày trong sách giáo trình [10]. Tuy về nguyên tắc
việc xác định các tham số trong mỗi thành phần của từng mô hình là giống nhau,
nhưng về cách tính cụ thể có nhiều cách xử lý riêng. Do vậy ở đây sẽ trình bày tóm
tắt cho một số dạng chính sau.
1) Xác định các thông số pi khi chỉ có thành phần AR(p) tham ra vào mô hình
Đây là dạng hàm tương tự như hồi quy tuyến tính nhiều biến, Tuy nhiên các
hệ số pi phải thỏa mản tiêu chuẩn hội tụ |p1 + pi ++ pi|<1 và được xác giải thông
qua hệ phương trình Yule_Walker [6]
Ck = p1Ck-l + p2Ck-2 + + ppCk-p; với k =1, 2,,p (2.14)
Ck là các mô men tương quan (Covarian - hiệp phương sai) giữa biến phụ
thuộc (yt) và các biến độc lập (yk-p).
Ck-l , Ck-2 , Ck-p là các mô men tương quan giữa các biến độc lập với nhau.
Đây là hệ p phương trình với các ẩn số là p1, p2 ,,pp và có thể giải được
theo nhiều phương pháp khác nhau như: dùng phương pháp ma trận nghịch đảo,
phương pháp khử Gauss, phương pháp Cramer
2) Xác định các thông số qi khi chỉ có thành phần MA(q) tham ra vào mô hình
Các thông số qi thoả mãn một hệ phương trình tương tự như hệ
Yule_Walker, được suy ra từ quan hệ:
22
2
2
1
11
0
k
...1
...
q
qkpkk
k
qqq
qqqqq
r
(2.15)
Trong đó : rk là hệ số tự tương quan của chuỗi yt,
34
γk và γ0 là các mô men tương quan bậc k và bậc 0 (k=0) của chuỗi
Cho k = 1,2,...,q ta được một hệ phương trình phi tuyến. Phương pháp để giải
hệ phương trình này được trình bày trong sách giáo trình [6,10].
3) Xác định các thông số pi, qi khi cả 2 thành phần AR(p) và MA(q) tham ra vào mô
hình
- Đối với thành phần AR(p): để giải quyết độc lập các giá trị pi theo công thức truy
hồi Durbin hay hệ phương trình Yule_Walker thì hệ thức (2.15) được viết bắt đầu
từ k > q tức là khi ấy các bi = 0. Do đó ta có hệ :
Cp+1 = p1Cp + p2Cp-2 + + ppCp+1-p
Cp+2 = p1Cp+1 + p2Cp+ + + ppCp+2-p
Cp+p = p1Cp+p-1 + p2Cp+p-2 + + ppCp
(2.16)
Hệ phương trình (2.26) cũng là hệ tuyến tính bậc nhất, do vậy có thể giải ra
tìm các nghiệm p1, p2,, pp.
- Đối với thành phần MA(q): để tìm các hệ số qi ta cũng xuất phát từ quan hệ (2.15),
nhưng khác với mô hình MA(q), ở đây γk và γ0 (hay Ck và C0) không phải là của yt
mà là của thành phần ngẫu nhiên at. Nghĩa là γk=γka, Ck=Cka. Quan hệ giữa γka và γky
hay Cka và Cky có dạng [6]:
ip
h
kih
p
i
ky
p
i
ika daaCaC
010
2
(2.17)
Với dk = Ck+i +Ck-i
Cky chính là mô men tương quan bậc k của y và có thể tính theo công thưc:
Cky = [(yt – ytb)(yt+k – ytb)] ; t = 1,2,...,n ; k = 1,2,,m ; (m < n)
Sau khi có Cky, sẽ tính được Cka theo (2.17), thay vào (2.15) được một hệ phương
trình, giải hệ này ta được các hệ số qi.
4) Xác định các thông số động thái Ui,t-b đối với các biến nhập Xi,t-b và các thông số
pi, qi trong mô hình động thái ARIMAX
35
Đối với các thông số Ui,t-b tương ứng với các biến nhập Xi,t-b được xác định
tương tự như mô hình hồi quy tuyến tính nhiều biến. Trong đó, biến phụ thuộc là
chuỗi Yt, biến độc lập là các chuỗi Xi,t-b (i = 1,2m là ký hiệu biến nhập, t-b là độ
chễ thời gian so với biến phụ thuộc). Phương pháp xác định các tham số hồi quy
tuyến tính nhiều biến (Ui,t-b) được trình bày trong giáo trình [10].
Đối với việc xác định các thông số pi, qi trong mô hình động thái ARIMAX
được xác định tương tự như trong mô hình ARIMA, tuy nhiên chuỗi tham gia vào
mô hình ARIMA trong trường hợp này là phần sai số của mô hình hồi quy tuyến
tính giữa Yt, biến độc lập Xi,t-b.
2.2.4. Kiểm định mô hình
Việc kiểm định để đánh giá sự phù hợp của mô hình là vô cùng quan trọng,
một mô hình được xem là phù hợp khi và chỉ khi thỏa mản các tiêu chuẩn như sau:
1) các giá trị sai số của mô hình phải độc lập nhau; 2) giá trị của mỗi thông số cần
phải đủ lớn để mô hình có ý nghĩa thống kê; 3) các thông số của mô hình phải độc
lập nhau. Dưới đây sẽ đi vào từng vấn đề cụ thể:
1) Kiểm tra tính độc lập của chuỗi sai số
Việc kiểm tra tính độc lập của chuỗi sai số chính là xem xét hàm tự tương
quan của các sai số đó. Nếu sai số có tính độc lập thì giữa chúng không có tương
quan, hay nói cách khác các hệ số tương quan giữa chúng phải thỏa mản biểu thức
giới hạn tin cậy (công thức 2.8 mục 2.2.1). Nếu không đạt tiêu chuẩn này thì giữa
chúng có mối quan hệ với nhau, như vậy mô hình được chọn chưa phù hợp với các
chuỗi số liệu được xem xét. Khi đó cần dựa vào hàm ACF và PACF của chuỗi sai
số để điều chỉnh lại dạng của mô hình.
2) Kiểm định độ lớn của các thông số
Thực chất các mô hình trong luận văn này là sử dụng phương pháp thống kê,
vì vậy các thông số được xác định trong mô hình đều phải kiểm định ý nghĩa thống
kê. Theo lý thuyết thống kê, việc kiểm định độ lớn các thông số trong mô hình
thường sử dụng công thức sau [6]:
ni
i
ni
s
n
t
(2.18)
36
Trong đó:tni là chỉ tiêu kiểm định độ lớn của thông số thứ i, ni là giá trị của
thông số thứ i, sni là độ lệch chuẩn của thông số thứ i.
- Nếu tni tα/2,n-1, (với = 0.05), thông số đó sẽ được giữ lại;
- Nếu tni < tα/2,n-1, thông số đó sẽ bị loại bỏ và quá trình tính toán để xác
định giá trị của thông số sẽ được thực hiện lại theo các thông số được giữ lại.
3) Kiểm định tính độc lập giữa các thông số
Quá trình phân tích trong bước kiểm định độ lớn của các thông số, đã loại đi
nhiều trường hợp. Tuy nhiên, vẫn có thể có một số thông số có sự tương quan mật
thiết với nhau. Để xem xét tính độc lập giữa các thông số, sẽ xem xét ma trận tương
quan của các thông số. Theo [6], nếu giữa các thông số có tương quan cao (r>0,8-
0,9) thì sẽ loại bỏ một trong 2 thông số tạo nên hệ số tương quan lớn này và một mô
hình có số thông số ít hơn sẽ được chọn.
2.2.5. Phần mềm thống kê SAS đối với mô hình ARIMA và ARIMAX
Việc tính toán đối với mô hình ARIMA và ARIMAX là bài toán khá phức tạp,
các bước kiểm nghiệm lặp lại nhiều lần. Để hỗ trợ trong quá trình tính toán đảm bảo
độ chính xác, nhiều phần mềm thống kê như: STATISTICA, NCSS, SYSTAT,
SASđều có chức năng tính toán các thông số trong mô hình ARIMA và ARIMAX.
Trong luận văn này đã sử dụng phần mềm thống kê SAS trong việc tính toán các
đặc trưng thống kê trong mô hình.
SAS (Statistical Analysis System) là hệ thống phần mềm thống kê do viện
nghiên cứu phần mềm thống kê của Mỹ xây dựng và đã được ở rất nhiều quốc gia
trên thế giới sử dụng. Trong phần mềm SAS có đầy đủ các chức năng thống kê như
lưu trữ, quản lý, kiểm tra, phân tích dử liệu. Đối với mô hình ARIMA và ARIMAX
trong phần mềm SAS bao gồm một hệ thống các câu lệnh được thực hiện trên mã
nguồn mở nhằm tính toán các tham số trong mô hình. Chi tiết về cách sử dụng phần
mềm này được trình bày trong [37,5].
2.3. Các nguồn số liệu đƣợc sử dụng
Trên cơ sở các mục tiêu và nội dung nghiên cứu luận văn, đã sử dụng các
loại số liệu sau: 1) Số liệu quan trắc lượng mưa tại 9 trạm đại diện cho khu vực
37
đồng bằng Bắc Bộ; 2) số liệu về các chỉ số khí hậu; 3) Số liệu vết đen Mặt Trời.
Trong đó, chuỗi số liệu từ năm 1951 đến 2008 để phát triển mô hình, từ năm 2009
đến 2013 để kiểm chứng mô hình. Chi tiết về các nguồn số liệu và cách xử lý đối
với từng loại số liệu sẽ được trình bày dưới đây:
2.3.1. Số liệu quan trắc mưa từ các trạm khí tượng
Bộ số liệu tổng lượng mưa tháng được thu thập tại Trung tâm Tư liệu Khí
tượng Thủy văn, hầu hết chuỗi dữ liệu được thu thập từ năm 1961 -2013, (riêng
trạm Hà Nội, từ năm 1951 -2013). Nhìn chung, nguồn số liệu này có độ tin cậy cao
do đã được kiểm tra, kiểm soát. Bộ số liệu này không chỉ có vai trò làm nhân tố dự
báo mà còn được sử dụng làm yếu tố dự báo trong mô hình ARIMA và ARIMAX.
2.3.2. Số liệu về các chỉ số khí hậu
Bộ số liệu các chỉ số khí hậu tham gia làm nhân tố dự báo bao gồm: chỉ số
dao đông nam (Southern Oscillation Index, SOI), dị thường nhiệt độ mặt nước biển
(ASST) trên các vùng NINO1.2, NINO3, NINO4, NINO3.4, vị trí của nhóm nhân tố
này được trình bày trên hình 2.5. Các nhân tố này phản ánh khá đầy đủ hoạt động
của hiện tượng ENSO trên khu vực xích đạo Thái Bình Dương. Nước ta nằm trong
khu vực chịu ảnh hưởng của hiện tượng này, do vậy diễn biến của 5 đặc trương khí
tượng hải dương này sẽ có quan hệ với diễn biến của lượng mưa trên các vùng lãnh
thổ nước ta. Bộ số liệu này có từ năm 1951 đến nay được cập nhật thường xuyên
qua mạng Internet thông qua website: edoc.
Hình 2.5 [2] Vị trí nhóm nhân tố ENSO
Tahiti
Darwin
SOI
38
Trong các chỉ số khí hậu thì nhóm nhân tố ENSO có quan hệ rất chặt chẽ với
nhau và gần như chúng có chung chu kỳ (hình2.6), do vậy, tại một bước trễ thời
gian, các nhân tố dự báo có thể chỉ có từ một đến 2 nhân tố thỏa mãn ý nghĩa thống
kê tham gia vào mô hình dự báo.
Hình 2.6. [9] Mật độ phổ của chỉ số SOI và nhiệt độ bề mặt nước biển
ở các vùng Nino.
2.3.3. Số liệu về số vết đen mặt trời (Sunspot Number)
Theo nhiều nghiên cứu [40] vết đen mặt trời có ảnh hưởng đến nhiều hiện
tượng vật lý khí quyển, trong đó ảnh hưởng rõ rệt nhất đến lượng mưa, hạn hán và
lũ lụt theo các chu kỳ và các thời gian khác nhau, số vết đen mặt trời trong một năm
nào đó thường ảnh hưởng đến lượng mưa ở nhiều năm sau đó. Với ý nghĩa quan
trọng như vậy, luận văn đã sử dụng số liệu này như là một trong những nhân tố dự
báo lượng mưa. Nguồn số liệu vết đen mặt trời được công bố bởi Trung tâm phân
tích ảnh hưởng mặt trời (Solar Influences Data analysis Center), thuộc đài quan sát
hoàng gia của Bỉ,. Đây là bộ số liệu được quan trắc hàng ngày từ năm 1818 đến nay
và được cập nhật thường xuyên trên trang webside:
data/dailyssn.php.
2.3.4. Xử lý số liệu
Trên cơ sở các nguồn số liệu được thu thập trong thời kỳ từ 1951 đến 2012,
tiến hành xử lý, tính toán theo các bước sau:
39
Bước 1: Tính toán bộ số liệu về tổng lượng mưa vụ đông xuân (từ tháng 10
năm trước đến tháng 1 năm sau) đối với 9 trạm khí tượng đại diện cho khu vực
đồng bằng Bắc Bộ. Đối với bộ số liệu nhân tố dự báo, nhóm các nhân tố ENSO là
các đặc trưng theo tháng ở dạng chuẩn sai nên được dữ nguyên, số liệu về số vết
đen Mặt trời được tính trung bình theo tháng.
Bước 2: Chuẩn hóa tập nhân tố dự báo
Do các nhân tố dự báo không có cùng thứ nguyên, ngoại tr
Các file đính kèm theo tài liệu này:
- luanvanthacsi_chuaphanloai_301_6075_1870186.pdf