乐趣区

关于算法:R语言多重比较方法

原文链接:http://tecdat.cn/?p=9578


假设检验的基本原理是小概率原理,即咱们认为小概率事件在一次试验中实际上不可能产生。

多重比拟的问题

当同一钻研问题下进行屡次假设检验时,不再合乎小概率原理所说的“一次试验”。如果在该钻研问题下只有有测验是阳性的,就对该问题下阳性论断的话,对该问题的测验的犯一类谬误的概率就会增大。如果同一问题下进行 n 次测验,每次的测验水准为 α(每次假阳性概率为 α),则 n 次测验至多呈现一次假阳性的概率会比 α 大。假如每次测验独立的条件下该概率可减少至


常见的多重比拟情景包含:

  • 多组间比拟
  • 多个次要指标
  • 临床试验中期中剖析
  • 亚组剖析

管制多重比拟舛误(Familywise error rate):Bonferroni改正

Bonferroni 法失去的改正 P 值 =P×n
Bonferroni 法非常简单,它的毛病在于十分激进(大略是各种办法中最激进的了),尤其当 n 很大时,通过 Bonferroni 法改正后总的一类谬误可能会远远小于既定 α。
 

管制 谬误发现率:Benjamini & Hochberg 法

简称 BH 法。首先将各 P 值从小到大排序,生成程序数
排第 k 的改正 P 值 =P×n/k
另外要保障改正后的各测验的 P 值大小程序不发生变化。

怎么做 测验

R 内置了一些办法来调整一系列 p 值,以管制多重比拟舛误(Familywise error rate)或管制谬误发现率。

Holm、Hochberg、Hommel 和 Bonferroni 办法管制了多重比拟舛误(Familywise error rate)。这些办法试图限度谬误发现的概率(I 型谬误,在没有实际效果时谬误地回绝无效假设),因而都是绝对较激进的。

办法 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

 

绘制(图表)

plot(X, Y,
        type="l",

调整后的 p 值与原始 p 值在 0 到 0.1 之间的一系列 5 个 p 值的绘图。请留神,Holm 和 Hochberg 的值与 Hommel 雷同,因而被 Hommel 暗藏。虚线示意一对一的线。


最受欢迎的见解

1.Matlab 马尔可夫链蒙特卡罗法(MCMC)预计随机稳定率(SV,Stochastic Volatility)模型

2.基于 R 语言的疾病制图中自适应核密度估计的阈值抉择办法

3.WinBUGS 对多元随机稳定率模型:贝叶斯预计与模型比拟

4. R 语言回归中的 hosmer-lemeshow 拟合优度测验

5.matlab 实现 MCMC 的马尔可夫切换 ARMA – GARCH 模型预计

6. R 语言区间数据回归剖析

7. R 语言 WALD 测验 VS 似然比测验

8.python 用线性回归预测股票价格

9. R 语言如何在生存剖析与 Cox 回归中计算 IDI,NRI 指标

退出移动版