1、概率統(tǒng)計(jì)的數(shù)值實(shí)驗(yàn)MATLAB在概率統(tǒng)計(jì)教學(xué)中的應(yīng)用,崔明濤 2012年10月11日,引言,而MATLAB 軟件具有簡(jiǎn)單易學(xué)、易操作和繪圖功能強(qiáng)等特點(diǎn)， 利用MATLAB 軟件的圖形可視功能將概率統(tǒng)計(jì)的內(nèi)容用圖形表示出來，通過圖形讓學(xué)生加深理解，以達(dá)到事半功倍的效果。,概率論與數(shù)理統(tǒng)計(jì)知識(shí)比較抽象，邏輯性較強(qiáng)。因此，建議讓學(xué)生結(jié)合理論和公式推導(dǎo)，進(jìn)行數(shù)值試驗(yàn)和相關(guān)調(diào)查，直觀地感受數(shù)學(xué)概念和理論，從而提高學(xué)生解決實(shí)際問題的信心和能力。,概率論,1.rand(m,n)：生成mn的隨機(jī)矩陣，每個(gè)元素都在(0,1) 間，生成方式為均勻分布。2.randn(m,n)：生成mn的隨機(jī)矩陣，每個(gè)元素都在(0

2、,1) 間，生成方式為正態(tài)分布。3.randperm(m)：生成一個(gè)1m的隨機(jī)整數(shù)排列。4.perms(1:n)：生成一個(gè)1n的全排列，共n!個(gè)。5.取整函數(shù)系列：（1）fix(x)：截尾法取整；（2）floor(x)：退一法取整（不超過x的最大整數(shù)）；（3）ceil(x)：進(jìn)一法取整（= floor(x)+1）；（4）round(x)：四舍五入法取整。6.unique(a)：合并a中相同的項(xiàng)。7.prod(x)：向量x的所有分量元素的積。,一、MATLAB常用的與隨機(jī)數(shù)產(chǎn)生相關(guān)的函數(shù)：,示例： rand(1) %生成一個(gè)(0,1)間的隨機(jī)數(shù) ans = 0.8147 rand(2,2) %生

3、成一個(gè)22階(0,1)間的隨機(jī)數(shù)矩陣ans = 0.9134 0.0975 0.6324 0.2785 randperm(5) %生成一個(gè)15的隨機(jī)整數(shù)排列ans = 4 1 5 2 3 a=1 2 4 2 3 3 2;unique(a)ans = 1 2 3 4,例1 隨機(jī)投擲均勻硬幣，觀察國徽朝上與國徽 朝下的頻率。,解 n= 3000100000000;m=0; for i=1:n t=randperm(2);%生成一個(gè)12的隨機(jī)整數(shù)排列 x=t-1;%生成一個(gè)01的隨機(jī)整數(shù)排列 y=x(1); if y=0; m=m+1; end end p1=m/n p2=1-p1,可見當(dāng) 時(shí)，,解

4、 記事件 為第i個(gè)人拿到自已槍，事件 為第i個(gè)人沒拿到自己槍，易知： 又記 為沒有一個(gè)人拿到自己槍的概率。有乘法公式可知：,例2 某班有n個(gè)人，每人各有一支槍，這些槍外形一樣。某次夜間緊急集合，若每人隨機(jī)地取走一支槍，問沒有一個(gè)人拿到自己槍的概率是多少？,于是 所以 特別地，當(dāng)n較大時(shí)， 。 因此，可隨機(jī)模擬出沒有人拿到自己槍的頻率，根據(jù)頻率的穩(wěn)定性，近似當(dāng)做概率，然后去估計(jì)自然常數(shù)e。算法如下：,1、產(chǎn)生n個(gè)隨機(jī)數(shù)的隨機(jī)序列；2、檢驗(yàn)隨機(jī)列與自然列是否至少有一個(gè)配對(duì)；3、對(duì)沒有一個(gè)配對(duì)的序列進(jìn)行累積 p；4、重復(fù)1、2、3步 m 次；5、估計(jì) 。,具體程序及相關(guān)結(jié)果為（注：自然常數(shù) e 2.

5、7183）：, m=40000; n=50; p=0; for j=1:m k=0; sui=randperm(n); for i=1:n if sui(i)=i k=k+1; else k=k; end end if k=0 p=p+1; else p=p; end end e=m/p e = 2.7313,設(shè)針與平行線的夾角為 ，針的中心與最近直線的距離為 。針與平行線相交的充要條件是 ，則所求概率為故可得 的近似計(jì)算公式 ，其中n為隨機(jī)試驗(yàn)次數(shù)，m為針與平行線相交的次數(shù)。,例3 Buffon投針實(shí)驗(yàn) 在畫有許多間距為 的等距平行線的白紙上，隨機(jī)投擲一根長為 的均勻直針，求針與平行線相交的

6、概率，并計(jì)算 的近似值。,解 clear,clf n=10000000;l=0.5;m=0;d=1; for i=1:n x=l/2*sin(rand(1)*pi);y=rand(1)*d/2; if x=y m=m+1; end end p1=m/n pai=2*n*l/(m*d),例4 在100個(gè)人的團(tuán)體中，不考慮年齡差異，研究是否有兩個(gè)以上的人生日相同。假設(shè)每人的生日在一年365天中的任意一天是等可能的，那么隨機(jī)找n個(gè)人(不超過365人)。 (1)求這n個(gè)人生日各不相同的概率是多少？從而求這n個(gè)人中至少有兩個(gè)人生日相同這一隨機(jī)事件發(fā)生的概率是多少？ (2)近似計(jì)算在30名學(xué)生的一個(gè)班中至

7、少有兩個(gè)人生日相同的概率是多少？,解: (1) clear,clf for n=1:100 p0(n)=prod(365:-1:365-n+1)/365n; p1(n)=1-p0(n); end p1=ones(1,100)-p0; n=1:100; plot(n,p0,n,p1,-) xlabel(人數(shù)),ylabel(概率) legend(生日各不相同的概率,至少兩人生日相同的概率) axis(0 100 -0.1 1.199),grid on,p1(30)=0.7063, p1(60)= 0.9941,分析：在30名學(xué)生中至少兩人生日相同的概率為70.63。下面進(jìn)行計(jì)算機(jī)仿真。 隨機(jī)產(chǎn)生

8、30個(gè)正整數(shù)，代表一個(gè)班30名學(xué)生的生日，然后觀察是否有兩人以上生日相同。當(dāng)30個(gè)人中有兩人生日相同時(shí)，輸出“1”，否則輸出“0”。如此重復(fù)觀察100次，計(jì)算出這一事件發(fā)生的頻率 。,(2) clear,clf n=0; for m=1:100 %做100次隨機(jī)試驗(yàn) y=0; x=1+fix(365*rand(1,30); %產(chǎn)生30個(gè)隨機(jī)數(shù) for i=1:29 %用二重循環(huán)尋找30個(gè)隨機(jī)數(shù) 中是否有相同數(shù) for j=i+1:30 if x(i)=x(j) y=1;break; end end end n=n+y; %累計(jì)有兩人生日相同的試驗(yàn)次數(shù) end f=n/m %計(jì)算頻率,f =0.

9、6900 f =0.7900 f =0.6700 f =0.7300 f =0.7500 f =0.6900 f =0.7200 f =0.6700 f =0.6800 ,重復(fù)觀察，數(shù)據(jù)如下：,例5 Galton釘板模型和二項(xiàng)分布 Galton釘板試驗(yàn)是由英國生物統(tǒng)計(jì)學(xué)家和人類學(xué)家Galton設(shè)計(jì)的。故而得名。 通過模擬Calton釘板試驗(yàn)，觀察和體會(huì)二項(xiàng)分布概率分布列的意義、形象地理解De Moivre -Laplace中心極限定理。,共15層小釘,高爾頓釘板試驗(yàn),小球最后落入的格數(shù) ?,記小球向右落下的次數(shù)為 則,記小球向左落下的次數(shù)為 則,符號(hào)函數(shù),大于0返回1,小于0返回-1,等于0返

10、回0,高爾頓( Francis Galton,1822-1911) 英國人類學(xué)家和氣象學(xué)家,O,記,則,近似,高爾頓釘板試驗(yàn),共15層小釘,模擬Galton釘板試驗(yàn)的步驟： (1) 確定釘子的位置：將釘子的橫、縱坐標(biāo)存儲(chǔ)在兩個(gè)矩陣X和Y中。 (2) 在Galton釘板試驗(yàn)中，小球每碰到釘子下落時(shí)都具有兩種可能性，設(shè)向右的概率為p，向左的概率為q1-p，這里p=0.5，表示向左向右的機(jī)會(huì)是相同的。 模擬過程如下：首先產(chǎn)生一均勻隨機(jī)數(shù)u，這只需調(diào)用隨機(jī)數(shù)發(fā)生器指令rand(m,n)。 rand(m,n)指令：用來產(chǎn)生mn個(gè)(0,1)區(qū)間中的隨機(jī)數(shù)，并將這些隨機(jī)數(shù)存于一個(gè)mn矩陣中，每次調(diào)用rand

11、(m,n)的結(jié)果都會(huì)不同。如果想保持結(jié)果一致，可與rand(seed,s)配合使用，這里s是一個(gè)正整數(shù)，例如 rand(seed,1),u=rand(1,6) u = 0.5129 0.4605 0.3504 0.0950 0.4337 0.7092 而且再次運(yùn)行該指令時(shí)結(jié)果保持不變。除非重設(shè)種子seed的值，如 rand(seed,2),u=rand(1,6) u = 0.0258 0.9210 0.7008 0.1901 0.8673 0.4185 這樣結(jié)果才會(huì)產(chǎn)生變化。,將0,1區(qū)間分成兩段，區(qū)間0,p)和p,1。如果隨機(jī)數(shù)u屬于0,p)，讓小球向右落下；若u屬于p,1 ，讓小球向左落下

12、。將這一過程重復(fù)n次，并用直線連接小球落下時(shí)所經(jīng)過的點(diǎn)，這樣就模擬了小球從頂端隨機(jī)地落人某一格子的過程。 (3) 模擬小球堆積的形狀。輸入扔球次數(shù)m(例如m50、 100、500等等)，計(jì)算落在第i個(gè)格子的小球數(shù)在總球數(shù)m中所占的比例，這樣當(dāng)模擬結(jié)束時(shí)，就得到了頻率 用頻率反映小球的堆積形狀。 (4)用如下動(dòng)畫指令制作動(dòng)畫： movien(n)：創(chuàng)建動(dòng)畫矩陣；制作動(dòng)畫矩陣數(shù)據(jù)； Getframe：拷貝動(dòng)畫矩陣； movie(Mat, m)：播放動(dòng)畫矩陣m次。 M文件如下：,解： clear,clf,m=100;n=5;y0=2;%設(shè)置參數(shù) ballnum=zeros(1,n+1); p=0.5

13、;q=1-p; for i=n+1:-1:1 %創(chuàng)建釘子的坐標(biāo)x，y x(i,1)=0.5*(n-i+1); y(i,1)=(n-i+1)+y0; for j=2:i x(i,j)=x(i,1)+(j-1)*1; y(i,j)=y(i,1); end end mm=moviein(m); %動(dòng)畫開始，模擬小球下落路徑 for i=1:m s=rand(1,n); %產(chǎn)生n個(gè)隨機(jī)數(shù) xi=x(1,1);yi=y(1,1);k=1;l=1; %小球遇到第一個(gè)釘子 for j=1:n plot(x(1:n,:),y(1:n,:),o,x(n+1,:),y(n+1,:),.-),%畫釘子的位置 axi

14、s(-2 n+2 0 y0+n+1),hold on,k=k+1; %小球下落一格 if s(j)p l=l+0;%小球左移 else l=l+1;%小球右移 end xt=x(k,l);yt=y(k,l);%小球下落點(diǎn)的坐標(biāo) h=plot(xi,xt,yi,yt);axis(-2 n+2 0 y0+n+1) %畫小球運(yùn)動(dòng)軌跡 xi=xt;yi=yt; end ballnum(l)=ballnum(l)+1; %計(jì)數(shù) ballnum1=3*ballnum./m; bar(0:n,ballnum1),axis(-2 n+2 0 y0+n+1) %畫各格子的頻率 mm(i)=getframe; %

15、存儲(chǔ)動(dòng)畫數(shù)據(jù) hold off end movie(mm,1) %播放動(dòng)畫一次,概率密度函數(shù)（pdf），求隨機(jī)變量X在x點(diǎn)處的概率密度值 累積分布函數(shù)（cdf），求隨機(jī)變量X在x點(diǎn)處的分布函數(shù)值 逆累積分布函數(shù)（inv），求隨機(jī)變量X在概率點(diǎn) 處的分 布函數(shù)反函數(shù)值 均值與方差計(jì)算函數(shù)（stat），求給定分布的隨機(jī)變量X的 數(shù)學(xué)期望E（X）和方差var（X）。 隨機(jī)數(shù)生成函數(shù)（rnd），模擬生成指定分布的樣本數(shù)據(jù)。,二、MATLAB為常見自然概率分布提供了下列5類函數(shù)：,具體函數(shù)的命名規(guī)則是：函數(shù)名分布類型名稱+函數(shù)類型名稱(pdf、cdf、inv、stat、rnd) 其中，分布類型名稱如下：

16、 分布類型 MATLAB名稱,正態(tài)分布 norm 指數(shù)分布 exp 均勻分布 unif 分布 beta 分布 gam 對(duì)數(shù)正態(tài)分布 logn rayleigh分布 rayl weibull 分布 weib 二項(xiàng)分布 bino Poisson分布 poiss 幾何分布 geo 超幾何分布 hyge 離散均勻分布 unid 負(fù)二項(xiàng)分布 nbin,例如，normpdf、normcdf、norminv、normstat和normrnd分別是正態(tài)分布的概率密度、累積分布、逆累積分布、數(shù)字特征和隨機(jī)數(shù)生成函數(shù)。,關(guān)于這5類函數(shù)的語法，請(qǐng)?jiān)斠娪嘘P(guān)書籍。 快捷的學(xué)習(xí)可借助MATLAB的系統(tǒng)幫助，通過指令doc

17、獲得具體函數(shù)的詳細(xì)信息，語法是 doc ,例6 到某服務(wù)機(jī)構(gòu)辦事總是要排隊(duì)等待的。設(shè)等待時(shí)間T是服從指數(shù)分布的隨機(jī)變量(單位：分鐘)，概率密度為設(shè)某人一個(gè)月內(nèi)要到此辦事10次，若等待時(shí)間超過15分鐘，他就離去。求： (1)恰好有兩次離去的概率； (2)最多有兩次離去的概率； (3)至少有兩次離去的概率； (4)離去的次數(shù)占多數(shù)的概率。,解 首先求任一次離去的概率，依題意 設(shè)10次中離去的次數(shù)為X，則 。, p=1-expcdf(15,10) %任一次離去的概率 p1=binopdf(2,10,p)%恰有兩次離去的概率 q=binopdf(0:2,10,p);p2=sum(q)%最多有兩次離去的

18、概率 q=binopdf(0:1,10,p);p3=1-sum(q)%最少有兩次離去的概率 q=binopdf(0:5,10,p);p4=1-sum(q)%離去的次數(shù)占多數(shù)的概率,p = 0.2231 p1 = 0.2972 p2 = 0.6073 p3 = 0.6899 p4 = 0.0112,例7 某一急救中心在長度為t的時(shí)間間隔內(nèi)收到的緊急呼救次數(shù)服從參數(shù)為t2的泊松分布，而與時(shí)間間隔的起點(diǎn)無關(guān)(時(shí)間以小時(shí)計(jì))，求： (1) 在某一天中午12時(shí)至下午3時(shí)沒有收到緊急呼救的概率； (2) 某一天中午12時(shí)至下午5時(shí)至少收到1次緊急呼救的概率。,（1） P1=poisscdf(0,3/2)

19、P1 =0.2231 或者 P1=poisspdf(0,3/2) P1 =0.2231 中午12時(shí)到下午3時(shí)沒有收到緊急呼救的概率為0.2231。 （2） P2=1-poisscdf(0,5/2) P2 = 0.9179 中午12時(shí)至下午5時(shí)至少收到1次緊急呼救的概率為0.9179。,解 本題計(jì)算需調(diào)用函數(shù)poisscdf，其格式為poisscdf(x,)， 返回 的值。,例8 某廠研發(fā)了一種新產(chǎn)品，現(xiàn)要設(shè)計(jì)它的包裝箱，要求每箱至少裝100件產(chǎn)品，且開箱驗(yàn)貨時(shí)，每箱至少裝有100件合格產(chǎn)品的概率不應(yīng)小于0.9，假設(shè)隨機(jī)裝箱時(shí)每箱中的不合格產(chǎn)品數(shù)服從參數(shù)為3的泊松分布。 問：要設(shè)計(jì)的這種包裝箱，

20、每箱至少應(yīng)裝多少件產(chǎn)品才能滿足要求？,解 設(shè)每箱至少裝100+m件產(chǎn)品，X表示每箱中的不合格品數(shù)，則X 服從參數(shù)為3的泊松分布，即 ，依題意，即要求按下面的不等式確定m。, clear;clf,m=0;p=0; while p0.9。即設(shè)計(jì)的包裝箱每箱至少應(yīng)裝106件產(chǎn)品。,例9 某種重大疾病的醫(yī)療險(xiǎn)種，每份每年需交保險(xiǎn)費(fèi)100元，若在這一年中，投保人得了這種疾病，則每份可以得到索賠額10000元，假設(shè)該地區(qū)這種疾病的患病率為0.0002，現(xiàn)該險(xiǎn)種共有10000份保單，問：(1)保險(xiǎn)公司虧本的概率是多少?(2)保險(xiǎn)公司獲利不少于80萬元的概率是多少?,解 設(shè) 表示這一年中發(fā)生索賠的份數(shù)，依題意

21、， 的統(tǒng)計(jì)規(guī)律可用二項(xiàng)分布 來描述。由二項(xiàng)分布與泊松分布的近似計(jì)算關(guān)系有 近似服從參數(shù)為2的泊松分布。 當(dāng)索賠份數(shù)超過100份時(shí)，則保險(xiǎn)公司發(fā)生虧本，虧本的概率為 當(dāng)索賠份數(shù)不超過20份時(shí)，則保險(xiǎn)公司獲利就不少于80萬元，其概率為, p=poisspdf(0:19,2);%計(jì)算出20個(gè)泊松分布概率值 或 p=binopdf(0:19,10000,0.0002); %按二項(xiàng)分布計(jì)算 p2=sum(p) %求出保險(xiǎn)公司獲利不少于80萬元的概率 p2 = 1.0000, p=poisspdf(0:100,2);%計(jì)算101個(gè)泊松分布概率值或 p=binopdf(0:100,10000,0.0002)

22、; %按二項(xiàng)分布計(jì)算 p1=1-sum(p) %求出保險(xiǎn)公司虧本的概率 p1 = 0.0000,例10 設(shè) ，求 ， 。,本題計(jì)算正態(tài)分布的累積概率值，調(diào)用函數(shù)normcdf，其格式為normcdf(x,)，返回 的值。 解： p1=normcdf(6,4,3)-normcdf(3,4,3) p1 = 0.3781 p2=1-normcdf(3,4,3) p2 = 0.6306,例11 繪制正態(tài)分布的密度函數(shù)、分布函數(shù)曲線，并求均值與方差。,解： clear mu=2.5;sigma=0.6; x=(mu-4*sigma):0.005:(mu+4*sigma); y=normpdf(x,mu,

23、sigma); f=normcdf(x,mu,sigma); plot(x,y,-g,x,f,:b) M,V=normstat(mu,sigma) legend(pdf,cdf,-1),M=2.5000 V=0.3600,從圖中可以看出，正態(tài)密度曲線是關(guān)于x對(duì)稱的鐘形曲線(兩側(cè)在處各有一個(gè)拐點(diǎn))，正態(tài)累積分布曲線當(dāng)x時(shí)F(x)0.5。,例12 觀察正態(tài)分布參數(shù)對(duì)密度曲線的影響。,解： clear mu1=2.5;mu2=3;sigma1=0.5;sigma2=0.6; x=(mu2-4*sigma2):0.01:(mu2+4*sigma2); y1=normpdf(x,mu1,sigma1);

24、 %考察均值的影響 y2=normpdf(x,mu2,sigma1); y3=normpdf(x,mu1,sigma1); %考察方差的影響 y4=normpdf(x,mu1,sigma2); subplot(1,2,1) %考察結(jié)果的可視化 plot(x,y1,-g,x,y2,-b) xlabel(fontsize1212,1=2 ) legend(1,2) subplot(1,2,2) plot(x,y3,-g,x,y4,-b) xlabel(fontsize121=2,12 ) legend(1,2),例13 正態(tài)分布參數(shù)和對(duì)變量x取值規(guī)律的約束3準(zhǔn)則。,解： clear,clf %（標(biāo)

25、準(zhǔn)）正態(tài)分布密度曲線下的面積 X=linspace(-5,5,100); Y=normpdf(X,0,1); yy=normpdf(-3,-2,-1,0,1,2,3,0,1); plot(X,Y,k-,0,0,0,yy(4),c-.) hold on plot(-2,-2,0,yy(2),m:,2,2,0,yy(6),m:,-2,-0.5,yy(6),yy(6),m:,0.5,2,yy(6),yy(6),m:) plot(-1,-1,0,yy(3),g:,1,1,0,yy(5),g:,-1,-0.5,yy(5), yy(5),g:,0.5,1,yy(5),yy(5),g:) plot(-3,-

26、3,0,yy(1),b:,3,3,0,yy(7),b:,-3,-0.5,yy(7), yy(7),b:,0.5,3,yy(7),yy(7),b:),hold off text(-0.5,yy(6)+0.005,fontsize1495.44%) text(-0.5,yy(5)+0.005,fontsize1468.26%) text(-0.5,yy(7)+0.005,fontsize1499.74%) text(-3.2,-0.03,fontsize10-3) text(-2.2,-0.03,fontsize10-2) text(-1.2,-0.03,fontsize10-) text(-0.

27、05,-0.03,fontsize10) text(0.8,-0.03,fontsize10+) text(1.8,-0.03,fontsize10+2) text(2.8,-0.03,fontsize10+3),例14 標(biāo)準(zhǔn)正態(tài)分布分位數(shù)的概念圖示。,解 %分位數(shù)示意圖（標(biāo)準(zhǔn)正態(tài)分布，=0.05） clear,clf data=normrnd(0,1,300,1); xalpha1=norminv(0.05,0,1); xalpha2=norminv(0.95,0,1); xalpha3=norminv(0.025,0,1); xalpha4=norminv(0.975,0,1); subp

28、lot(3,1,1) capaplot(data,-inf,xalpha1);axis(-3,3,0,0.45) subplot(3,1,2) capaplot(data,xalpha2,inf);axis(-3,3,0,0.45) subplot(3,1,3) capaplot(data,-inf,xalpha3);axis(-3,3,0,0.45) hold on capaplot(data,xalpha4,inf);axis(-3,3,0,0.45) hold off xalpha1 xalpha2 xalpha3 xalpha4,xalpha1 = -1.6449 xalpha2 =

29、1.6449 xalpha3 = -1.9600 xalpha4 = 1.9600,數(shù)理統(tǒng)計(jì)基礎(chǔ),Matlab統(tǒng)計(jì)工具箱中常見的統(tǒng)計(jì)命令,1、基本統(tǒng)計(jì)量 對(duì)于隨機(jī)變量x，計(jì)算其基本統(tǒng)計(jì)量的命令如下： 均值：mean(x) 標(biāo)準(zhǔn)差：std(x) 中位數(shù)：median(x) 方差：var(x) 偏度：skewness(x) 峰度：kurtosis(x) 2、頻數(shù)直方圖的描繪 A、給出數(shù)組data的頻數(shù)表的命令為： N,X=hist(data,k) 此命令將區(qū)間min(data),max(data)分為k個(gè)小區(qū)間(缺省為10)，返回?cái)?shù)組data落在每一個(gè)小區(qū)間的頻數(shù)N和每一個(gè)小區(qū)間的中點(diǎn)X。 B、描

30、繪數(shù)組data的頻數(shù)直方圖的命令為： hist(data,k),3、參數(shù)估計(jì) A、對(duì)于正態(tài)總體，點(diǎn)估計(jì)和區(qū)間估計(jì)可同時(shí)由以下命令獲得： muhat,sigmahat,muci,sigmaci=normfit(x,alpha) 此命令在顯著性水平alpha下估計(jì)x的參數(shù)（alpha缺省值為5%），返回值muhat是均值的點(diǎn)估計(jì)值，sigmahat是標(biāo)準(zhǔn)差的點(diǎn)估計(jì)值，muci是均值的區(qū)間估計(jì)，sigmaci是標(biāo)準(zhǔn)差的區(qū)間估計(jì)。 B、對(duì)其他分布總體，兩種處理辦法：一是取容量充分大的樣本，按中心極限定理，它近似服從正態(tài)分布，仍可用上面估計(jì)公式計(jì)算；二是使用特定分布總體的估計(jì)命令，常用的命令如： muh

31、at,muci=expfit(x,alpha) lambdahat, lambdaci=poissfit(x,alpha) phat, pci=weibfit(x,alpha),4、正態(tài)總體假設(shè)檢驗(yàn) A、單總體均值的z檢驗(yàn)： h,sig,ci=ztest(x,m,sigma,alpha,tail) 檢驗(yàn)數(shù)據(jù)x關(guān)于總體均值的某一假設(shè)是否成立，其中sigma為已知方差，alpha為顯著性水平，究竟檢驗(yàn)什么假設(shè)取決于tail的取值： tail=0，檢驗(yàn)假設(shè)“x的均值等于m” tail=1，檢驗(yàn)假設(shè)“x的均值大于m” tail=-1，檢驗(yàn)假設(shè)“x的均值小于m” tail的缺省值為0， alpha的缺省

32、值為5%。 返回值h為一個(gè)布爾值，h=1表示可拒絕原假設(shè)， h=0表示不可拒絕原假設(shè)，sig為假設(shè)成立的概率，ci為均值的1- alpha置信區(qū)間。 B、單總體均值的t檢驗(yàn)： h,sig,ci=ttest(x,m,alpha,tail) C、雙總體均值的t檢驗(yàn)： h,sig,ci=ttest2(x,y,alpha,tail),5、非參數(shù)檢驗(yàn)：總體分布的檢驗(yàn) Matlab統(tǒng)計(jì)工具箱提供了兩個(gè)對(duì)總體分布進(jìn)行檢驗(yàn)的命令： A、 h=normplot(x) 此命令顯示數(shù)據(jù)矩陣x的正態(tài)概率圖，如果數(shù)據(jù)來自于正態(tài)分布，則圖形顯示出直線形態(tài)，而其他概率分布函數(shù)顯示出曲線形態(tài)。 B、h=weibplot(x)

33、 此命令顯示數(shù)據(jù)矩陣x的Weibull概率圖，如果數(shù)據(jù)來自于Weibull分布，則圖形顯示出直線形態(tài)，而其他概率分布函數(shù)顯示出曲線形態(tài)。,例15 一道工序用自動(dòng)化車床連續(xù)加工某種零件，由于刀具損壞等會(huì)出現(xiàn)故障。故障是完全隨機(jī)的，并假定生產(chǎn)任一零件時(shí)出現(xiàn)故障機(jī)會(huì)均相同，工作人員是通過檢查零件來確定工序是否出現(xiàn)故障的?，F(xiàn)積累有100次故障紀(jì)錄，故障出現(xiàn)時(shí)該刀具完成的零件數(shù)如下：459 362 624 542 509 584 433 748 815 505 612 452 434 982 640 742 565 706 593 680 926 653 164 487 734 608 428 1153

34、 593 844 527 552 513 781 474 388 824 538 862 659 775 859 755 49 697 515 628 954 771 609 402 960 885 610 292 837 473 677 358 638 699 634 555 570 84 416 606 1062 484 120 447 654 564 339 280 246 687 539 790 581 621 724 531 512 577 496 468 499 544 645 764 558 378 765 666 763 217 715 310 851試觀察該刀具出現(xiàn)故障時(shí)完成

35、的零件數(shù)屬于哪種分布？, %數(shù)據(jù)輸入 x1=459 362 624 542 509 584 433 748 815 505; x2=612 452 434 982 640 742 565 706 593 680; x3=926 653 164 487 734 608 428 1153 593 844; x4=527 552 513 781 474 388 824 538 862 659; x5=775 859 755 49 697 515 628 954 771 609; x6=402 960 885 610 292 837 473 677 358 638; x7=699 634 555 57

36、0 84 416 606 1062 484 120; x8=447 654 564 339 280 246 687 539 790 581; x9=621 724 531 512 577 496 468 499 544 645; x10=764 558 378 765 666 763 217 715 310 851; x=x1 x2 x3 x4 x5 x6 x7 x8 x9 x10; %作頻數(shù)直方圖 hist(x,10) N,X=hist(x,10) %分布的正態(tài)性檢驗(yàn) normplot(x),N = 3 3 7 14 24 22 14 8 3 2 X = 1.0e+003 * 0.1042

37、0.2146 0.3250 0.4354 0.5458 0.6562 0.7666 0.8770 0.9874 1.0978, %參數(shù)估計(jì) muhat,sigmahat,muci,sigmaci=normfit(x) muhat = 594 sigmahat = 204.1301 muci = 553.4962 634.5038 sigmaci = 179.2276 237.1329 刀具壽命服從正態(tài)分布，均值估計(jì)值為594，方差估計(jì)值為204.1301，均值的95%置信區(qū)間為553.4962，634.5038，方差的95%置信區(qū)間為179.2276，237.1329, %假設(shè)檢驗(yàn) h,sig

38、,ci=ttest(x,594)%已知刀具壽命服從正態(tài)分布，方差未知的情況下，檢驗(yàn)壽命均值是否等于594。 h = 0 sig = 1 ci = 553.4962 634.5038 檢驗(yàn)結(jié)果：布爾變量h=0，表示不可拒絕原假設(shè)，說明假設(shè)壽命均值等于594是合理的。 95%置信區(qū)間為553.4962，634.5038完全包括594，估計(jì)精度較高。 sig = 1遠(yuǎn)超過0.05，不可拒絕原假設(shè) 所以可以認(rèn)為刀具平均壽命為594（件）,例16 用模擬試驗(yàn)的方法直觀地驗(yàn)證教材6.2抽樣分布定理一的結(jié)論。 假定變量 ，用隨機(jī)數(shù)生成的方法模擬對(duì) 的500次簡(jiǎn)單隨機(jī)抽樣，每個(gè)樣本的容量為16。利用這5001

39、6個(gè)樣本數(shù)據(jù)直觀地驗(yàn)證樣本均值 的抽樣分布為均值等于60、方差等于2516的正態(tài)分布，即,解 %1、用隨機(jī)數(shù)生成的方法模擬簡(jiǎn)單隨機(jī)抽樣 x=;%生成一個(gè)存放樣本數(shù)據(jù)的空表（維數(shù)可變的動(dòng)態(tài)矩陣） for byk=1:500 %循環(huán)控制，循環(huán)執(zhí)行下面的指令500次，本例中相當(dāng)于500次抽樣 xx=normrnd(60,5,16,1);%生成一個(gè)來自N(60,25)的容量為16 的樣本（列向量） x=x,xx; %將樣本數(shù)據(jù)逐列存入數(shù)表x，可從matlab的變量瀏覽器（workspace）中觀察這個(gè)數(shù)表 end %2、計(jì)算每個(gè)樣本的樣本均值（1500） xmean=mean(x);%可從變量瀏覽器中

40、觀察這500個(gè)數(shù)據(jù) %3、繪制500個(gè)樣本均值數(shù)據(jù)的直方圖k=ceil(1.87*(length(x)-1)(2/5);%確定分組數(shù) h=histfit(xmean,k);%繪制附正態(tài)參考曲線的數(shù)據(jù)直方圖 set(h(1),FaceColor,c,EdgeColor,w)%修飾，設(shè)置直方圖線條顏色與填充色 %4、用這500個(gè)樣本均值數(shù)據(jù)驗(yàn)證樣本均值的均值和方差 M=mean(xmean) %求（1500）樣本的樣本均值的均值 V=var(xmean)%求（1500）樣本的樣本均值的方差,M = 59.9879 V = 1.4129 M = 60.0117 V = 1.3900 M = 59.9

41、749 V = 1.5158 M = 59.9929 V = 1.5757 M = 59.8809 V = 1.6855 ,例17 觀察：用binornd模擬5000次投球過程，觀察小球堆積的情況。, clear;clf, n=5;p=0.5;m=5000;x=0:1:n rand(seed,3) R=binornd(n,p,1,m);%模擬服從二項(xiàng)分布的隨機(jī)數(shù)，相當(dāng)于模擬 投球m次 for I=1:n+1 %開始計(jì)數(shù) k=; k=find(R=(I-1);%find是一個(gè)有用的指令，本語句的作用是找出R中等于(I-1)元素下標(biāo)，并賦予向量k中 h(I)=length(k)/m;%計(jì)算落于編號(hào)

42、(I-1)的格子中的小球頻率 end bar(x,h),axis(-1 6 0 1)%畫頻率圖 title(fontsize18fontname華文新魏5000次投球小球堆積的頻率圖), f=binopdf(x,n,p), bar(x,f), axis(-1 6 0 1) title(fontsize18fontname華文新魏B(5,0.5)理論分布圖),例18 利用隨機(jī)數(shù)樣本驗(yàn)證中心極限定理。,獨(dú)立同分布的隨機(jī)變量的和的極限分布服從正態(tài)分布，通過產(chǎn)生容量為n的poiss分布和exp分布的樣本，研究其和的漸近分布。 算法如下： 產(chǎn)生容量為n的獨(dú)立同分布的隨機(jī)數(shù)樣本，得其均值和標(biāo)準(zhǔn)差； 將隨機(jī)

43、數(shù)樣本和標(biāo)準(zhǔn)化； 重復(fù)、； 驗(yàn)讓所得標(biāo)準(zhǔn)化的隨機(jī)數(shù)樣本和是否服從標(biāo)準(zhǔn)正態(tài)分布。 具體程序如下：, clear n=2000; means=0; s=0; y=; lamda=4; a=lamda; for i=1:n r=poissrnd(a,n,1);%可換成r=exprnd(a,n,1)； means=mean(r);%計(jì)算樣本均值 s=std(r);%計(jì)算樣本標(biāo)準(zhǔn)差 y(i)=sqrt(n).*(means-a)./sqrt(s); end normplot(y);%分布的正態(tài)性檢驗(yàn) title(poiss分布，中心極限定理),例19 在同一坐標(biāo)軸上畫box圖，并對(duì)兩個(gè)班的成績(jī)進(jìn)行初步的

44、分析比較。 兩個(gè)教學(xué)班各30名同學(xué)，在數(shù)學(xué)課程上，A班用新教學(xué)方法組織教學(xué)，B班用傳統(tǒng)方法組織教學(xué)，現(xiàn)得期末考試成績(jī)?nèi)缦隆?A：82，92，77，62，70，36，80，100，74，64，63，56，72，78，68，65，7270，58，92，79,92，65，56，85，73，61，71，42，89 B：57，67，64，54，77，65，71，58，59，69，67，84，63，95，81，46，49, 60, 64，66，74，55，58，63，65，68，76，72，48，72,解 clearx=82,92,77,62,70,36,80,100,74,64,63,56,72,78,

45、68,65,72,70,58,92,79,92,65,56,85,73,61,71,42,89;57,67,64,54,77,65,71,58,59,69,67,84,63,95,81,46,49,60,64,66,74,55,58,63,65,68,76,72,48,72; boxplot(x),從圖中直觀地看出，兩個(gè)班成績(jī)的分布是正態(tài)（對(duì)稱）的，A班成績(jī)較為分散(方差大)，B班成績(jī)則較集中(方差小)。A班成績(jī)明顯高于B班(均值比較并且A班25低分段上限接近B班中值線，A班中值線接近B班25高分段下限)。A班的平均成績(jī)約為70分(中值)，B班約為65分(中值)。A班有一名同學(xué)的成績(jī)過低(離群

46、)，而B班成績(jī)優(yōu)秀的只有一人(離群)。,需要注意的是，從圖中我們不能得出新教學(xué)方法一定優(yōu)于傳統(tǒng)教學(xué)方法的結(jié)論，因?yàn)槲覀儾⒉恢纼蓚€(gè)班級(jí)原有的數(shù)學(xué)基礎(chǔ)是怎樣的。,三、MATLAB也為常用的三大統(tǒng)計(jì)分布提供了相應(yīng)的pdf、cdf、inv、stat、rnd類函數(shù)，具體分布類型函數(shù)名稱如下： 分布類型 MATLAB名稱,分布 chi2 t分布 t F分布 f 非中心 分布 ncx2 非中心t分布 nct 非中心F分布 ncf,例20 分布的密度函數(shù)曲線。,解： %繪制不同自由度的卡方分布概率密度曲線 clear,clf X=linspace(0,20,100); Y1=chi2pdf(X,1);%自由

47、度等于1 Y2=chi2pdf(X,3);%自由度等于3 Y3=chi2pdf(X,6);%自由度等于6 plot(X,Y1,-g,X,Y2,-b,X,Y3,-k) title(fontsize18fontname華文新魏不同自由度的chi2分布概率密度曲線的比較) text(0.6,0.6,fontsize12df:n=1) text(2.6,0.2,fontsize12df:n=3) text(8.6,0.09,fontsize12df:n=6) legend(df:n=1,df:n=3,df:n=6),解： %繪制t分布概率密度曲線 clear,clf X=linspace(-4,4,100); Y0=normp