1. 程式人生 > >R語言入門:使用函式sample進行抽樣

R語言入門:使用函式sample進行抽樣

在醫學統計學或者流行病學裡的現場調查、樣本選擇經常會提到一個詞:隨機抽樣。隨機抽樣是為了保證各比較組之間均衡性的一個很重要的方法。那麼今天介紹的第一個函式就是用於抽樣的函式sample

> x=1:10

> sample(x=x)

 [1]  3  5  9  6 10  7  2  1  8  4

第一行程式碼表示給x向量賦值1~10,第二行程式碼表示對x向量進行隨機抽樣。結果輸出為每次抽樣抽得的結果,可以看出該抽樣為無放回抽樣------最多抽n次,n為x向量中元素的個數。

如果想指定在該向量中抽取元素的個數,需要加一個引數size

> x=1:1000

> sample(x=x,size=20)

 [1]  66 891 606 924 871 374 879 573 284 305 914 792 398 497 721 897 324 437

[19] 901  33

這是在1~1000的正整數中抽樣,其中size指定抽樣的次數,抽了20次,結果如上所示。

這些都是無放回抽樣。所謂無放回抽樣,也就是說某個元素一旦被選擇,該總體中就不會再有該元素。如果是有放回抽樣,則需新增一個引數repalce=T

> x=1:10

> sample(x=x,size=5,replace=T)

[1] 4 7 2 4 8

“replace”就是重複的意思。即可以重複對元素進行抽樣,也就是所謂的有放回抽樣。我們看上面的結果,元素4在5次隨機抽樣的過程中被抽取了兩次。

R語言程式碼有一個特性就是“對位性”,也許我的詞不專業,但是它的意思就是:如果我們輸入程式碼的位置與某個函式中引數的位置一一對應的話,我們可以不寫該函式的引數,如:

> x=1:10

> sample(x,20,T)

 [1] 1 2 2 1 5 5 5 9 9 5 2 9 8 3 4 8 8 8 1 1

在上述程式碼中我們省略了引數xsizerepalce,但是仍然可以運算並且表示對x向量有放回隨機抽取20次。我們之所以儘量在每次編寫程式碼時帶上引數是因為我覺得這個習慣比較好,而且看起來清楚明白。另外,省略引數的前提是你非常熟悉某個函式引數的位置,否則一旦沒有“對位”,那麼結果肯定是錯誤的。而且很多函式有較多的引數,想記住它們的位置是困難的。而如果帶上引數,那麼即使位置不對應,也一樣可以運算:

> x=1:10

> sample(size=20,replace=T,x=x)

 [1]  4  9  2  6  4  5  4  7 10  5  2  2  3  4  2  4  6  8  7  8

這種優點顯而易見,不僅清楚,而且無需對應。另外我們也可以看出,有放回抽樣的話size可以無窮大,而無放回抽樣size的大小就取決於總體的容量了。

對於擲骰子,投硬幣(這可能是介紹抽樣必介紹的內容),都屬於有放回抽樣。

這裡要說明,對於sample函式,引數x可以是數值,也可以是字元,實際上引數x代表任意一個向量:

> a=c("A","B")

> sample(x=a,size=10,replace=T)

 [1] "B" "A" "A" "A" "B" "A" "A" "B" "A" "A"

上述程式碼可以理解為擲硬幣,拋了10次,其中正面(A)與反面(B)出現的次數。

上述抽樣過程,每一個元素被抽取的概率相等,稱為隨機抽樣。

有時候我們的抽取元素的概率未必相等(如常見的二項分佈概率問題),此時我們需要新增一個引數prob,也就是“probability”(概率)的縮寫。假設一名醫生給患者做某手術成功的概率是80%,那麼現在他給20例病人做手術,可能有哪幾次是成功的呢?程式碼如下:

> x=c("S","F")

> sample(x,size=20,replace=T,prob=c(0.8,0.2))

 [1] "F" "S" "S" "S" "S" "S" "S" "S" "S" "S" "S" "S" "F" "S" "S" "F" "S" "S"

[19] "F" "S"

其中“S”代表成功,“F”代表失敗。

> x=c(1,3,5,7)

> sample(x,size=20,replace=T,prob=c(0.1,0.2,0.3,0.9))

 [1] 3 5 7 3 7 3 7 5 3 7 7 7 1 5 7 5 7 7 3 7

這些程式碼告訴我們,對每一個元素都可以給定一個概率,且每個概率是獨立的,即在引數prob中,不一定所有元素的概率加起來等於1,它只代表某元素被抽取的概率而已。

對於sample函式,其引數x可以是R中任意一個物件(如上述對字元的抽樣)。其中還有一個功能相同的函式sample.int,“int”即“intger”的縮寫,也就是“整數”。它的引數n必須是正整數:

> x=-10.5:7.5

> sample(x=x,size=3);sample.int(n=x,size=3)

[1] -5.5 -7.5  0.5

Error in sample.int(x, size = 3) : invalid first argument

第一行程式碼生成了-10.5到7.5的等差數列,結果輸出的第一行是sample的結果;第二行是sample.int的結果,提示錯誤:“第一個自變數無效”,因為它不是正整數。其餘的用法與sample是一樣的。

接下來介紹第二個函式max/min,以及pmax/pmin

> x=sample(1:1000,size=20)

> x

 [1] 514 164 107 354 719 317 204 372 714 522 900 635 441  7  230 653 420 804

[19] 817  36

> max(x);min(x)

[1] 900

[1] 7

從1:1000中隨機抽取了20個數,並得出它們的最大值和最小值。如果有缺失值,那麼需要一個引數na.rm=T,即移除缺失值,這在前面的mean函式中已經介紹過了:

> x=c(x,NA)

> mean(x);mean(x,na.rm=T);max(x);max(x,na.rm=T);min(x);min(x,na.rm=T)

[1] NA

[1] 446.5

[1] NA

[1] 900

[1] NA

[1] 7

第一行程式碼我們給原來的向量x添加了一個缺失值“NA”,第二行程式碼及結果分別顯示了在計算x的均值、最大值、最小值時是否帶引數na.rm=T的情況。

> x=sample(1:1000,size=20)

> y=sample(1:1000,size=20)

> x;y

 [1] 596 333 401 746 970 121 665 831 853 741 855  99 331 275 286 370 558 373

[19]  55 786

 [1] 291 360 964 173 497 553 155  34 488  29 661 736 591 602 548 527 450 416

[19] 596 965

這行程式碼我們對1:1000分別進行20次無放回隨機性抽樣得到向量x和y。

> pmax(x,700)

 [1] 700 700 700 746 970 700 700 831 853 741 855 700 700 700 700 700 700 700

[19] 700 786

> pmin(y,700)

 [1] 291 360 700 173 497 553 155  34 488  29 661 700 591 602 548 527 450 416

[19] 596 700

第一個程式碼表示用700和x向量中的每一個元素進行比較,將每次比較結果較大者輸出。

第二個程式碼表示用700和y向量中的每一個元素進行比較,將每次比較結果較小者輸出。

這就是pmax/pmin的用法,p在這裡表示“parallel”,可以理解為對向量做平行比較。

> pmax(x,y)

 [1] 596 360 964 746 970 553 665 831 853 741 855 736 591 602 548 527 558 416

[19] 596 965

該結果顯示了將x向量與y向量中的每個元素平行比較後得出的較大值的結果。如第一個元素的比較,x為596,y為291,所以結果輸出為596;第二個元素比較,x向量為333,y向量為360,因此輸出結果為360。後面的結果也是一樣。

因為x向量與y向量是等長的,所以我們可以對x與y每一個元素進行平行比較。因為R中向量運算時具有“迴圈原則”,所以在上面的例子中我們也可以用x或y向量與一個元素“700”比較(相當於包含一個元素“700”的向量重複了20次)。

如果迴圈完整,那麼較短向量中的元素依次與較長向量中的元素進行比較

> x=c(1,3,5,7)

> y=c(4,6)

> pmax(x,y)

[1] 4 6 5 7

結果解釋:第一個是1與4比較的最大值,第二個是3與6比較的最大值,第三個是5與4比較的最大值,最後一個是7與6比較的最大值。如果迴圈不完整,自然會出現警告,但輸出結果依然是各元素依次平行比較的結果:

> x=c(1,3,5,7,9)

> y=c(4,6)

> pmax(x,y)

[1] 4 6 5 7 9

Warning message:

In pmax(x, y) : an argument will be fractionally recycled

警告的內容是:自變數將部分迴圈。 

對於函式pmax/pmin,也有引數na.rm=T,用法是一樣的。同時有對應的pmax.int/pmn.int,也就是比較向量內的元素必須是整數。

接下來給大家介紹今天的最後一個函式length,即檢視某一物件的長度。

> x=c(1,3,5,7,9,11)

> y=c("a","b","c")

> z=matrix(1:12,nrow=3)

> a=cbind(x,y)

> length(x);length(y);length(z);length(a)

[1] 6

[1] 3

[1] 12

[1] 12

結果告訴我們向量x的長度是6,向量y的長度是3,矩陣z的長度是12(3×4)。以列結合的物件a的長度為12(6×2)。

對於資料框,其長度為列的長度:

> x=1:5

> y=6:10

> z=7:11

> d=data.frame(x,y,z)

> length(d)

[1] 3

函式length也可以“反向使用”,即對於已知向量,先給定其長度然後輸出結果:

> x=1:20

> length(x)=12

> x

 [1]  1  2  3  4  5  6  7  8  9 10 11 12

X向量的真實長度為20,但我們輸出長度為12的x向量。