數(shù)理統(tǒng)計直觀教學(xué)的實驗設(shè)計與R程序?qū)崿F(xiàn)

呂書龍， 劉文麗， 梁飛豹， 葉福玲

(福州大學(xué) 數(shù)學(xué)與計算機科學(xué)學(xué)院， 福建 福州　350116)

?

數(shù)理統(tǒng)計直觀教學(xué)的實驗設(shè)計與R程序?qū)崿F(xiàn)

呂書龍， 劉文麗， 梁飛豹， 葉福玲

(福州大學(xué) 數(shù)學(xué)與計算機科學(xué)學(xué)院， 福建 福州350116)

針對數(shù)理統(tǒng)計中關(guān)于格列汶科定理、正態(tài)抽樣定理、點估計和區(qū)間估計在教學(xué)中的常見疑惑和問題進行討論與分析，通過實驗設(shè)計和R程序，以圖形方式加以直觀解決。該實驗設(shè)計一方面促進了學(xué)生對于數(shù)理統(tǒng)計理論與方法的理解、提升了教學(xué)效果，另一方面也為R軟件在教學(xué)上的應(yīng)用提供一種思路。

實驗設(shè)計； R軟件； 數(shù)理統(tǒng)計； 直觀教學(xué)

數(shù)理統(tǒng)計獨特的思維方式、抽象的理論、多學(xué)科知識的融合和豐富的應(yīng)用，使得學(xué)生在學(xué)習(xí)這門課程時存在著一定的困難[1]。數(shù)理統(tǒng)計教學(xué)的重點和難點在于如何把統(tǒng)計思想和統(tǒng)計方法闡述清楚，并將其應(yīng)用到實際問題中，這一點在文獻[2-3]中有充分闡述。關(guān)于統(tǒng)計思想和方法論的教育，有很多學(xué)者發(fā)表了富有建設(shè)性的觀點，例如在教學(xué)中充分利用統(tǒng)計軟件[4]和仿真技術(shù)[5]，強化實驗教學(xué)[6]和設(shè)計性實驗[7-10]等。

本文借助R軟件強大的隨機模擬和繪圖功能[11]，針對數(shù)理統(tǒng)計教學(xué)中的疑難點，特別是對抽象的統(tǒng)計思想及原理，提出以實驗設(shè)計和圖形展示為主的直觀教學(xué)[12]模式，并以格列汶科定理、正態(tài)抽樣定理、點估計和區(qū)間估計的教學(xué)為例予以說明，目的是促進學(xué)生對課程內(nèi)容的理解、提升教學(xué)效果，同時也嘗試培養(yǎng)學(xué)生的動手能力和創(chuàng)新應(yīng)用能力。

1　格列汶科定理

設(shè)總體X的分布函數(shù)為F(x)，經(jīng)驗分布函數(shù)為Fn(x)，則格列汶科定理可以描述成：

ns=c(30,50,100,200,400,500,800,1000)

y=seq(-4,4,by=0.001)

par(mfrow=c(2,4))

for (i in 1:8)

{plot(y, pnorm(y), type=′l′)

x=rnorm(ns[i])

lines(ecdf(x),cex=0.1)

text(-2.5,0.9,paste(′n=′,ns[i],sep=″′))

}

由圖1可粗略地看出，經(jīng)驗分布函數(shù)要較好地逼近分布函數(shù)，n不能太小。當n為100和200時，從逼近程度看還是可以接受的；若想得到更好的逼近效果，n應(yīng)不小于200。從教學(xué)角度講，直觀地展示定理結(jié)論是一方面，而了解定理的適用條件也是必要的。因此建議把尋找滿足某個逼近標準的最小樣本容量n的問題作為課外探索題，這樣更能增強學(xué)生對定理的理解程度。

圖1　經(jīng)驗分布函數(shù)與分布函數(shù)比較

該定理的特別之處在于它研究的是|Fn(x)-F(x)|變量的偏差上界，雖然理論上已經(jīng)對它的分布做了證明并能給出相應(yīng)的概率，但還是沒給出n多大時才能較可靠地使用這個定理。不妨通過抽樣模擬來刻畫這個過程：以標準正態(tài)分布為例，讓樣本容量n從10變到500，對每個n，模擬100次得到100個樣本，并計算每個樣本的經(jīng)驗分布函數(shù)與分布函數(shù)偏差的最大值，再計算最大值序列的極差，最后畫極差圖(見圖2)，R程序如下：

getmax=function(size,n)

{a=numeric(size)

for(i in 1:size){

x=rnorm(n); tmp=ks.test(x,′pnorm′)

a[i]=tmp$statistic} #利用ks檢驗得到最大偏差

return(a)

}

n=seq(10,500,by=10);nlen=length(ns)

mat=matrix(0,nrow=2,ncol=nlen)

for(i in 1:nlen) mat[,i]=range(getmax(100, n[i]))

plot(n, mat[1,], ylim=c(0,max(mat[2,])+0.01), cex=0.5)

points(n, mat[2,], cex=0.5)

for(i in 1:nlen) lines(c(n[i], n[i]), mat[,i])

從圖2可直觀看出，極差變化先隨n的變大快速變窄，當n≈200時，極差變化趨于平緩，之后基本保持不變?？梢妌≥200是個不錯的選擇，這也進一步肯定了圖1的結(jié)論。

圖2　最大偏差的極差圖

2　正態(tài)抽樣定理的一個直觀描述

n=100;times=1000;mean=1;sd=2;N1=99; N2=100;sn=numeric(times);tmp=(n-1)/(sd^2)

for(i in 1:times){

set.seed(i); x=rnorm(n,mean,sd);sn[i]=tmp*var(x)}

plot(ecdf(sn),verticals=TRUE,do.point=FALSE)

nx=seq(70,140,by=1)

lines(nx,pchisq(nx,N1),lty=2,lwd=2)

lines(nx,pchisq(nx,N2),lty=4)

legend(130,0.8,c(′經(jīng)驗分布函數(shù)′,′chisq(99)′, ′chisq(100)′),lty=c(1,2,4), cex=0.75, lwd=c(1,2,1)

圖3　(n-1)S2/σ2與χ2(99)，χ2(100)的比較

3　矩估計和極大似然估計的比較

矩估計方法簡單易用，但也有明顯的缺點，例如估計量表示不唯一、樣本信息利用不充分等。極大似然估計能克服矩估計的缺點，在實際中應(yīng)用更廣泛。相信學(xué)生并不喜歡這樣的比較，而若能配置一個實例進行直觀展示，將會收到良好的教學(xué)效果。實驗設(shè)計為：從總體U(0,θ) (不妨設(shè)θ=10)中抽取容量n=20的一個樣本，計算θ的矩估計和極大似然估計，重復(fù)40次后繪制估計值(見圖4)，R程序如下：

#est1和est2 分別記錄每次的矩估計和極大似然估計值

n=20; times=40;est1=rep(0,m);est2=rep(0,m)

for(i in 1:times)

{ x=runif(n,0,10); est1[i]=2*mean(x);est2[i]=max(x)}

plot(1:m,est1,type=″l″, ylim=c(min(c(est1,est2))-0.1,max(c(est1,est2))+0.1) )

points(1:m,est1); lines(1:m,est2)

points(1:m,est2,pch=16);abline(h=10)#畫出真值直線

圖4　均勻分布參數(shù)兩種估計的比較

4　參數(shù)的置信區(qū)間

給定置信度1-α，從正態(tài)總體中抽取50個隨機數(shù)，計算均值的置信區(qū)間；重復(fù)該過程N次(N=100)得到N個置信區(qū)間并繪制成空心端點垂線圖，然后繪制真值水平線，將沒被水平線穿過的區(qū)間用實心端點表示(見圖5)，R程序如下：

Confidence = function(n=50,times=100,mu=0,sd=1,alpha=0.05)

{interval=matrix(0,nrow=times, ncol=2)

for(i in 1:times){

x=rnorm(n,mu,sd);mx=mean(x)

u_half=qnorm(1-alpha/2,mu,sd)*sd/sqrt(n)

interval[i,]=c(mx-u_half,mx+u_half)}

miny=min(interval)-0.01;maxy=max(interval)+0.01

plot(1:times,seq(miny,maxy,length=times),type=′n′)

abline(h=mu)

for(i in 1:times){

if(prod(interval[i,]-mu)>0) points(c(i,i),interval[i,],pch=16)

else points(c(i,i),interval[i,],pch=1)

lines(c(i,i),interval[i,])}

}

Cofidence(50,100,0,1,0.05)

在置信度為95%的前提下，不包含參數(shù)真值的區(qū)間數(shù)不超過5次，從圖5可看出：該模擬有3次不包含真值，符合這個要求。顯然Nα只是不包含真值這個事件的頻數(shù)，受樣本的隨機波動影響，不具有嚴格的相等性約束。

圖5　置信度的一次模擬圖形

5　關(guān)于最短雙側(cè)置信區(qū)間的討論

ns=c(2,4,6,8,10,15,20,25,30)

alpha=0.05

par(mfrow=c(3,3))

for(n in ns)

{ mat=matrix(0,nrow=200,ncol=2)

standar=c(qchisq(alpha/2,n),qchisq(1-alpha/2,n))

mat[,1]=seq(standar[1]/20,standar[1]*1.2,length=200)

for(i in 1:200)

{ p=pchisq(mat[i,1],n)

mat[i,2]=qchisq(p+1-alpha,n)

}

plot(mat[,2]-mat[,1],type=′l′,xlab=paste(′n=′,n,sep=′′),ylab=′區(qū)間長度′)

abline(h=standar[2]-standar[1])

}

上述過程直觀地展示了標準置信區(qū)間的漸進收斂性，是對區(qū)間估計教學(xué)的直觀補充，更是對習(xí)慣性使用標準置信區(qū)間的一種直觀解釋。此外，還可以將尋找最短置信區(qū)間留作課后探索性實驗題，讓學(xué)生探討置信區(qū)間的常規(guī)求解方法。

從圖6還看到：當自由度n≥20時，標準置信區(qū)間長度已經(jīng)很接近最短區(qū)間長度。此處n≥20可作為區(qū)間估計實際應(yīng)用的一種參考；其次也說明此時χ2(n)的分布已經(jīng)很接近對稱型分布了。實際上，從χ2(n)的構(gòu)造定義知其滿足獨立同分布中心極限定理的條件，于是容易推出χ2(n)的漸進分布為正態(tài)分布N(n,2n)，而正態(tài)分布就是對稱型分布，這就不難解釋上述的結(jié)論。以下給出不同自由度下χ2(n)的密度函數(shù)比較(見圖7)，R程序如下：

ns=seq(10,40,by=10)

par(mfrow=c(1,4))

for(n in ns)

{

x=seq(0.01,n*2,by=0.1)

plot(x,dchisq(x,n),type=′l′,xlab=paste(′n=′,n,sep=″))

}

圖6　χ2(n)不同自由度下置信區(qū)間的比較

圖7　χ2(n)的密度函數(shù)

6　結(jié)語

在數(shù)理統(tǒng)計教學(xué)中引入R統(tǒng)計軟件是一種趨勢，利用R強大的隨機模擬和繪圖功能，通過設(shè)計合理的實驗和編寫程序，可將數(shù)理統(tǒng)計的概念、理論和方法演繹得更生動、更直觀。這種直觀的處理方式可增強學(xué)生對課程知識的理解，也能提升課程內(nèi)容的趣味性和實用性，更能激發(fā)學(xué)生的探索興趣和創(chuàng)造性。現(xiàn)在已有越來越多的學(xué)生主動關(guān)注利用R軟件和統(tǒng)計方法去解決實際的問題。

References)

[1] 王巖.工科專業(yè)數(shù)理統(tǒng)計實驗教學(xué)的實踐探索[J].數(shù)學(xué)教育學(xué)報,2007,16(3):95-98.

[2] 劉超,吳喜之.統(tǒng)計教學(xué)面臨的挑戰(zhàn)[J].統(tǒng)計研究,2012,29(2):105-108.

[3] 史書良.統(tǒng)計思想教育重于統(tǒng)計方法教育[J].中國統(tǒng)計,2008(2):56-57.

[4] 吉祖勤,蔡長安.NS2仿真技術(shù)在網(wǎng)絡(luò)實驗教學(xué)中的應(yīng)用[J].實驗技術(shù)與管理,2011，28(12):96-99.

[5] 關(guān)彥輝.R軟件在《概率統(tǒng)計》教學(xué)中的應(yīng)用[J].現(xiàn)代計算機:專業(yè)版,2009(12):87-90.

[6] 顧光同,張香云,徐光輝.統(tǒng)計實驗寓于概率統(tǒng)計教學(xué)的探索與實踐[J].統(tǒng)計與決策,2007(21):165-167.

[7] 陳懷俠,蔡火操,黃建林，等.設(shè)計性實驗教學(xué)的實踐與思考[J].實驗技術(shù)與管理2006,23(11):105-107.

[8] 郝小江,繆志農(nóng),黃昆.基于DSP的數(shù)字信號處理實驗設(shè)計[J].2012,29(2):44-47.

[9] Robert C P, Casella G. Monte Carlo Statistical Methods[M].New York:Springer,2004.

[10] Efron B, Rogosa D, Tibshiran R. Resampling Methods of Estimation[M]//International Encyclopedia of the Social & Behavioral Sciences.2001:13216-13220.

[11] 薛毅,陳立萍.統(tǒng)計建模與R軟件[M].北京：清華大學(xué)出版社,2007.

[12] 呂書龍,劉文麗,涂淑珍,等.概率論直觀教學(xué)與R統(tǒng)計軟件實現(xiàn)的探索[J].龍巖學(xué)院學(xué)報,2015,33(5):98-101.

Experimental design and R program realization of intuitive teaching of mathematical statistics

Lü Shulong, Liu Wenli , Liang Feibao, Ye Fuling

(College of Mathematics and Computer Science, Fuzhou University, Fuzhou 350116. China)

This paper discusses the doubts and problems of the Glivenko-Cantelli theorem, sampling theory, point estimation and interval estimation in mathematical statistics teaching, and puts forward the experimental design and R programs to graphically solve them directly. It not only promotes the students’ understanding of mathematical statistics and improves the teaching effects, but also enhances the application of R software in the teaching fields.

experimental design; R software; mathematical statistics; intuitive teaching

10.16791/j.cnki.sjg.2016.10.036

2016-03-29

福建省本科高校教育教學(xué)改革研究項目(JAS151395)；福州大學(xué)研究生優(yōu)質(zhì)課程建設(shè)項目(52004634,52004612)；福州大學(xué)高等教育教學(xué)改革工程(52001024,52001069)

呂書龍(1977—)，男，福建閩侯，碩士，副教授，主要研究方向為概率統(tǒng)計、統(tǒng)計計算和統(tǒng)計應(yīng)用.E-mail：wujispace@126.com

O212.1；TP311

A

1002-4956(2016)10-0142-05