Recursos para el análisis de series temporales interrumpidas en R

12

Soy bastante nuevo en R. He intentado leer sobre análisis de series de tiempo y ya he terminado

  1. Análisis de series de tiempo de Shumway y Stoffer y sus aplicaciones 3ra Edición ,
  2. Excelente pronóstico de Hyndman : principios y práctica
  3. El uso de R de Avril Coghlan para el análisis de series temporales
  4. A. Ian McLeod et al Análisis de series de tiempo con R
  5. Análisis de series de tiempo aplicadas del Dr. Marcel Dettling

Editar: no estoy seguro de cómo manejar esto, pero encontré un recurso útil fuera de Cross Validated. Quería incluirlo aquí en caso de que alguien se tope con esta pregunta.

Análisis de regresión segmentada de estudios de series de tiempo interrumpido en la investigación del uso de medicamentos.

Tengo una serie temporal univariada de la cantidad de artículos consumidos (datos de conteo) medidos diariamente durante 7 años. Se aplicó una intervención a la población de estudio aproximadamente a la mitad de la serie temporal. No se espera que esta intervención produzca un efecto inmediato y el momento del inicio del efecto es esencialmente incognoscible.

Usando el forecastpaquete de Hyndman, he ajustado un modelo ARIMA a los datos previos a la intervención auto.arima(). Pero no estoy seguro de cómo usar este ajuste para responder si ha habido un cambio estadísticamente significativo en la tendencia y cuantificar la cantidad.

# for simplification I will aggregate to monthly counts
# I can later generalize any teachings the community supplies
count <- c(2464, 2683, 2426, 2258, 1950, 1548, 1108,  991, 1616, 1809, 1688, 2168, 2226, 2379, 2211, 1925, 1998, 1740, 1305,  924, 1487, 1792, 1485, 1701, 1962, 2896, 2862, 2051, 1776, 1358, 1110,  939, 1446, 1550, 1809, 2370, 2401, 2641, 2301, 1902, 2056, 1798, 1198,  994, 1507, 1604, 1761, 2080, 2069, 2279, 2290, 1758, 1850, 1598, 1032,  916, 1428, 1708, 2067, 2626, 2194, 2046, 1905, 1712, 1672, 1473, 1052,  874, 1358, 1694, 1875, 2220, 2141, 2129, 1920, 1595, 1445, 1308, 1039,  828, 1724, 2045, 1715, 1840)
# for explanatory purposes
# month <- rep(month.name, 7)
# year <- 1999:2005
ts <- ts(count, start(1999, 1))
train_month <- window(ts, start=c(1999,1), end = c(2001,1))
require(forecast)
arima_train <- auto.arima(train_month)
fit_month <- Arima(train_month, order = c(2,0,0), seasonal = c(1,1,0), lambda = 0)
plot(forecast(fit_month, 36)); lines(ts, col="red")

¿Hay algún recurso que se ocupe específicamente del análisis de series temporales interrumpidas en R? Encontré esto relacionado con ITS en SPSS pero no he podido traducir esto a R.

dais.johns
fuente
¿Desea hacer una inferencia sobre si la intervención tuvo un efecto estadísticamente significativo o desea modelar la intervención para obtener mejores pronósticos ? ¿Y podría hacer que los datos estén disponibles?
Stephan Kolassa
@StephanKolassa ¡Ciertamente! Mi objetivo es hacer inferencia. Proporcionaré datos ficticios en una Edición para ilustrar mejor mi punto.
dais.johns
@StephanKolassa Datos proporcionados a lo mejor de mis habilidades.
dais.johns
Investigaciones previas sugieren que el efecto de la intervención está en la escala de +/- 5% de cambio.
dais.johns
@StephanKolassa Proporcionó datos utilizables reales
dais.johns

Respuestas:

4

Esto se conoce como análisis de punto de cambio. El paquete R changepointpuede hacer esto por usted: consulte la documentación aquí (incluidas las referencias a la literatura): http://www.lancs.ac.uk/~killick/Pub/KillickEckley2011.pdf

Brent Kerby
fuente
Gracias. Estoy investigando esto. Por lo que puedo decir, esto calcularía los posibles puntos de cambio en la serie, pero no analizará la diferencia de tendencia. Pido disculpas si esta suposición es incorrecta. No he podido revisar el paquete de manera superficial.
dais.johns
Después de identificar el punto de cambio, puede dividir los datos en dos series de tiempo (antes y después del punto de cambio) y estimar los parámetros de las dos series de tiempo por separado. Un par de sugerencias más: como sus datos tienen una fuerte tendencia estacional, esto debe eliminarse antes del análisis del punto de cambio; y si va a usar un modelo ARIMA, entonces la diferenciación también debe realizarse antes del análisis del punto de cambio (o, alternativamente, deberá usar un procedimiento más especializado).
Brent Kerby
Gracias por sus sugerencias que intentaré implementar y marcaré como "respondidas" si esto resuelve el problema.
dais.johns
2

Sugeriría un modelo jerárquico de medidas repetidas. Este método debería proporcionar resultados sólidos, ya que cada individuo actuará como su propio control. Intente revisar este enlace de UCLA.

kblansit
fuente
0

Para un enfoque bayesiano, puede usar mcppara ajustar un modelo de Poisson o Binomial (porque tiene recuentos de períodos de intervalo fijo) con autorregresión aplicada a los residuos (en el espacio logarítmico). Luego compare un modelo de dos segmentos con un modelo de un segmento utilizando validación cruzada.

Antes de comenzar, tenga en cuenta que para este conjunto de datos, este modelo no se ajusta bien y la validación cruzada parece inestable. Por lo tanto, me abstendría de usar lo siguiente en escenarios de alto riesgo, pero ilustra un enfoque general:

# Fit the change point model
library(mcp)
model_full = list(
  count ~ 1 + ar(1),  # intercept and AR(1)
  ~ 1  # New intercept
)
fit_full = mcp(model_full, data = df, family = poisson(), par_x = "year")


# Fit the null model
model_null = list(
  count ~ 1 + ar(1)  # just a stable AR(1)
)
fit_null = mcp(model_null, data = df, family = poisson(), par_x = "year")

# Compare predictive performance using LOO cross-validation
fit_full$loo = loo(fit_full)
fit_null$loo = loo(fit_null)
loo::loo_compare(fit_full$loo, fit_null$loo)

Para el presente conjunto de datos, esto da como resultado

       elpd_diff se_diff
model2    0.0       0.0 
model1 -459.1      64.3 

Es decir, una elpd_diff/se_diffrelación de alrededor de 7 a favor del modelo nulo (sin cambios). Las posibles mejoras incluyen:

  • modelando la tendencia periódica usando sin()o cos().
  • adición de información previa acerca de la ubicación probable del cambio, por ejemplo, prior = list(cp_1 = dnorm(1999.8, 0.5).

Obtenga más información sobre el modelado de autorregresión, la comparación de modelos y la configuración previa del mcpsitio web . Divulgación: soy el desarrollador de mcp.

Jonas Lindeløv
fuente