R Language - сортировка данных по диапазонам; в среднем; игнорировать выбросы - PullRequest
5 голосов
/ 30 января 2011

Я анализирую данные с ветряной турбины, обычно это то, что я бы делал в Excel, но количество данных требует чего-то сверхмощного. Я никогда раньше не использовал R, поэтому я просто ищу несколько указателей.

Данные состоят из 2 столбцов WindSpeed ​​ и Power , так что до сих пор я пришел к импорту данных из файла CSV и нанес их на график с разбивкой.

Далее я хотел бы отсортировать данные по диапазонам; например, все данные, где WindSpeed ​​ находится между x и y, а затем найти среднее значение мощности, генерируемой для каждого диапазона, и построить график кривой.

Из этого среднего я хочу пересчитать среднее на основе данных, которые попадают в одно из двух стандартных отклонений среднего (в основном игнорируя выбросы).

Любые указатели приветствуются.

Для тех, кто заинтересован, я пытаюсь создать график, аналогичный this . Это довольно стандартный тип графиков, но, как я уже сказал, для сдвига количества данных требуется нечто более тяжелое, чем Excel.

Ответы [ 5 ]

5 голосов
/ 30 января 2011

Поскольку вы больше не находитесь в Excel, почему бы не использовать современную статистическую методологию, которая не требует грубого объединения данных и специальных методов для устранения выбросов: локально сглаженная регрессия, как это реализовано в loess.1002 * Используя небольшую модификацию выборочных данных csgillespie:

w_sp <- sample(seq(0, 100, 0.01), 1000)
power <- 1/(1+exp(-(w_sp -40)/5)) + rnorm(1000, sd = 0.1)

plot(w_sp, power)

x_grid <- seq(0, 100, length = 100)
lines(x_grid, predict(loess(power ~ w_sp), x_grid), col = "red", lwd = 3)
2 голосов
/ 30 января 2011

Добавьте эту версию, похожую по мотивации на @ hadley's, в микс, используя аддитивную модель с адаптивным сглаживателем, используя пакет mgcv:

Сначала фиктивные данные, используемые @ hadley

w_sp <- sample(seq(0, 100, 0.01), 1000)
power <- 1/(1+exp(-(w_sp -40)/5)) + rnorm(1000, sd = 0.1)
df <- data.frame(power = power, w_sp = w_sp)

Подгонка аддитивной модели с использованием gam() с использованием адаптивного сглаживания и выбора сглаживания с помощью REML

require(mgcv)
mod <- gam(power ~ s(w_sp, bs = "ad", k = 20), data = df, method = "REML")
summary(mod)

Прогнозирование из нашей модели и получение стандартных ошибок подгонки, используйте последнюю для получения приблизительного значения 95% доверительный интервал

x_grid <- with(df, data.frame(w_sp = seq(min(w_sp), max(w_sp), length = 100)))
pred <- predict(mod, x_grid, se.fit = TRUE)
x_grid <- within(x_grid, fit <- pred$fit)
x_grid <- within(x_grid, upr <- fit + 2 * pred$se.fit)
x_grid <- within(x_grid, lwr <- fit - 2 * pred$se.fit)

Изобразите все и Лесс подходит для сравнения

plot(power ~ w_sp, data = df, col = "grey")
lines(fit ~ w_sp, data = x_grid, col = "red", lwd = 3)
## upper and lower confidence intervals ~95%
lines(upr ~ w_sp, data = x_grid, col = "red", lwd = 2, lty = "dashed")
lines(lwr ~ w_sp, data = x_grid, col = "red", lwd = 2, lty = "dashed")
## add loess fit from @hadley's answer
lines(x_grid$w_sp, predict(loess(power ~ w_sp, data = df), x_grid), col = "blue",
      lwd = 3)

adaptive smooth and loess fits

2 голосов
/ 30 января 2011

Сначала мы создадим несколько примеров данных, чтобы конкретизировать проблему:

w_sp = sample(seq(0, 100, 0.01), 1000)
power = 1/(1+exp(-(rnorm(1000, mean=w_sp, sd=5) -40)/5))

Предположим, мы хотим поместить значения power между [0,5), [5,10) и т. Д. Тогда

bin_incr = 5
bins = seq(0, 95, bin_incr)
y_mean = sapply(bins, function(x) mean(power[w_sp >= x & w_sp < (x+bin_incr)]))

Мы создали средние значения между интересующими диапазонами. Обратите внимание: если вы хотите получить медианные значения, просто измените значение mean на median. Все, что осталось сделать, это построить их:

plot(w_sp, power)
points(seq(2.5, 97.5, 5), y_mean, col=3, pch=16)

Чтобы получить среднее на основе данных, которые попадают в два стандартных отклонения от среднего, нам нужно создать чуть более сложную функцию:

noOutliers = function(x, power, w_sp, bin_incr) {
  d = power[w_sp >= x & w_sp < (x + bin_incr)]
  m_d = mean(d)
  d_trim = mean(d[d > (m_d - 2*sd(d)) & (d < m_d + 2*sd(d))])
  return(mean(d_trim))
}

y_no_outliers = sapply(bins, noOutliers, power, w_sp, bin_incr)
1 голос
/ 30 января 2011

Вот некоторые примеры подогнанных кривых (анализ Вейбулла) для коммерческих турбин:

http://www.inl.gov/wind/software/

http://www.irec.cmerp.net/papers/WOE/Paper%20ID%20161.pdf

http://www.icaen.uiowa.edu/~ie_155/Lecture/Power_Curve.pdf

0 голосов
/ 30 января 2011

Я бы также порекомендовал поиграть с собственным ggplot2 Хэдли.Его сайт - отличный ресурс: http://had.co.nz/ggplot2/.

    # If you haven't already installed ggplot2:
    install.pacakges("ggplot2", dependencies = T)

    # Load the ggplot2 package
    require(ggplot2)

    # csgillespie's example data
    w_sp <- sample(seq(0, 100, 0.01), 1000)
    power <- 1/(1+exp(-(w_sp -40)/5)) + rnorm(1000, sd = 0.1)

    # Bind the two variables into a data frame, which ggplot prefers
    wind <- data.frame(w_sp = w_sp, power = power)

    # Take a look at how the first few rows look, just for fun
    head(wind)


    # Create a simple plot
    ggplot(data = wind, aes(x = w_sp, y = power)) + geom_point() + geom_smooth()

    # Create a slightly more complicated plot as an example of how to fine tune
    # plots in ggplot
    p1 <- ggplot(data = wind, aes(x = w_sp, y = power))
    p2 <- p1 + geom_point(colour = "darkblue", size = 1, shape = "dot") 
    p3 <- p2 + geom_smooth(method = "loess", se = TRUE, colour = "purple")
    p3 + scale_x_continuous(name = "mph") + 
             scale_y_continuous(name = "power") +
             opts(title = "Wind speed and power")
...