Erreur dans rollapply: indice en dehors des limites

Question

Je voudrais d'abord décrire mon problème: Ce que je veux faire est de calculer le nombre de pics sur les prix dans une fenêtre de 24 heures, alors que je possède des données demi-heure.

J'ai vu tous les posts de Stackoverflow comme e.g. celui-ci: Rollapply for time series

(S'il y a les plus pertinents, s'il vous plaît laissez-moi savoir;))

Comme je ne peux pas et devrait probablement pas télécharger mes données, voici un exemple minimal: je simule un échantillon aléatoire variable, le convertit en un objet xts, et utilise une fonction définie par l'utilisateur pour détecter les "pointes" (bien sûr assez ridicule dans ce cas, mais illustre l'erreur).

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) 

La fonction fonctionne comme elle le devrait: elle prend les 10 valeurs précédentes et calcule la moyenne. Puis, j'ai défini une fonction propre pour trouver des pics: Un pic est un maximum local (supérieur à m points autour de lui) ET est au moins aussi grand que la moyenne des séries temporelles + h. Cela conduit à:

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 
} 

Et fonctionne très bien: Retour à l'exemple:

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") 

Cependant, en utilisant la fonction rollapply() conduit à:

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 

j'ai pensé, eh bien, peut-être que l'erreur se produit car elle peut s'exécuter dans un index négatif pour les premiers points, à cause du paramètre m. Malheureusement, mettre m à zéro ne change pas l'erreur.

J'ai essayé de tracer cette erreur aussi, mais je ne trouve pas la source. Quelqu'un peut-il m'aider ici?

Edit: Une image de pointes: Spikes on the australian Electricity Market. find_peaks(20,50) determines the red points to be spikes, find_peaks(0,50) additionally finds the blue ones to be spikes (therefore, the second parameter h is important, because the blue points are clearly not what we want to analyse when we talk about spikes).

Answer 1

Je ne suis toujours pas tout à fait sûr de ce qu'il est que vous êtes après. En supposant que donné une fenêtre de données que vous voulez identifier si son centre est plus grand que le reste de la fenêtre en même temps comme étant supérieure à la moyenne de la fenêtre + h alors vous pourriez faire ce qui suit:

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) 

Bien qu'il me semble que si le point central est plus grand que ses voisins, il est nécessairement plus grand que la moyenne locale aussi cette partie de la fonction ne serait pas nécessaire si h >= 0. Si vous voulez utiliser la moyenne globale de la série temporelle, remplacez simplement le calcul de meandat par la moyenne globale pré-calculée passée en argument à peakfinder.