Python實(shí)現(xiàn)蒙特卡洛算法小實(shí)驗(yàn)過程詳解

蒙特卡洛算法思想

蒙特卡洛（Monte Carlo）法是一類隨機(jī)算法的統(tǒng)稱，提出者是大名鼎鼎的數(shù)學(xué)家馮·諾伊曼，他在20世紀(jì)40年代中期用馳名世界的賭城—摩納哥的蒙特卡洛來命名這種方法。

通俗的解釋一下蒙特卡洛算法的思想。假如籃子里有1000個(gè)蘋果，讓你每次閉著眼睛拿1個(gè)，挑出最大的。于是你閉著眼睛隨機(jī)拿了一個(gè)，然后再隨機(jī)拿一個(gè)與第一個(gè)比，留下大的，再隨機(jī)拿一個(gè)，與前次留下的比較，又可以留下大的……你每拿一次，留下的蘋果至少是當(dāng)前最大的，循環(huán)往復(fù)這樣，拿的次數(shù)越多，挑出最大蘋果的可能性也就越大，但除非你把1000個(gè)蘋果都挑一遍，否則你無法肯定最終挑出來的就是最大的一個(gè)。也就是說，蒙特卡洛算法是樣本越多，越能找到最佳的解決辦法，但只是盡量找最好的，不保證一定是最好的。

與它形成對比的是拉斯維加斯算法思想。假如有一把鎖，有1000把鑰匙進(jìn)行選擇，但只有1把是對的。于是你每次隨機(jī)拿1把鑰匙去試，打不開就再換1把。你試的次數(shù)越多，打開最優(yōu)解的機(jī)會就越大，但在打開之前，那些錯(cuò)的鑰匙都是沒有用的。所以拉斯維加斯算法就是盡量找最好的解決辦法，但是不保證能找到。假設(shè)試了999次后沒有任何一把鑰匙能打開鎖，真正的鑰匙是第1000把，但是樣本并沒有第1000次選擇，那么拉斯維加斯算法就不能找到打開鎖的鑰匙。

蒙特卡洛和拉斯維加斯本身是兩座著名賭城，因?yàn)橘€博中體現(xiàn)了許多隨機(jī)算法，所以借此命名。它們只是概括了隨機(jī)算法的特性，算法本身可能復(fù)雜，也可能簡單，在這兩類隨機(jī)算法之間的選擇，往往受到問題的局限。如果問題要求在有限采樣內(nèi)，必須給出一個(gè)解，但不要求是最優(yōu)解，那就要用蒙特卡羅算法。反之，如果問題要求必須給出最優(yōu)解，但對采樣沒有限制，那就要用拉斯維加斯算法。

蒙特卡洛算法實(shí)驗(yàn)

這么看來蒙特卡洛方法的理論支撐其實(shí)是概率論或統(tǒng)計(jì)學(xué)中的大數(shù)定律?；驹砗唵蚊枋鍪窍却罅磕M，然后計(jì)算一個(gè)事件發(fā)生的次數(shù)，再通過這個(gè)發(fā)生次數(shù)除以總模擬次數(shù)，得到想要的結(jié)果。下面我們以三個(gè)經(jīng)典的小實(shí)驗(yàn)來學(xué)習(xí)下蒙特卡洛算法思想。

1.計(jì)算圓周率pi（π）值

實(shí)驗(yàn)原理：在正方形內(nèi)部有一個(gè)相切的圓，圓面積/正方形面積之比是(PixRxR)/(2Rx2R)= Pi/4。在這個(gè)正方形內(nèi)隨機(jī)產(chǎn)生n個(gè)點(diǎn)，假設(shè)點(diǎn)落在圓內(nèi)的概率為P，那么P=圓面積/正方形面積，則P= Pi/4。如何計(jì)算點(diǎn)落在圓內(nèi)的概率P？可以計(jì)算點(diǎn)與中心點(diǎn)的距離，判斷是否落在圓的內(nèi)部，若這些點(diǎn)均勻分布，用M表示落到圓內(nèi)投點(diǎn)數(shù) ， N表示總的投點(diǎn)數(shù)，則圓周率Pi=4P=4xM/N。

實(shí)驗(yàn)步驟：

（1）將圓心設(shè)在原點(diǎn)（0，0），以R為半徑形成圓，則圓面積為PixRxR

（2）將該圓外接正方形, 坐標(biāo)為（-R，-R）（R，-R）（R， R）（-R，R），則該正方形面積為R*R

（3）隨即取點(diǎn)（X，Y），使得-R <=X<=R并且-R <=Y<=R，即點(diǎn)在正方形內(nèi)

（4）通過公式 XxX+YxY<= RxR判斷點(diǎn)是否在圓周內(nèi)（直角三角形邊長公式）。

（5）設(shè)所有點(diǎn)（也就是實(shí)驗(yàn)次數(shù)）的個(gè)數(shù)為N，落在圓內(nèi)的點(diǎn)（滿足步驟4的點(diǎn)）的個(gè)數(shù)為M，則P=M/N，于是Pi=4xM/N。

（6）運(yùn)行結(jié)果為3.143052

def cal_pai_mc(n=1000000):
 r = 1.0
 a, b = (0.0, 0.0)
 x_neg, x_pos = a - r, a + r
 y_neg, y_pos = b - r, b + r
 m = 0
 for i in range(0, n+1):
 x = random.uniform(x_neg, x_pos)
 y = random.uniform(y_neg, y_pos)
 if x**2 + y**2 <= 1.0:
 m += 1
 return (m / float(n)) * 4

2.計(jì)算函數(shù)定積分值

實(shí)驗(yàn)原理：若要求函數(shù)f(x)從a到b的定積分，我們可以用一個(gè)比較容易算得面積的矩型包圍在函數(shù)的積分區(qū)間上（假設(shè)其面積為Area），定積分值其實(shí)就是求曲線下方的面積。隨機(jī)地向這個(gè)矩形框里面投點(diǎn)，統(tǒng)計(jì)落在函數(shù)f(x)下方的點(diǎn)數(shù)量占所有點(diǎn)數(shù)量的比例為P，那么就可以據(jù)此估算出函數(shù)f(x)從a到b的定積分為Area×P。此處我們將a和b設(shè)為0和1，函數(shù)f(x)=x2。

運(yùn)行結(jié)果為0.333749

def cal_integral_mc(n = 1000000):
 x_min, x_max = 0.0, 1.0
 y_min, y_max = 0.0, 1.0
 m = 0
 for i in range(0, n+1):
 x = random.uniform(x_min, x_max)
 y = random.uniform(y_min, y_max)
 # x*x > y 表示該點(diǎn)位于曲線的下面。
 if x*x > y:
 m += 1
 #所求的積分值即為曲線下方的面積與正方形面積的比
 return m / float(n)

3.計(jì)算函數(shù)極值，可避免陷入局部極值

實(shí)驗(yàn)原理：極值是“極大值” 和 “極小值”的統(tǒng)稱。如果一個(gè)函數(shù)在某點(diǎn)的一個(gè)鄰域內(nèi)處處都有確定的值，函數(shù)在該點(diǎn)的值大于或等于在該點(diǎn)附近任何其他點(diǎn)的函數(shù)值，則稱函數(shù)在該點(diǎn)的值為函數(shù)的“極大值”。如果函數(shù)在該點(diǎn)的值小于或等于在該點(diǎn)附近任何其他點(diǎn)的函數(shù)值，則稱函數(shù)在該點(diǎn) 的值為函數(shù)的“極小值”。此處在區(qū)間[-2,2]上隨機(jī)生成一個(gè)數(shù)，求出其對應(yīng)的y，找出其中最大值認(rèn)為是函數(shù)在[-2,2]上的極大值。

運(yùn)行結(jié)果發(fā)現(xiàn)極大值185.1204262706596, 極大值點(diǎn)為1.5144491499169481

def cal_extremum_mc(n = 1000000):
 y_max = 0.0
 x_min, x_max = -2.0, 2.0
 y = lambda x:200*np.sin(x)*np.exp(-0.05*x)#匿名函數(shù)
 for i in range(0, n+1):
 x0 = random.uniform(x_min, x_max)
 if y(x0) > y_max:
 y_max = y(x0)
 x_max = x0
 return y_max, x_max

以上三個(gè)例子也稱為基于蒙特卡洛的投點(diǎn)法，由此得出的值并不是一個(gè)精確值，而是一個(gè)近似值。當(dāng)投點(diǎn)的數(shù)量越來越大時(shí)，這個(gè)近似值也越接近真實(shí)值。

以上就是本文的全部內(nèi)容，希望對大家的學(xué)習(xí)有所幫助，也希望大家多多支持億速云。