2011-02-11 12 views
14

Tengo una serie temporal para la cual quiero interpolar inteligentemente los valores perdidos. El valor en un momento particular está influenciado por una tendencia de varios días, así como su posición en el ciclo diario.Interpolar valores perdidos en una serie temporal con un ciclo estacional

Aquí se muestra un ejemplo en el que el décimo observación no se encuentra en myzoo

start <- as.POSIXct("2010-01-01") 
freq <- as.difftime(6, units = "hours") 
dayvals <- (1:4)*10 
timevals <- c(3, 1, 2, 4) 
index <- seq(from = start, by = freq, length.out = 16) 
obs <- (rep(dayvals, each = 4) + rep(timevals, times = 4)) 
myzoo <- zoo(obs, index) 
myzoo[10] <- NA 

Si tuviera que implementar esto, me gustaría usar algún tipo de media ponderada de los tiempos estrechos en los días próximos, o añadir una valor del día para una línea de función adaptada a la tendencia más amplia, pero espero que ya exista algún paquete o funciones que se apliquen a esta situación?

EDITAR: Modifiqué el código ligeramente para aclarar mi problema. Existen na.* métodos que se interpolan desde los vecinos más cercanos, pero en este caso no reconocen que el valor faltante se encuentra en el momento que es el valor más bajo del día. Quizás la solución sea redirigir los datos a formato ancho y luego interpolar, pero no me gustaría ignorar por completo los valores contiguos del mismo día. Vale la pena señalar que diff(myzoo, lag = 4) devuelve un vector de 10's. La solución puede estar en alguna combinación de reshape, na.spline y diff.inv, pero no puedo entenderlo.

Aquí hay tres enfoques que no funcionan: enter image description here

Edit2. Imagen producida usando el siguiente código.

myzoo <- zoo(obs, index) 
myzoo[10] <- NA # knock out the missing point 
plot(myzoo, type="o", pch=16) # plot solid line 
points(na.approx(myzoo)[10], col = "red") 
points(na.locf(myzoo)[10], col = "blue") 
points(na.spline(myzoo)[10], col = "green") 
myzoo[10] <- 31 # replace the missing point 
lines(myzoo, type = "o", lty=3, pch=16) # dashed line over the gap 
legend(x = "topleft", 
     legend = c("na.spline", "na.locf", "na.approx"), 
     col=c("green","blue","red"), pch = 1) 
+0

Este código no se ejecuta. índice y obs no están definidos. 'na.approx',' na.spline', 'na.locf' y otras funciones' na. * 'en el paquete del zoo pueden completar los valores' NA'. –

+0

Gracias, pegado el bloque correcto. –

+0

Por favor, muestre el código que usó para crear la trama y explique qué significa "no funciona". –

Respuesta

17

Prueba esto:

x <- ts(myzoo,f=4) 
fit <- ts(rowSums(tsSmooth(StructTS(x))[,-2])) 
tsp(fit) <- tsp(x) 
plot(x) 
lines(fit,col=2) 

La idea es utilizar un modelo estructural básico para la serie de tiempo, que maneja el valor que falta fino usando un filtro de Kalman. Luego se usa un Kalman suave para estimar cada punto de la serie temporal, incluido el omitido.

Tuve que convertir tu objeto de zoológico en un objeto ts con frecuencia 4 para usar StructTS. Es posible que desee volver a cambiar los valores ajustados al zoológico.

+0

Gracias, eso lo hace casi exactamente. Pero aquí hay algo extraño: la trama muestra que el primer punto de "ajuste" está bastante lejos (por .85), y la suma de (ajuste x)^2 es ~ 0.96. Pero, si reemplaza x con 'x <- ts (rev (myzoo), f = 4)', el ajuste se vuelve perfecto. ¿Alguna idea de lo que pasó? –

+0

La función 'zoo :: na.StructTS' realiza líneas 2-3 más fácilmente:' fit2 <- na.StructTS (x) 'crea una serie idéntica a' x' con NA llena a través del filtro estacional de Kalman (30.66, mismo valor que 'fit' en esta respuesta). –

2

En este caso, creo que desea una corrección de estacionalidad en el modelo ARIMA. No hay fecha suficiente aquí para ajustarse al modelo estacional, pero esto debería comenzar.

library(zoo) 
start <- as.POSIXct("2010-01-01") 
freq <- as.difftime(6, units = "hours") 
dayvals <- (1:4)*10 
timevals <- c(3, 1, 2, 4) 
index <- seq(from = start, by = freq, length.out = 16) 
obs <- (rep(dayvals, each = 4) + rep(timevals, times = 4)) 
myzoo <- myzoo.orig <- zoo(obs, index) 
myzoo[10] <- NA 

myzoo.fixed <- na.locf(myzoo) 

myarima.resid <- arima(myzoo.fixed, order = c(3, 0, 3), seasonal = list(order = c(0, 0, 0), period = 4))$residuals 
myzoo.reallyfixed <- myzoo.fixed 
myzoo.reallyfixed[10] <- myzoo.fixed[10] + myarima.resid[10] 

plot(myzoo.reallyfixed) 
points(myzoo.orig) 

En mis pruebas, el ARMA (3, 3) está muy cerca, pero eso es solo suerte. Con una serie temporal más larga, debería poder calibrar la corrección estacional para proporcionarle buenas predicciones. Sería útil tener una buena información previa sobre qué mecanismos subyacentes para la señal y la corrección estacional pueden mejorar el rendimiento de la muestra.

+0

agregó una imagen. trazar no es ningún problema: 'points (na.locf (myzoo) [10], col =" blue ")' –

+0

@jonw - ¡Oh! Entendí mal. Pensé que el problema era obtener un punto. El problema es obtener una buena estimación con esta "estacionalidad", que en realidad es una estacionalidad diaria. Debería haber intentado trazar (simplemente intenté '? Points.zoo'). –

2

forecast::na.interp es un buen enfoque. Del documentation

Utiliza la interpolación lineal para series no estacionales y una descomposición stl periódica con series estacionales para reemplazar valores perdidos.

library(forecast) 
fit <- na.interp(myzoo) 
fit[10] # 32.5, vs. 31.0 actual and 32.0 from Rob Hyndman's answer 

This paper evalúa varios métodos de interpolación en contra de series de tiempo real, y concluye que na.interp es a la vez precisa y eficiente:

De las implementaciones de I ensayados en este trabajo, na.interp de la previsión paquete y na.StructTS del paquete del zoológico mostraron los mejores resultados globales.

La función na.interp tampoco es mucho más lenta que na.approx [el método más rápido], por lo que la descomposición del loess no parece ser muy exigente en términos de tiempo de cálculo.

También digno de mención que Rob Hyndman escribió el paquete forecast, e incluyó na.interp después de proporcionar su respuesta a esta pregunta. Es probable que na.interp sea una mejora con respecto a este enfoque, aunque en este caso empeoró (probablemente debido a la especificación del período en StructTS, donde na.interp lo resuelve).

Cuestiones relacionadas