R語言多重比較示例：Bonferroni校訂法和Benjamini & Hochberg法

時間 2021-04-11

標籤 python spa code 排序事件 rem get table ast 欄目 Python 简体版

原文原文鏈接

原文連接：http://tecdat.cn/?p=21825

=python

假設檢驗的基本原理是小几率原理，即咱們認爲小几率事件在一次試驗中實際上不可能發生。spa

=code

多重比較的問題

當同一研究問題下進行屢次假設檢驗時，再也不符合小几率原理所說的「一次試驗」。若是在該研究問題下只要有檢驗是陽性的，就對該問題下陽性結論的話，對該問題的檢驗的犯一類錯誤的機率就會增大。若是同一問題下進行n次檢驗，每次的檢驗水準爲α（每次假陽性機率爲α），則n次檢驗至少出現一次假陽性的機率會比α大。假設每次檢驗獨立的條件下該機率可增長至排序

常見的多重比較情景包括：事件

多組間比較
多個主要指標
臨牀試驗中期中分析
亞組分析

=rem

控制多重比較謬誤（Familywise error rate)：Bonferroni矯正

Bonferroni法獲得的矯正P值=P×n
Bonferroni法很是簡單，它的缺點在於很是保守（大概是各類方法中最保守的了），尤爲當n很大時，通過Bonferroni法矯正後總的一類錯誤可能會遠遠小於既定α。get

控制錯誤發現率：Benjamini & Hochberg法

簡稱BH法。首先將各P值從小到大排序，生成順序數
排第k的矯正P值=P×n/k
另外要保證矯正後的各檢驗的P值大小順序不發生變化。it

怎麼作檢驗

R內置了一些方法來調整一系列p值，以控制多重比較謬誤（Familywise error rate)或控制錯誤發現率。table

Holm、Hochberg、Hommel和Bonferroni方法控制了多重比較謬誤（Familywise error rate)。這些方法試圖限制錯誤發現的機率（I型錯誤，在沒有實際效果時錯誤地拒絕無效假設），所以都是相對較保守的。ast

方法BH（Benjamini-Hochberg，與R中的FDR相同）和BY（Benjamini & Yekutieli）控制錯誤發現率，這些方法試圖控制錯誤發現的指望比例。

請注意，這些方法只須要調整p值和要比較的p值的數量。這與Tukey或Dunnett等方法不一樣，Tukey和Dunnett也須要基礎數據的變異性。Tukey和Dunnett被認爲是多重比較謬誤（Familywise error rate)方法。

要了解這些不一樣調整的保守程度，請參閱本文下面的兩個圖。

關於使用哪一種p值調整度量沒有明確的建議。通常來講，你應該選擇一種你的研究領域熟悉的方法。此外，可能有一些邏輯容許你選擇如何平衡犯I型錯誤和犯II型錯誤的機率。例如，在一項初步研究中，你可能但願保留儘量多的顯著值，來避免在將來的研究中排除潛在的顯著因素。另外一方面，在危及生命而且治療費用昂貴的醫學研究中，得出一種治療方法優於另外一種治療方法的結論以前，你應該有很高的把握。

具備25個p值的多重比較示例

====================

### --------------------------------------------------------------
### 多重比較示例
### --------------------------------------------------------------

Data = read.table(Input,header=TRUE)

按p值排序數據

Data = Data[order(Data$Raw.p),]

檢查數據是否按預期的方式排序

執行p值調整並添加到數據框

Data$Bonferroni =
      p.adjust(Data$Raw.p,
               method = "bonferroni")

Data$BH =
      p.adjust(Data$Raw.p,
               method = "BH")

Data$Holm =
      p.adjust(Data$ Raw.p,
               method = "holm")

Data$Hochberg =
      p.adjust(Data$ Raw.p,
               method = "hochberg")

Data$Hommel =
      p.adjust(Data$ Raw.p,
               method = "hommel")

Data$BY =
      p.adjust(Data$ Raw.p,
               method = "BY")

Data

繪製圖表

=====

plot(X, Y,
        xlab="原始的p值",
        ylab="矯正後的P值"
        lty=1,
        lwd=2

調整後的p值與原始的p值的圖爲一系列的25個p值。虛線表示一對一的線。

5個p值的多重比較示例

================

### --------------------------------------------------------------
### 多重比較示例，假設示例
### --------------------------------------------------------------
Data = read.table(Input,header=TRUE)

執行p值調整並添加到數據幀

Data$Bonferroni =
      p.adjust(Data$Raw.p,
               method = "bonferroni")

Data$BH =
      signif(p.adjust(Data$Raw.p,
               method = "BH"),
             4)

Data$Holm =
      p.adjust(Data$ Raw.p,
               method = "holm")

Data$Hochberg =
      p.adjust(Data$ Raw.p,
               method = "hochberg")

Data$Hommel =
      p.adjust(Data$ Raw.p,
               method = "hommel")

Data$BY =
      signif(p.adjust(Data$ Raw.p,
               method = "BY"),
             4)

Data

繪製(圖表)

=======