李斯儒 譚靜
DOI:10.19850/j.cnki.2096-4706.2021.08.004
摘? 要:眾所周知,隨機(jī)游動(dòng)作為一種特殊的馬爾科夫鏈在諸多領(lǐng)域有著廣泛的應(yīng)用,而研究系統(tǒng)的閾值狀態(tài)的首達(dá)概率則是重中之重。文章首先對(duì)對(duì)稱一維隨機(jī)游動(dòng)某狀態(tài)的首達(dá)概率、首達(dá)時(shí)進(jìn)行計(jì)算分析;然后采用了對(duì)遞推公式求母函數(shù)的方法求解對(duì)稱一維隨機(jī)游動(dòng)的首次返回概率,最后通過蒙特卡洛方法求其首次返回概率的模擬值與理論值對(duì)照,通過大樣本容量的計(jì)算證實(shí)理論解的精確性。
關(guān)鍵詞:隨機(jī)游動(dòng);首達(dá)時(shí)間;首達(dá)概率;馬爾科夫鏈
中圖分類號(hào):O211.1? ? 文獻(xiàn)標(biāo)識(shí)碼:A? ? 文章編號(hào):2096-4706(2021)08-0013-04
Research and Analysis on the First Arrival Probability in Random Walk
LI Siru,TAN Jing
(Nanhang Jincheng College,Nanjing? 211156,China)
Abstract:As we all know,the random walk is a special Markov chain that has wide applications in many fields,and the research on the first arrival probability of the threshold state of the system is the most important thing. First,the paper calculates and analyzes the first arrival probability and first arrival time of a certain state of a symmetric one-dimensional random walk;then uses the method of finding the generating function of the recurrence formula to solve the first return probability of the symmetric one-dimensional random walk,and finally uses Monte Carlo method compares the simulated value of the probability of its first return with the theoretical value,and confirms the accuracy of the theoretical solution through the calculation of a large sample volume.
Keywords:random walk;first arrival time;first arrival probability;Markov chain
0? 引? 言
在馬爾科夫鏈中有一類非常特殊的隨機(jī)過程被稱為隨機(jī)游走[1](Random Walk,RW)。隨機(jī)游走是一種數(shù)學(xué)統(tǒng)計(jì)模型,由一連串的軌跡組成,其中的每一步都是隨機(jī)的。隨機(jī)游走的理論最早在1905年由Karl Pearson提出[2],就像一個(gè)人喝醉了之后酒后漫步一樣,一定時(shí)間中所記錄的軌跡就是隨機(jī)游動(dòng)。
一維隨機(jī)游走作為一種特殊的馬爾科夫鏈,可以有如下定義:假設(shè)游走者在一維坐標(biāo)軸上從起始點(diǎn)x處出發(fā),與時(shí)間無關(guān)以概率p向左移動(dòng)一個(gè)單位(以1-p的概率向右移動(dòng)一個(gè)單位,很多時(shí)候取p=0.5)。不妨記第游走者的第k次游動(dòng)為Xk,n時(shí)刻游走著的位置是Ln,于是有:
Ln=x+X1+X2+…+Xk+…+Xn
這里X1,X2,…,Xn要滿足P(Xi=1)=p=1-P(Xi=-1)。
再列舉兩個(gè)常用的隨機(jī)游走問題:假設(shè)你是一位賭徒,你在賭場(chǎng)賭博,與時(shí)間無關(guān),你每一場(chǎng)贏的概率為p,輸?shù)糍€局的概率為1-p。你一開始身上帶了k元錢,贏一場(chǎng)加一元錢,輸一場(chǎng)少一元錢,在不能借錢賭博的情況下問你輸光所有賭錢的概率是多少,還有平均賭博幾場(chǎng)會(huì)輸光所有本錢。或是假設(shè)一個(gè)只會(huì)向前或者向后移動(dòng)的醉鬼站在一個(gè)一邊是懸崖的道路上,他再向前跨一步就會(huì)跌入懸崖,每秒醉鬼向前一步或者是向后一步的概率都是50%,問醉鬼掉下懸崖的概率是多少。
1? 一維隨機(jī)游走的首次返回問題求解
1.1? 游動(dòng)點(diǎn)返回原點(diǎn)的理論概率
在這個(gè)特殊的一維對(duì)稱隨機(jī)游動(dòng)中,假設(shè)游動(dòng)點(diǎn)從0點(diǎn)出發(fā)經(jīng)過m步返回0點(diǎn),顯然m必須是偶數(shù)才有可能返回,所以我們可以只計(jì)算游動(dòng)點(diǎn)經(jīng)過2n步之后返回0點(diǎn)的概率[1]。計(jì)算2n步之后首次返回0點(diǎn)的概率雖然有些難度,但我們很容易就能寫出2n步之后游動(dòng)點(diǎn)停在0點(diǎn)的概率:
我們以此為計(jì)算的起始點(diǎn),來分析2n步之后游動(dòng)點(diǎn)首次回到0點(diǎn)的概率[3]。
2? 數(shù)軸上的隨機(jī)游動(dòng)仿真與蒙特卡洛方法
2.1? 蒙特卡洛方法
蒙特卡洛方法是一種一般使用偽隨機(jī)數(shù),通過模擬的方式抽取系統(tǒng)狀態(tài)來解決很多計(jì)算問題的方法,與之對(duì)應(yīng)的是確定性算法。
事實(shí)上很早以前布豐通過投針實(shí)驗(yàn)[7]來計(jì)算π就算是一種蒙特卡洛方法,目前在計(jì)算復(fù)雜的積分時(shí)經(jīng)常會(huì)先選擇一個(gè)能完全包含積分曲線的長(zhǎng)方形區(qū)域,再在區(qū)域內(nèi)大量平均撒點(diǎn)將落在積分區(qū)域內(nèi)點(diǎn)數(shù)比上全部的撒點(diǎn)次數(shù)就可以得到曲線與x軸包圍的面積與長(zhǎng)方形區(qū)域的面積之比,從而求得積分?jǐn)?shù)值解。
2.2? 蒙特卡洛模擬首達(dá)概率
在MATLAB環(huán)境中我用蒙特卡洛方法模擬了1.1中對(duì)稱隨機(jī)游走中游動(dòng)點(diǎn)從0點(diǎn)出發(fā)經(jīng)過2k次游動(dòng)是否首次返回0點(diǎn)的問題,MATLAB蒙特卡洛實(shí)驗(yàn)代碼為:
1n1=10000000;? ? ? ? ? ? ? ? %每次中層嵌套的蒙特卡洛要模擬一千萬次
i0=1;? ? ? ? ? ? ? ? ? ? ? ? %初始化i0=1
A=zeros(2,200);? ? ? ? ? ? ? ? %創(chuàng)造一個(gè)2行200列的矩陣
neighbour=[-1 1];? ? ? ? ? ? ? ?%游動(dòng)點(diǎn)游動(dòng)矩陣,數(shù)軸坐標(biāo)+1或者是-1
for n=2:2:400? ? ? ? ? ? %最外層循環(huán)確定在一千萬次模擬中游動(dòng)的步數(shù)(n為偶數(shù))
z=0;? ? ? ? ? ? ? %初始化z,z為一千萬次實(shí)驗(yàn)中成功經(jīng)過n步首次返回原點(diǎn)的次數(shù)
for i1=1:n1? ? ? ? ? ? ? ? ?%中層循環(huán)蒙特卡洛實(shí)驗(yàn)開始
x=0;? ? ? ? ? ? ? ? ? ? ?%x代表坐標(biāo)
y=0;? ? ? ? ? ? ? ? ? ? ? %y表示在本次游動(dòng)中一共返回0點(diǎn)的次數(shù)
for i2=1:n
r=floor(1+2*rand());? ? ?%偽隨機(jī)數(shù)來判定x是向左游動(dòng)一個(gè)單位還是向右移動(dòng)一個(gè)單位
x=x+neighbour(1,r);
if x==0? ? ? ? ? ? ? ?%游動(dòng)過程中x每返回一次原點(diǎn)y自增1
y=y+1;
else end
end? ? ? ? ? ? ? ? ? ? ? %內(nèi)層循環(huán)結(jié)束,一次n步的隨機(jī)游走完成
if (x==0)&(y==1)? ? ? ? ? ? ? ?%一次蒙特卡洛模擬結(jié)束,當(dāng)且僅當(dāng)%
z=z+1;? ? ? ? ? ? ? ? ? ? ? %x停在0點(diǎn)且只有這次停在0點(diǎn)時(shí)z自增1
else end
end? ? ? ? ? ? ? ? ? ? ? ? ? %記錄完一千萬次實(shí)驗(yàn)中有幾次符合要求
A(:,i0)=[n;z];? ? ? ? ? ? ? ? ? ? ?%將實(shí)驗(yàn)數(shù)據(jù)放置在矩陣A中
i0=i0+1;
end? ? ? ? ? ? ? ? ? ? ? ? ? ? ?%外層循環(huán)結(jié)束蒙特卡洛實(shí)驗(yàn)完成
xlswrite('A.xlsx', A);? ? ? ? ? ? ? ?%將實(shí)驗(yàn)數(shù)據(jù)導(dǎo)出到Excel中防止丟失
B=zeros(2,200);
B=A;
i0=1;
for n=2:2:400
w=(1./(n-1)).*nchoosek(n,n/2).*((1/2)^n);? ?%將之前計(jì)算的概率分布的理論值付給w
B(:,i0)=[n;log10(w)];? ?%由于某些概率過低將理論值取以10為底數(shù)的對(duì)數(shù)儲(chǔ)存
i0=i0+1; y1=A(2,:);
end
y1=y1./10000000;? ? ?%求出蒙特卡洛法從0點(diǎn)出發(fā)經(jīng)過n次游動(dòng)首次返回0點(diǎn)的概率
for i1=1:200
q=y1(:,i1);
y1(:,i1)=log10(q);? ? %同理將蒙特卡洛方法求得的概率的數(shù)值解去對(duì)數(shù)儲(chǔ)存
end
y2=B(2,:);
x0=A(1,:);
plot(x0,y2)? ? ? ? ? ? ? ? ? ?%繪制出橫坐標(biāo)為游走步數(shù)縱坐標(biāo)為從0點(diǎn)出發(fā)經(jīng)過
hold on;? ? ? ? ? ? ? % n次游動(dòng)首次返回0點(diǎn)的概率的對(duì)數(shù)的理論值和模擬實(shí)驗(yàn)值
plot(x0,y1)? ? ? %理論值為光滑曲線,蒙特卡洛實(shí)驗(yàn)值為震蕩曲線
運(yùn)行結(jié)果如圖1、圖2所示。
3? 結(jié)? 論
綜合全文,本文使用遞推和母函數(shù)性質(zhì)的方法著重研究了一維對(duì)稱隨機(jī)游動(dòng)中粒子的首次返回原點(diǎn)的概率分布,給出MATLAB蒙特卡洛方法的具體代碼,經(jīng)過多次試驗(yàn)將理論精確解與模擬數(shù)值解進(jìn)行比對(duì)。
在一維對(duì)稱隨機(jī)游動(dòng)的狀態(tài)首達(dá)概率的理論計(jì)算時(shí),由于使用的遞推公式較為繁瑣導(dǎo)致后面不得不使用數(shù)學(xué)歸納法證明f(x)的通項(xiàng)公式。
另外在進(jìn)行一維對(duì)稱隨機(jī)游動(dòng)蒙特卡洛方法的仿真中,若實(shí)驗(yàn)次數(shù)稍少就會(huì)出現(xiàn)理論解與實(shí)際解不匹配,震蕩嚴(yán)重的情況。本文的代碼運(yùn)行用時(shí)過長(zhǎng)超過8小時(shí),今后應(yīng)采取效率更高的算法進(jìn)行仿真。
參考文獻(xiàn):
[1] 陸大絟.隨機(jī)過程及其應(yīng)用 [M].北京:清華大學(xué)出版社,1986.
[2] 何聲武.隨機(jī)過程導(dǎo)論 [M].上海:華東師范大學(xué)出版社,1989.
[3] 李林海,常建平.隨機(jī)信號(hào)分析 [M].北京:科學(xué)出版社,2006.
[4] 陳懷琛,吳大正,高西全.MATLAB及在電子信息課程中的應(yīng)用:第4版 [M].北京:電子工業(yè)出版社,2013.
[5] 何書元.概率引論 [M].北京:高等教育出版社,2011.
[6] 陳紀(jì)修.數(shù)學(xué)分析:第2版 [M].北京:高等教育出版社,2004.
[7] SPITZER F. Principles of RANDOM WALK [M].New York:Springer,2001.
作者簡(jiǎn)介:李斯儒(1999—),男,漢族,江蘇揚(yáng)州人,本科在讀,研究方向:信息與信號(hào)處理;通訊作者:譚靜(1978—),女,漢族,江蘇南通人,副教授,碩士,研究方向:雷達(dá)信號(hào)處理。
收稿日期:2021-03-02