Nội dung

Mô phỏng Monte Carlo (Monte Carlo Simulation)

1. Giới thiệu

Đã bao giờ bạn tự hỏi làm thế nào để giải quyết một bài toán quá phức tạp để có thể tính toán bằng các công thức toán học thông thường? Câu trả lời thường là: hãy để sự ngẫu nhiên làm việc đó cho bạn. Đó chính là triết lý đằng sau Mô phỏng Monte Carlo (Monte Carlo Simulation).

Tên gọi “Monte Carlo” được nhà toán học Stanislaw Ulam và John von Neumann đặt vào những năm 1940, lấy cảm hứng từ sòng bạc nổi tiếng ở Monaco, nơi chú của Ulam thường xuyên lui tới. Ý tưởng cốt lõi rất đơn giản: sử dụng lấy mẫu ngẫu nhiên (random sampling) lặp đi lặp lại để xấp xỉ kết quả của những bài toán khó hoặc không thể giải quyết bằng phương pháp giải tích (analytic).

Từ vật lý hạt nhân, đồ họa máy tính, cho đến tài chính định lượng, Monte Carlo đã trở thành một công cụ vô giá. Trong bài viết này, chúng ta sẽ tìm hiểu cách Monte Carlo hoạt động và cách nó được ứng dụng để đánh giá rủi ro tài chính.

2. Ý tưởng cốt lõi: Luật số lớn

Sức mạnh của Monte Carlo bắt nguồn từ một định lý nền tảng trong xác suất: Luật số lớn (Law of Large Numbers).

Giả sử bạn muốn tính kỳ vọng $E[f(X)]$ của một hàm số $f$ áp dụng lên biến ngẫu nhiên $X$. Thay vì vật lộn với các tích phân phức tạp, bạn có thể làm theo cách sau:

  1. Sinh ra một số lượng lớn $N$ các giá trị ngẫu nhiên $X_1, X_2, \dots, X_N$ độc lập và có cùng phân phối với $X$.
  2. Tính trung bình cộng của $f(X_i)$.
$$ \frac{1}{N} \sum_{i=1}^N f(X_i) \approx E[f(X)] $$

Ví dụ kinh điển: Ước lượng số $\pi$ Hãy tưởng tượng bạn ném phi tiêu ngẫu nhiên vào một bảng vuông có cạnh bằng $2$. Bên trong bảng vuông này vẽ một hình tròn bán kính $R=1$.

  • Diện tích hình vuông là $2 \times 2 = 4$.
  • Diện tích hình tròn là $\pi \times 1^2 = \pi$.

Nếu bạn ném phi tiêu hoàn toàn ngẫu nhiên, xác suất một mũi phi tiêu rơi trúng vào bên trong hình tròn sẽ là tỷ lệ diện tích của hình tròn so với hình vuông: $p = \frac{\pi}{4}$. Bằng cách ném một số lượng lớn phi tiêu ($N$), đếm số phi tiêu nằm trong hình tròn ($M$), ta có thể xấp xỉ $p \approx \frac{M}{N}$, từ đó suy ra $\pi \approx 4 \times \frac{M}{N}$. Chúng ta đã tính được một hằng số toán học phức tạp chỉ bằng cách… ném phi tiêu!

(Mẹo: Trong các bài toán tính tích phân hoặc hằng số như tính $\pi$, để tăng tốc độ hội tụ thay vì dùng số giả ngẫu nhiên thuần tuý, người ta ưu tiên dùng Chuỗi Sobol (Quasi-Monte Carlo) để các điểm phi tiêu phân bố đều đặn lấp kín mặt phẳng và không bị vón cục).

Hãy thử ngay bên dưới — bấm Bắt đầu để xem phi tiêu rơi (mô phỏng này đang sử dụng chuỗi Sobol ở chế độ chạy ngầm) và $\pi$ hội tụ dần:

Phi tiêu
0
Trong vòng tròn
0
π ≈
π thật = 3.14159265…

3. Sinh số ngẫu nhiên từ con số không

Trong thế giới thực, máy tính không thể thực sự ném phi tiêu. Chúng chỉ có thể sinh ra các số giả ngẫu nhiên (pseudo-random numbers). Nền tảng của mọi mô phỏng Monte Carlo là khả năng sinh ra các số ngẫu nhiên phân phối đều đặn trên khoảng $[0, 1]$, ký hiệu là $U \sim \text{Uniform}(0, 1)$. Phương pháp phổ biến để làm việc này là Linear Congruential Generator, một thuật toán toán học đơn giản nhưng tạo ra chuỗi số nhìn có vẻ ngẫu nhiên.

Tuy nhiên, trong tài chính và thống kê, chúng ta hiếm khi chỉ cần phân phối đều. Chúng ta cần phân phối Chuẩn (Normal), phân phối Student-t, v.v. Làm sao để biến đổi từ Uniform sang các phân phối khác?

Phương pháp biến đổi ngược (Inverse CDF Method)

Nếu bạn biết hàm phân phối tích lũy (CDF) $F(x)$ của một biến ngẫu nhiên, bạn có thể sinh ra nó bằng cách lấy $X = F^{-1}(U)$, với $U \sim \text{Uniform}(0, 1)$. Tuy nhiên, với phân phối Chuẩn, hàm $F^{-1}$ không có dạng biểu diễn đơn giản, khiến phương pháp này kém hiệu quả.

Phép biến đổi Box-Muller (Box-Muller Transform)

Để sinh ra các biến ngẫu nhiên phân phối Chuẩn $N(0,1)$, phép biến đổi Box-Muller là một giải pháp cực kỳ thanh lịch. Ý tưởng là thay vì sinh từng biến một trên trục tọa độ Descartes, ta làm việc trong hệ tọa độ cực.

$$ Z_1 = \sqrt{-2 \ln(U_1)} \cos(2\pi U_2) $$$$ Z_2 = \sqrt{-2 \ln(U_1)} \sin(2\pi U_2) $$

Trực giác: $\sqrt{-2 \ln(U_1)}$ đóng vai trò như bán kính (có bình phương tuân theo phân phối Exponential, tương đương với phân phối Chi-bình phương bậc tự do 2 của $Z_1^2 + Z_2^2$), và $2\pi U_2$ là góc quay phân phối đều từ $0$ đến $2\pi$.

Sinh biến ngẫu nhiên Student-t

Trong tài chính, lợi suất cổ phiếu hiếm khi hoàn toàn tuân theo phân phối Chuẩn. Chúng thường có “đuôi dày” (heavy tails) — nghĩa là các sự kiện cực đoan (như khủng hoảng thị trường) xảy ra thường xuyên hơn so với dự đoán của phân phối Chuẩn.

$$ T = \frac{Z}{\sqrt{V / \nu}} $$

sẽ tuân theo phân phối Student-t với $\nu$ bậc tự do.

Biến $\chi^2_\nu$ bản thân nó có thể được sinh ra bằng cách cộng tổng bình phương của $\nu$ biến phân phối chuẩn tắc độc lập (nếu $\nu$ nguyên), hoặc sử dụng các phép biến đổi từ phân phối Gamma.

4. Monte Carlo trong ước lượng rủi ro tài chính

Một trong những ứng dụng quan trọng nhất của Monte Carlo là trong quản trị rủi ro, cụ thể là để tính toán Giá trị chịu rủi ro (Value at Risk - VaR) và Thiếu hụt kỳ vọng (Expected Shortfall - ES).

Giả sử bạn đang nắm giữ một danh mục đầu tư và muốn biết: “Trong 20 ngày giao dịch tới, với độ tự tin 95%, mức lỗ tối đa của tôi là bao nhiêu?” Đó chính là VaR. Nếu tình huống tệ nhất xảy ra (vượt qua mốc VaR), mức lỗ trung bình sẽ là bao nhiêu? Đó là ES.

Để ước lượng các giá trị này, chúng ta thực hiện các bước sau:

  1. Mô hình hóa lợi suất: Sử dụng log-return (lợi suất logarit) vì chúng có tính chất cộng dồn. Lợi suất trong $h$ ngày ($h=20$) bằng tổng lợi suất của $h$ ngày đơn lẻ.
  2. Mô phỏng đường dẫn (Paths): Sinh ra $N$ kịch bản (ví dụ $N = 100,000$). Trong mỗi kịch bản, ta mô phỏng một chuỗi $h$ biến ngẫu nhiên đại diện cho lợi suất từng ngày. Chúng ta có thể dùng phân phối Chuẩn hoặc tốt hơn là phân phối Student-t (để bắt được hiện tượng đuôi dày).
  3. Tính toán phân phối tổn thất: Tổng hợp lại, ta được một tập hợp $N$ giá trị lợi suất giả định cho 20 ngày. Nếu nhân với giá trị danh mục, ta có phân phối của tổn thất.
  4. Ước lượng VaR và ES:
    • VaR(95%): Chính là bách phân vị (percentile) thứ 5 từ dưới lên của phân phối lợi suất giả định.
    • ES(95%): Là trung bình của tất cả các kịch bản tổn thất tồi tệ hơn mức VaR(95%).

Bằng cách sử dụng Student-t thay vì phân phối Chuẩn trong bước 2, mô hình Monte Carlo sẽ sinh ra nhiều giá trị cực đoan hơn ở phần đuôi, dẫn đến ước lượng VaR và ES bảo thủ và an toàn hơn, phản ánh đúng thực tế khắc nghiệt của thị trường tài chính. Sau đó, người ta thường dùng các phương pháp kiểm định ngược (backtesting) với cửa sổ cuộn (rolling-window) và kiểm định Kupiec để đánh giá xem mô hình VaR có thực sự nắm bắt đúng rủi ro hay không.

5. Ví dụ minh hoạ bằng Python

Dưới đây là một đoạn mã Python ngắn gọn thể hiện việc sử dụng Box-Muller để sinh số ngẫu nhiên và áp dụng Monte Carlo để ước lượng VaR và ES cho chân trời 20 ngày.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
import numpy as np

# 1. Thuật toán Box-Muller sinh phân phối chuẩn
def box_muller_normal(n):
    u1 = np.random.uniform(0, 1, n)
    u2 = np.random.uniform(0, 1, n)
    z1 = np.sqrt(-2 * np.log(u1)) * np.cos(2 * np.pi * u2)
    return z1

# Thông số giả định của thị trường
mu_daily = 0.0005     # Lợi suất kỳ vọng hằng ngày
sigma_daily = 0.01    # Độ lệch chuẩn hằng ngày
h = 20                # Số ngày (20-day horizon)
N = 10000             # Số đường dẫn mô phỏng (paths)

# 2. Mô phỏng N kịch bản lợi suất cho 20 ngày
# Mỗi kịch bản là tổng của h ngày lợi suất
# Ta cần tổng cộng N * h biến chuẩn
z_random = box_muller_normal(N * h).reshape((N, h))

# Tính log-return cho từng ngày rồi cộng dồn lại theo từng path
daily_returns = mu_daily + sigma_daily * z_random
cumulative_returns = np.sum(daily_returns, axis=1)

# 3. Tính toán rủi ro (VaR và ES tại mức 95%)
alpha = 0.05
# VaR là bách phân vị thứ 5 của phân phối lợi suất
var_95 = np.percentile(cumulative_returns, alpha * 100)

# ES là trung bình của các lợi suất nhỏ hơn hoặc bằng VaR
es_95 = np.mean(cumulative_returns[cumulative_returns <= var_95])

print(f"Ước lượng qua {N} kịch bản Monte Carlo (chân trời {h} ngày):")
print(f"VaR 95%: {-var_95:.2%}")
print(f"ES 95% : {-es_95:.2%}")

(Lưu ý: VaR và ES thường được báo cáo dưới dạng số dương đại diện cho tỷ lệ tổn thất).

6. Hội tụ và độ chính xác

Một câu hỏi hiển nhiên là: liệu chạy $10,000$ kịch bản đã đủ chưa? Tại sao không phải là $100,000$ hay $1,000,000$?

Độ chính xác của phương pháp Monte Carlo bị chi phối bởi phương sai của kết quả ước lượng. Sai số chuẩn (standard error) của ước lượng Monte Carlo tỷ lệ thuận với $1/\sqrt{N}$. Điều này có nghĩa là, để giảm một nửa sai số, bạn phải tăng số lượng mô phỏng lên gấp 4 lần.

Tốc độ hội tụ $O(1/\sqrt{N})$ này tuy khá chậm, nhưng lợi thế lớn nhất của nó là không phụ thuộc vào số chiều của bài toán (điều mà các phương pháp tích phân số trị truyền thống rất sợ). Trong thực tế, để mô phỏng tính toán VaR và ES ổn định, người ta thường sử dụng hàng chục ngàn đến hàng trăm ngàn kịch bản ($N = 100,000$).

Để tiết kiệm chi phí tính toán, các kỹ thuật giảm phương sai (variance reduction techniques) như biến đối nghịch (antithetic variates) hay lấy mẫu phân tầng (stratified sampling) thường được kết hợp để tăng độ chính xác mà không cần tăng $N$ lên quá cao.

7. Tổng kết

Mô phỏng Monte Carlo là một sự kết hợp tuyệt đẹp giữa lý thuyết xác suất và sức mạnh tính toán hiện đại. Bằng cách dùng máy tính “tung xúc xắc” hàng trăm nghìn lần thông qua các phép biến đổi toán học khéo léo như Box-Muller, chúng ta có thể dựng lên bức tranh toàn cảnh về tương lai của một danh mục đầu tư. Từ những con số vô hồn, Monte Carlo giúp chúng ta nhìn thấy hình hài của rủi ro, cho dù đó là VaR hay ES, đặc biệt là khi kết hợp với các phân phối đuôi dày như Student-t để tôn trọng sự khắc nghiệt của thị trường.