隨機(jī)樣本的生成

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

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

u=unifrnd(0,1,1,n1);i=1:n1;%生成1個(gè)服從標(biāo)準(zhǔn)正態(tài)分布的隨機(jī)變量

y(j)=sum(u(i))-6;end%求生成的服從標(biāo)準(zhǔn)正態(tài)分布的隨機(jī)變量的樣本均值和方差k=1:n2;yjz=sum(y(k))/n2yfc=sum((y(k)-yjz).^2)/(n2-1)%根據(jù)生成的隨機(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%與標(biāo)準(zhǔn)正態(tài)分布密度函數(shù)比較x=-5:0.1:5;yx=normpdf(x,0,1);plot(x,yx,'r')axis([-5500.7])離散隨機(jī)變量的生成：若離散隨機(jī)變量可以取值，則有可以得到生成離散隨機(jī)變量的步驟：1）生成區(qū)間上服從均勻分布的隨機(jī)變量；2）如果，則令，這里，定義且。例4（5.6.5）：設(shè)服從區(qū)間上的均勻分布，我們可以用來生成服從參數(shù)為的二項(xiàng)分布隨機(jī)變量。程序演示：clearallclcn=1000;%生成均勻隨機(jī)變量u=unifrnd(0,1,1,n);%生成服從參數(shù)為（4,5/8）的二項(xiàng)分布的隨機(jī)變量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）的二項(xiàng)分布的樣本均值和方差k=1:n;yjz=sum(y(k))/nyfc=sum((y(k)-yjz).^2)/(n-1)5.6.2間接法舍選法：（無法通過直接變換生成隨機(jī)變量）先看一個(gè)簡單的例子（5.6.7）：設(shè)待生成的目標(biāo)變量服從參數(shù)為的貝塔分布.。如果，均為整數(shù)，則可用直接變換法生成；如果，不是整數(shù)，則直接法失效。不妨令

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

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

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

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