如何重复 Grubbs 测试并标记异常值

Question

我想对一组数据重复应用 Grubbs 检验，直到它不再发现异常值。我希望标记而不是删除异常值，以便我可以将数据绘制为直方图，异常值具有不同的颜色。我使用了异常值包中的 grubbs.test 来手动识别异常值，但无法弄清楚如何循环遍历它们并成功标记它们。我的目标输出如下：

X   Outlier
152.36  Yes
130.38  Yes
101.54  No
96.26   No
88.03   No
85.66   No
83.62   No
76.53   No
74.36   No
73.87   No
73.36   No
73.35   No
68.26   No
65.25   No
63.68   No
63.05   No
57.53   No

Answer 1

Sam Dickson 的 answer 很棒，但是如果您到达一个点，即除了两个值之外的所有值都被标记为异常值，或者如果您一开始只使用三个值 (grubbs.test() won'如果输入向量中只有两个值，则返回一个 p 值）。

我为这种意外情况在 while 循环中添加了一个断点，如果发生这种情况，它也会发出警告。此外，当您从少于两个输入值开始时，它会引发信息错误。

grubbs.flag <- function(x) {
  outliers <- NULL
  test <- x
  grubbs.result <- grubbs.test(test)
  pv <- grubbs.result$p.value
  # throw an error if there are too few values for the Grubb's test
  if (length(test) < 3 ) stop("Grubb's test requires > 2 input values")
  while(pv < 0.05) {
    outliers <- c(outliers,as.numeric(strsplit(grubbs.result$alternative," ")[[1]][3]))
    test <- x[!x %in% outliers]
    # stop if all but two values are flagged as outliers
    if (length(test) < 3 ) {
      warning("All but two values flagged as outliers")
      break
    }
    grubbs.result <- grubbs.test(test)
    pv <- grubbs.result$p.value
  }
  return(data.frame(X=x,Outlier=(x %in% outliers)))
}

当然，如果您一开始只有三个数据点，那么可能进行异常值测试没有多大意义，但我不了解您的业务。

Answer 2

看起来你需要一个简短的函数来做你想做的事：

library(outliers)
library(ggplot2)

X <- c(152.36,130.38,101.54,96.26,88.03,85.66,83.62,76.53,
       74.36,73.87,73.36,73.35,68.26,65.25,63.68,63.05,57.53)

grubbs.flag <- function(x) {
  outliers <- NULL
  test <- x
  grubbs.result <- grubbs.test(test)
  pv <- grubbs.result$p.value
  while(pv < 0.05) {
    outliers <- c(outliers,as.numeric(strsplit(grubbs.result$alternative," ")[[1]][3]))
    test <- x[!x %in% outliers]
    grubbs.result <- grubbs.test(test)
    pv <- grubbs.result$p.value
  }
  return(data.frame(X=x,Outlier=(x %in% outliers)))
}

这是输出：

grubbs.flag(X)
         X Outlier
1   152.36    TRUE
2   130.38    TRUE
3   101.54   FALSE
4    96.26   FALSE
5    88.03   FALSE
6    85.66   FALSE
7    83.62   FALSE
8    76.53   FALSE
9    74.36   FALSE
10   73.87   FALSE
11   73.36   FALSE
12   73.35   FALSE
13   68.26   FALSE
14   65.25   FALSE
15   63.68   FALSE
16   63.05   FALSE
17   57.53   FALSE

如果你想要一个不同颜色的直方图，你可以使用以下：

ggplot(grubbs.flag(X),aes(x=X,color=Outlier,fill=Outlier))+
  geom_histogram(binwidth=diff(range(X))/30)+
  theme_bw()