0

我首先想描述我的问题:我想要做的是计算 24 小时窗口内价格峰值的数量,而我拥有半小时的数据。

我看过所有 Stackoverflow 的帖子,比如这个: Rollapply for time series

(如果有更多相关的,请告诉我;))

因为我不能也可能也不应该上传我的数据,这里有一个最小的例子:我模拟一个随机变量,将它转换为一个 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参数的原因而对第一个点运行 int 负索引。可悲的是,设置m为零不会改变错误。

我也试图追踪这个错误,但没有找到源头。有谁可以帮我离开这里吗?

编辑:尖峰图片:澳大利亚电力市场的尖峰。find_peaks(20,50) 确定红色点为尖峰,find_peaks(0,50) 额外发现蓝色点为尖峰(因此,第二个参数 h 很重要,因为蓝色点显然不是我们要分析的当我们谈论尖峰时)。

4

1 回答 1

0

我仍然不完全确定你追求的是什么。假设给定一个数据窗口,您想确定其中心是否大于窗口的其余部分,同时大于窗口的平均值,+ 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

于 2016-02-04T13:53:56.670 回答