第五章 蒙特卡羅方法與隨機(jī)數(shù)

第五章蒙特卡羅方法與隨機(jī)數(shù)

第五章蒙特卡羅方法與隨機(jī)數(shù)

蒙特卡羅（MonteCarlo）方法，又稱(chēng)隨機(jī)抽樣法或統(tǒng)計(jì)實(shí)驗(yàn)法，是基于對(duì)概率模型的觀察或抽樣試驗(yàn)得到問(wèn)題近似解的方法。本章討論：

1、蒙特卡羅方法的基本思想；

2、隨機(jī)數(shù)的產(chǎn)生和檢驗(yàn)。

蒙特卡羅（Monte-Carlo）方法

蒙特卡羅方法以概率統(tǒng)計(jì)理論為其主要理論基礎(chǔ)，以隨機(jī)抽樣為主要手段。先建立一個(gè)概率模型或隨機(jī)過(guò)程，使它的參數(shù)等于問(wèn)題的解；然后通過(guò)對(duì)模型或過(guò)程的觀察或抽樣試驗(yàn)來(lái)計(jì)算所求隨機(jī)參數(shù)計(jì)算所求隨機(jī)參數(shù)的統(tǒng)計(jì)特征，最后給出所求解的近似值。

若射擊彈著點(diǎn)到靶心的距離r，則得分為Y，Y與r的關(guān)系為g(r)，即

r用統(tǒng)計(jì)試驗(yàn)的方法求其數(shù)學(xué)期望

N次射擊的彈著點(diǎn)股票價(jià)格的預(yù)測(cè)研究

（ITO）過(guò)程

維納過(guò)程（Wienerprocess）

服從正態(tài)分布

σ＝0.03，Δt=0.01,股票的初始值S=30元，用變換

每周期開(kāi)始價(jià)隨機(jī)樣本v1∈N(0,1)隨機(jī)樣本v2股價(jià)變化ΔS每周期結(jié)束價(jià)30.000-0.12675-0.0026-0.0780829.92229.9220.7216090.0228480.68366430.60630.606-1.08988-0.0315-0.9639729.64229.642-0.60426-0.01693-0.5017729.14029.1401.9316120.0591481.72357530.86330.863-0.12521-0.00256-0.078930.78530.785-1.266-0.03678-1.1322629.65229.6521.5679780.0482391.43040631.08331.083-0.26258-0.00668-0.2075630.87530.8750.0898150.0038940.12024230.99530.9951.0504160.0327121.01393632.009（1）構(gòu)造或描述概率過(guò)程。（2）實(shí)現(xiàn)抽樣。（3）建立估計(jì)量。實(shí)施蒙特卡羅仿真的主要步驟：蒙特卡羅仿真的優(yōu)點(diǎn)：（1）收斂速度與問(wèn)題的維數(shù)無(wú)關(guān)。（2）受問(wèn)題的條件限制的影響小。（3）程序結(jié)構(gòu)簡(jiǎn)單、清晰，便于編制和調(diào)試。（4）對(duì)于一些物理問(wèn)題，具有其它數(shù)值方法不能替代的作用。隨機(jī)數(shù)的生成和檢驗(yàn)

隨機(jī)數(shù)就是按隨機(jī)方法而生成的數(shù)碼。即0，1，2……9這十個(gè)數(shù)字出現(xiàn)的機(jī)會(huì)等概率，排列的順序隨機(jī)。均勻分布隨機(jī)數(shù)隨機(jī)數(shù)的產(chǎn)生方法0123468932765140675454324678568087586545均勻分布隨機(jī)數(shù)f(x)F(x)xx1010111

均勻分布概率密度及分布函數(shù)

隨機(jī)抽取等概且服從均勻分布的小數(shù)（從0到1，但不包括1）稱(chēng)為均勻隨機(jī)數(shù)

產(chǎn)生均勻分布隨機(jī)數(shù)方法

手工方法2.隨機(jī)數(shù)表

物理方法4.數(shù)學(xué)方法

手工方法擲骰子、抽簽、發(fā)紙牌或從“攪拌均勻的容器”中模編號(hào)球的方法等，但效率很低。隨機(jī)數(shù)表物理方法對(duì)隨機(jī)數(shù)表法的一種改進(jìn)辦法是在計(jì)算機(jī)上裝一臺(tái)物理隨機(jī)數(shù)發(fā)生器，把具有隨機(jī)性質(zhì)的物理過(guò)程（如以放射性物質(zhì)為隨機(jī)源），在計(jì)算機(jī)上直接轉(zhuǎn)換成隨機(jī)數(shù)字。這種方法雖然產(chǎn)生了隨機(jī)數(shù)，但無(wú)可追溯性。由此得到的隨機(jī)數(shù)一去不返，不能重復(fù)出現(xiàn)，因此無(wú)法再用原來(lái)的隨機(jī)數(shù)試算或檢查。又因?yàn)閷?duì)設(shè)備的要求較高，而實(shí)用價(jià)值較低

數(shù)學(xué)方法利用數(shù)學(xué)遞推公式，在計(jì)算機(jī)上產(chǎn)生隨機(jī)數(shù)，是目前使用較廣、發(fā)展較快的一種方法。由于這種方法產(chǎn)生的數(shù)，只能近似地具備隨機(jī)性質(zhì)，因此稱(chēng)為偽隨機(jī)數(shù)（pseudorandomnumber）。產(chǎn)生偽隨機(jī)數(shù)的方法

（一）歷史上曾經(jīng)采用過(guò)的方法1.平方取中法(mid-squaremethod)

平方取中法計(jì)算簡(jiǎn)單，但有明顯得缺點(diǎn)，一是無(wú)法說(shuō)明用什么樣的種子數(shù)可保證有足夠長(zhǎng)的周期，二是容易退化到某一常數(shù)或零，而一旦有一個(gè)零，則以后的數(shù)全為零以21為初始種子數(shù)，用平方取中法求隨機(jī)數(shù)。xixi2xi+1ui21449364986300441193686494096008100640036000900004493649863000.440.930.640.090.080.060.030.000.002.乘法取中法(mid-productmethod）

設(shè)初始種子為12和34，用乘法取中法求隨機(jī)數(shù)。xixi?xi+1xi+1ui1234403644585519472100408136014491584255231901045007600280014000200004036445855194721000.400.360.440.580.550.190.040.070.020.010.000.003.常數(shù)乘子法（constantmultipliermethod）

常數(shù)乘子法的周期較長(zhǎng)，均勻分布特性較好，但仍有平方取中和乘法取中的缺點(diǎn)——容易退化。同時(shí)，成功與否取決于所選常數(shù)。xiki?xixi+1ui341870870.87874785780.78784290290.29291595590.59593245240.24241320320.32321760760.76764180180.1818990990.99995445440.44442420420.42422310310.31311705700.07703850850.85854675670.67673685680.68683740740.7474407070.07341870870.87874785780.78選擇常數(shù)k=55，初始種子數(shù)為34，用常數(shù)乘子法求隨機(jī)數(shù)。

（二）線(xiàn)性同余法1.混合線(xiàn)性同余法

2.乘法線(xiàn)性同余法3.二次同余法

4.加法同余法

用線(xiàn)性同余、混合線(xiàn)性同余、二次同余及加法同余，計(jì)算隨機(jī)數(shù)。其中：a=21，b=59，c=11；模m=128；種子數(shù)是：x0

=34,x-1

=51。參數(shù)m，a，c滿(mǎn)足條件：①m與c互素；②如果q是m的一個(gè)素因子，則q也是a-l的因子；③如果m能被4整除，則a-1也能被4整除。通常取m=2b,c為奇數(shù)，而a-1可被4整除，可得到滿(mǎn)周期。隨機(jī)數(shù)的檢驗(yàn)

由于計(jì)算機(jī)產(chǎn)生的隨機(jī)數(shù)，是按某一迭代公式計(jì)算出來(lái)的，因此實(shí)際上由此而得的數(shù)列是確定性而不是隨機(jī)的。為此，通常要對(duì)隨機(jī)數(shù)進(jìn)行測(cè)試，觀察它與均勻分布U(0,1)的類(lèi)似程度。隨機(jī)數(shù)檢驗(yàn)可分為:參數(shù)檢驗(yàn)均勻性檢驗(yàn)獨(dú)立性檢驗(yàn)

檢驗(yàn)隨機(jī)數(shù)的均值、二階矩及方差等估計(jì)量與均勻分布U(0,1)的理論值有無(wú)顯著差別

參數(shù)檢驗(yàn)用參數(shù)檢驗(yàn)由某隨機(jī)數(shù)發(fā)生器生成的100個(gè)隨機(jī)數(shù)的均勻性。

0.38200.10070.59650.89910.88460.95850.01450.40740.86320.13860.24500.04550.03240.16410.21960.01710.28500.34310.55360.35740.37180.35560.91030.46600.42620.30390.97570.80670.99120.25630.95170.05340.70500.81650.97250.46630.30020.75020.35150.77570.07430.19840.06410.35830.48700.51120.37350.98590.04070.23070.00500.92610.10030.25670.77570.67960.80910.72430.08510.13230.75620.62650.17370.40480.55230.71150.55520.18120.97030.68690.52880.79670.80570.26220.17800.86680.11480.05950.76160.73840.98630.92560.90390.54500.50080.67500.48980.14580.03800.79630.67160.73170.58450.15220.89220.37780.20050.20580.33400.3251

檢驗(yàn)隨機(jī)數(shù)在各子區(qū)間的頻率與理論值N/n有無(wú)顯著差異。常用的檢驗(yàn)方法有:頻率檢驗(yàn)（擬合優(yōu)度檢驗(yàn)）

累積頻率檢驗(yàn)（柯?tīng)柲缏宸颍姑谞栔Z夫檢驗(yàn)）

均勻性檢驗(yàn)頻率檢驗(yàn)是按以下方法，對(duì)隨機(jī)數(shù)獨(dú)立且服從U(0,1)分布的假設(shè)作顯著性檢驗(yàn)。具體步驟是：①將區(qū)間[0,1]分成k個(gè)不相交的等長(zhǎng)子區(qū)間②計(jì)算經(jīng)驗(yàn)頻數(shù)③計(jì)算統(tǒng)計(jì)量④判斷統(tǒng)計(jì)量χ2是否小于臨界值在5％的顯著性水平上，用χ

2檢驗(yàn)所給的100個(gè)隨機(jī)數(shù)的均勻性。將區(qū)間[0,1]分成10個(gè)等間距的子區(qū)間，各子區(qū)間的隨機(jī)數(shù)個(gè)數(shù)ni

如表4.2.5:

經(jīng)驗(yàn)頻數(shù)表子區(qū)間i12345678910

ni121291379512912據(jù)表的經(jīng)驗(yàn)頻數(shù)計(jì)算統(tǒng)計(jì)量χ2

χ2＝6.2<16.92

累積頻率檢驗(yàn)又稱(chēng)柯?tīng)柲缏宸颍姑谞栔Z夫檢驗(yàn)，是連續(xù)分布的擬合性檢驗(yàn)（簡(jiǎn)稱(chēng)K-S檢驗(yàn)）。按以下步驟對(duì)樣本的經(jīng)驗(yàn)分布函數(shù)Fn(x)與總體分布函數(shù)F(x)進(jìn)行顯著性檢驗(yàn)。

①將隨機(jī)數(shù)按大小順序排列成u(1)，u(2)，…，u(n)，并計(jì)算隨機(jī)數(shù)的經(jīng)驗(yàn)分布函數(shù)Fn(x)②計(jì)算

③計(jì)算統(tǒng)計(jì)量④根據(jù)顯著性水平α及樣本大小n查K-S表檢索臨界值⑤比較統(tǒng)計(jì)量D及臨界值

【例4.2.7】在5％的顯著性水平，用K-S法檢驗(yàn)例4.2.5前10個(gè)隨機(jī)數(shù)的均勻性。解：將隨機(jī)數(shù)按大小順序排列，得相關(guān)數(shù)據(jù)如表

Dn+=0.161,D-n=0.263，

Dn=max{0.161,0.263}=0.263

D0.05(10)=0.41

i12345678910u(i)0.0150.1010.1390.3820.4070.5970.8630.8850.8990.959i/n0.10.20.30.40.50.60.70.80.91.0i/n-u(i)0.0850.0990.1610.0180.0930.003-0.163-0.0850.0010.041u(i)-(i-1)/n0.0150.001-0.0610.0820.0070.0970.2630.1850.0990.059

獨(dú)立性檢驗(yàn)就是檢查一個(gè)序列的隨機(jī)數(shù)之間是否存在相關(guān)性。常用檢驗(yàn)的方法有：自相關(guān)檢驗(yàn)

組合規(guī)律檢驗(yàn)

連貫性檢驗(yàn)

獨(dú)立性檢驗(yàn)

利用自相關(guān)系數(shù)值的大小來(lái)判斷隨機(jī)數(shù)的獨(dú)立性，因此也稱(chēng)相關(guān)系數(shù)檢驗(yàn)。具體步驟是：

①計(jì)算j階自相關(guān)系數(shù)的估計(jì)值

②構(gòu)造統(tǒng)計(jì)量vj（j=1，2，…，m）

③根據(jù)顯著性水平α，比較zα/2和統(tǒng)計(jì)量判斷待檢隨機(jī)數(shù)是否具有獨(dú)立性。

>zα/則拒絕相關(guān)系數(shù)ρj=0的假設(shè)

自相關(guān)檢驗(yàn)

通常用撲克檢驗(yàn)（pokertest）檢驗(yàn)樣本中重復(fù)出現(xiàn)數(shù)字的頻率，從而檢查其組合規(guī)律性。

組合規(guī)律檢驗(yàn)

若Ak為隨機(jī)數(shù)中有且僅有k個(gè)數(shù)字相同的事件，Ak出現(xiàn)的頻數(shù)為mk。則三種可能事件發(fā)生的理論頻數(shù)為μk=nPk。構(gòu)造統(tǒng)計(jì)量

由于統(tǒng)計(jì)量χ2漸近服從χ2(2)分布。因此可利用χ2分布檢驗(yàn)隨機(jī)數(shù)列的組合規(guī)律性來(lái)判斷其獨(dú)立性。?=-=3122)(kkkkmmmc有1000個(gè)用三位小數(shù)表示的隨機(jī)數(shù)，經(jīng)過(guò)統(tǒng)計(jì)知有680個(gè)三位數(shù)都不相同，289個(gè)恰有兩位數(shù)字相同，31個(gè)三位數(shù)字都相同。試檢驗(yàn)這些數(shù)據(jù)是否相互獨(dú)立。解：三種不同組合情況的有關(guān)數(shù)據(jù)如下：據(jù)此計(jì)算統(tǒng)計(jì)量,若顯著性水平α＝0.05，查分布表得，因?yàn)椋?7.66>5.991，拒絕符合組合規(guī)律，表明這1000個(gè)隨機(jī)數(shù)不獨(dú)立。類(lèi)型概率觀測(cè)頻數(shù)mk

理論頻數(shù)μk

A10.7268072016002.22A20.272892703611.34A30.01311044144.1連貫性檢驗(yàn)通過(guò)考察數(shù)字的排列來(lái)檢驗(yàn)隨機(jī)數(shù)獨(dú)立性假設(shè)，即檢驗(yàn)隨機(jī)數(shù)列先后出現(xiàn)的隨機(jī)性，重點(diǎn)是檢驗(yàn)其連貫現(xiàn)象是否異常，從而判斷這些數(shù)據(jù)是否獨(dú)立。連貫檢驗(yàn)可分為：1、升降連檢驗(yàn)2、正負(fù)連檢驗(yàn)。升降連檢驗(yàn)由u1，u2，…，un，構(gòu)造數(shù)列{vi}，vi＝ui-ui-1(i=2,3,…,n)，根據(jù)數(shù)列{vi}的符號(hào)，分成升降兩類(lèi)連，vi連續(xù)為正的連稱(chēng)為升連，vi連續(xù)為的連稱(chēng)為降連。用連來(lái)表示隨機(jī)數(shù)增減及其長(zhǎng)度的變化規(guī)律的方法稱(chēng)為升降連檢驗(yàn)。ui0.7440.9060.2640.2310.2570.0070.1590.3630.7140.070Vi0.162-0.642-0.0330.026-0.2500.1520.2040.351-0.644+--+-+++-ui0.3690.2680.8980.13620.0860.4730.7080.2450.4380.132vi0.299-0.1010.630-0.762-0.0500.3870.236-0.4630.193-0.306+-+--++-+-ui0.4420.2270.49390.5810.6240.5980.7450.8060.6180.703vi0.310-0.2150.2670.0870.043-0.0270.1480.061-0.1880.085+-+++-++-+ui0.2740.5580.47430.1630.4350.9830.1970.6360.6400.530vi-0.4290.284-0.084-0.3110.2720.548-0.7860.4390.004-0.110-+--++-++-ui0.4390.4090.94630.9740.7160.0640.7910.7110.5520.496vi-0.091-0.0300.5380.027-0.258-0.6520.727-0.080-0.159-0.056--++--+---{vi}及升降連升連、降連共32個(gè),其中長(zhǎng)度為1的有19個(gè),長(zhǎng)度為2的有9個(gè),長(zhǎng)度為3的有4個(gè)T=32，n=50，統(tǒng)計(jì)量

又zα/2＝z0.025=1.96，可見(jiàn)|Z|<z0.025，所以不能拒絕獨(dú)立性假設(shè)。正負(fù)連檢驗(yàn)正負(fù)連是將數(shù)值在均值以上的標(biāo)以“＋”，表示正連；在均值以下的標(biāo)以“－”，表示負(fù)連。然后，用類(lèi)似升降連檢驗(yàn)的方法進(jìn)行檢驗(yàn)。

假設(shè)n1和n2分別是在均值以上和以下的數(shù)據(jù)個(gè)數(shù)，T為連的總數(shù)，其均值和方差分別是

，

（5.3.11）當(dāng)n1或n2大于20時(shí)，T近似服從正態(tài)分布，統(tǒng)計(jì)量

（5.3.12）服從標(biāo)準(zhǔn)正態(tài)分布N(0,1)。同樣地，給定顯著性水平α，若|Z|≤zα/2，則不拒絕獨(dú)立性的假設(shè)。

顯著性水平α＝5％，用正負(fù)連法檢驗(yàn)上例所給的隨機(jī)數(shù)列的獨(dú)立性。ui0.7440.9060.2640.2310.2570.0070.1580.3630.7140.070vi0.2580.420-0.222-0.255-0.229-0.479-0.327-0.1230.228-0.416符號(hào)++------++ui0.3690.2680.8980.1360.0860.4730.7080.2450.4380.132vi-0.117-0.2180.412-0.350-0.400-0.0130.222-0.241-0.048-0.354符號(hào)--+-