2009-12-14 17 views
7

Tengo una serie temporal estacionaria a la que quiero ajustar un modelo lineal con un término autorregresivo para corregir la correlación serial, es decir, usando la fórmula At = c1 * Bt + c2 * Ct + ut, donde ut = r * ut-1 + etGLM con término autorregresivo para corregir la correlación serial

(ut es un AR (1) plazo para corregir la correlación serial en los términos de error)

¿alguien sabe qué usar en R para modelar ¿esta?

Gracias Karl

Respuesta

8

El paquete GLMMarp se ajustará a estos modelos. Si solo desea un modelo lineal con errores gaussianos, puede hacerlo con la función arima() donde las covariables se especifican a través del argumento xreg.

+0

El paquete GLMMarp se eliminó del repositorio CRAN. ¿Conoces algún otro paquete que logre esto? – TMS

+0

Todos los paquetes que conozco están listados en http://cran.r-project.org/web/views/TimeSeries.html –

2

¿Cuál es su función de enlace?

La forma de describirlo suena como una regresión lineal básica con errores autocorrelacionados. En ese caso, una opción es usar lm para obtener una estimación consistente de sus coeficientes y usar Newey-West HAC standard errors.

No estoy seguro de la mejor respuesta para GLM en general.

3

Hay varias maneras de hacer esto en R. He aquí dos ejemplos utilizando el "Seatbelts" time series dataset in the datasets package que viene con R.

La función arima() viene en el paquete: las estadísticas que se incluyen con R. La función toma un argumento de el formulario order=c(p, d, q) donde puede especificar el orden del componente autoregresivo, integrado y promedio móvil. En su pregunta, sugiere que desea crear un modelo AR (1) para corregir la autocorrelación de primer orden en los errores y eso es todo. Podemos hacer eso con el siguiente comando:

arima(Seatbelts[,"drivers"], order=c(1,0,0), 
     xreg=Seatbelts[,c("kms", "PetrolPrice", "law")]) 

El valor para la orden especifica que queremos un (1) modelo AR. El componente xreg debe ser una serie de otras X que queremos agregar como parte de una regresión. La salida se ve un poco como la salida de summary.lm() puesta de lado.

Otro proceso alternativo podría ser más familiar para la forma en que usted ha ajustado los modelos de regresión es usar gls() en el nlme package. El siguiente código convierte el objeto de series de tiempo del cinturón de seguridad en una trama de datos y luego extrae y añade una nueva columna (t) que es sólo un contador en el objeto de series de tiempo Ordenada:

Seatbelts.df <- data.frame(Seatbelts) 
Seatbelts.df$t <- 1:(dim(Seatbelts.df)[1]) 

Las dos líneas anteriores sólo son poniendo los datos en forma. Como la función arima() está diseñada para series de tiempo, puede leer objetos de series de tiempo más fácilmente. Para ajustar el modelo con nlme usted entonces ejecute:

library(nlme) 
m <- gls(drivers ~ kms + PetrolPrice + law, 
     data=Seatbelts.df, 
     correlation=corARMA(p=1, q=0, form=~t)) 
summary(m) 

La línea que comienza con "correlación" es la forma en que se pasa en la estructura de correlación ARMA de GLS. Los resultados no serán exactamente los mismos porque arima() usa la máxima probabilidad de estimar modelos y gls() usa la máxima verosimilitud restringida de forma predeterminada. Si agrega method="ML" a la llamada al gls(), obtendrá las estimaciones idénticas que obtuvo con la función ARIMA anterior.