劉二鋼+馬建強(qiáng)

摘要：R語言在統(tǒng)計(jì)分析中的應(yīng)用已經(jīng)越來越廣泛。文章通過一些實(shí)例，簡要地介紹了R語言在統(tǒng)計(jì)分析中的一些技巧和應(yīng)用，希望能夠?yàn)閺V大讀者學(xué)習(xí)R語言提供一些啟迪和借鑒作用。

關(guān)鍵詞：R語言；統(tǒng)計(jì)分析；統(tǒng)計(jì)函數(shù)；檢驗(yàn)

中圖分類號：TP312 文獻(xiàn)標(biāo)識碼：A 文章編號：1009-3044（2017）01-0251-03

Abstract： The R language in the application of statistical analysis has been more and more widely.This paper introduces the R language techniques and applications in statistical analysis through some examples， and wants to be able to provide some enlightenment and reference for the readers to study the R language.

Key words： R language； statistical analysis； statistical function； test

R語言主要誕生于上個世紀(jì)九十年代，最初是S語言的一種實(shí)現(xiàn)。R是由來自于新西蘭奧克蘭大學(xué)的Ross Ihaka和Robert Gentleman開發(fā)（也因此稱為R），現(xiàn)在由“R開發(fā)核心團(tuán)隊(duì)”負(fù)責(zé)開發(fā)[1]。它的主要特點(diǎn)是具備一系列連貫而又完整的數(shù)據(jù)分析中間工具，擁有一整套數(shù)組和矩陣操作運(yùn)算符能有效地處理保存數(shù)據(jù)，其圖形統(tǒng)計(jì)功能可以對數(shù)據(jù)直接進(jìn)行分析和顯示，加上R是一種面向?qū)ο蟮目删幊陶Z言，和其他編程語言、數(shù)據(jù)庫之間有很好的接口[2]。

R語言能夠提供強(qiáng)大的統(tǒng)計(jì)分析、數(shù)據(jù)處理以及圖形可視化功能。與其說R是一種軟件，不如說R是一種統(tǒng)計(jì)應(yīng)用的環(huán)境。另外R語言還是一個完全免費(fèi)的統(tǒng)計(jì)軟件，許多程序員提供了大量豐富的程序包，并且擁有強(qiáng)大的社區(qū)資源，有許多不同領(lǐng)域的專家、學(xué)者在此處交流與探討。正因如此，R語言在國外的統(tǒng)計(jì)教學(xué)及科研工作中已經(jīng)得到了廣泛應(yīng)用，但是在國內(nèi)，R語言還比較陌生，研究人員相對較少。不過可喜的是，已經(jīng)有不少學(xué)者已經(jīng)在關(guān)注R語言，并為推廣和使用R語言而努力。

本文就是在此環(huán)境之下，介紹R語言的一些簡單統(tǒng)計(jì)功能，希望能夠?yàn)閺V大學(xué)者和教學(xué)人員提供一些借鑒之處。

1 R語言中的對象與數(shù)據(jù)類型

R語言把操作對象分為多種類型，如向量、矩陣、因子（factor）、清單（list）及數(shù)據(jù)框（data frame）等。部分R函數(shù)會根據(jù)類型的不同而作不同的處理。比如函數(shù)plot會隨著不同的數(shù)據(jù)類型而得到不同的結(jié)果：輸入數(shù)字向量，會得到點(diǎn)散圖；輸入因子，會得到盒形圖；輸入回歸模型分析則得到圖像分析。這種處理方式在R中被稱為“對象導(dǎo)向式程序編寫”[3]。

1.1 標(biāo)量與向量

標(biāo)量是由簡單的數(shù)字或者字符串所組成的，如果是字符串，則在左右兩端需加上雙引號。比如3代表數(shù)字，而“3”則代表一個符號。向量則是一組一維標(biāo)量的集合。如果一個向量同時(shí)包含了數(shù)字和字符，則R語言會統(tǒng)一按照文字處理。例如：

> a=c（1，2，3，4，"f"）

1.2 矩陣

矩陣就是二維的向量，輸入矩陣和輸入向量類似。例如：

> matrix（c（1：9），nrow=3，ncol=3，byrow="T"）

[，1] [，2] [，3]

[1，] 1 2 3

[2，] 4 5 6

[3，] 7 8 9

這里輸入的3行3列的一個矩陣，byrow="T"表示是按照行進(jìn)行排列矩陣。如果是一行的矩陣，就可以按照向量進(jìn)行處理了。

1.3 清單

清單可以將多種類型的數(shù)據(jù)放在一起。如果將向量看作是一組標(biāo)量的排列，則清單可以理解為一組不同形態(tài)對象的排列。例如：

> b=list（c（"one"，2，3，"four"）， matrix（c（1，2，3，4，5，6，7，8，9）， nrow=3， ncol=3）， c（3，9，15））

想輸出第一組內(nèi)容，可以使用：

> b[[1]]

[1] "one" "2" "3" "four"

如果輸入：> b[[2]][1：3，2]，則顯示的是矩陣中第二列的內(nèi)容。

2 R的統(tǒng)計(jì)功能

R語言中內(nèi)置了許多統(tǒng)計(jì)函數(shù)，用這些函數(shù)可以方便地實(shí)現(xiàn)大量的統(tǒng)計(jì)功能，方便地解決在統(tǒng)計(jì)中概率計(jì)算、臨界值、分位數(shù)等問題。表1是R語言中常見的統(tǒng)計(jì)函數(shù)名稱及相關(guān)參數(shù)。

為了便于方便理解，R將統(tǒng)計(jì)分布設(shè)為四個基本項(xiàng)目，即概率密度函數(shù)、累積分布函數(shù)、分位數(shù)和偽隨機(jī)數(shù)[3]。每個統(tǒng)計(jì)函數(shù)所對應(yīng)的項(xiàng)目函數(shù)，其命名規(guī)則是用d、p、q、r做為名稱前綴。比如幾何分布函數(shù)，R語言對應(yīng)名稱為geom，則可以用dgeom、pgeom、qgeom、rgeom去表示幾何分布的這四個基本項(xiàng)目。

下面舉一個例子表明R語言對于統(tǒng)計(jì)函數(shù)的處理方法，簡單的例句及釋義如下：

> pnorm（5） #計(jì)算標(biāo)準(zhǔn)正太隨機(jī)變量Z<5的方差

> pnorm（5，3，1.2） #計(jì)算正太隨機(jī)變量平均數(shù)是3，標(biāo)準(zhǔn)方差是2，分布函數(shù)<5的概率

> pnorm（5，3，1.2）-pnorm（4，3，1.2） #求分布函數(shù)4～5之間的概率

> qnorm（0.85，3，1.2） #求85%的分位數(shù)

> dnorm（4，3，1.2） #求值為4時(shí)的概率密度

> rnorm（20，3，1.2） #產(chǎn)生20個隨機(jī)數(shù)

如上可以看出，對于相關(guān)的函數(shù)，只需要在前面加上相應(yīng)的前綴，就可以進(jìn)行計(jì)算了，非常方便。表1中列舉的其它函數(shù)也可以以此類推。

3 R與概率曲線

3.1 正態(tài)分布

正態(tài)分布是統(tǒng)計(jì)學(xué)中最重要、最常見的分布之一，現(xiàn)實(shí)生活中的很多隨機(jī)現(xiàn)象都可以用正態(tài)分布或類似正態(tài)分布進(jìn)行描述。正態(tài)分布對于統(tǒng)計(jì)學(xué)中的三大抽樣分布具有引導(dǎo)作用。對于正態(tài)分布，可以用下面語句畫出整體曲線：

> curve（dnorm（x）， xlim=c（-3，3），col="red"， lwd=3， lty=1，main="標(biāo)準(zhǔn)正態(tài)分布"）

其中，curve是繪畫曲線的常用命令，是根據(jù)（x，dnorm（x））的坐標(biāo)關(guān)系作出的。橫軸的取值范圍為-3，3，曲線顏色為紅色，線條寬度3，main的內(nèi)容為圖形標(biāo)題。如果想改為累積分布函數(shù)，則按照上面第2點(diǎn)所敘述，只需要將dnorm函數(shù)改為pnorm即可，語句如下：

> curve（pnorm（x）， xlim=c（-3，3），col="red"， lwd=3， lty=1，main="標(biāo)準(zhǔn)正態(tài)累積分布"）

做出的圖形分別如圖1和圖2所示。

這里再舉一個例子，對于正態(tài)分布進(jìn)行簡單模擬。

假設(shè)正態(tài)分布樣本數(shù)目為5000，進(jìn)行模擬后畫出正態(tài)分布曲線和樣本直方圖，所用語句如下：

> n=5000 #設(shè)定樣本數(shù)目

> x=rnorm（n） #產(chǎn)生偽隨機(jī)樣本

> curve（dnorm（x）， xlim=c（-3，3），col="red"， lwd=3， lty=1） #畫出標(biāo)準(zhǔn)正態(tài)分布曲線

> hist（x，probability=T，add=T） #畫出直方圖

> lines（density（x），col="red"， lwd=3， lty=2） #在圖上添加x的模擬密度曲線

> title（main="正態(tài)分布曲線與直方圖"） #添加圖表標(biāo)題

所做圖形如圖3所示。由圖可以看出，所產(chǎn)生的模擬曲線與標(biāo)準(zhǔn)正態(tài)分布曲線是比較吻合的，如果樣本空間越大，則吻合度越高。

3.2 二項(xiàng)分布

二項(xiàng)分布與貝努力實(shí)驗(yàn)有著極大的聯(lián)系。一個實(shí)驗(yàn)中有兩個實(shí)驗(yàn)：成功（記為1）與失?。ㄓ洖?），出現(xiàn)的概率分別為p和1-p，則一次實(shí)驗(yàn) （稱為貝努力實(shí)驗(yàn)） 成功的次數(shù)服從一個參數(shù)為p的貝努力分布。如果將實(shí)驗(yàn)次數(shù)重復(fù)n次，則實(shí)驗(yàn)成功的次數(shù)服從參數(shù)為（n，p）的二項(xiàng)分布。所以當(dāng)二項(xiàng)分布實(shí)驗(yàn)為1時(shí)，就是貝努力分布。

現(xiàn)在假設(shè)二項(xiàng)分布參數(shù)為f（n，p） n=30，p=0.5的計(jì)算機(jī)模擬。由貝努力實(shí)驗(yàn)與二項(xiàng)分布的關(guān)系可以將此過程設(shè)計(jì)為：

> n=30 #設(shè)置f（n，p）中n的值為30

> x=c（1：5000） #將實(shí)驗(yàn)過程模擬為5000次，并將x保存為一個長度為5000的向量記錄實(shí)驗(yàn)結(jié)果

> for（i in 1：5000）{x[i]=sum（sample（c（0，1），n，replace=T））} #設(shè)置循環(huán)過程，在等概率抽樣條件下，將每次試驗(yàn)中1出現(xiàn)的次數(shù)賦值給相應(yīng)的X向量。

> hist（x，col="light blue"， xlim=c（min（x），max（x））， probability=T， main="二項(xiàng)分布模擬圖"） #畫出X的頻率直方圖，并添加圖表標(biāo)題

> lines（density（x，bw=1），col="red"，lwd=3） #添加模擬分布密度曲線

畫出的圖形如圖4所示。由此可以看出頻率直方圖及模擬的密度曲線與標(biāo)準(zhǔn)的正態(tài)分布曲線還是比較相似的[4]。如果根據(jù)中心極限定理，當(dāng)隨機(jī)變量足夠大時(shí)，正態(tài)分布與二項(xiàng)分布是非常接近的。

4 R語言與統(tǒng)計(jì)檢驗(yàn)

在數(shù)理統(tǒng)計(jì)分析中，只能由估計(jì)量估計(jì)總體的參數(shù)，總體參數(shù)始終是不可知的，只能通過統(tǒng)計(jì)檢驗(yàn)，由統(tǒng)計(jì)量推斷總體參數(shù)[5]。一般在統(tǒng)計(jì)推斷過程中，對參數(shù)先提出假設(shè)，然后根據(jù)假設(shè)進(jìn)行假設(shè)檢驗(yàn)。在R語言中包含了多種參數(shù)假設(shè)辦法。下面舉一個檢驗(yàn)的例子。

某種元件的壽命X（以小時(shí)計(jì)算）服從正態(tài)分布Ν（μ，σ2），現(xiàn)測得16只元件的壽命如下：

149 260 485 210 239 280 191 212

224 379 179 264 222 362 198 250

假設(shè)顯著性水平a=0.05，問是否可以認(rèn)定元件的平均壽命不大于225小時(shí)？

分析：這是一個關(guān)于均值的檢驗(yàn)問題，可以采用單邊檢測方法，需要檢驗(yàn)

由于總體方差未知，故采用t-檢驗(yàn)方法。檢驗(yàn)統(tǒng)計(jì)量公式為：

在a=0.05前提條件下，樣本數(shù)量n=16，臨界值t1-a（n-1）=1.7531，則拒絕域{t> 1.7531}，根據(jù)公式計(jì)算X平均值=256.5，樣本方差s=86.02，所以t=（256.5-225）/（86.02/）=1.4648，因此樣本離差t=1.4648小于臨界值1.7531。因此接受原假設(shè)，可以認(rèn)為元件的平均壽命不大于225小時(shí)。

而在R語言中，計(jì)算上述問題只需要下面兩個步驟就可以了。

> X=c（149，260，485，210，239，280，191，212，224，379，179，264，222，362，198，250）

> t.test（X，mu=225，alternative='greater'）

得出結(jié)果為：

One Sample t-test

data： X

t = 1.4648， df = 15， p-value = 0.0818

alternative hypothesis： true mean is greater than 225

95 percent confidence interval：

218.8023 Inf

sample estimates：

mean of x

256.5

在這里沒有給出結(jié)論，但是可以得出t值和p值，從p而判斷是否拒絕原假設(shè)。如果從p值判斷，得出p值為0.0818>0.05，不能拒絕原假設(shè)。則接受H0，因此在a=0.05情況下，可以認(rèn)為平均壽命不大于225小時(shí)。t值的結(jié)果與判斷方式和前面所述一樣。另外參數(shù)alternative可以為“l(fā)ess”，“greater”和“two side”，可以反映出假設(shè)是單邊假設(shè)還是雙邊假設(shè)，默認(rèn)是雙邊假設(shè)。

5 結(jié)束語

R語言語法簡單，容易編寫。通過文本的敘述，相信讀者對R語言已經(jīng)有了一定的了解。本文所敘述的只是R語言的一小部分內(nèi)容，讀者可以通過查看相關(guān)參考資料對R語言做進(jìn)一步的掌握。另外由于R語言是開源軟件，自由免費(fèi)，非常適合做統(tǒng)計(jì)分析和教學(xué)。在目前注重知識產(chǎn)權(quán)的大背景下，掌握和利用好R語言對于統(tǒng)計(jì)分析具有重大的現(xiàn)實(shí)作用。

參考文獻(xiàn)：

[1] smalllleopard. R語言為Hadoop集群數(shù)據(jù)統(tǒng)計(jì)分析帶來革命性變化[EB/OL]. （2012-11-10）. http：//blog.sina.com.cn/s/blog_6cfc336b01018wvt.html.

[2] 陳希. 基于R語言數(shù)據(jù)挖掘的社交網(wǎng)絡(luò)客戶細(xì)分研究[D]. 北京： 北京郵電大學(xué)， 2011： 34-40

[3] 陳毅恒， 梁沛霖. R軟件操作入門[M]. 北京： 中國統(tǒng)計(jì)出版社， 2006.

[4] 葉文春. 淺談R語言在統(tǒng)計(jì)學(xué)中的應(yīng)用[J]. 中共貴州省委黨校學(xué)報(bào)， 2008（4）： 123-125.

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