BỘ GIÁO DỤC VÀ ĐÀO TẠO
TRƯỜNG ĐẠI HỌC SƯ PHẠM HÀ NỘI 2
—————————————-

HOÀNG THỊ ANH SANG

PHƯƠNG PHÁP XÍCH MARKOV
MONTE CARLO VÀ ỨNG DỤNG

LUẬN VĂN THẠC SĨ TOÁN HỌC

HÀ NỘI, 2017


BỘ GIÁO DỤC VÀ ĐÀO TẠO
TRƯỜNG ĐẠI HỌC SƯ PHẠM HÀ NỘI 2
—————————————-

HOÀNG THỊ ANH SANG

PHƯƠNG PHÁP XÍCH MARKOV
MONTE CARLO VÀ ỨNG DỤNG
Chuyên ngành: Toán ứng dụng
Mã số: 60 46 01 12

LUẬN VĂN THẠC SĨ TOÁN HỌC

Người hướng dẫn khoa học: TS. NGÔ HOÀNG LONG

HÀ NỘI, 2017


Học viên cao học

Hoàng Thị Anh Sang

ii


Mục lục
Lời cảm ơn

i

Lời cam đoan

ii

Danh mục các kí hiệu, chữ viết tắt

v

Danh sách bảng

vi

Danh sách hình vẽ

vii

MỞ ĐẦU


Mô phỏng Xích Markov . . . . . . . . . . . . . . . . . . . . . .

7

1.2.3

Xích Markov tối giản và phi tuần hoàn . . . . . . . . . . . . . .

12

1.2.4

Phân phối dừng . . . . . . . . . . . . . . . . . . . . . . . . . . .

13

1.2.5

Xích Markov khả nghịch . . . . . . . . . . . . . . . . . . . . . .

15

2 Phương pháp Monte-Carlo
2.1

16

Số ngẫu nhiên và phương pháp sinh số ngẫu nhiên . . . . . . . . . . . .

16

Phương pháp Monte-Carlo . . . . . . . . . . . . . . . . . . . . . . . . .

24

2.2.1

Phương pháp lấy mẫu phân tầng . . . . . . . . . . . . . . . . .

24

2.2.2

Chọn mẫu quan trọng . . . . . . . . . . . . . . . . . . . . . . .

27

2.2.3

Các số ngẫu nhiên chung . . . . . . . . . . . . . . . . . . . . . .

28

3 Xích Markov Monte Carlo

32

3.1

Phương pháp Xích Markov Monte Carlo . . . . . . . . . . . . . . . . .


Tài liệu tham khảo

60

PHỤ LỤC: Code mô phỏng bằng Python

61

iv


Danh mục các kí hiệu, chữ viết tắt

E[X]

Kì vọng của biến X

Var[X]

Phương sai của biến X

ψ(x)

Hàm khởi tạo

φ(x)

Hàm cập nhật

dTV


Danh sách hình vẽ
3.1

Một cấu hình thực hiện được (lựa chọn bất kì theo độ đo xác suất µG ). . . .

33

3.2

Đồ thị G với 4 đỉnh. . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

57

3.3

Đồ thị G với 5 đỉnh. . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

58

vii


MỞ ĐẦU

1. Lý do chọn đề tài
Trong nhiều ứng dụng thực tế, bài toán thường được quy về việc xác định kì vọng
của hàm của một đại lượng ngẫu nhiên theo phân phối nhất định. Khi đó, người ta
thường tìm cách xây dựng dãy biến ngẫu nhiên độc lập có phân phối này trên máy tính
và sử dụng luật số lớn để xấp xỉ kì vọng cần tính. Tuy nhiên, trong nhiều trường hợp,


5. Phương pháp nghiên cứu
❼ Nghiên cứu lý thuyết.
❼ Nghiên cứu thực nghiệm mô phỏng trên máy tính.

6. Dự kiến đóng góp mới
Luận văn hệ thống hoá và làm rõ cơ sở lý thuyết của phương pháp xích Markov
Monte Carlo và ứng dụng phương pháp này trong một vấn đề thực tiễn.

NỘI DUNG
Luận văn tốt nghiệp sẽ được chia làm ba chương cộng với phần Mở đầu, Kết luận
và Tài liệu tham khảo. Nội dung cụ thể trong Chương 1, Chương 2, Chương 3 của luận
văn dự kiến sẽ được phân bổ như sau:
Chương 1: Xích Markov
1.1. Một số khái niệm cơ bản của lý thuyết xác suất
ix


1.2. Xích Markov
Chương 2: Phương pháp Monte-Carlo
2.1. Số ngẫu nhiên và phương pháp sinh số ngẫu nhiên
2.2. Phương pháp Monte-Carlo
Chương 3: Xích Markov Monte Carlo
3.1. Phương pháp xích Markov Monte Carlo
3.2. Một số phương pháp cải tiến
3.3. Mô phỏng

x



P(A ∩ B)
.
P(B)

Giải thích trực quan của P(A|B) là xem xét khả năng xảy ra của biến cố A khi ta
đã biết biến cố B đã xảy ra.
Hai biến cố A và B được gọi là độc lập nếu P(A ∩ B) = P(A)P(B).
Tổng quát hơn, các biến cố A1 , A2 , . . . , Ak được gọi là độc lập nếu với mỗi l ≤ k và

1


i1 , . . . , il ∈ {1, 2, . . . , k} bất kì với i1 < i2 < · · · < il ta có
l

P(Ai1 ∩ Ai2 ∩ · · · ∩ Ail ) =

P(Ain ).
n=1

Chú ý rằng nếu P(B) > 0, A và B độc lập thì ta có P(A|B) = P(A), nghĩa là sự xảy
ra của biến cố B không làm ảnh hưởng đến khả năng xảy ra của biến cố A.
Biến ngẫu nhiên là đại lượng mà giá trị của nó phụ thuộc vào kết quả ngẫu nhiên
của phép thử. Thông thường, một biến ngẫu nhiên là một giá trị thực, trong trường
hợp này nó là một hàm X : Ω → R. Tuy nhiên, xem xét các biến ngẫu nhiên trong
một ý nghĩa tổng quát hơn, đó là những hàm X : Ω → S, S là tập bất kì.
Biến cố A được xác định dưới dạng biến ngẫu nhiên X nếu ta có thể có hoặc
không có A xảy ra từ giá trị của X. Ví dụ các biến cố được xác định dưới dạng biến
ngẫu nhiên X là
A = {X ≤ 4.7} = {ω ∈ Ω : X(ω) ≤ 4.7}

phần tiếp theo là một lớp đặc biệt của quá trình ngẫu nhiên.
Ta xét hai loại biến ngẫu nhiên giá trị thực: biến ngẫu nhiên rời rạc và biến ngẫu
nhiên liên tục. Biến ngẫu nhiên rời rạc thì các giá trị của chúng là tập con hữu hạn
hoặc đếm được của R.
Biến ngẫu nhiên liên tục X là một biến ngẫu nhiên mà tồn tại một hàm mật
độ fX : R → [0, ∞) sao cho
x

fX (x)dx = FX (x) = P(X ≤ x),
−∞

với mọi x ∈ R. Một ví dụ rất nổi tiếng của biến ngẫu nhiên liên tục X sinh ra bằng
cách cho X bởi hàm mật độ Gaussian
fX (x) = √

1

2 /2σ 2

2πσ 2

e−(x−µ)

,

với các thông số µ và σ > 0. Tuy nhiên, các biến ngẫu nhiên liên tục đó có phân phối
đều trên [0, 1], khi đó ta có hàm mật độ

fX (x) =


3


khoảng I có độ dài a trong đoạn [0, 1], ta có: P(X ∈ I) = a.
Kì vọng (hoặc giá trị kì vọng, hoặc trung bình) E[X] của một biến ngẫu nhiên
giá trị thực X là, trong một nghĩa nào đó, giá trị "trung bình" mong đợi từ x. Nếu X
là một biến ngẫu nhiên liên tục với hàm mật độ fX (x), khi đó kì vọng được xác định
như sau

∞

E[X] =

xfX (x)dx.
−∞

Trong trường hợp X có phân phối đều trên [0, 1] thì
1

1
xdx = .
2

E[X] =
0

Trong trường hợp X là biến ngẫu nhiên có giá trị nguyên không âm, kì vọng được xác
định như sau

∞

Nếu c là một hằng số thì
E[cX] = cE[X].
4

(1.3)


Đối với các phương sai ta có
Var[cX] = c2 Var[X].
Khi X1 , X2 , . . . , Xn độc lập
Var[X1 + X2 + · · · + Xn ] = Var[X1 ] + Var[X2 ] + · · · + Var[Xn ].
Định lý 1.1. (Bất đẳng thức Chebyshev) Cho X là một biến ngẫu nhiên với kì
vọng µ và phương sai σ 2 . Cho a > 0 bất kì, ta có xác suất P[|X − µ| ≥ a] của độ lệch
trung bình nhỏ nhất, thỏa mãn
P(|X − µ| ≥ a) ≤

σ2
.
a2

Định lý 1.2. (Luật Số Lớn) Cho X1 , X2 , ... là các biến ngẫu nhiên độc lập khả tích
và có cùng phân phối với trung bình µ. Cho Mn biểu thị giá trị trung bình của n biến
Xi đầu tiên, có nghĩa là:
Mn =

1
(X1 + X2 + · · · + Xn ).
n

Khi đó, Mn hội tụ hầu chắc chắn tới µ.

1, 2, . . . , k}. Một quá trình ngẫu nhiên (X0 , X1 , ...) với không gian trạng thái hữu hạn
5


S = {s1 , s2 , . . . , sk } là một xích Markov (thuần nhất) với ma trận chuyển P ,
nếu với tất cả n, tất cả i, j ∈ {1, 2, . . . , k} và với mọi i0 , i1 , . . . , in−1 ∈ {1, 2, . . . , k} ta
có
P(Xn+1 = sj |X0 = si0 , X1 = si1 , . . . , Xn−1 = sin−1 , Xn = si )
= P(Xn+1 = sj |Xn = si )
= Pi,j .
Mỗi phần tử của ma trận chuyển P được gọi là xác suất chuyển. Xác suất chuyển
Pi,j là xác suất hệ sẽ ở trạng thái sj vào "ngày mai" với điều kiện hệ đang trong trạng
thái si vào "ngày hôm nay".
Mỗi ma trận chuyển thỏa mãn
Pi,j ≥ 0 với mọi i, j ∈ {1, 2, . . . , k},
và

(1.4)

k

Pi,j = 1 với mọi i ∈ {1, 2, . . . , k}.

(1.5)

j=1

Tính chất (1.5) có nghĩa là
P(Xn+1 = s1 |Xn = si ) + P(Xn+1 = s2 |Xn = si ) + · · ·
+ P(Xn+1 = sk |Xn = si ) = 1.

6


Định lý 1.5. Cho xích Markov (X0 , X1 , ...) với không gian trạng thái S = {s1 , s2 , . . . , sk },
phân phối ban đầu µ(0) và ma trận chuyển P , với n bất kì thì phân phối µ(n) tại thời
điểm n thỏa mãn
µ(n) = µ(0) P n .

(1.6)

Định nghĩa 1.6. Cho P (1) , P (2) , ... là một dãy các ma trận cỡ k × k, mỗi ma trận
n

thỏa mãn các điều kiện: Pi,j ≥ 0 với mọi i, j ∈ {1, 2, . . . , k} và

Pi,j = 1 với i ∈
j=1

{1, 2, . . . , k}. Một quá trình ngẫu nhiên (X0 , X1 , ...) với không gian trạng thái hữu
hạn S = {s1 , s2 , ..., sk } được gọi là một xích Markov không thuần nhất với các
ma trận chuyển P (1) , P (2) , ... nếu với mọi n, với mọi i, j ∈ {1, 2, . . . , k} và tất cả
i0 , i1 , . . . , in−1 ∈ {1, 2, . . . , k} ta có
P(Xn+1 = sj |X0 = si0 , X1 = si1 , . . . , Xn−1 = sin−1 , Xn = si )
= P(Xn+1 = sj |Xn = si )
(n+1)

= Pi,j

.


Hàm khởi tạo ψ : [0, 1] → S là một hàm từ khoảng đơn vị tới không gian trạng
thái S mà chúng ta sử dụng để tạo ra các giá trị X0 ban đầu. Giả sử
(i) ψ là hàm hằng từng khúc (tức là [0, 1] có thể chia thành hữu hạn những khoảng
con sao cho ψ là hằng số trên mỗi khoảng con),
(ii) với mỗi s ∈ S, tổng độ dài của các khoảng con thoả mãn ψ(x) = s bằng µ(0) (s).
Một cách diễn đạt khác của (ii) đó là
1

I{ψ(x)=s} dx = µ(0) (s),

(1.7)

0

với mỗi s ∈ S; ở đây I{ψ(x)=s} được gọi là hàm chỉ tiêu của {ψ(x) = s}, nghĩa là

I{ψ(x)=s} =



1 nếu ψ(x) = s

0 ngược lại.

Với điều kiện ta có một hàm ψ, ta có thể tạo ra X0 từ một số ngẫu nhiên đầu tiên U0
bằng cách đặt X0 = ψ(U0 ). Điều này cho ta phân phối chính xác của X0 , vì với s ∈ S
bất kì ta có
1

I{ψ(x)=s} dx = µ(0) (s).



..


.

với x ∈ 0, µ(0) (s1 )
với x ∈ µ(0) (s1 ), µ(0) (s1 ) + µ(0) (s2 )

i−1


si với x ∈







..


.







µ(0) (sj ) −

I{ψ(x)=si } dx =
j=1

0

µ(0) (sj ) = µ(0) (si ).
j=1

với i = 1, 2, . . . , k. Điều này có nghĩa là ψ xác định như (1.8) là một hàm khởi tạo hợp
lí cho xích Markov (X0 , X1 , ...).
Ta đã biết cách tạo ra X0 và bây giờ ta nghiên cứu cách tạo ra Xn+1 từ Xn với mọi
n, sau đó ta có thể dùng quy trình này để sinh ra xích (X0 , X1 , ...). Ta sẽ dùng số ngẫu
nhiên Un+1 và hàm cập nhật φ : S × [0, 1] → S, đầu vào là trạng thái s ∈ S và một
số nằm giữa 0 và 1 sẽ tạo ra một trạng thái khác s ∈ S là đầu ra. Tương tự như hàm
ψ thì hàm φ cần thỏa mãn:
(i) với si cố định, hàm φ(si , x) hằng từng khúc theo biến x,
(ii) với mỗi si , sj ∈ S cố định, tổng độ dài các khoảng mà trên đó φ(si , x) = sj là
bằng Pi,j .
Giống như hàm khởi tạo thì tính chất (ii) có thể được viết lại là
1

I{φ(si ,x)=sj } dx = Pi,j .
0

9

(1.9)








s2





.


..

với x ∈ 0, Pi,1
với x ∈ Pi,1 , Pi,1 + Pi,2
(1.11)

j

j−1


sj với x ∈



Đó là một hàm cập nhật hợp lí, chú ý với mỗi si , sj ∈ S ta có
1

j

j−1

Pi,l −

I{φ(si ,x)=sj } dx =
l=1

0

Pi,l
l=1

= Pi,j .
Vì vậy, ta có công thức hoàn thiện cho mô phỏng của một xích Markov: Trước hết xây
dựng hàm khởi tạo hợp lí và hàm cập nhật hợp lí ψ và φ, và khi đó

10





X0




Ta chỉ ra làm thế nào phương pháp trên có thể khái quát với mô phỏng xích Markov
không thuần nhất. Cho X0 , X1 , ... là một xích Markov không thuần nhất với không gian
trạng thái S = {s1 , . . . , sk }, phân phối ban đầu µ(0) và các ma trận chuyển P (0) , P (1) , ....
Ta có được hàm khởi tạo ψ và bắt đầu với giá trị X0 như trong trường hợp thần nhất.
Việc cập nhật thực hiện giống như trong trường hợp thuần nhất, từ xích không thuần
nhất ta cần một số hàm cập nhật khác nhau: φ(0) , φ(1) , ... và từ đó ta có:
1
(n)

I{φ(n) (si ,x)=sj } dx = Pi,j ,
0

với mỗi n và mỗi si , sj ∈ S. Ta có hàm khái quát của (1.11)


(n)


s1 với x ∈ 0, Pi,1





(n)
(n)
(n)








k−1

(n)


Pi,l , 1 .
sk với x ∈
l=1

11

(1.12)


Xích Markov không thuần nhất được mô phỏng bằng cách thiết lập



X0 = ψ(U0 )







được từ tất cả xích khác. Cụ thể, xét xích Markov (X0 , X1 , ...) với không gian trạng
thái S = {s1 , s2 , . . . , sk } và ma trận chuyển P . Ta nói trạng thái si thông với trạng
thái sj , kí hiệu si → sj , nếu tồn tại n ≥ 0 sao cho
P(Xm+n = sj |Xm = si ) > 0.
Nếu si → sj và sj → si thì ta nói trạng thái si và sj là liên thông, và viết si ↔ sj .
Định nghĩa 1.8. Xích Markov (X0 , X1 , ...) với không gian trạng thái S = {s1 , s2 , . . . , sk }
và ma trận chuyển P được gọi là tối giản nếu mọi si , sj ∈ S ta có si ↔ sj . Trái lại thì
xích Markov đó được gọi là không tối giản.
Một cách phát biểu khác của định nghĩa xích tối giản nếu si , sj ∈ S bất kì ta có
thể tìm được n sao cho (P n )i,j > 0.
Ta xét khái niệm phi tuần hoàn. Cho một tập hữu hạn hoặc không hữu hạn các số
nguyên dương {a1 , a2 , ...}, ta viết ƯCLN{a1 , a2 , ...} là ước chung lớn nhất của a1 , a2 , ....
Chu kì d(si ) của trạng thái si ∈ S được xác định như sau
d(si ) = ƯCLN{n ≥ 1 : (P n )i,i > 0}.
Chu kì của si là ước chung lớn nhất của tập thời gian mà xích trở lại tới si , ta bắt
đầu với X0 = si , nếu d(si ) = 1 thì ta nói trạng thái si là phi tuần hoàn.
Định nghĩa 1.9. Xích Markov gọi là phi tuần hoàn nếu tất cả các trạng thái của
xích Markov là phi tuần hoàn. Trái lại xích Markov đó được gọi là tuần hoàn.
12


Định lý 1.10. Giả sử ta có một xích Markov không tuần hoàn (X0 , X1 , ...) với không
gian trạng thái S = {s1 , s2 , ..., sk } và ma trận chuyển P . Khi đó tồn tại N < ∞
sao cho:
(P n )i,i > 0,
với i ∈ {1, 2, . . . , k} và với mọi n ≥ N .
Hệ quả 1.11. Cho (X0 , X1 , ...) là một xích Markov tối giản và phi tuần hoàn với không
gian trạng thái S = {s1 , s2 , . . . , sk } và ma trận chuyển P . Khi đó tồn tại M < ∞ sao
cho
(P n )i,j > 0,

hơn là cho xích Markov).
Ta xét 3 vấn đề: sự tồn tại của phân phối dừng, tính duy nhất của phân phối
dừng và sự hội tụ tới dừng bắt đầu từ phân phối ban đầu bất kì.
13