Nhiễu trắng (tiếng Anh: White noise) có vai trò hết sức quan trọng trong các mô hình tài chính, kinh tế, vận tải, y học… đặc biệt là các mô hình dự báo giá trong tài chính. Show
Hình minh họa (Nguồn: Medium) Nhiễu trắngXem thêm: Anaconda python là gì ? Cài đặt anaconda Xem thêm: renown Khái niệm Nhiễu trắng hay tiếng ồn trắng hay tạp âm trắng trong tiếng Anh gọi là: White noise. Tiếng ồn trắng là sự kết hợp của tất cả các tần số âm thanh, vì vậy, tại một thời điểm nào đó, sự xuất hiện của một tần số là bất kì, hoàn toàn ngẫu nhiên, có thể cao hoặc thấp, không thể biết trước được. Sự biến thiên hoàn toàn mang tính ngẫu nhiên và không có các phần tử mang tính hệ thống nào khi xuất hiện trong các quá trình ngẫu nhiên tài chính thay vì gọi là tiếng ồn trắng, chúng ta thường gọi là nhiễu trắng. Ứng dụng của nhiễu trắng trong tài chínhTiếng ồn trắng có vai trò hết sức quan trọng trong các mô hình tài chính, kinh tế, vận tải, y học… đặc biệt là các mô hình dự báo giá trong tài chính. Dựa vào tính chất ngẫu nhiên của nhiễu trắng mà nhiều mô hình toán học đã được xây dựng để ước lượng các tham số trong các mô hình dự báo giá tài sản như dự báo giá chứng khoán. Mục đích chính trong việc mô hình thống kê là để tách ra thông tin cơ bản của các quá trình nhiều nhất có thể và để lại phần dư xấp xỉ một nhiễu trắng. Nhiễu trắng là một trong những quá trình ngẫu nhiên quan trọng trong nhiều lĩnh vực. Một cách toán học, nó có phổ xác định, giống như ánh sáng trắng chúng ta nhìn thấy bằng mắt. Các nhà phân tích tài chính đã áp dụng nhiễu trắng để mô hình các thị trường chứng khoán. Thực tế nhiễu trắng hay tiếng ồn trắng đã được sử dụng trong xung lực, tổng hợp âm thanh, nghệ thuật, điều trị giấc ngủ… Nhiễu trắng là dạng cơ bản của quá trình ngẫu nhiên cung cấp nền tảng cho hầu hết tất cả các mô hình thống kê ứng dụng được dùng trong khoa học tự nhiên và khoa học xã hội. Các mô hình chuỗi thời gian dự báo giá tài sản, giá chứng khoán/ cổ phiếu: – Mô hình một bước thời gian; – Mô hình với thời gian rời rạc; – Mô hình tự hồi qui tích hợp trung bình trượt ARIMA; – Mô hình ARCH (Autoregressive conditional heteroskedasticity) – mô hình phương sai có điều kiện của sai số thay đổi tự hồi qui và mô hình GARCH (Mở rộng Autoregressive conditioned heteroskedasticity); – Mô hình Black Scholes và ứng dụng trong tài chính của định lí Black S – Ocone. (Tài liệu tham khảo: Tiếng ồn trắng và ứng dụng trong tài chính, ThS. Vũ Thị Hương Sắc, Tạp chí Công thương, 2017) Dự báo chuỗi thời gian là một lớp mô hình quan trọng trong thống kê, kinh tế lượng và machine learning. Sở dĩ chúng ta gọi lớp mô hình này là chuỗi thời gian (time series) là vì mô hình được áp dụng trên các chuỗi đặc thù có yếu tố thời gian. Một mô hình chuỗi thời gian thường dự báo dựa trên giả định rằng các qui luật trong quá khứ sẽ lặp lại ở tương lai. Do đó xây dựng mô hình chuỗi thời gian là chúng ta đang mô hình hóa mối quan hệ trong quá khứ giữa biến độc lập (biến đầu vào) và biến phụ thuộc (biến mục tiêu). Dựa vào mối quan hệ này để dự đoán giá trị trong tương lai của biến phụ thuộc. Do là dữ liệu chịu ảnh hưởng bởi tính chất thời gian nên chuỗi thời gian thường xuất hiện những qui luật đặc trưng như : yếu tố chu kỳ, mùa vụ và yếu tố xu hướng. Đây là những đặc trưng thường thấy và xuất hiện ở hầu hết các chuỗi thời gian.
Hình 1: Đồ thị về chuỗi nhiệt độ trung bình theo tháng thể hiện yếu tố mùa vụ.
Hình 2: Đồ thị về yếu tố xu hướng trong chuỗi thời gian của chuỗi giá. Các dự báo chuỗi thời gian có tính ứng dụng cao và được sử dụng rất nhiều lĩnh vực như tài chính ngân hàng, chứng khoán, bảo hiểm, thương mại điện tử, marketing, quản lý chính sách. Bên dưới là một số ứng dụng của dự báo chuỗi thời gian:
Vai trò của chuỗi thời gian rất quan trọng đối với nền kinh tế và hoạt động của doanh nghiệp nên trong machine learning và thống kê có những ngành học nghiên cứu chuyên sâu về chuỗi thời gian như kinh tế lượng, định giá tài sản tài chính. Khác với các mô hình dự báo thông thường trong machine learning, các mô hình trong dự báo chuỗi thời gian trong kinh tế lượng có những đặc trưng rất riêng. Đòi hỏi phải tuân thủ nghiêm ngặt các điều kiện về chuỗi dừng, nhiễu trắng và tự tương quan. Có rất nhiều lớp mô hình chuỗi thời gian khác nhau và mỗi một lớp mô hình sẽ có một tiêu chuẩn áp dụng cụ thể. Chúng ta có thể liệt kê một số mô hình phổ biến:
2. Mô hình ARIMAHiện tại cả R và python đều support xây dựng các mô hình chuỗi thời gian ARIMA, SARIMA, ARIMAX, GARCH,…. Trên R chúng ta có thể sử dụng các packages như 1 2 3 4 và 1 2 3 5 để xây dựng các mô hình này khá dễ dàng. Đối với thống kê và các mô hình chuỗi thời gian R đang support tốt hơn python. Một lý do đó là các nhà thống kê và kinh tế lượng ưa chuộng sử dụng R hơn. Hướng dẫn xây dựng mô hình ARIMA trên R các bạn có thể xem tại ARIMA turtorial và GARCH time series model. Bài hôm nay tôi sẽ hướng dẫn các bạn sử dụng python để xây dựng các model ARIMA. Nhưng trước tiên chúng ta cần tìm hiểu lý thuyết của mô hình trước khi đi vào phần thực hành ở mục 3 và 4. 2.1. Lý thuyết mô hình ARIMALý thuyết: Chúng ta biết rằng hầu hết các chuỗi thời gian đều có sự tương quan giữa giá trị trong quá khứ đến giá trị hiện tại. Mức độ tương quan càng lớn khi chuỗi càng gần thời điểm hiện tại. Chính vì thể mô hình ARIMA sẽ tìm cách đưa vào các biến trễ nhằm tạo ra một mô hình dự báo fitting tốt hơn giá trị của chuỗi. ARIMA model là viết tắt của cụm từ Autoregressive Intergrated Moving Average. Mô hình sẽ biểu diễn phương trình hồi qui tuyến tính đa biến (multiple linear regression) của các biến đầu vào (còn gọi là biến phụ thuộc trong thống kê) là 2 thành phần chính:
\[\text{AR}(p) = \phi_{0}+\phi_{1}x_{t-1} + \phi_{2}x_{t-2} + \dots + \phi_{p} x_{t-p}\]
\[\begin{equation} \left\{ \begin{array}{l l} \text{E}(\epsilon_t) & = 0 & (1) \\ \sigma(\epsilon_t) & = \alpha & (2)\\ \rho(\epsilon_t, \epsilon_{t-s}) & = 0, \forall s <= t & (3) \end{array} \right. \end{equation}\] Vế (1) có nghĩa rằng kì vọng của chuỗi bằng 0 để đảm bảo chuỗi dừng không có sự thay đổi về trung bình theo thời gian. Vế (2) là phương sai của chuỗi không đổi. Do kì vọng và phương sai không đổi nên chúng ta gọi phân phối của nhiễu trắng là phân phối xác định (identical distribution) và được kí hiệu là $\epsilon_{t}∼\text{WN}(0,\sigma^2)$. Nhiễu trắng là một thành phần ngẫu nhiên thể hiện cho yếu tố không thể dự báo của model và không có tính qui luật. Qúa trình trung bình trượt được biểu diễn theo nhiễu trắng như sau: \[\text{MA}(q) = \mu+\sum_{i=1}^{q} \theta_i\epsilon_{t-i}\] Qúa trình này có thể được biểu diễn theo dịch chuyển trễ - backshift operator $B$ như sau: \[\text{MA}(q) = \mu+(1+\theta_1 B + \dots + \theta_q B^{q})\epsilon_{t}\] Như vậy bạn đọc đã hình dung ra moving average là gì rồi chứ? Về mặt ý tưởng thì đó chính là quá trình hồi qui tuyến tính của giá trị hiện tại theo các giá trị hiện tại và quá khứ của sai số nhiễu trắng (white noise error term) đại diện cho các yếu tố shock ngẫu nhiên, những sự thay đổi không lường trước và giải thích bởi mô hình.
Thông thường chuỗi sẽ dừng sau quá trình đồng tích hợp $\text{I}(0)$ hoặc $\text{I}(1)$. Rất ít chuỗi chúng ta phải lấy tới sai phân bậc 2. Một số trường hợp chúng ta sẽ cần biến đổi logarit hoặc căn bậc 2 để tạo thành chuỗi dừng. Phương trình hồi qui ARIMA(p, d, q) có thể được biểu diễn dưới dạng: \[\Delta x_{t} = \phi_{1} \Delta x_{t-1}+\phi_{2} \Delta x_{t-2}+...+\phi_{p}\Delta x_{t-p}+ \theta_{1}\epsilon_{t-1}+\theta_{2}\epsilon_{t-2}+...+\theta_{q}\epsilon_{t-q}\] Trong đó $\Delta x_t$ là giá trị sai phân bậc $d$ và $\epsilon_t$ là các chuỗi nhiễu trắng. Như vậy về tổng quát thì ARIMA là mô hình kết hợp của 2 quá trình tự hồi qui và trung bình trượt. Dữ liệu trong quá khứ sẽ được sử dụng để dự báo dữ liệu trong tương lai. Trước khi huấn luyện mô hình, cần chuyển hóa chuỗi sang chuỗi dừng bằng cách lấy sai phân bậc 1 hoặc logarit. Ngoài ra mô hình cũng cần tuân thủ điều kiện ngặt về sai số không có hiện tượng tự tương quan và phần dư là nhiễu trắng. Đó là lý thuyết của kinh tế lượng. Còn theo trường phái machine learning thì tôi chỉ cần quan tâm đến làm sao để lựa chọn một mô hình có sai số dự báo là nhỏ nhất. Tiếp theo chúng ta sẽ sử dụng package vnquant, một package được tôi viết để hỗ trợ cộng đồng khai thác dữ liệu chứng khoán thuận tiện hơn. 3. Ứng dụng vnquant trong thu thập dữ liệu.3.1. Thu thập dữ liệuHiện tại package vnquant đã cho phép chúng ta thu thập được hầu hết mã chứng khoán trên thị trường chứng khoán Việt Nam, ở phiên bản R (package VNDS) còn thu thập được báo cáo tài chính và dòng tiền của doanh nghiệp. Ngoài ra vnquant còn hỗ trợ vẽ biểu đồ nến theo thời gian và kết hợp giữa giá và khối lượng giao dịch. Đây là một package mà mình nghĩ là rất hữu ích đối với các bạn làm trong lĩnh vực quant tại Việt Nam. Để cài đặt package này bạn làm như hướng dẫn trong phần read me nhé. Bây giờ chúng ta sẽ cùng lấy dữ liệu chỉ số VNINDEX 30 thông qua package vnquant nào. 1 2 3 4 5 6 7 8 9 10 from vnquant.DataLoader import DataLoader loader = DataLoader(symbols="VN30", data = loader.download()
data.head()1 2 2019-12-12 04:22:07,057 : INFO : NumExpr defaulting to 2 threads. 2019-12-12 04:22:07,076 : INFO : data VN30 from 2019-01-01 to 2019-12-09 have already cloned! Attributes change_perc1 change_perc2 open high low close avg volume_match volume_reconcile volume Symbols VN30 VN30 VN30 VN30 VN30 VN30 VN30 VN30 VN30 VN30 date 2019-01-02 0.00 0.000000 855.01 864.64 852.79 855.66 855.66 37199310.0 3099070.0 40298380.0 2019-01-03 -16.87 0.019716 855.20 856.21 834.81 838.79 838.79 55933762.0 5880342.0 61814104.0 2019-01-04 1.38 0.001645 836.60 842.05 821.83 840.17 840.17 41706643.0 4266803.0 45973446.0 2019-01-07 11.24 0.013378 843.95 858.69 843.95 851.41 851.41 35905579.0 5233979.0 41139558.0 2019-01-08 -5.98 0.007024 851.20 852.92 841.86 845.43 845.43 33310330.0 5386790.0 38697120.0 Sử dụng hàm visualization trên chính vnquant để visualize dữ liệu lịch sử giá và khối lượng giao dịch. 1 2 3 4 5 from vnquant import Plot Plot._vnquant_candle_stick(data = data,
Nhận xét, năm 2019 là một năm thị trường chứng khoán có nhiều thăng trầm biến động khi chỉ số có lúc vượt ngưỡng 950 nhưng có những giai đoạn hạ xuống thấp hơn 850 điểm. Tuy nhiên chưa thể phục hồi lại ngưỡng đỉnh cao trên 1000 điểm của năm 2018. 3.2. Khởi tạo chuỗi lợi suất và khảo sát dữ liệuĐể thuận tiện cho việc xây dựng mô hình ARIMA ta sẽ chuyển chuỗi giá close sang chuỗi dừng bằng cách lấy lợi suất theo công thức sai phân bậc 1 của logarit như bên dưới: \[r_t = \log(\frac{x_t}{x_{t-1}})\] Mục tiêu của mô hình sẽ là dự báo chuỗi $r_t$. Từ chuỗi $r_t$ ta có thể dễ dàng biến đổi ngược lại thành giá đóng cửa của chuỗi VnIndex 30. Hàm 1 2 3 6 sẽ giúp ta lấy trễ bậc 1 của chuỗi giá close. Công thức tính toán lợi suất như sau: 1 2 3 import numpy as np Tính chuỗi returnr_t = np.log(data['close']/data['close'].shift(1)).values[:, 0] Fill giá trị 1 2 3 7 bằng trung bình chuỗi. 1 2 3 mean = np.nanmean(r_t) r_t[0]=mean r_t[:5] from vnquant.DataLoader import DataLoader loader = DataLoader(symbols="VN30", data = loader.download()
data.head()0 from vnquant.DataLoader import DataLoader loader = DataLoader(symbols="VN30", data = loader.download()
data.head()1 Các biểu đồ:
Để hiểu rõ hơn về xu hướng, biến động, ta hãy cùng vẽ biểu đồ chuỗi lợi suất $r_t$. from vnquant.DataLoader import DataLoader loader = DataLoader(symbols="VN30", data = loader.download()
data.head()2 from vnquant.DataLoader import DataLoader loader = DataLoader(symbols="VN30", data = loader.download()
data.head()3 Nhận xét: Biểu đồ chuỗi lợi suất cho thấy nó là một biến động ngẫu nhiên dạng nhiễu trắng, có trung bình gần như bằng 0 và phương sai không đổi.
Ta có thể vẽ biểu đồ biểu diễn chuỗi $r_t$ dựa trên chuỗi $r_{t-1}$ để xem chúng có quan hệ tuyến tính hay ngẫu nhiên. from vnquant.DataLoader import DataLoader loader = DataLoader(symbols="VN30", data = loader.download()
data.head()4 from vnquant.DataLoader import DataLoader loader = DataLoader(symbols="VN30", data = loader.download()
data.head()5 Đồ thị cho thấy 2 chuỗi $r_t$ và $r_{t-1}$ không có mối quan hệ tương quan. Biểu đồ của chúng là một tợp hợp các điểm không tuân theo một trend cụ thể.
Ta hãy cũng khảo sát qua biểu đồ phân phối xác suất của chuỗi lợi suất. from vnquant.DataLoader import DataLoader loader = DataLoader(symbols="VN30", data = loader.download()
data.head()2 from vnquant.DataLoader import DataLoader loader = DataLoader(symbols="VN30", data = loader.download()
data.head()7 Từ biểu đồ phân phối lợi suất ta nhận thấy chuỗi phân phối lợi suất dường như có dạng phân phối chuẩn và có kì vọng bằng 0. Chúng ta có thể kiểm định phân phối chuẩn thông qua biểu đồ qqplot (quantiles quantiles plot).
\[P(r_t < r_i| r_t \sim N(\mu, \sigma^2)) = 0.25\] thì nó sẽ nằm trong khoảng ngũ phân vị thứ 2. Biểu đồ qqplot sẽ biểu diễn các điểm có tọa độ $(x, y)$ sao cho giá trị $x$ chính là khoảng ngũ phân vị theo phân phối chuẩn và giá trị $y$ chính là giá trị thực nghiệm. Từ biểu đồ ta kết luận 2 chuỗi có phân phối tương tự nhau nếu như các điểm trên đồ thì nằm trên một đường thẳng. Khi đó $r_t$ có thể được coi như là một phân phối chuẩn. Cách kiểm định phân phối chuẩn dựa trên biểu đồ qqplot được sử dụng khá phổ biến nhưng nó có một nhược điểm là không cung cấp một tiêu chuẩn xác định của việc chấp nhận/bác bỏ giả thuyết. Để vẽ biểu đồ qqplot ta có thể làm theo 2 cách khác nhau. Cách 1: Sử dụng trực tiếp hàm 1 2 3 1 2 3 from vnquant.DataLoader import DataLoader loader = DataLoader(symbols="VN30", data = loader.download()
data.head()9 Cách 2: Tính ra chuỗi lý thuyết và vẽ đồ thị Chúng ta cũng có thể vẽ biểu đồ qqplot bằng cách tính ra trực tiếp giá trị phân phối thực nhiệm (Theoretical Quantiles) và vẽ biểu đồ lợi suất thực tế đã được sắp xếp theo thứ tự tăng dần (Sample Quantiles). 1 2 0 1 2 1 Như vậy từ đồ thị ta có thể khẳng định chuỗi lợi suất có phân phối chuẩn vì đồ thị của phân phối lý thuyết và phân phối thực nghiệm nằm trên cùng một đường thẳng. 3.3. Kiểm tra tính dừng.Một trong những điều kiện tiền đề khi hồi qui các mô hình chuỗi thời gian đó là chuỗi phải dừng. Để kiểm định tính dừng chúng ta có thể sử dụng kiểm định Argument Dickey Fuller hay còn gọi là kiểm định nghiệm đơn vị. Giả sử ta có một quá trình tự hồi qui $\text{AR}(1)$ đối với chuỗi $y_t$ được xác định như sau: \(y_t = \alpha + \phi y_{t-1} + \epsilon_t\) Phương trình trên tương ứng với trường hợp chuỗi $y_t$ có hệ số chặn và có xu hướng (tổng quát nhất).
Với $\epsilon_t \sim WN(0, \sigma^2)$ là một chuỗi sai số phân phối nhiễu trắng. Khi khai triển $y_t$ một cách liên tục theo giá trị trễ ta có: \[\begin{eqnarray} y_t& = & \phi y_{t-1} + \epsilon_t \\ & = & \phi(\phi y_{t-2} + \epsilon_{t-1}) + \epsilon_t \\ & = &\dots \\ & = & \phi(\phi( \dots (\phi y_{0} + \epsilon_{0})+ \dots) + \epsilon_{t-1}) + \epsilon_t \\ & = & \phi^{t}y_0 + \phi^{t-1}\epsilon_{1} + \dots + \phi \epsilon_{t-1} + \epsilon_t \end{eqnarray}\] Do đó \(\mathbf{E}(y_t) = \phi^{t}\mathbf{E}(y_0)+\mathbf{E}(\phi^{t-1}\epsilon_{1} + \dots + \phi \epsilon_{t-1} + \epsilon_t) = \phi^{t}\mathbf{E}(y_0)\) Dấu bằng thứ 2 xảy ra là do $\text{E}(\epsilon_t)=0,\forall t\in N$ Mặt khác: \[\begin{equation} \begin{cases} \phi > 1: \lim\limits_{t \rightarrow \infty} \phi^{t} \mathbf{E}(y_0) = \infty \\ \phi = 1: \lim\limits_{t \rightarrow \infty} \phi^{t} \mathbf{E}(y_0) = \mathbf{E}(y_0) \\ \phi < 1: \lim\limits_{t \rightarrow \infty} \phi^{t} \mathbf{E}(y_0) = 0 \end{cases} \end{equation}\] Như vậy tính chất dừng của chuỗi $y_t$ sẽ phụ thuộc vào phương trình đặc trưng $\theta(y)=1−\theta L=0$ có nghiệm đơn vị hay không. Nếu phương trình đặc trưng có nghiệm đơn vị (unit root), chuỗi $y_t$ sẽ không dừng. Trái lại ta có thể khẳng định $y_t$ là chuỗi dừng như đồ thị mô tả bên dưới. Đồ thị mô tả một khả năng của nghiệm đơn vị. Đường màu đỏ hiển thị sự sụt giảm của output và đường phục hồi nếu chuỗi thời gian có nghiệm đơn vị. Màu xanh hiển thị sự phục hồi nếu không có nghiệm đơn vị và chuỗi là chuỗi dừng có xu hướng (trend-stationary). (Nguồn Unit Root: Simple Definition, Unit Root Tests). Nhắc lại lý thuyết về kiểm định. Một kiểm định thống kê sẽ bao gồm 2 cặp giả thuyết kiểm định đó là: Gỉa thuyết null (null hypothesis), còn được gọi là 1 2 3 9, kí hiệu là $H_0$. Đây là giả thuyết bị nghi ngờ xảy ra và được sử dụng để kiểm chứng những tính chất liên quan đến mẫu mà chúng ta chưa biết rằng mẫu sở hữu trên thực tế. Lấy ví dụ, một giáo viên nói rằng các học sinh của trường thi đạt học đạt điểm toán trung bình là 7. Tuy nhiên đây mới chỉ là một giả thuyết cần kiểm chứng. Chúng ta có thể kiểm tra bằng cách thu thập 30 mẫu các học sinh của trường và kiểm chứng giả thuyết của cô giáo bằng cách kiểm định xem trung bình điểm thi toán của 30 học sinh trên có bằng 7 hay không?
Để kiếm tra phương trình đặc trưng của chuỗi có nghiệm đơn vị hay không chúng ta sử dụng kiểm định ADF. Giả thuyết null được đặt ra đó là phương trình đặc trưng có nghiệm đơn vị. Trong trường hợp p-value < 0.05 thì ta sẽ loại bỏ giả thuyết null, chấp nhận giả thuyết thay thế. Khi đó ta có thể khẳng định rằng chuỗi không có nghiệm đơn vị và có tính chất dừng. \[\begin{equation} \begin{cases} H_0: \phi = 1, \Rightarrow \text{unit root, non-stationary} \\ H_1: |\phi| < 1, \Rightarrow \text{non-unit root, stationary} \end{cases} \end{equation}\] Giá trị ngưỡng kiểm định: \(DF = \frac{\hat\phi - 1}{SE(\hat\phi)}\) Chúng ta sẽ so sánh giá trị ngưỡng kiểm định này với giá trị tới hạn của phân phối Dickey - Fuller để đưa ra kết luận về chấp nhận hoặc bác bỏ giả thuyết $H_0$. Trên python đã hỗ trợ kiểm định ADF thông qua package import numpy as np Tính chuỗi returnr_t = np.log(data['close']/data['close'].shift(1)).values[:, 0] 2. Ta sẽ kiểm định ADF cho chuỗi lợi suất. from vnquant.DataLoader import DataLoader loader = DataLoader(symbols="VN30", data = loader.download()
data.head()4 1 2 3 1 2 4 1 2 5 Gía trị p-value < 0.05, kết luận chúng ta sẽ bác bỏ giả thuyết $H_0$. Phương trình đặc trưng không có nghiệm đơn vị. Do đó chuỗi lợi suất có tính chất dừng. 4. Xây dựng mô hìnhỞ mục 2.1. ta đã trình bày về mặt lý thuyết của mô hình ARIMA đã khá đầy đủ. Tuy nhiên lý thuyết chỉ là lý thuyết. Chúng ta sẽ phải biến lý thuyết thành thực tiễn thông qua việc thực hành. Phần này chúng ta sẽ cùng thử nghiệm xây dựng một mô hình ARIMA trên python sử dụng dữ liệu được lấy từ package vnquant. Mục tiêu của mô hình là dự báo lợi suất theo ngày của chuỗi VnIndex 30 từ dữ liệu được khai thác trong giai đoạn từ import numpy as np Tính chuỗi returnr_t = np.log(data['close']/data['close'].shift(1)).values[:, 0] 3 đến import numpy as np Tính chuỗi returnr_t = np.log(data['close']/data['close'].shift(1)).values[:, 0] 4 theo phương pháp ARIMA. Đầu tiên cần phải lựa chọn bậc phù hợp cho mô hình ARIMA. 4.1. Lựa chọn tham số ARIMA(p, d, q)Tự tương quan (ACF - AutoCorrelation Function): Tự tương quan là một khái niệm quan trọng trong chuỗi thời gian. Hầu hết các chuỗi thời gian sẽ có sự tương quan với giá trị trễ của nó và các giá trị càng gần nhau thì tương quan càng mạnh hoặc các giá trị cùng thuộc 1 chu kì của chuỗi thì sẽ có tương quan cao (chẳng hạn như cùng tháng trong chu kì năm hay cùng quí trong chu kì năm). Chính vì vậy hệ số này mới có tên là tự tương quan. Hệ số tự tương quan được viết tắt là ACF và thường dùng để tìm ra độ trễ của quá trình trung bình trượt $\text{MV}(q)$ để xây dựng các mô hình như ARIMA, GARCH, ARIMAX,… và kiểm tra yếu tố mùa vụ. Hệ số tự tương quan bậc $s$ được xác định như sau: \[\rho(s,t) = \frac{cov(x_s,x_t)}{\sqrt{\sigma_s\sigma_t}}\] giá trị $\rho(s,t)$ đo lường khả năng dự báo của biến $x_t$ nếu chỉ sử dụng biến $x_s$. Trong trường hợp 2 đại lượng có tương quan hoàn hảo tức $\rho(s,t)=±1$ ta có thể biểu diễn $x_t=\beta_0+\beta_1 x_s$. Hệ số của $\beta_1$ sẽ ảnh hưởng lên chiều của hệ số tương quan. Theo đó $\rho(s, t)=1$ khi $\beta_1>0$ và $\rho(s,t)=−1$ khi $\beta_1<0$. Chúng ta có thể vẽ biểu đồ các hệ số tự tương quan ACF theo các bậc liên tiếp thông qua hàm plot_acf của import numpy as np Tính chuỗi returnr_t = np.log(data['close']/data['close'].shift(1)).values[:, 0] 2 như bên dưới: 1 2 0 1 2 7 Trục hoành là độ trễ, trục tung là giá trị của hệ số tự tương quan tương ứng với độ trễ. Dải màu hồng chính là khoảng tin cậy 95% để giá trị hệ số tự tương quan bằng 0. Nếu tại một độ trễ nhỏ nhất mà đoạn thẳng (vuông góc với trục hoành) mà độ dài đại diện cho giá trị của hệ số tự tương quan nằm ngoài khoảng tin cậy thì đó chính là độ trễ phù hợp lớn nhất mà ta nên lựa chọn cho quá trình trung bình trượt $\text{MV}(q)$. Nhìn chung bậc $q$ không nên quá lớn. Thông thường tôi chỉ chọn tối đa là 5. Đối với bài toán này toàn bộ các hệ số tự tương quan với bậc nhỏ hơn hoặc bằng 5 đều có giá trị nằm trong khoảng tin cậy 95% của 0. Do đó chúng ta có thể linh hoạt lựa chọn bậc $q=5$ là vị trí mà hệ số tự tương quan lớn nhất. Tự tương quan riêng phần (PACF - Partitial AutoCorrelation Function): Về cơ bản tương quan riêng phần cũng là chỉ số đo lường hệ số tương quan như ACF. Tuy nhiên vẫn có sự khác biệt đó là hệ số tương quan này loại bỏ ảnh hưởng của các chuỗi độ trễ trung gian (là các chuỗi trễ $x_{t-1}, \dots, x_{t-k+1}$ nằm giữa $x_{t}$ và $x_{t-k}$). Một phương trình hồi qui tuyến tính giữa chuỗi hiện tại với các chuỗi độ trễ trung gian được xây dựng nhằm đánh giá ảnh hưởng của các chuỗi độ trễ lên chuỗi hiện tại. Sau đó, để tính hệ số tương quan riêng phần chúng ta sẽ loại bỏ ảnh hưởng của các độ trễ trung gian khỏi chuỗi hiện tại bằng cách trừ đi giá trị ước lượng từ phương trình hồi qui. Lấy ví dụ: Để tính tự tương quan riêng phần PACF bậc $k$ của chuỗi $x_t$. Đầu tiên ta sẽ hồi qui tuyến tính $x_t$ theo các chuỗi trễ của nó là $x_{t-1},\dots,x_{t-k}$. Khi đó ta thu được phương trình hồi qui tuyến tính tổng quát bậc $k$ là: \[x_t = \epsilon_t + \alpha_0 + \alpha_1 x_{t-1} + \dots + \alpha_k x_{t-k}\] $\epsilon_t$ là thành phần đại diện cho sai số. Gía trị ước lượng của mô hình đối với $x_t$ chính là: \[P_{t, k}(x_t) = \alpha_0 + \alpha_1 x_{t-1} + \dots + \alpha_k x_{t-k}\] Hệ số tự tương quan tuyến tính sau đó sẽ chính bằng: \[\phi_{k} = corr(x_t-P_{t,k}(x_t), x_{t-k}-P_{t,k}(x_{t-k}))\] Trong đó $corr()$ là hàm tính hệ số tương quan. Đó là tất cả về PACF. Khá dễ hiểu phải không nào? PACF sẽ có tác dụng tìm ra hệ số bậc tự do $p$ của quá trình tự hồi qui $\text{AR}(p)$. Tương tự như ACF, thông qua một biểu đồ PACF về giá trị các hệ số tương quan riêng phần tương ứng với các độ trễ khác nhau, chúng ta sẽ tìm ra được các bậc tự do $p$ phù hợp. Đó chính là vị trí mà giá trị của hệ số tương quan riêng phần nằm ngoài ngưỡng tin cậy 95% của giả thuyết hệ số tương quan riêng phần bằng 0. 1 2 0 1 2 9 Tương tự như ACF, bậc của PACF cũng thường nhỏ hơn 5. Như vậy ta cũng có thể lựa chọn bậc tự do của PACF là một giá trị nào đó từ 1 đến 5. Kết hợp giữa bậc của $p$ và $q$ và giá trị của $d = 0$ do chuỗi $r_t$ đã là một chuỗi dừng ta có thể thu được một số kịch bản
4.2. Chỉ số AIC - Akaike Information CriteriaNhư vậy chúng ta đã thu được 5 kịch bản mô hình khác nhau và cần chọn ra một mô hình phù hợp nhất. Một trong những tiêu chí thường được sử dụng để lựa chọn mô hình đó là chỉ số AIC (Akaike Information Criteria). Tiêu chí thông tin này là một công cụ ước tính lỗi dự báo và do đó đánh giá chất lượng tương đối của các mô hình thống kê trên một tập hợp dữ liệu nhất định. Gỉa sử có một tập hợp các mô hình được xây dựng trên cùng một bộ dữ liệu, AIC ước tính chất lượng của từng mô hình trong mối liên quan đến từng mô hình khác. Do đó, AIC cung cấp một phương tiện để lựa chọn mô hình. AIC được hình thành dựa trên lý thuyết thông tin (information theory). Khi một mô hình thống kê được sử dụng để dự báo, kết quả sẽ gần như không bao giờ chính xác hoàn toàn. Vì vậy một số thông tin sẽ bị mất do không thể dự báo từ mô hình. AIC ước tính lượng thông tin tương đối bị mất bởi một mô hình nhất định: mô hình mất càng ít thông tin thì chất lượng của mô hình đó càng cao. Giả sử rằng chúng ta có một mô hình thống kê tương ứng với một bộ dữ liệu. Gọi $k$ là số lượng tham số ước tính trong mô hình. Đặt $\hat{L}$ là giá trị tối đa của hàm hợp lý (maximum likelihood function) của mô hình. Khi đó, giá trị AIC của mô hình được tính như sau: \[AIC = 2k-2\ln(\hat{L})\] Tóm lại rằng giá trị của AIC càng nhỏ thì mô hình của chúng ta càng phù hợp. 4.3. Đọc hiểu kết quả một mô hình ARIMA.Mô hình ARIMA có thể được xây dựng khá dễ dàng trên python thông qua package import numpy as np Tính chuỗi returnr_t = np.log(data['close']/data['close'].shift(1)).values[:, 0] 2. Điều mà chúng ta cần thực hiện chỉ là khai báo bậc của mô hình ARIMA. Giả sử cần xây dựng một mô hình ARIMA(2, 0, 0) ta thực hiện như sau: 1 2 3 4 5 2019-12-12 04:22:07,057 : INFO : NumExpr defaulting to 2 threads. 2019-12-12 04:22:07,076 : INFO : data VN30 from 2019-01-01 to 2019-12-09 have already cloned! 1 2019-12-12 04:22:07,057 : INFO : NumExpr defaulting to 2 threads. 2019-12-12 04:22:07,076 : INFO : data VN30 from 2019-01-01 to 2019-12-09 have already cloned! 2 2019-12-12 04:22:07,057 : INFO : NumExpr defaulting to 2 threads. 2019-12-12 04:22:07,076 : INFO : data VN30 from 2019-01-01 to 2019-12-09 have already cloned! 3 Bảng trên chính là summary kết quả từ mô hình ARIMA.
import numpy as np Tính chuỗi returnr_t = np.log(data['close']/data['close'].shift(1)).values[:, 0] 8 là độ lệch chuẩn của hệ số ước lượng. Từ giá trị ước lượng và độ lệch chuẩn ta có thể tính toán ra được khoảng tin cậy. Cận trên và dưới của khoảng tin cậy là các cột import numpy as np Tính chuỗi returnr_t = np.log(data['close']/data['close'].shift(1)).values[:, 0] 9 và 1 2 3
1 2 3 1 chính là giá trị ngưỡng tới hạn được suy ra từ phân phối chuẩn hóa. Gía trị ngưỡng tới hạn được tính bằng: \(z = \frac{\beta-\mu}{u_{\alpha/2}\sigma}\) với $\beta$ là giá trị ước lượng, $\mu, \sigma$ là các giá trị trung bình, độ lệch chuẩn của hệ số ước lượng. $u_{\alpha/2}$ là mức phân vị 97.5% của phân phối chuẩn hóa. Do ở đây ta đang kiểm định các giá trị $\beta$ bằng 0 nên trung bình được giả định bằng 0. Tức là $\mu=0$. Nhìn chung là ta không cần quan tâm đến cột này lắm bởi cột 1 2 3 2 đã thực hiện chức năng thay thế nó trong việc kiểm định các hệ số ước lượng có ý nghĩa thống kê hay không. Cột 1 2 3 2 là xác suất để giá trị $P(\mid X\mid > 0 \mid X \sim N(0, \sigma^2))$. Đây chính là giá trị P-value của cặp giả thuyết dấu bằng. Gía trị của P-value < 0.05 sẽ cho thấy hệ số ước lượng lớn hơn 0 là có ý nghĩa thống kê. Các chỉ số ở góc trên bên phải lần lượt là:
Như vậy xét trên khía cạnh mô hình thống kê thì tất cả các hệ số ước lượng đều có ý nghĩa thống kê với mức ý nghĩa 95% ngoại trừ giá trị ước lượng của hệ số tự do. Ta có thể thấy các mô hình trong thống kê và kinh tế lượng tuy đơn giản (chỉ là hồi qui tuyến tính) nhưng rất chú trọng tới các yếu tố như:
Trên đây là những điểm sơ đẳng nhất mà tôi rút ra từ kinh nghiệm của mình. Ngoài ra còn rất nhiều những sự khác biệt nữa giữa machine learning và thống kê, kinh tế lượng mà làm nhiều chúng ta sẽ tự đúc kết ra. Quay trở lại việc lựa chọn mô hình tốt nhất trong lớp các mô hình ARIMA, đơn giản chúng ta có thể căn cứ trên AIC như sau: 2019-12-12 04:22:07,057 : INFO : NumExpr defaulting to 2 threads. 2019-12-12 04:22:07,076 : INFO : data VN30 from 2019-01-01 to 2019-12-09 have already cloned! 4 2019-12-12 04:22:07,057 : INFO : NumExpr defaulting to 2 threads. 2019-12-12 04:22:07,076 : INFO : data VN30 from 2019-01-01 to 2019-12-09 have already cloned! 5 1 2 0 2019-12-12 04:22:07,057 : INFO : NumExpr defaulting to 2 threads. 2019-12-12 04:22:07,076 : INFO : data VN30 from 2019-01-01 to 2019-12-09 have already cloned! 7 Ta nhận thấy mô hình ARIMA(2,0,2) là phù hợp nhất với bộ dữ liệu lợi suất vì nó tương ứng với chỉ số AIC là nhỏ nhất. 4.4. Phương pháp Auto ARIMAChúng ta thấy rằng việc lựa chọn mô hình tốt nhất chỉ đơn thuần dựa trên chỉ số AIC, khá đơn giản. Do đó chúng ta hoàn toàn có thể tự động thực hiện qui trình này. Trên python đã hỗ trợ tìm kiếm mô hình ARIMA phù hợp thông qua package auto arima. Chúng hoạt động như một grid search mà tham số chúng ta truyền vào chỉ là các hệ số giới hạn trên của các bậc $(p, d, q)$. Mọi việc còn lại hãy để thuật toán tự giải quyết. Nếu bạn làm quen với R thì cũng hỗ trợ chức năng auto ARIMA tương tự như vậy. Để sử dụng phương pháp auto arima chúng ta phải dùng tới package pyramid. Cài đặt pyramid như bên dưới: ` pip install pyramid-arima ` Xây dựng phương trình hồi qui theo phương pháp Auto ARIMA 1 2 3 4 5 6 7 8 9 10 2019-12-12 04:22:07,057 : INFO : NumExpr defaulting to 2 threads. 2019-12-12 04:22:07,076 : INFO : data VN30 from 2019-01-01 to 2019-12-09 have already cloned! 9 1 2 4 1 2 3 4 5 1 Ở đây chúng ta sẽ khai báo điểm bắt đầu và kết thúc của $p, q$, bậc của $d$ và phương pháp điều chỉnh mô hình là stepwise. Chiến lược stepwise: Là một chiến lược được sử dụng khá phổ biến trong việc lựa chọn biến cho mô hình hồi qui. Chiến lược này sẽ trải qua nhiều bước. Mỗi một bước tìm cách thêm vào hoặc bớt đi một biến giải thích nhằm tạo ra một mô hình mới sao cho sai số giảm so với mô hình gốc. Như vậy càng qua nhiều bước mô hình sẽ càng chuẩn xác hơn và sau khi kết thúc chiến lược ta sẽ thu được một mô hình cuối cùng là mô hình tốt nhất. Tiêu chuẩn để đo lường sai số và ra quyết định lựa chọn thường là AIC. Trong trường hợp mô hình có yếu tố mùa vụ thì ta sẽ cần thiết lập 1 2 3 4 và kết hợp thêm chu kì của mùa vụ. Chẳng hạn chu kì là 12 tháng thì có thể khai báo 1 2 3 5. Khi đó mô hình ARIMA sẽ trở thành mô hình SARIMA (seasonal ARIMA). Kết quả mô hình tốt nhất thu được là ARIMA(0, 0, 0): ARMA Model Results Dep. Variable: y No. Observations: 234 Model: ARMA(0, 0) Log Likelihood 822.143 Method: css S.D. of innovations 0.007 Date: Thu, 12 Dec 2019 AIC -1640.286 Time: 04:39:00 BIC -1633.376 Sample: 0 HQIC -1637.500 coef std err z P>|z| [0.025 0.975] const 0.0001 0.000 0.255 0.799 -0.001 0.001 4.5. Kiếm tra yếu tố mùa vụTrong một số chuỗi thời gian thường xuất hiện yếu tố mùa vụ. Việc tìm ra chu kì và qui luật mùa vụ sẽ giúp cho mô hình dự báo chuẩn xác hơn. Yếu tố mùa vụ cũng không phải là một trong những yếu tố quá khó nhận biết. Chúng ta có thể dễ dàng phát hiện ra chúng thông qua đồ thị của chuỗi. Chẳng hạn bên dưới là dữ liệu sản xuất công nghiệp điện và khí đốt tại Hoa Kỳ từ năm 1985 đến năm 2019, với tần suất theo tháng. Chúng ta hãy cùng xem biểu diễn đồ thị của chuỗi. from vnquant.DataLoader import DataLoader loader = DataLoader(symbols="VN30", data = loader.download()
data.head()4 1 2 3 4 5 3 from vnquant.DataLoader import DataLoader loader = DataLoader(symbols="VN30", data = loader.download()
data.head()0 1 2 3 4 5 5 Ta nhận thấy chuỗi có chu kì là 1 năm. Nhu cầu tiêu thụ điện và gas tăng vào những tháng mùa đông do nhu cầu sưởi ấm tăng cao. Ngoài ra chúng ta có thể sử dụng một phép phân rã mùa vụ (seasonal decompose) để trích lọc ra các thành phần cấu thành nên chuỗi bao gồm: xu hướng (trend), mùa vụ (seasonal), phần dư (residual) như bên dưới: 1 2 0 1 2 3 4 5 7 Như vậy các thành phần đã được tách ra khá rõ ràng như thể hiện trong biểu đồ trên. Tiếp theo ta sẽ cùng hồi qui mô hình SARIMA. 4.6. Hồi qui mô hình SARIMAPhân chia tập train/test Đầu tiên để thuận tiện cho việc kiêm định mô hình dự báo chúng ta sẽ phân chia tập train/test sao cho năm 2019 sẽ được sử dụng làm dữ liệu test và dữ liệu còn lại được sử dụng để huấn luyện mô hình. 1 2 3 1 2 3 4 5 9 1 2 from vnquant import Plot Plot._vnquant_candle_stick(data = data,
1 Chúng ta sẽ cùng kiểm tra xem các đặc tính tự tương quan và tương quan riêng phần của chuỗi tiêu thụ điện và gas ra sao. Từ đó quyết định xem quá trình tự hồi qui và trung bình trượt của mô hình ARIMA nên nằm trong khoảng giá trị bao nhiêu và sử dụng phương pháp stepwise để tìm kiếm mô hình phù hợp nhất. 1 2 0 from vnquant import Plot Plot._vnquant_candle_stick(data = data,
3 Như vậy từ biểu đồ ta có thể lựa chọn bậc tự tương quan riêng phần PACF và tự tương quan ACF là các giá trị nhỏ hơn hoặc bằng 5. Do chuỗi có trend nên chúng ta sẽ lấy sai phân bậc 1 để tạo chuỗi dừng, hay nói cách khác bậc của intergration $d=1$. Ngoài ra chúng ta cần phải xác định thêm các bậc $(P, D, Q)$ của yếu tố mùa vụ được trích xuất từ chuỗi ban đầu. Để mô hình hiểu được chúng ta đang hồi qui trên mô hình SARIMA thì cần thiết lập tham số 1 2 3 6 và chu kì của mùa vụ 1 2 3 7. Chiến lược stepwise sẽ tự động tìm cho ta một mô hình tốt nhất dựa trên tham số đã thiết lập. 1 2 3 4 5 6 7 8 9 10 from vnquant import Plot Plot._vnquant_candle_stick(data = data,
5 from vnquant.DataLoader import DataLoader loader = DataLoader(symbols="VN30", data = loader.download()
data.head()4 from vnquant import Plot Plot._vnquant_candle_stick(data = data,
7Phương pháp stepwise đã giúp chúng ta tìm được mô hình SARIMA tốt nhất cho bài toán dự báo như bên dưới: Statespace Model Results Dep. Variable: y No. Observations: 408 Model: SARIMAX(1, 1, 1)x(2, 1, 2, 12) Log Likelihood -908.295 Date: Thu, 12 Dec 2019 AIC 1832.589 Time: 09:21:50 BIC 1864.420 Sample: 0 HQIC 1845.201 - 408 Covariance Type: opg coef std err z P>|z| [0.025 0.975] intercept -0.0022 0.001 -2.509 0.012 -0.004 -0.000 ar.L1 0.5137 0.041 12.601 0.000 0.434 0.594 ma.L1 -0.9769 0.013 -76.118 0.000 -1.002 -0.952 ar.S.L12 0.5542 0.128 4.319 0.000 0.303 0.806 ar.S.L24 -0.3529 0.053 -6.663 0.000 -0.457 -0.249 ma.S.L12 -1.2539 0.131 -9.593 0.000 -1.510 -0.998 ma.S.L24 0.5087 0.106 4.790 0.000 0.301 0.717 sigma2 5.5948 0.346 16.155 0.000 4.916 6.274 Ljung-Box (Q): 46.03 Jarque-Bera (JB): 12.84 Prob(Q): 0.24 Prob(JB): 0.00 Heteroskedasticity (H): 3.05 Skew: -0.02 Prob(H) (two-sided): 0.00 Kurtosis: 3.88 Warnings: Đó chính là mô hình SARIMA(p=1, d=1, q=1)(p=2, D=1, q=2, m=12). Mô hình cho kết quả khá tốt khi các hệ số hồi qui đều có ý nghĩa thống kê (toàn bộ cột 1 2 3 2 đều nhỏ hơn 0.05). 4.7. Dự báoSau khi đã tìm ra được mô hình ARIMA tốt nhất. Chúng ta sẽ dự báo cho khoảng thời gian tiếp theo. Dự báo cho chuỗi thời gian khá đặc thù và khác biệt so với các lớp mô hình dự báo khác vì giá trị time step liền trước sẽ được sử dụng để dự báo cho time step liền sau. Do đó đòi hỏi phải có một vòng lặp liên tiếp dự báo qua các bước thời gian. Rất may mắn là hàm 1 2 3 9 đã tự động giúp ta thực hiện việc đó. Ta sẽ chỉ phải xác định số lượng phiên tiếp theo muốn dự báo là bao nhiêu. Chẳng hạn bên dưới ta sẽ dự báo cho 10 tháng tới tương ứng trên tập mean = np.nanmean(r_t) r_t[0]=mean r_t[:5] 0 kèm theo giá trị ngưỡng tin cậy của nó. from vnquant import Plot Plot._vnquant_candle_stick(data = data,
8from vnquant import Plot Plot._vnquant_candle_stick(data = data,
9Chúng ta biết rằng một mô hình có thể fit với tập huấn luyện nhưng chưa chắc đã tốt khi dự báo. Chính vì thế cần kiểm tra chất lượng của mô hình trên tập dự báo. Trong mô hình phân loại chúng ta thường quan tâm đến tỷ lệ chính xác accuracy, trong trường hợp mẫu mất cân bằng thì precision, recall, f1 là những chỉ số đo lường độ chính xác khác được thay thế. Tuy nhiên với lớp mô hình dự báo thì sẽ sử dụng một tập hợp các tham số khác liên quan đến đo lường sai số giữa giá trị dự báo và giá trị thực tế. Đó là các chỉ số:
Các chỉ số này càng nhỏ thì chứng tỏ mô hình dự báo càng khớp với giá trị thực tế. Mô hình có thể được tính toán trên tập test như sau: 1 2 3 0 1 2 3 1 1 2 0 1 2 3 3 Giải thích ý nghĩa các thông số:
5. Tổng kếtNhư vậy tôi đã giới thiệu tới các bạn từ lý thuyết đến thực hành các mô hình ARIMA và SARIMA áp dụng trong các bài toán dự báo chuỗi thời gian. Đây là những lớp mô hình phổ biến và có độ chính xác cao không thua kém gì những mô hình deep learning như LSTM, CNN. Thậm chí một số trường hợp còn cho độ chính xác cao hơn. Tuy nhiên mọi sự so sánh đều chỉ là tương đối, tùy vào bài toán và tùy vào dữ liệu mà chúng ta sẽ cần phải lựa chọn phương pháp phù hợp. Hoặc thậm chí là thử nghiệm nhiều phương pháp khác nhau. Về nội dung của bài viết tôi hi vọng các bạn nắm được: |