R之組間差異的非參數檢驗

# t 檢驗
# 在研究中最常見的行為就是對兩個組進行比較。接受某種新藥治療的患者是否較使用某種現
# 有藥物的患者表現出了更大程度的改善？某種制造工藝是否較另外一種工藝制造出的不合格品
# 更少？兩種教學方法中哪一種更有效？如果你的結果變量是類別型的，那么可以直接使用7.3節
# 中闡述的方法。這里我們將關注結果變量為連續型的組間比較，并假設其呈正態分布。
# 為了闡明方法，我們將使用MASS包中的UScrime數據集。它包含了1960年美國47個州的刑
# 罰制度對犯罪率影響的信息。我們感興趣的結果變量為Prob（監禁的概率）、U1（14~24歲年齡
# 段城市男性失業率）和U2（35~39歲年齡段城市男性失業率）。類別型變量So（指示該州是否位
# 于南方的指示變量）將作為分組變量使用。數據的尺度已被原始作者縮放過
[python] view plain copy
print?

    <span style="font-family:SimHei;font-size:18px;">library(MASS)</span> 
獨立樣本的t檢驗
如果你在美國的南方犯罪，是否更有可能被判監禁？我們比較的對象是南方和非南方各州，
因變量為監禁的概率。一個針對兩組的獨立樣本t檢驗可以用于檢驗兩個總體的均值相等的假設
這里假設兩組數據是獨立的，并且是從正態總體中抽得。檢驗的調用格式為:

[python] view plain copy
print?

    <span style="font-family:SimHei;font-size:18px;">    
    t.test(y~X,data)  
      
    其中的y是一個數值型變量，x是一個二分變量  
      
    t.test(y1,y2)  
      
    其中的y1和y2為數值型向量（即各組的結果變量）?？蛇x參數data的取值為一個包含了這些  
    變量的矩陣或數據框,里的t檢驗默認假定方差不相等，并使  
    用Welsh的修正自由度。你可以添加一個參數var.equal=TRUE以假定方差相等，并使用合并方  
    差估計。默認的備擇假設是雙側的（即均值不相等，但大小的方向不確定）。你可以添加一個參  
    數alternative="less"或alternative="greater"來進行有方向的檢驗。</span> 
我們使用了一個假設方差不等的雙側檢驗，比較了南方（group 1）和非南
方（group 0）各州的監禁概率
[python] view plain copy
print?

    <span style="font-family:SimHei;font-size:18px;">t.test(Prob~So,data = UScrime)</span>  

[python] view plain copy
print?

    <span style="font-family:SimHei;font-size:18px;">#   
    # > t.test(Prob~So,data = UScrime)  
    #   
    # Welch Two Sample t-test  
    #   
    # data:  Prob by So  
    # t = -3.8954, df = 24.925, p-value = 0.0006506  
    # alternative hypothesis: true difference in means is not equal to 0  
    # 95 percent confidence interval:  
    #   -0.03852569 -0.01187439  
    # sample estimates:  
    #   mean in group 0 mean in group 1   
    # 0.03851265      0.06371269 </span>  

你可以拒絕南方各州和非南方各州擁有相同監禁概率的假設（p < .001）。
非獨立樣本的t檢驗
再舉個例子，你可能會問：較年輕（14~24歲）男性的失業率是否比年長（35~39歲）男性的
失業率更高？在這種情況下，這兩組數據并不獨立。你不能說亞拉巴馬州的年輕男性和年長男性
的失業率之間沒有關系。在兩組的觀測之間相關時，你獲得的是一個非獨立組設計（dependent
groups design）。前—后測設計（pre-post design）或重復測量設計（repeated measures design）同樣
也會產生非獨立的組。
非獨立樣本的t檢驗假定組間的差異呈正態分布,對于本例，檢驗的調用的格式為:
[python] view plain copy
print?

    <span style="font-family:SimHei;font-size:18px;">t.text(y1,y2,pairred=TRUE)  
      
    其中的y1和y2為兩個非獨立組的數值向量  
      
    sapply(UScrime[c("U1","U2")],function(x){c(mean=mean(x),sd=sd(x))})  
      
    with(UScrime,t.test(U1,U2,paired = TRUE))  
      
    #   
    # > sapply(UScrime[c("U1","U2")],function(x){c(mean=mean(x),sd=sd(x))})  
    # U1       U2  
    # mean 95.46809 33.97872  
    # sd   18.02878  8.44545  
    # > with(UScrime,t.test(U1,U2,paired = TRUE))  
    #   
    # Paired t-test  
    #   
    # data:  U1 and U2  
    # t = 32.407, df = 46, p-value < 2.2e-16  
    # alternative hypothesis: true difference in means is not equal to 0  
    # 95 percent confidence interval:  
    #   57.67003 65.30870  
    # sample estimates:  
    #   mean of the differences   
    # 61.48936   
    </span>  

差異的均值（61.5）足夠大，可以保證拒絕年長和年輕男性的平均失業率相同的假設。
年輕男性的失業率更高。事實上，若總體均值相等，獲取一個差異如此大的樣本的概率小于
0.000 000 000 000 000 22（即2.2e?16）

組間差異的非參數檢驗
如果數據無法滿足t檢驗或ANOVA的參數假設，可以轉而使用非參數方法
兩組的比較
若兩組數據獨立，可以使用Wilcoxon秩和檢驗（更廣為人知的名字是Mann–Whitney U檢驗）
來評估觀測是否是從相同的概率分布中抽得的（即，在一個總體中獲得更高得分的概率是否比另
一個總體要大）。調用格式為：wilcox.text(y~x,text)
其中的y是數值型變量，而x是一個二分變量:wilcox.test(y1,y2)
其中的y1和y2為各組的結果變量?？蛇x參數data的取值為一個包含了這些變量的矩陣或數據框。默
認進行一個雙側檢驗。你可以添加參數exact來進行精確檢驗，指定alternative="less"或
alternative="greater"進行有方向的檢驗。
如果你使用Mann–Whitney U檢驗回答上一節中關于監禁率的問題，將得到這些結果：
[python] view plain copy
print?

    <span style="font-family:SimHei;font-size:18px;">with(UScrime,by(Prob,So,median))  
    # > with(UScrime,by(Prob,So,median))  
    # So: 0  
    # [1] 0.038201  
    # ------------------------------------------------------------------   
    #   So: 1  
    # [1] 0.055552  
    # >   
      
    wilcox.test(Prob~So,data = UScrime)  
    #   
    # 你可以再次拒絕南方各州和非南方各州監禁率相同的假設（p < 0.001）  
    > wilcox.test(Prob~So,data = UScrime)  
    #   
    # Wilcoxon rank sum test  
    #   
    # data:  Prob by So  
    # W = 81, p-value = 8.488e-05  
    # alternative hypothesis: true location shift is not equal to 0  
    # </span>  


Wilcoxon符號秩檢驗是非獨立樣本t檢驗的一種非參數替代方法
它適用于兩組成對數據和
無法保證正態性假設的情境。調用格式與Mann–Whitney U檢驗完全相同，不過還可以添加參數
paired=TRUE。讓我們用它解答上一節中的失業率問題

[python] view plain copy
print?

    <span style="font-family:SimHei;font-size:18px;">sapply(UScrime[c("U1","U2")],median)  
    #> sapply(UScrime[c("U1","U2")],median)  
    # U1 U2   
    # 92 34   
      
    with(UScrime,wilcox.test(U1,U2,paired = TRUE))  
    # > with(UScrime,wilcox.test(U1,U2,paired = TRUE))  
    #   
    # Wilcoxon signed rank test with continuity correction  
    #   
    # data:  U1 and U2  
    # V = 1128, p-value = 2.464e-09  
    # alternative hypothesis: true location shift is not equal to 0  
    </span>  

你再次得到了與配對t檢驗相同的結論.
多于兩組的比較
在要比較的組數多于兩個時，必須轉而尋求其他方法?？紤]7.4節中的state.x77數據集。
它包含了美國各州的人口、收入、文盲率、預期壽命、謀殺率和高中畢業率數據。如果你想比較
美國四個地區（東北部、南部、中北部和西部）的文盲率，應該怎么做呢？這稱為單向設計（one-way
design），我們可以使用參數或非參數的方法來解決這個問題
如果無法滿足ANOVA設計的假設，那么可以使用非參數方法來評估組間的差異
如果各組獨立，則Kruskal—Wallis檢驗將是一種實用的方法
Kruskal–Wallis檢驗的調用格式為:kruskal.test(y~A,data)
其中的y是一個數值型結果變量，A是一個擁有兩個或更多水平的分組變量（grouping variable）。
（若有兩個水平，則它與Mann–Whitney U檢驗等價。
如果各組不獨立（如重復測量設計或隨機區組設計），那么Friedman檢驗會更合適
friedman.test(y~A|B,data)其中的y是數值型結果變量，A是一個分組變量，而B是一個用以認定匹配觀測的區組變量（blocking
variable）
讓我們利用Kruskal–Wallis檢驗回答文盲率的問題。首先，你必須將地區的名稱添加到數據
集中。這些信息包含在隨R基礎安裝分發的state.region數據集中:

[python] view plain copy
print?

    <span style="font-family:SimHei;font-size:18px;">states<-as.data.frame(cbind(state.region,state.x77))  
      
    kruskal.test(Illiteracy~state.region,data = states)  
      
    # > kruskal.test(Illiteracy~state.region,data = states)  
    #   
    # Kruskal-Wallis rank sum test  
    #   
    # data:  Illiteracy by state.region  
    # Kruskal-Wallis chi-squared = 22.672, df = 3, p-value = 4.726e-05</span>  


顯著性檢驗的結果意味著美國四個地區的文盲率各不相同（p <0.001）
雖然你可以拒絕不存在差異的原假設，但這個檢驗并沒有告訴你哪些地區顯著地與其他地區
不同。要回答這個問題，你可以使用Mann–Whitney U檢驗每次比較兩組數據。一種更為優雅的
方法是在控制犯第一類錯誤的概率（發現一個事實上并不存在的差異的概率）的前提下，執行可
以同步進行的多組比較，這樣可以直接完成所有組之間的成對比較。npmc包提供了所需要的非
參數多組比較程序

[python] view plain copy
print?

    <span style="font-family:SimHei;font-size:18px;">install.packages("npmc")  
    var<-state.x77[,c("Illiteracy")]  
    mydata<-as.data.frame(cbind(class,var))  
    rm(class,var)  
    library(npmc)  
    summary(npmc(mydata),type="BF")  
    </span>  

注意：npmc 包已經被棄用了?。?！
但是下面的截圖是以前沒有棄用時的圖