仿真法，自助法

• 仿真法/統(tǒng)計模擬方法(simulation based methods)
• 自助法(Bootstrap method)
• 三. 結(jié)語

在最后一章，小弟為大家介紹兩個統(tǒng)計方法，說老實話，這兩個方法可以說的是上古秘籍。統(tǒng)計的起源最早可追溯到1910年，那時候計算機(jī)還沒出生呢，當(dāng)中國還是清朝偽滿洲的時候，有這么一群數(shù)學(xué)家利用紙和筆在條件很艱苦的情況下，想出了一些統(tǒng)計方法，這些方法就是本章要介紹的。我們不能忘記過去數(shù)學(xué)家們的智慧和對如今大數(shù)據(jù)時代的貢獻(xiàn)。

仿真法/統(tǒng)計模擬方法(simulation based methods)

在18世紀(jì)，也就是中國乾隆時期，西方有這么個問題—蒲豐投針問題:
假設(shè)我們有一個以平行且等距木紋(黑線)鋪成的地板（如圖），

隨意拋一支長度比木紋之間距離小的針(紅線段)，求針和其中一條木紋相交的概率。
經(jīng)過蒲豐計算，如果條紋間距為a，針的長度為$l$,那么針與其中平行線中任意一條相交的概率為$\frac{2l}{\pi a}$.這里推導(dǎo)也是比較簡單的(相關(guān)推導(dǎo)網(wǎng)上均有)，但小弟在這還是寫幾筆:

設(shè)X表示針的中點到平行線的距離，X服從均勻分布(0,$\frac{a}{2}$),也就是(0<X<$\frac{a}{2}$)的概率是均勻的(數(shù)據(jù)分析模型第二章，連續(xù)均勻分布)
$$f(X)=\{\begin{array}{ll}\frac{1}{\frac{a}{2}?0}=\frac{2}{a},& 0<X<\frac{a}{2}\\ 0,& 其他\end{array}$$

設(shè)Y表示針與平行線的夾角。Y服從均勻分布(0,$\frac{\pi }{2}$)
$$f(Y)=\{\begin{array}{ll}\frac{1}{\frac{\pi }{2}?0}=\frac{2}{\pi },& 0<Y<\frac{\pi }{2}\\ 0,& 其他\end{array}$$

X和Y相互獨立,那么X和Y的聯(lián)合概率為:
$$f(X,Y)=\{\begin{array}{ll}\frac{2}{a}?\frac{2}{\pi }=\frac{4}{\pi a},& 0<X<\frac{a}{2},0<Y<\frac{\pi }{2}\\ 0,& 其他\end{array}$$

當(dāng)$X<\frac{l}{2}sinY$時，針與線相交。那么:
$$P\{X<\frac{l}{2}sinY\}={\int }_0^{\frac{\pi }{2}}{\int }_0^{\frac{1}{2}sin(y)}f(X,Y)dxdy={\int }_0^{\frac{\pi }{2}}{\int }_0^{\frac{1}{2}sin(y)}\frac{4}{\pi a}dxdy=\frac{2l}{\pi a}$$
（根據(jù)數(shù)據(jù)分析模型第二章知識點和內(nèi)容）
因為$X<\frac{l}{2}sinY$,所以${\int }_0^{\frac{1}{2}sin(y)}$,因為$0<Y<\frac{\pi }{2}$所以${\int }_0^{\frac{\pi }{2}}$. 先算${\int }_0^{\frac{1}{2}sin(y)}(?)dx$，后算${\int }_0^{\frac{\pi }{2}}(?)dy$。

有位叫蒙特卡洛(Monte Carlo)同學(xué)給了另外一個解決方法，那就是實驗，硬干它。經(jīng)過大量實驗，記錄多少次針與條紋相交的次數(shù)，從而得到那么針與其中平行線中任意一條相交的概率。我們也知道根據(jù)弱大數(shù)理論，如果實驗多次，那么它估計的概率其實也是很準(zhǔn)確的。

但當(dāng)時人家蒙特卡洛同學(xué)本意不是為了估計針與條紋相交的概率，人家統(tǒng)計這個概率是為了估計$\pi$值。利用實驗得到的概率具體值和蒲豐利用微積分計算的概率$\frac{2l}{\pi a}$,來反向估計出$\pi$值。所以蒲豐投針問題也是用來計算$\pi$值的。

以現(xiàn)在看，當(dāng)時人家蒙特卡洛同學(xué)通過實驗的方法來估計某一事件的概率，大家覺得這不就是小兒科么。但這就是早期統(tǒng)計的起源，也是最早的統(tǒng)計模擬方法。

有的時候，對于某件事，我們很難純用數(shù)學(xué)來計算它的概率，舉幾個例子:

1.如果$X_1,X_2～N(0,1)$，那么$P(\sqrt{|X_1|}+\sqrt{|X_1|}>2)$是多少？

2.如果$X_1～Poi({\lambda }_1),X_2～Poi({\lambda }_2)$,那么$P(X_1?X_2^{\frac{3}{2}}>0)$是多少？

這些問題，要是靠純數(shù)學(xué)計算推導(dǎo)是很難的，這個時候我們需要模擬這個事件，也就是做實驗來計算它的概率。

實驗概率/經(jīng)驗概率收斂(Convergence of Empirical Probabilities)
我們可以利用統(tǒng)計模擬方法來解決這些問題:

第一步:
反復(fù)m次，每一次從它的分布中產(chǎn)生$X_1,X_2..,X_p$個數(shù)據(jù)。

第二步:估計概率，利用實驗
$$P(X_1,X_2,..,X_p\in Z)\approx \frac{1}{m}\sum _{i=1}^{m}I(X_1,X_2,..,X_p\in Z)$$
$I(?)$為指標(biāo)方程，即滿足條件則為1，不滿足條件則為0

舉個例子:
$X～N(0,1)$,計算$P\{\sqrt{|X|}>1\}$的概率

我們先利用純數(shù)學(xué)來計算下:
設(shè):$Q=\sqrt{|X|}$,那么根據(jù)正態(tài)分布$N(0,1)$,得:
$$P(Q=q)=(\frac{2^3}{\pi }{)}^{\frac{1}{2}}qexp(?\frac{q^4}{2})$$
那么:
$$P(Q>1)={\int }_1^{\infty }P(q)dq=1?\frac{2}{\sqrt{\pi }}{\int }_0^{x}exp(?t^2)dt=0.3171$$

我們在通過實驗概率/經(jīng)驗概率收斂的方法來算下:

#MATLAB或者R X=rnorm(m,0,1) mean(sqrt(abs(X))>1) >m=100, 0.3300 >m=1000, 0.3310 >m=100000, 0.3176 >m=10^9, 0.317102

那么根據(jù)抽取$10^9$次數(shù)據(jù)(每次實驗抽一個數(shù)據(jù))，我們利用實驗估計$P\{\sqrt{|X|}>1\}$的概率大約為0.317102.

我們可以發(fā)現(xiàn)其實通過實驗來計算概率和純數(shù)學(xué)計算的結(jié)果幾乎一樣，當(dāng)我們的實驗次數(shù)足夠多的情況下。

當(dāng)然了這個方法也有缺點:
那就是實驗多少次才算足夠，當(dāng)我們計算出的概率結(jié)果比$\frac{1}{m}$還要小，說明我們的計算已經(jīng)被$\frac{1}{m}$所影響了，那樣計算出的概率結(jié)果是不準(zhǔn)確的。因為如果當(dāng)m足夠大時，我們的$\frac{1}{m}$會趨近于0.

但我們有解決這個缺點的辦法:
我們會先叫實驗的次數(shù)和每次實驗抽取多少數(shù)據(jù)均為m，如果最后計算的概率結(jié)果小于$\frac{1}{m}$,那就說明估計不準(zhǔn)確，對于較小值的概率結(jié)果，我們則需要較大的實驗次數(shù)，試驗次數(shù)至少是$\frac{1}{概率結(jié)果}$的10倍。

偽隨機(jī)數(shù)(Pseudo Random Numbers)
在做仿真法/統(tǒng)計模擬法時，我們會需要一些隨機(jī)數(shù)，例如:根據(jù)分布隨機(jī)生成數(shù)據(jù)，隨機(jī)生成數(shù)字。
雖然計算機(jī)可以隨機(jī)生成數(shù)字，但如果計算機(jī)是根據(jù)一些特定規(guī)則進(jìn)而隨機(jī)生成數(shù)字，這類隨機(jī)數(shù)被稱為偽隨機(jī)數(shù)。偽隨機(jī)數(shù)可以用計算機(jī)大量生成，是模擬研究中重要的一步，可以提高效率。

偽隨機(jī)數(shù)最初是利用均勻分布的思想，但如今小弟就不知道了，因為創(chuàng)造偽隨機(jī)數(shù)的方法多了，例如馬特塞特旋轉(zhuǎn)演算，線性同余法，乘同余法，后面的兩種方法均是離散數(shù)學(xué)范疇但還是用到了均勻分布還是比較簡單的。據(jù)小弟所知，大部分編程語言創(chuàng)造偽隨機(jī)數(shù)還是利用均勻分布的思想。

小弟在此就分享下利用均勻分布來創(chuàng)建偽隨機(jī)數(shù)，若有其他創(chuàng)造偽隨機(jī)數(shù)的方法，那基本上是大同小異。

均勻分布隨機(jī)變量產(chǎn)生其他分布的隨機(jī)變量
首先我們要有個均勻分布$U(0,1)$,選取變量∈(0,1)。那么我們可以利用均勻分布變量∈(0,1)來得到其他分布的變量。

第一步:給定其他分布或者概率密度函數(shù)$Z=F(Y)$

第二步:計算它的反函數(shù)$F^{?1}(Z)=Y$,那么我們都知道這個$Z$∈(0,1)。有的沒有反函數(shù)就分段。

第三步：生成隨機(jī)變量X從均勻分布$U(0,1)$

第四步: 帶入反函數(shù)$F^{?1}(X)=Y$那么我們便利用均勻分布的變量$X$，得到了其他分布的變量$Y$了。

我們既然可以利用均勻分布隨機(jī)變量得到其他分布的隨機(jī)變量。那么，如何叫它隨機(jī)呢？

隨機(jī)化:
隨機(jī)化主要是在于隨機(jī)選取均勻分布的變量上。
$$f(X_0,X_1,...,X_{n?2},X_{n?1})=X_n$$
舉個例子:
我們先初始給一個$X_0$的值，該值被稱為種子(seed)。$f(?)$被稱為轉(zhuǎn)換方程。這里轉(zhuǎn)換方程就很多了，各有各的優(yōu)點，例如有的方法可以使隨機(jī)數(shù)更均勻，有的影響隨機(jī)數(shù)列的周期長短，有的需要參數(shù)間互為質(zhì)數(shù)等等。
當(dāng)我們有了$X_0$,我們利用$f(?)$得到$X_1$,再利用$X_0,X_1$帶入$f(?)$中得到$X_2$。但無論怎樣我們的$X_0,X_1,...,X_{n?2},X_{n?1}\in (0,1)$.

這樣我們的實現(xiàn)了在均勻分布上隨機(jī)選取變量，從而隨機(jī)產(chǎn)生其他分布上的變量。

偽隨機(jī)數(shù)缺點:
好的種子比較關(guān)鍵。另外如果使用相同的種子，那么它們產(chǎn)生的偽隨機(jī)數(shù)是相同的。

自助法(Bootstrap method)

自助法也經(jīng)常被用在統(tǒng)計模擬方法中，自助法和統(tǒng)計模擬法均為大家族多次采樣方法中的一員。我們在前兩章其實已經(jīng)了解了它們家族的另外一員，那就是CV(交叉驗證)。

對于交叉驗證，我們估計對于未來數(shù)據(jù)(測驗數(shù)據(jù))預(yù)測的差值。

對于自助法, 它是用來看看我們的估計是否穩(wěn)定(例如標(biāo)準(zhǔn)差，置信區(qū)間等等)

在第三章中，我們經(jīng)過多次取樣，我們會得到樣本分布，即如下圖:

Population, 則是我們的總體數(shù)據(jù)，多次取樣(samples)，取相同數(shù)量的數(shù)據(jù)，再利用每次取樣得來的樣本數(shù)據(jù)估計總體參數(shù)。

舉個例子:
我們目前有一組很感興趣的樣本，那么我們可以利用上圖的自助法分布來估計這組由感興趣樣本所估計出的參數(shù)$\bar{\theta }$偏差和方差

假如我們已經(jīng)反復(fù)取樣了M次，分別估計它們的參數(shù)，偏有了M個${\bar{\theta }}^{(i)}$，那么我們可以估計出該感興趣樣本參數(shù)的$\bar{\theta }$的偏差和方差

估計$\bar{\theta }$的偏差:
$$bias=\frac{1}{M}\sum _{i=1}^M(\bar{\theta }?{\bar{\theta }}^{(i)})$$

估計$\bar{\theta }$的方差:
$$Var=(\frac{1}{M?1})\sum _{i=1}^M({\bar{\theta }}^{(i)}?\bar{\theta }{)}^2$$

除此之外，經(jīng)過多次取樣后，我們自然也會得到樣本的分布。就像我們之前第三章和第四章內(nèi)容所講，經(jīng)過利用樣本的估計參數(shù)后(就是一些數(shù)學(xué)轉(zhuǎn)換)，我們自然也會得到估計參數(shù)的分布。

舉個例子：
我們抽一次樣本得$y=(y_1,y_2,..,y_n)$,該樣本有n個數(shù)據(jù)，并作為我們感興趣的樣本。

我們利用該感興趣的樣本來估計參數(shù)$\bar{\theta }$

利用自助法邏輯：
1.抽取M次數(shù)樣本，每個樣本均也需要有n個數(shù)據(jù)
2.每個樣本分別估計參數(shù)記作${\bar{\theta }}^{(i)},i=1,2,3..,M$

那么我們便可以做很多事情了，估計我們興趣樣本的偏差，方差。當(dāng)然了，因為我們有了${\bar{\theta }}^{(i)},i=1,2,3..,M$的分布，我自然也可以計算這個分布的置信區(qū)間。

大家不要被小弟上述例子從而局限了想象力，自助法的本質(zhì)就是多次取樣。還記得第八章中的一幅圖么，即如下:

我們不看綠線，就看圖中橘色虛線和線段，即平均擬合模型和95%平均擬合模型。這里95%其實就是95%${\bar{\theta }}^{(i)},i=1,2,3..,M$的分布置信區(qū)間，然后把估計的參數(shù)依次帶入到模型里去，然后在把你想預(yù)測的數(shù)據(jù)代入進(jìn)去，這樣你的預(yù)測就會有個范圍即95%。至于橘色線段即平均擬合，則是算的預(yù)測值的均值罷了。

在第6章，我們也接觸過簡單線性回歸對于預(yù)測的置信區(qū)間，還記得么$(Y?{\bar{\beta }}_0?{\bar{\beta }}_1{x}_0)/[\sqrt{\frac{n+1}{n}+\frac{(x_0?\bar{x}{)}^2}{({\sum }_{i=1}^n{x}_i^2?n{\bar{x}}^2)}}\sqrt{\frac{RSS}{n?2}}]～t_{n?2}$，$x_0$則是測試數(shù)據(jù)。對沒錯，有了分布我們便可以得到關(guān)于簡單線性回歸的預(yù)測95%置信區(qū)間。你可以用純數(shù)學(xué)的方式來找到關(guān)于某個模型的分布，但現(xiàn)實很殘酷，有些模型過于復(fù)雜，利用純數(shù)學(xué)的方式我們是很難做的，旦旦一個簡單線性回歸的分布都如此復(fù)雜，就不要說比簡單新型回歸更復(fù)雜的模型了。還是以上圖為例，該圖其實是多項式分布模型，綠色線為真多項式分布，并且有5項，要想用數(shù)學(xué)找到該模型的分布是很困難的。所以才會有蒙特塔羅的思想，就是多做實驗，多抽樣，用簡單粗暴的方式。

總之一句話，多做實驗，多動手，多取樣，再看分布，有時候你覺得能用數(shù)學(xué)推，到最后你會發(fā)現(xiàn)現(xiàn)實很殘酷。

自助法的優(yōu)點:
簡單粗暴，很容易做，也很容易理解。
估計準(zhǔn)確度會隨著樣本的增加而增加
自助法經(jīng)常用在計算偏差，方差，預(yù)測差，置信區(qū)間。

自助法缺點:
算起來太慢了。
對于某些問題需要轉(zhuǎn)換，例如拉索

排列檢驗(Permutation tests)
在最后，小弟介紹一個另一個統(tǒng)計方法，即排列檢驗，經(jīng)常用在樣本分布上面。其實大家早就知道了，它就是p值。

簡單回顧下
$$H_0:\beta =0,H_1=0$$
計算p值，看是否拒絕原假設(shè)。

那么如何利用p值用在多次抽樣上面呢？
我們用例子來解釋:

我們要用排列檢驗計算p值，來看是否$\bar{\beta }=1.4310$

現(xiàn)在我們從這些數(shù)據(jù)(總體)里取樣。取個m次，分別估計參數(shù)${\bar{\beta }}^{(i)}$.
此時我們計算p值:
$$p\approx \frac{1}{m}\sum _{i=1}^{m}I(|{\bar{\beta }}^{(i)}|>|\bar{\beta }|)$$
$I(?)$符合條件為1，不符合條件為0

我們來張圖：

經(jīng)過m=10000次的抽樣，我們估計10000個${\bar{\beta }}^{(i)}$即我們的x軸.計算頻率，則有了y軸。我們利用所有數(shù)據(jù)估計的$\bar{\beta }$=1.4310.可以看到$P(|{\bar{\beta }}^{(i)}|>|\bar{\beta }|)$即上圖±1.4310的兩個外端概率是很小的約為0.0013(p值太小)拒絕原假設(shè)。我們其實能夠清楚看到，也側(cè)面說明了，$\bar{\beta }=1.4310$。

排列檢驗與我們之間的假設(shè)檢驗區(qū)別在于，排列檢驗不需要基于原假設(shè)，但假設(shè)檢驗需要基于原假設(shè)去計算。

在排列檢驗里，其實它更注重的是多次取樣后的數(shù)據(jù)分布。然后根據(jù)條件，計算最小概率即p值來進(jìn)行觀察原假設(shè)或者備擇假設(shè)。

排列檢驗的優(yōu)點:
對于非正態(tài)分布的模型計算p值更準(zhǔn)確
好寫代碼，很好應(yīng)用

排列檢驗的缺點
結(jié)果很容易受到采樣的次數(shù)影響。
可能會運行很慢

三. 結(jié)語

有公式推導(dǎo)錯誤的或理論有謬誤的，請大家指出，我好及時更正，感謝。

感謝大家來閱讀小弟我分享的一些關(guān)于基礎(chǔ)統(tǒng)計的知識。如果將來有時間，小弟會陸續(xù)分享更多知識，例如高等數(shù)據(jù)分析，機(jī)器學(xué)習(xí)，計算機(jī)視覺，神經(jīng)網(wǎng)絡(luò)，優(yōu)化等幾門課的知識分享。也希望小弟可以堅持分享下來吧。

總結(jié)

以上是生活随笔為你收集整理的数据分析模型 第十一章的全部內(nèi)容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網(wǎng)站內(nèi)容還不錯，歡迎將生活随笔推薦給好友。