Ошибка в rollapply: индекс за пределами

Question

Сначала я хотел бы описать мою проблему: я хочу рассчитать количество скачков цен в 24-часовом окне, в то время как у меня есть полчасовые данные.

Я видел все сообщения Stackru, например, вот этот: Rollapply для временных рядов

(Если есть более подходящие, пожалуйста, дайте мне знать;))

Поскольку я не могу и, вероятно, также не должен загружать свои данные, вот минимальный пример: я имитирую случайную переменную, преобразую ее в объект xts и использую пользовательскую функцию для обнаружения "всплесков" (конечно, довольно смешно в этом случае, но иллюстрирует ошибку).

library(xts)
##########Simulate y as a random variable
y <- rnorm(n=100)
##########Add a date variable so i can convert it to a xts object later on
yDate <- as.Date(1:100)
##########bind both variables together and convert to a xts object
z <- cbind(yDate,y)
z <- xts(x=z, order.by=yDate)
##########use the rollapply function on the xts object:
x <- rollapply(z, width=10, FUN=mean)

Функция работает так, как она должна: она принимает 10 предыдущих значений и вычисляет среднее значение.

Затем я определил собственную функцию для поиска пиков: пик - это локальный максимум (выше, чем m точек вокруг него), И он, по крайней мере, такой же большой, как среднее для временных рядов +h. Это ведет к:

find_peaks <- function (x, m,h){
  shape <- diff(sign(diff(x, na.pad = FALSE)))
  pks <- sapply(which(shape < 0), FUN = function(i){
    z <- i - m + 1
    z <- ifelse(z > 0, z, 1)
    w <- i + m + 1
    w <- ifelse(w < length(x), w, length(x))
    if(all(x[c(z : i, (i + 2) : w)] <= x[i + 1])&x[i+1]>mean(x)+h) return(i + 1) else return(numeric(0))
  })
  pks <- unlist(pks)
  pks
}

И работает нормально: Вернемся к примеру:

plot(yDate,y)
#Is supposed to find the points which are higher than 3 points around them
#and higher than the average:
#Does so, so works.
points(yDate[find_peaks(y,3,0)],y[find_peaks(y,3,0)],col="red")

Однако, используя rollapply() Функция приводит к:

x <- rollapply(z,width = 10,FUN=function(x) find_peaks(x,3,0))
#Error in `[.xts`(x, c(z:i, (i + 2):w)) : subscript out of bounds 

Сначала я подумал, что, может быть, ошибка возникает потому, что для первых точек она может работать с отрицательным индексом из-за m параметр. К сожалению, установка m ноль не меняет ошибку.

Я тоже пытался отследить эту ошибку, но не нашел источника. Может кто-нибудь помочь мне здесь?

Редактировать: изображение шипов: шипы на австралийском рынке электроэнергии. find_peaks (20,50) определяет красные точки, которые будут шипами, find_peaks(0,50) дополнительно находит синие точки как шипы (поэтому второй параметр h важен, потому что синие точки явно не то, что мы хотим анализировать когда мы говорим о шипах).

Answer 1

Я до сих пор не совсем уверен, что вы ищете. Предполагая, что для данного окна данных вы хотите определить, больше ли его центр, чем остальная часть окна, в то же время, что он больше среднего значения окна. + h тогда вы можете сделать следующее:

peakfinder = function(x,h = 0){
  xdat = as.numeric(x)
  meandat = mean(xdat)
  center = xdat[ceiling(length(xdat)/2)]
  ifelse(all(center >= xdat) & center >= (meandat + h),center,NA)
}

y <- rnorm(n=100)
z = xts(y, order.by = as.Date(1:100))
plot(z)
points(rollapply(z,width = 7, FUN = peakfinder, align = "center"), col = "red", pch = 19)

Хотя мне кажется, что если центральная точка больше, чем ее соседи, она обязательно будет больше, чем локальное среднее значение, поэтому эта часть функции не будет нужна, если h >= 0, Если вы хотите использовать глобальное среднее временного ряда, просто подставьте вычисление meandat с предварительно рассчитанным глобальным средним значением, переданным в качестве аргумента peakfinder,