隨機樣本的生成

5.6隨機樣本的生成主要問題：我們已經(jīng)介紹了許多描述隨機變量行為的方法：（變換、分布、矩的計算以及極限定理等等）實際應(yīng)用中，通常是先采集隨機變量的觀測值，然后將這些隨機變量用于描述實際現(xiàn)象并進行建模。通常考慮：給定隨機變量，利用分布函數(shù)的性質(zhì)刻畫隨機變量的行為。如何生成服從給定的分布的隨機數(shù)？這里主要研究：給定分布，生成相應(yīng)的隨機樣本。我們來看一個例子（5.6.1）：設(shè)某類電子元件的使用壽命服從參數(shù)為的指數(shù)分布，我們考慮個元件中至少有個使用壽命大于等于小時的概率。

假設(shè)各元件之間獨立，則可將對個元件的檢測看作伯努利實驗（獨立重復(fù)試驗），有明顯地，的計算量有些繁重，不過這里服從指數(shù)分布，因此能表示為初等函數(shù)：計算的一個模擬算法是生成具有相應(yīng)分布的隨機變量，并用弱大數(shù)定律證實模擬的可靠性，令，獨立同分布，則有：樣本空間足夠大時，樣本均值依概率收斂于分布的期望。計算步驟：對任意的1）生成獨立同分布且服從參數(shù)為的指數(shù)分布的隨機變量；2）如果至少有個使用壽命，則令；否則令。于是由且可知，主要思路：首先,是研究怎樣生成我們需要的隨機變量。其次,要用大數(shù)定律證明模擬算法所得近似結(jié)果的可靠性。問題轉(zhuǎn)化：假定能夠生成獨立同分布的均勻隨機變量。事實上，有許多可以生成偽隨機數(shù)的算法。注：偽隨機數(shù)，即不是真正的隨機數(shù)，是一種具有非常長周期的能通過數(shù)理統(tǒng)計中均勻性檢驗的數(shù)列，是均勻隨機數(shù)的一種可行近似。如何由獨立同分布均勻隨機變量得到滿足其他分布的隨機變量？主要方法：直接法——變換間接法——舍選法5.6.1直接法理論基礎(chǔ)：概率積分變換（定理2.1.10）設(shè)隨機變量有連續(xù)累積分布函數(shù)，令，則服從上的均勻分布，即這里我們主要是利用上述變換的逆，在已知一系列上均勻分布的情況下得到滿足其他已知分布的隨機變量。5.6.1直接法不難看出，若函數(shù)可以表示成初等函數(shù)，隨機變量服從區(qū)間上的均勻分布時，變換后的隨機變量滿足某指定分布，則可以使用直接法生成隨機變量。換句話說，即要求方程有初等函數(shù)解。直接法的關(guān)鍵在于求。例1（5.6.3）：設(shè)隨機變量的累積分布函數(shù)為，服從區(qū)間上的均勻分布，則隨機變量服從分布，若服從參數(shù)為的指數(shù)分布，則也服從參數(shù)為的指數(shù)分布。程序演示：clearallclc%生成均勻隨機變量n=10000;u=unifrnd(0,1,1,n);%求均勻隨機變量樣本均值和方差i=1:n;jz=sum(u(i))/nfc=sum((u(i)-jz).^2)/(n-1)%生成服從參數(shù)為2的指數(shù)分布的隨機變量y(i)=-2*log(1-u(i));%求生成的服從參數(shù)為2的指數(shù)分布的隨機變量的樣本均值和方差yjz=sum(y(i))/nyfc=sum((y(i)-yjz).^2)/(n-1)%根據(jù)生成的隨機變量作頻率直方圖a=min(y);b=max(y);m=200;de=(b-a)/m;[r,xout]=hist(y,[a:de:b]);f=r./(n*de);bar(xout,f)holdon%與參數(shù)為2的指數(shù)分布密度函數(shù)比較x=0:0.1:20;yx=exppdf(x,2);plot(x,yx,'r')axis([02000.5])鑒于指數(shù)分布與其他分布的聯(lián)系，我們可以進一步生成許多不同類型的隨機變量。例如：根據(jù)指數(shù)分布與伽馬分布的關(guān)系及伽馬分布的可加性容易得到：若則不存在初等函數(shù)解時，我們尋找新類型的變換：例2（5.6.4）：Box-Muller算法是一對獨立的隨機變量。所以這里盡管沒有生成單個隨機變量的快速變換，我們卻可以同時生成兩個隨機變量。首先生成服從區(qū)間上均勻分布的兩個隨機變量和，然后令則程序演示：clearallclc%生成兩個均勻隨機變量n=10000;u=unifrnd(0,1,1,n);v=unifrnd(0,1,1,n);%求均勻隨機變量樣本均值和方差i=1:n;ujz=sum(u(i))/nufc=sum((u(i)-ujz).^2)/(n-1)vjz=sum(v(i))/nvfc=sum((v(i)-vjz).^2)/(n-1)%生成兩個服從n（0,1）分布的隨機變量forj=1:n

x(j)=sqrt(-2*log(u(j)))*cos(2*pi*v(j));

y(j)=sqrt(-2*log(u(j)))*sin(2*pi*v(j));endx;y;%求生成的兩個服從n（0,1）分布的隨機變量的樣本均值和方差xjz=sum(y(i))/nxfc=sum((y(i)-xjz).^2)/(n-1)yjz=sum(y(i))/nyfc=sum((y(i)-yjz).^2)/(n-1)%根據(jù)生成的隨機變量作頻率直方圖與n（0,1）分布的密度函數(shù)比較subplot(1,2,1);a=min(x);b=max(x);m=200;de=(b-a)/m;[r,xout]=hist(x,[a:de:b]);f=r./(n*de);bar(xout,f)holdonz1=-5:0.1:5;xz1=normpdf(z1,0,1);plot(z1,xz1,'r')axis([-5500.5])subplot(1,2,2);a=min(y);b=max(y);m=200;de=(b-a)/m;[r,xout]=hist(y,[a:de:b]);f=r./(n*de);bar(xout,f)holdonz2=-5:0.1:5;yz2=normpdf(z2,0,1);plot(z2,yz2,'r')axis([-5500.5])例3：標準正態(tài)隨機數(shù)的生成設(shè)服從區(qū)間上的均勻分布，我們可以利用中心極限定理來生成標準正態(tài)分布的隨機數(shù)。這里用到了林德貝格-勒維中心極限定理：設(shè)是獨立同分布的隨機變量序列，且，。記對任意實數(shù)，有程序演示：clearallclcn1=12;n2=1000;forj=1:n2%生成12個均勻隨機變量

u=unifrnd(0,1,1,n1);i=1:n1;%生成1個服從標準正態(tài)分布的隨機變量

y(j)=sum(u(i))-6;end%求生成的服從標準正態(tài)分布的隨機變量的樣本均值和方差k=1:n2;yjz=sum(y(k))/n2yfc=sum((y(k)-yjz).^2)/(n2-1)%根據(jù)生成的隨機變量作頻率直方圖a=min(y);b=max(y);m=200;de=(b-a)/m;[r,xout]=hist(y,[a:de:b]);f=r./(n2*de);bar(xout,f)holdon%與標準正態(tài)分布密度函數(shù)比較x=-5:0.1:5;yx=normpdf(x,0,1);plot(x,yx,'r')axis([-5500.7])離散隨機變量的生成：若離散隨機變量可以取值，則有可以得到生成離散隨機變量的步驟：1）生成區(qū)間上服從均勻分布的隨機變量；2）如果，則令，這里，定義且。例4（5.6.5）：設(shè)服從區(qū)間上的均勻分布，我們可以用來生成服從參數(shù)為的二項分布隨機變量。程序演示：clearallclcn=1000;%生成均勻隨機變量u=unifrnd(0,1,1,n);%生成服從參數(shù)為（4,5/8）的二項分布的隨機變量fori=1:nif0<u(i)<=0.020

y(i)=0

elseif0.020<u(i)<=0.152

y(i)=1

elseif0.152<u(i)<=0.481

y(i)=2

elseif0.481<u(i)<=0.847

y(i)=3

elseif0.847<u(i)<=1

y(i)=4endend%求生成的服從參數(shù)為（4,5/8）的二項分布的樣本均值和方差k=1:n;yjz=sum(y(k))/nyfc=sum((y(k)-yjz).^2)/(n-1)5.6.2間接法舍選法：（無法通過直接變換生成隨機變量）先看一個簡單的例子（5.6.7）：設(shè)待生成的目標變量服從參數(shù)為的貝塔分布.。如果，均為整數(shù)，則可用直接變換法生成；如果，不是整數(shù)，則直接法失效。不妨令

，。在下圖中我們將貝塔概率密度函數(shù)置于長為、寬為的矩形中?，F(xiàn)在我們計算：設(shè)是一對獨立隨機變量且都服從區(qū)間上的均勻分布，則圖中陰影部分概率為：于是，我們可以根據(jù)均勻分布概率計算貝塔概率，這就表明我們可以利用均勻隨機變量生成貝塔隨機變量。令，則有，因此于是，我們可以得到生成參數(shù)為的貝塔隨機變量的下列算法：1）生成區(qū)間上服從均勻分布的隨機變量；2）如果，則令；否則返回步驟1）。只要，該算法就能生成參數(shù)為的貝塔隨機變量。注意：上面的算法并不完整，因為我們不知道要用多少對才能生成。如果定義隨機變量則由可知，服從參數(shù)為的幾何分布，因此，生成一個通常需要對，所以越小算法越好，所以的最優(yōu)值為。算法改進觀察發(fā)現(xiàn)，上面的算法中并沒有用到區(qū)域，由此我們可以換用與貝塔隨機變量相近的其他變量以改進算法。設(shè)，可以計算。

則步驟2）可以推廣為比較。

比值越大，“看起來越像”概率密度函數(shù)

為的隨機變量，的可能性也越大。5.6.3舍選法：設(shè)，其中有相同的支撐集且按下列步驟可生成隨機變量：1）生成區(qū)間上服從均勻分布的隨