2016-03-13 3 views
2

Я работаю над функцией, которая избавится от выбросов в заданном наборе данных на основе правила 3 ​​сигмы. Мой код представлен ниже. «data» - это набор данных для обработки.Функция: sapply in apply, удаление выбросов

rm.outlier <- function(data){ 

    apply(data, 2, function(var) { 
     sigma3.plus <- mean(var) + 3 * sd(var) 
     sigma3.min <- mean(var) - 3 * sd(var) 
     sapply(var, function(y) { 
     if (y > sigma3.plus){ 
      y <- sigma3.plus 
     } else if (y < sigma3.min){ 
      y <- sigma3.min 
     } else {y <- y} 
     }) 
    }) 
    as.data.frame(data) 
} 

Для того, чтобы проверить, если функция работает, я написал короткий тест:

set.seed(123) 
a <- data.frame("var1" = rnorm(10000, 0, 1)) 
b <- a 
sum(a$var1 > mean(a$var1) + 3 * sd(a$var1)) # number of outliers in a 

В результате, я получаю:

[1] 12

Таким образом, переменная var1 в кадре данных a имеет 12 выбросов. Далее я попытаюсь применить свою функцию на этом объекте:

a2 <- rm.outlier(a) 
sum(b$var1 - a2$var1) 

К сожалению, это дает 0, что ясно указывает на то, что что-то не работает. Я уже разработал, что реализация sapply верна, поэтому в моем применении должна быть ошибка. Любая помощь будет оценена по достоинству.

ответ

1

Кажется, вы просто забыли присвоить свои результаты функции apply новому кадру данных. (Сравните третью строку с кодом)

rm.outlier <- function(data){ 

    # Assign the result to a new dataframe 
    data_new <- apply(data, 2, function(var) { 
    sigma3.plus <- mean(var) + 3 * sd(var) 
    sigma3.min <- mean(var) - 3 * sd(var) 
    sapply(var, function(y) { 
     if (y > sigma3.plus){ 
     y <- sigma3.plus 
     } else if (y < sigma3.min){ 
     y <- sigma3.min 
     } else {y <- y} 
    }) 
    }) 

    # Print the new dataframe 
    as.data.frame(data_new) 
} 

set.seed(123) 
a <- data.frame("var1" = rnorm(10000, 0, 1)) 
sum(a$var1 > mean(a$var1) + 3 * sd(a$var1)) # number of too big outliers 
# 15 
sum(a$var1 < mean(a$var1) - 3 * sd(a$var1)) # number of too small outliers 
# 13 
# Overall 28 outliers 

# Check the function for the number of outliers 
a2 <- rm.outlier(a) 
sum(a2$var1 == a$var1) - length(a$var1) 
3

Если время выполнения важно для вас, вы можете рассмотреть другой подход. Вы можете процитировать эту фильтрацию, например. используя pmin и pmax, которые являются одинаково читаемыми и в 15 раз быстрее. Если вам нравится это немного более сложным, вы могли бы использовать findInterval и получить еще больше скорости:

rm.outlier2 <- function(x) { 
    ## calculate -3/3 * sigma borders 
    s <- mean(x) + c(-3, 3) * sd(x) 
    pmin(pmax(x, s[1]), s[2]) 
} 

rm.outlier3 <- function(x) { 
    ## calculate -3/3 * sigma borders 
    s <- mean(x) + c(-3, 3) * sd(x) 
    ## sorts x into intervals 0 == left of s[1], 2 == right of s[2], 1 
    ## between both s 
    i <- findInterval(x, s) 
    ## which values are left/right of the interval 
    j <- which(i != 1L) 
    ## add a value between s to directly use output of findInterval for subsetting 
    s2 <- c(s[1], 0, s[2]) 
    ## replace all values that are left/right of the interval 
    x[j] <- s2[i[j] + 1L] 
    x 
} 

Бенчмаркинг материал:

## slightly modified OP version 
rm.outlier <- function(x) { 
    sigma3 <- mean(x) + c(-3,3) * sd(x) 
    sapply(x, function(y) { 
    if (y > sigma3[2]){ 
     y <- sigma3[2] 
    } else if (y < sigma3[1]){ 
     y <- sigma3[1] 
    } else {y <- y} 
    }) 
} 

set.seed(123) 
a <- rnorm(10000, 0, 1) 

# check output 
all.equal(rm.outlier(a), rm.outlier2(a)) 
all.equal(rm.outlier2(a), rm.outlier3(a)) 

library("rbenchmark") 

benchmark(rm.outlier(a), rm.outlier2(a), rm.outlier3(a), 
      order = "relative", 
      columns = c("test", "replications", "elapsed", "relative")) 
#   test replications elapsed relative 
#3 rm.outlier3(a)   100 0.028 1.000 
#2 rm.outlier2(a)   100 0.102 3.643 
#1 rm.outlier(a)   100 1.825 65.179 
+0

Ого, это действительно ускоряли вещи, спасибо большое! – kaksat

Смежные вопросы