蒙特卡羅算法ppt課件

1、Monte Carlo Simulation Methods(蒙特卡羅模擬方法)主要內(nèi)容: 一. M-C方法概述. 二. 隨機(jī)數(shù)的生成. 三. 模擬訓(xùn)練. 四. 實(shí)驗(yàn)標(biāo)題.成信院數(shù)學(xué)與信息科學(xué)系 李勝坤 蒙特卡洛蒙特卡洛(Monte Carlo)方法，或稱(chēng)計(jì)算機(jī)方法，或稱(chēng)計(jì)算機(jī)隨機(jī)模擬方法，是一種基于隨機(jī)模擬方法，是一種基于“隨機(jī)數(shù)的計(jì)隨機(jī)數(shù)的計(jì)算方法。這一方法源于美國(guó)在第二次世界大算方法。這一方法源于美國(guó)在第二次世界大戰(zhàn)中研制原子彈的戰(zhàn)中研制原子彈的“曼哈頓方案。該方案曼哈頓方案。該方案的主持人之一、數(shù)學(xué)家馮的主持人之一、數(shù)學(xué)家馮諾伊曼用著名世諾伊曼用著名世界的賭城界的賭城摩納哥的摩納哥的M

2、onte Carlo來(lái)命名來(lái)命名這種方法，為它蒙上了一層奧秘顏色。這種方法，為它蒙上了一層奧秘顏色。一.M-C方法概述根本思想很早以前就被人們所發(fā)現(xiàn)和利用。17世紀(jì)，人們就知道用事件發(fā)生的“頻率來(lái)決議事件的“概率。19世紀(jì)人們用投針實(shí)驗(yàn)的方法來(lái)決議。高速計(jì)算機(jī)的出現(xiàn)，使得用數(shù)學(xué)方法在計(jì)算機(jī)上大量模擬這樣的實(shí)驗(yàn)成為能夠。從Buffon(蒲豐)投針問(wèn)題談起 220, /2 0, sin,:sin.llaXAX隨機(jī)投針可以理解成針的中心點(diǎn)與最近的平行線(xiàn)的距離X是均勻地分布在區(qū)間 上的r.v.，針與平行線(xiàn)的夾角 是均勻地分布在區(qū)間 上的r.v.，且X與 相互獨(dú)立，于是針與平行線(xiàn)相交的充要條件為 即相交

3、 2sin0022(sin)2lllpP Xdxdaa 于是有：2lap試驗(yàn)者時(shí)間(年)針長(zhǎng)投針次數(shù)相交次數(shù)的估計(jì)值Wolf18500.80500025323.15956Smith18550.60320412183.15665Fox18840.7510304893.15951Lazzarini19250.83340818083.14159292數(shù)值積分問(wèn)題niiixfnxfExxf110)(1)(d)(選取(0, 1)中隨機(jī)數(shù)序列x1， x2， x3， xn。那么誤差約 ，它并不能和一些高級(jí)的數(shù)值積分算法比較， 但對(duì)多維情況，MC方法卻很有吸引力。n1 niiiizyxfnzyxzyxf110

4、1010),(1ddd),(我們可產(chǎn)生一系列隨機(jī)數(shù)可簡(jiǎn)單取3個(gè)隨機(jī)數(shù)構(gòu)成一個(gè)隨機(jī)點(diǎn)，即,.,654321),.,(),(654321相應(yīng)地，niibaxfnxxfab1)(1d)(1 niiidcbayxfnyxyxfabcd1),(1dd),()(1普通地，A)in points randomn over of (average) of measure(d)(fAxfAMonte Carlo數(shù)值積分的優(yōu)點(diǎn)與普通的數(shù)值積分方法比較，Monte Carlo方法具有以下優(yōu)點(diǎn)：1. Monte Carlo一般的數(shù)值方法很難推廣到高維積分的情形，而方法很容易推廣到高維情形2/1/22. ()() dO

5、 nO n一般的數(shù)值積分方法的收斂階為 ，而由中心極限定理可以保證 Monte Carlo 方法的收斂階為 。此收斂階與維數(shù)無(wú)關(guān)，且在高維時(shí)明顯優(yōu)于一般的數(shù)值方法。M-C的根本思緒1.針對(duì)實(shí)踐問(wèn)題建立一個(gè)簡(jiǎn)單且便于實(shí)現(xiàn)的概率統(tǒng)計(jì)模型，使所求的量或解恰好是該模型某個(gè)目的的概率分布或者數(shù)字特征。2.對(duì)模型中的隨機(jī)變量建立抽樣方法，在計(jì)算機(jī)上進(jìn)展模擬測(cè)試，抽取足夠多的隨機(jī)數(shù)，對(duì)有關(guān)事件進(jìn)展統(tǒng)計(jì)3.對(duì)模擬實(shí)驗(yàn)結(jié)果加以分析，給出所求解的估計(jì)及其精度(方差)的估計(jì)4.必要時(shí)，還應(yīng)改良模型以降低估計(jì)方差和減少實(shí)驗(yàn)費(fèi)用，提高模擬計(jì)算的效率回想幾種延續(xù)型分布1.均勻分布U(a,b)其概率密度函數(shù)為其他,0,1)

6、(bxaabxf有),(), c (,badabdcdxcP其中 均勻性特點(diǎn):均勻分布隨機(jī)變量X落在(a,b)內(nèi)恣意子區(qū)間的概率只與子區(qū)間的長(zhǎng)度有關(guān),而與子區(qū)間的位置無(wú)關(guān). 可假設(shè)有這種特性的隨機(jī)變量服從均勻分布.均勻分布隨機(jī)變量X的取值具有均勻性2. 正態(tài)分布正態(tài)分布隨機(jī)變量X的概率密度函數(shù)是Rxxxf,)(21exp21)(2正態(tài)分布由兩個(gè)參數(shù) 和 獨(dú)一確定.其中 是X的均值(數(shù)學(xué)期望): =E(X),它確定了概率曲線(xiàn)的中心位置,而 是X的規(guī)范差: ,它確定了概率曲線(xiàn)的寬窄程度.)(XD 在許多實(shí)踐問(wèn)題中,有一類(lèi)隨機(jī)變量可以表示成為許多相互獨(dú)立的隨機(jī)變量之和,而其中每個(gè)隨機(jī)變量對(duì)總和只起微

7、小的影響,這類(lèi)隨機(jī)變量往往服從或近似服從正態(tài)分布.在實(shí)踐運(yùn)用中,假設(shè)我們分析到一個(gè)隨機(jī)變量遭到較多獨(dú)立的微小要素的疊加影響,就可以用正態(tài)分布來(lái)模擬這個(gè)變量.如:工廠(chǎng)產(chǎn)品的丈量尺寸,農(nóng)作物的收獲量,某地域成年人的身高,體重等可看成服從正態(tài)分布的隨機(jī)變量.3. 指數(shù)分布指數(shù)分布隨機(jī)變量X的概率密度為0, 00,)(xxexfx指數(shù)分布常用來(lái)描畫(huà)壽命問(wèn)題.二.隨機(jī)數(shù)的生成1.蒙特卡羅模擬的關(guān)鍵是生成優(yōu)良的隨機(jī)數(shù)。2.在計(jì)算機(jī)實(shí)現(xiàn)中,我們是經(jīng)過(guò)確定性的算法生成 隨機(jī)數(shù),所以這樣生成的序列在本質(zhì)上不是隨機(jī) 的,只是很好的模擬了隨機(jī)數(shù)的性質(zhì)(如可以經(jīng)過(guò) 統(tǒng)計(jì)檢驗(yàn))。我們通常稱(chēng)之為偽隨機(jī)數(shù)(pseudo-r

8、andom numbers)。3.在模擬中，我們需求產(chǎn)生各種概率分布的隨機(jī)數(shù)，而大多數(shù)概率分布的隨機(jī)數(shù)產(chǎn)生均基于均勻分布U(0,1)的隨機(jī)數(shù)。U(0,1)隨機(jī)數(shù)的生成乘同余法：1101 mod , , /iiiiixaxmuaxxmxm其中 均為整數(shù), 可以任意選取。 稱(chēng)為種子,a 是乘因子,m是模數(shù)0 x一個(gè)簡(jiǎn)單的例子1i+116 mod11, u/11 6,11iiixxxam（）0 1 ,x 1,6,3,7,9,10,5,8,4當(dāng)時(shí) 得到序列:,1,6,3.,2.003,1 ,1,3,9.3,2,2,1,3,9,5,42,6,7,10,8,6.axax如果令 得到序列:如果令 得到序列：

9、一個(gè)簡(jiǎn)單的例子(續(xù))上面的例子中，第一個(gè)隨機(jī)數(shù)生成器的周期長(zhǎng)度是 10，而后兩個(gè)生成器的周期長(zhǎng)度只需它的一半。我們自然希望生成器的周期越長(zhǎng)越好，這樣我們得到的分布就更接近于真實(shí)的均勻分布。0 (max在給定 的情況下，生成器的周期與和初值種子)選擇有關(guān)。線(xiàn)性同余生成器(混合同余法)(Linear Congruential Generator )111 () mod /iiiixaxcmuxm一般形式：02. 1 mmaxm線(xiàn)性同余器可以達(dá)到的最長(zhǎng)周期為 ，我們可以通過(guò)適當(dāng)?shù)倪x擇 和，使無(wú)論選取怎樣的初值 都可以達(dá)到最大周期(一般選取 為質(zhì)數(shù)) c是非負(fù)整數(shù).經(jīng)過(guò)適中選取參數(shù)c可以改善 隨機(jī)數(shù)的

10、統(tǒng)計(jì)性質(zhì)(獨(dú)立性,均勻性).常用的線(xiàn)性同余生成器Modulus mMultiplier aReference231-1=214748364716807Lewis, Goodman, and Miller39373LEcuyer742938285Fishman and Moore950706376Fishman and Moore1226874159Fishman and Moore214748339940692LEcuyer214748356340014LEcuyer復(fù)雜一些的生成器Multiple recursive generator1122102(.) mo,.)d/iiiki kikk

11、ixa xa xa xmuxxxxm需(要選取種子12-1-1 (,.) 1ijiii kkkxmxxxmm每個(gè)有種選擇，所以向量 可以取個(gè)不同的值，所以這樣的隨機(jī)數(shù)生成器的最大周期可以達(dá)到 ，大大提高了簡(jiǎn)單同余生成器的周期。算法實(shí)現(xiàn)許多程序文語(yǔ)中都自帶生成隨機(jī)數(shù)的方法，如 c 中的 random() 函數(shù)，Matlab中的rand()函數(shù)等。但這些生成器生成的隨機(jī)數(shù)效果很不一樣，比如 c 中的函數(shù)生成的隨機(jī)數(shù)性質(zhì)就比較差，假設(shè)用 c ，最好本人再編一個(gè)程序。Matlab 中的 rand() 函數(shù)，經(jīng)過(guò)了很多優(yōu)化。可以產(chǎn)生性質(zhì)很好的隨機(jī)數(shù)，可以直接利用。從U(0,1)到其它概率分布的隨機(jī)數(shù)1.

12、離散型隨機(jī)數(shù)的模擬2.延續(xù)型隨機(jī)數(shù)的模擬3.正態(tài)隨機(jī)數(shù)的模擬1.離散型隨機(jī)數(shù)的模擬設(shè)隨機(jī)變量 的分布律為令將 作為區(qū)間(0,1)的分點(diǎn).假設(shè)隨機(jī)變量 ,有X), 2 , 1(ipxXPii),2, 1(,)(, 0)0(1npnPPnii) 1 , 0( UU),2, 1( ,)1()()()1(npnPnPnPUnPPn)(nP令那么有據(jù)此,可得產(chǎn)生 的隨機(jī)數(shù)的詳細(xì)過(guò)程為:每產(chǎn)生一個(gè)(0,1)區(qū)間上均勻分布隨機(jī)數(shù) ,假設(shè) 那么令 取值 .nx)() 1(nxXnPUnPXnnpxXPX)() 1(nPUnPU例1:離散型隨機(jī)變量X有如下分布律: X 0 1 2P(x) 0.3 0.3 0.4

13、設(shè) 是(0,1)上均勻分布的隨機(jī)數(shù),令那么 是具有X分布律的隨機(jī)數(shù).iiiiUUUx6 . 0, 26 . 03 . 0, 13 . 00, 0NUUU,21Nxxx,212.延續(xù)型隨機(jī)數(shù)的模擬a.逆變換方法(常用) (Inverse Transform Method)b.舍取方法 (Acceptance-Rejection Method)定理: 設(shè)隨機(jī)變量Y的分布函數(shù)F(y)是延續(xù)函數(shù), 而U是在(0,1)上均勻分布的隨機(jī)變量, 令, 那么Y與X有一樣的分布.)(1UFX1 1 ( ) (0,1) ( ) F xUuFu由 定 理， 要 產(chǎn) 生 來(lái) 自的 隨 機(jī) 數(shù) ， 只 要 先產(chǎn) 生 來(lái)

14、 自隨 機(jī) 數(shù)， 然 后 計(jì) 算 即可 。 具 體 步 驟 如 下 ：(1) (0,1)U上生 成均 勻 分 布 的 隨 機(jī) 數(shù)。-1(2) , ( ) ( ) XXUFFx計(jì)算則為來(lái)自 分布的隨機(jī)數(shù).11() ()()( )( )FUP XxP FUxP UF xF x 由 的 定 義 和 均 勻 分 布 的 分 布 函明數(shù) 可證：得 ：-11 ( , ) 0 ( ) 1, ( )() 01 (0,1) ( - ) ( , ) XU a bxaxaF xaxbbaxbFyaba yyUUab a UU a b例 ： 設(shè)，則其分布函數(shù)為 ，生成隨機(jī)數(shù)，則是來(lái)自的隨機(jī)數(shù)。/12 exp( ) (

15、 )1, 0 ( )log(1) log(1) ( ) 1 xXXF xexFyyXUUUU 例 ：設(shè)服從指數(shù)分布，則的分布函數(shù)為：通過(guò)計(jì)算得，則：服從指數(shù)分布 其中服從均勻分布又因?yàn)?和有著同樣的分布，所以也可以取： log( )XU 舍取方法舍取方法(Acceptance-Rejection )方法最早由 Von Neumann提出，如今曾經(jīng)廣泛運(yùn)用于各種隨機(jī)數(shù)的生成。根本思緒：經(jīng)過(guò)一個(gè)容易生成的概率分布 g 和一個(gè)取舍準(zhǔn)那么生成另一個(gè)與 g 相近的概率分布 f 。詳細(xì)步驟： ( ) ( ) ( )( ) 1, f xg xf xcg xcx 假設(shè)和均為集合 上的概率密度函數(shù)。且滿(mǎn)足： 。

16、易于從g(x)中生成樣本X.用如下步驟生成有Y:1. ( ) ;2. (0,1),U3. ()/() g xXUXUf Xcg X生成 的樣本生成U且 與 獨(dú)立；如果，則取Y=X,返回步驟1， 否則舍去X,返回步驟1.下面我們驗(yàn)證由上述步驟生成的隨機(jī)數(shù) Y 確實(shí)具有密度函數(shù) f(x) ()(|()/()( , ()/() ()/()( )( )1/( )( )()( )( )( )AAAP YAP XA Uf Xcg XP XAUf Xcg XP Uf Xcg Xf xg x dxccg xf xP YAcg x dxf x dxcg xx對(duì)于任何的 ，我們有：）分母即Y的概率密度函數(shù)為f(

17、).()/()1/ ( ) ( )(0,1) P Uf Xcg Xcf xcg xU由上面的證明過(guò)程我們可以看出每次接受的概率為 ）。也就是說(shuō)為了生成一個(gè) 的樣本，我們需要平均生成個(gè)和 的樣本。所以為了提高舍取法的效率，我們應(yīng)該使 c 的取值盡能夠的小，也就是使 f 和 g 的分布更為相近。3.規(guī)范正態(tài)分布隨機(jī)數(shù)的生成正態(tài)分布是概率統(tǒng)計(jì)中最重要的分布，在此我們著重討論如何生成規(guī)范正態(tài)分布隨機(jī)數(shù)。引理：121/21121/221212, (0,1) ( 2ln)cos(2) ( 2ln)sin(2) U UUZUUZUUZZ 設(shè) 是獨(dú)立同分布的變量，令：則證明：與 獨(dú)由隨機(jī)立，且均服從向量的變換

18、公標(biāo)準(zhǔn)式易正態(tài)分布。得，略。Box-Muller 算法121. ,U U生成兩個(gè)獨(dú)立的均勻分布隨機(jī)數(shù) 1/21121/221212 ( 2ln)cos(2) ( 2ln)sin(2) ZUUZUUZZ 2.令 則，為獨(dú)立的正態(tài)隨機(jī)數(shù)利用中心極限定理設(shè) 是n個(gè)相互獨(dú)立的在(0,1)上均勻分布的隨機(jī)變量,由中心極限定理知 漸近服從正態(tài)分布N(0,1).普通取n=10即可,假設(shè)取n=12,那么上式簡(jiǎn)化為 再由公式 即可得到正態(tài)分布 的隨機(jī)數(shù).nUUU,2112/ )2(1nnUXniiXY),(2N1216iiUXMatlab程序Function r=rnd-u(a,b)%產(chǎn)生在a,b間均勻分布的隨

19、機(jī)數(shù)r=a+(b-a)*rand;returnMatlab程序Function r=rnd-beta(lamda)%模擬指數(shù)分布%lamda表示指數(shù)分布的參數(shù)r=-log(rand)/lamda;returnMatlab程序Function y= rnd(mean, segema)%模擬均值為mean,方差為segma的正態(tài)分布r=rand(1,12);x=sum(r)-6;y=segma*x+mean;return三. 模擬訓(xùn)練例1. 模擬求近似圓周率在邊長(zhǎng)為1的正方形內(nèi)有一半徑為0.5的內(nèi)切圓.如今模擬產(chǎn)生在正方形內(nèi)均勻分布的點(diǎn)n個(gè).如果有m個(gè)在圓內(nèi),那么圓面積與正方形的面積比可近似為m/

20、n.即/4m/n 4m/nn=10000m=0;For i=1:n if rand2+rand2 =0,那么可以經(jīng)過(guò)模擬估算.構(gòu)造一個(gè)矩形包含曲邊梯形,dmax f(x).產(chǎn)生n(足夠大)個(gè)在矩形區(qū)域內(nèi)的點(diǎn),假設(shè)落在由函數(shù)f(x)構(gòu)成的曲邊梯形內(nèi)的點(diǎn)為m個(gè),那么所求定積分為 .dabnmxfba)()(n=106;a=0;b=1;d=1;m=0;for i=1:n x=a+rand*(b-a); y=d*rand; if y=x2 m=m+1; endends=m/n*(b-a)*dn=106;x=rand(1,n);y=x.2;s=sum(y)/n采用前面講的方法:niiixfnxfExxf

21、110)(1)(d)(例3. 渡口模型問(wèn)題描畫(huà):一個(gè)渡口的渡船營(yíng)運(yùn)者擁有一只甲板長(zhǎng)32米,可以并排停放兩列車(chē)輛的渡船.他在思索怎樣在甲板上安排過(guò)河車(chē)輛的位置,才干平安地運(yùn)過(guò)最多數(shù)量的車(chē)輛.分析:怎樣安排過(guò)河車(chē)輛,關(guān)懷一次可以運(yùn)多少輛各類(lèi)車(chē).預(yù)備任務(wù):察看數(shù)日,發(fā)現(xiàn)每次情況不盡一樣,得到以下數(shù)據(jù)和情況:(1)車(chē)輛隨機(jī)到達(dá),構(gòu)成一個(gè)等待上船的車(chē)列;(2)來(lái)到車(chē)輛,轎車(chē)約占40%,卡車(chē)約占55%,摩托車(chē)約占5%;(3)轎車(chē)車(chē)身長(zhǎng)3.55.5米,卡車(chē)車(chē)身長(zhǎng)為810米.問(wèn)題分析:這是一個(gè)機(jī)理較復(fù)雜的隨機(jī)問(wèn)題,是遵先到先效力的隨機(jī)排隊(duì)問(wèn)題.處理方法:采用模擬模型方法.因此需思索以下問(wèn)題: (1)應(yīng)該怎樣安

22、排摩托車(chē)? (2)下一輛到達(dá)的車(chē)是什么類(lèi)型? (3)怎樣描畫(huà)一輛車(chē)的車(chē)身長(zhǎng)度? (4)如何安排到達(dá)車(chē)輛參與甲板上兩列車(chē)隊(duì)中的哪一列中去?本實(shí)驗(yàn)主要模擬裝載車(chē)輛的情況,暫時(shí)不思索渡船的平安.模型建立:設(shè)到達(dá)的卡車(chē),轎車(chē)長(zhǎng)度分別為隨機(jī)變量 , .結(jié)合實(shí)踐,這里無(wú)妨設(shè)卡車(chē),轎車(chē)的車(chē)身長(zhǎng)度 , 均服從正態(tài)分布.由于卡車(chē)車(chē)身長(zhǎng)810米,所以卡車(chē)車(chē)長(zhǎng) 的均值為(8+10)/2=9米,由概率知識(shí)中的3原那么,其規(guī)范差為(9-8)/3=1/3,所以得到 .同理可得 .1L2L1L2L1L)91, 9(1NL)91, 5 . 4(2NL模擬程序設(shè)計(jì):由以上的分析,程序設(shè)計(jì)時(shí)應(yīng)劃分的主要模塊(函數(shù))如下:(1)確

23、定下一輛到達(dá)車(chē)輛的類(lèi)型:(2)根據(jù)車(chē)的類(lèi)型確定到達(dá)車(chē)輛的長(zhǎng)度;(3)根據(jù)一定的停放規(guī)那么,確定放在哪一列.function r =makeid%模擬下一輛到達(dá)車(chē)的類(lèi)型t=rand;if t=0.55 r=1;%到達(dá)卡車(chē)elseif tL(2)v if L(2)+newlen32v pos=2;v elsev full=1;v endvelsev if L(1)+newlen32v pos=1;v elsev full=1;v endvendvfunction sim_dukou v%渡口模型的模擬vn=input(輸入模擬次數(shù):);vif isempty(n)|(n500)v n=500;ve

24、ndvN=zeros(1,3);vfor i=1:nv isfull=0;v L=0,0;v while isfullv id=makeid;v N(id)=N(id)+1;v newlen=getlength(id);v isfull,pos=getiffull(L,newlen);v if isfullv L(pos)=L(pos)+newlen;v endv endvendvdisp(平均每次渡船上的車(chē)數(shù))vmean_n=N/n四.實(shí)驗(yàn)標(biāo)題趕上火車(chē)的概率問(wèn)題描畫(huà):如圖,一列火車(chē)從A站開(kāi)往B站,某人每天趕往B站上這趟火車(chē).他已了解到:(1)火車(chē)從A站到B站的運(yùn)轉(zhuǎn)時(shí)間是均值為30分鐘,規(guī)范差為2分鐘的隨機(jī)變量;(2)火車(chē)在下午大約1點(diǎn)分開(kāi)A站