汲守峰,劉 卉
(唐山學(xué)院 基礎(chǔ)教學(xué)部,河北 唐山 063000)
概率統(tǒng)計(jì)是研究自然界中隨機(jī)現(xiàn)象統(tǒng)計(jì)規(guī)律的一門(mén)數(shù)學(xué)方法,廣泛應(yīng)用在金融、經(jīng)濟(jì)、生物、醫(yī)學(xué)、運(yùn)籌管理和工程技術(shù)等領(lǐng)域。在概率統(tǒng)計(jì)過(guò)程中幾乎每個(gè)環(huán)節(jié)都離不開(kāi)統(tǒng)計(jì)軟件的輔助。目前的統(tǒng)計(jì)軟件主要包括Matlab,SAS,SPSS等商用專(zhuān)業(yè)軟件,但這些軟件除運(yùn)行環(huán)境封閉、下載安裝復(fù)雜、內(nèi)存占用較高外,還需要對(duì)不同的插件支付額外的專(zhuān)利費(fèi)用,而R軟件則不受這些條件的約束。由Ross Ihaka和Robert Gentleman開(kāi)發(fā)的面向?qū)ο蟮腞軟件,是一款免費(fèi)開(kāi)源且能夠自由有效地用于統(tǒng)計(jì)計(jì)算和繪圖的計(jì)算機(jī)軟件[1],它提供了廣泛的統(tǒng)計(jì)分析和繪圖技術(shù),其功能包括程序編輯與運(yùn)算、數(shù)據(jù)存儲(chǔ)與處理、數(shù)組運(yùn)算(其向量、矩陣運(yùn)算功能尤其強(qiáng)大),除此之外,用戶還可以根據(jù)自己的需要安裝現(xiàn)成的統(tǒng)計(jì)軟件包,它支持對(duì)包的代碼進(jìn)行修改和重新編寫(xiě)[2-3],因此,此軟件被統(tǒng)計(jì)學(xué)家、工程師和科學(xué)家廣泛使用。如張志成[4]使用R軟件對(duì)經(jīng)典的蒲豐投針實(shí)驗(yàn)進(jìn)行了隨機(jī)模擬,從理論和實(shí)踐中得到了圓周率π的近似值;李秀敏等[5]使用R軟件通過(guò)隨機(jī)模擬實(shí)驗(yàn)研究了統(tǒng)計(jì)學(xué)中幾種比較重要的抽樣分布問(wèn)題,驗(yàn)證了中心極限定理的正確性;熊炳忠[6]則是利用R軟件對(duì)概率分布、大數(shù)定律、中心極限定理與假設(shè)檢驗(yàn)等進(jìn)行了模擬驗(yàn)證和分析。本文應(yīng)用R軟件,對(duì)概率統(tǒng)計(jì)中中心極限定理的理論結(jié)果進(jìn)行數(shù)值模擬,借助R軟件強(qiáng)大的數(shù)據(jù)處理和計(jì)算功能簡(jiǎn)化假設(shè)檢驗(yàn)中龐雜的數(shù)據(jù)計(jì)算過(guò)程,并利用其數(shù)據(jù)解析和圖形繪制功能進(jìn)行線性回歸分析,由此簡(jiǎn)化抽象復(fù)雜的問(wèn)題,提高運(yùn)算效率,增加檢驗(yàn)結(jié)果的可視化程度。
中心極限定理在概率統(tǒng)計(jì)中具有十分重要的地位,主要研究獨(dú)立的隨機(jī)變量序列之和的分布近似服從正態(tài)分布的有關(guān)問(wèn)題。中心極限定理主要包含三大定理,其他定理的證明以及某些統(tǒng)計(jì)推斷都是建立在三大定理的理論之上。中心極限定理描述的內(nèi)容較為抽象,當(dāng)涉及的統(tǒng)計(jì)模型數(shù)據(jù)較多時(shí)也不太容易計(jì)算和驗(yàn)證,可借助R軟件對(duì)概率模型進(jìn)行編程模擬,然后輸出統(tǒng)計(jì)制圖和數(shù)據(jù)計(jì)算結(jié)果,增加結(jié)果的可視化程度。下面通過(guò)對(duì)服從二項(xiàng)分布和指數(shù)分布的隨機(jī)變量序列之和的直方圖與正態(tài)分布密度曲線對(duì)比,驗(yàn)證棣莫弗-拉普拉斯定理和萊維定理理論結(jié)果的正確性。
棣莫弗-拉普拉斯中心極限定理是關(guān)于二項(xiàng)分布漸近趨于正態(tài)分布的極限定理,也稱(chēng)二項(xiàng)分布的中心極限定理。假設(shè)隨機(jī)變量X~B(n,p),依據(jù)棣莫弗-拉普拉斯中心極限定理,隨著n→∞,X的分布將依概率收斂于正態(tài)分布N(np,np(1-p))。除用嚴(yán)格的數(shù)學(xué)證明外,可應(yīng)用R軟件對(duì)其進(jìn)行統(tǒng)計(jì)模擬并加以驗(yàn)證。
>layout(matrix(c(1,2,3,4),ncol=2,byrow=T))
sim<-function(m=20,n=50,p=0.2)
{y<-rbinom(m,n,p)
x=(y-n*p)/sqrt(n*p*(1-p))
hist(x,prob=T,breaks=30,main=paste("n=",n,"p=",p,"m=",m))
curve(dnorm(x),add=T)}
sim()
sim(200)
sim(2000)
sim(20000)
輸出統(tǒng)計(jì)制圖結(jié)果,圖1為隨機(jī)產(chǎn)生的四組來(lái)自于B(50,0.2)的隨機(jī)變量序列統(tǒng)計(jì)直方圖(柱狀圖)與對(duì)應(yīng)的正態(tài)分布密度曲線(曲線圖),比較發(fā)現(xiàn),隨著產(chǎn)生的個(gè)數(shù)m越大,兩者近似效果越好,直觀地解釋和驗(yàn)證了棣莫弗-拉普拉斯中心極限定理。
>layout(matrix(c(1:4),ncol=2,byrow=T))
lambda<-0.05
for(n in c(5,15,30,50)){
mu<-n/lambda
sumx<-numeric(1000)
sdsumx<-sqrt(n)/lambda
for(i in 1:1000){
sumx[i]<-sum(rexp(n,rate=0.05))}
hist(sumx,prob=T,main=paste("hist ogram.sumx,n=",n),col=gray(.5),lwd=2)
real<-dnorm(seq(mu-3*sdsumx,mu+3*sdsumx,0.01),mu,sdsumx)
lines(seq(mu-3*sdsumx,mu+3*sdsumx,0.01),real,lty=1,col=2,lwd=2)
box()}
由圖2可知,當(dāng)n>15時(shí),樣本值和的直方圖與正態(tài)分布概率密度曲線近似精度較高,直觀解釋了中心極限定理中獨(dú)立隨機(jī)變量和隨變量個(gè)數(shù)增加而趨近于正態(tài)分布的結(jié)論。
假設(shè)檢驗(yàn)往往會(huì)涉及大量的計(jì)算,有些計(jì)算簡(jiǎn)單但需要多次重復(fù),而有些計(jì)算需要一些特殊的技巧或查閱相關(guān)分布表,如F分布、t分布、χ2分布等分位數(shù)和分位點(diǎn)的計(jì)算。實(shí)踐中很多統(tǒng)計(jì)計(jì)算都需借助計(jì)算機(jī)軟件,否則會(huì)使統(tǒng)計(jì)工作難以高效開(kāi)展。下面應(yīng)用R軟件對(duì)大學(xué)生的期末考試成績(jī)進(jìn)行統(tǒng)計(jì)分析。
選取兩個(gè)年級(jí)GM18和GM19《概率統(tǒng)計(jì)》期末考試成績(jī),并假設(shè)GM18~N(mu1,sigma1^2),GM19~N(mu2,sigma2^2),應(yīng)用R軟件計(jì)算兩個(gè)年級(jí)置信水平為90%的平均成績(jī)差的置信區(qū)間:
>GM18<-c(81,69,76,62,67,60,77,63,71,76,70,42,76,86,88,74,100,37,81,31,60,68,68,82,79,74,81,98,80,5,20,60,36,57,47,68,61,56,48,43,63,63,78,68,42,69,38,81,74,83,76,83,96,44,69,54,85,50,64,78,84,66,70,52,92,87,93,37,58,36,12,26,4,42,36,36,61,51,63,4)
GM19<-c(51,84,68,84,93,61,99,75,74,62,54,66,86,56,61,57,77,78,73,92,18,59,48,87,68,36,68,76,53,69,67,65,49,50,57,81,52,78,48,94,45,59,68,82,27,77,77,96,59,73,58,63,55,86,77,62,67,89,47,54,73,77,65,52,48,81,70,72,81,77,80,65,67,60,64,34,61,20,46,51,59,57,82,51,62)
t.test(GM18,GM19,conf.level=0.9)
Welch Two Sample t-test
data:GM18 and GM19
t=-1.1433,df=145.35,p-value=
0.2548
alternative hypothesis: true difference in means is not equal to 0
90 percent confidence interval:
-8.492335 1.554100
sample estimates:
mean of x mean of y
61.82500 65.29412
結(jié)果顯示,在方差不等的情況下,均值差的置信水平為90%的置信區(qū)間為(-8.492 335,1.554 100),0被包含在區(qū)間內(nèi)部,由此認(rèn)為兩個(gè)年級(jí)的“成績(jī)相同”。但兩個(gè)年級(jí)成績(jī)均值分別為61.825 00和65.294 12,二者是否存在顯著性差異?可通過(guò)R軟件進(jìn)一步驗(yàn)證。
>var.test(GM18,GM19)
F test to compare two variances
data:GM18 and GM19
F=1.8175,num df=79,denom df=84,p-value=0.007404
alternative hypothesis: true ratio of variances is not equal to 1
95 percent confidence interval:
1.175249 2.819119
sample estimates:
ratio of variances
1.81751
可以看出,兩個(gè)正態(tài)總體方差以95%的置信水平落入?yún)^(qū)間(1.175 249,2.819 119)內(nèi),區(qū)間端點(diǎn)明顯大于1,因此可認(rèn)為兩個(gè)年級(jí)的成績(jī)存在顯著性差異,且p=0.007 404(<0.05),也支持該結(jié)論的正確性。
在刑事科學(xué)技術(shù)中,辦案人員往往根據(jù)現(xiàn)場(chǎng)遺留的蛛絲馬跡尋找案件的突破口。例如利用從現(xiàn)場(chǎng)提取的足跡長(zhǎng)度推算出犯罪嫌疑人身高的近似值?,F(xiàn)隨機(jī)抽取10個(gè)樣本,并測(cè)得以下數(shù)據(jù)(見(jiàn)表1),應(yīng)用R軟件對(duì)其進(jìn)行回歸分析。
表1 樣本的足跡長(zhǎng)度與身高 cm
利用R軟件做散點(diǎn)圖(圖3),觀察兩者的大致關(guān)系。
圖3 足跡長(zhǎng)度與身高關(guān)系的散點(diǎn)圖
>x<-c(21.6,22.3,23.6,24.6,25.3,25.8,26.7,27.1,28.2,28.4)
y<-c(156.3,160.8,165.3,170.1,172.6,173.6,179.1,179.2,185.2,186.6)
plot(x,y,ylab="f(x)",type="b",col=2,lwd=2)
應(yīng)用R軟件內(nèi)嵌函數(shù)做出線性回歸直線,與散點(diǎn)圖進(jìn)行比較(圖4)。
>cb<-lm(formula=y~x)
summary(cb)
summary(cb)$coefficients[,1]
nihe<-predict(cb)
plot(x,y,ylab="f(x)",type="p",col=1,lwd=2)
lines(x,nihe,col=1,lwd=2)
legend("topleft",c("nihe","sandian"),lty=1:2,col=1:1,lwd=2)
圖4 足跡長(zhǎng)度與身高的散點(diǎn)圖與擬合曲線
圖3反映出足跡長(zhǎng)度與身高存在較為明顯的線性關(guān)系,與圖4對(duì)比不難看出,回歸直線與散點(diǎn)圖存在較為顯著的近似關(guān)系,這為辦案人員根據(jù)足跡長(zhǎng)度預(yù)測(cè)犯罪嫌疑人身高提供了科學(xué)依據(jù)。
應(yīng)用R軟件對(duì)概率統(tǒng)計(jì)問(wèn)題進(jìn)行了輔助研究。在數(shù)據(jù)處理中通過(guò)R軟件的合理應(yīng)用,省去了復(fù)雜的計(jì)算和繁瑣的推導(dǎo)過(guò)程,增強(qiáng)了數(shù)據(jù)的處理速度和統(tǒng)計(jì)制圖的繪制能力,有效提升了工作效率。