Tópicos En Series De Tiempo

   EMBED

Share

Preview only show first 6 pages with water mark for full document please download

Transcript

´ picos en Series de Tiempo To Autor: Grisel M. Britos Directora: Dra. Silvia Mar´ıa Ojeda C´ordoba, Diciembre de 2012 Agradecimientos A mi familia por la paciencia, la comprensi´on y tanto que me han dado en todos estos a˜ nos de la carrera. A Mart´ın Ocampo por las sonrisas, la compa˜ n´ıa y el amor que me di´o y sigue dando. A mi directora, por la confianza, la dedicaci´on y sobre todo por las palabras de apoyo y aliento. A los compa˜ neros y amigos de la Facultad: gracias por hacer m´as llevaderas las largas horas de estudio. A los amigos de la vida, por alegrarla. A aquellos profesores que en la Facultad me ense˜ naron con paciencia y entrega. De nuevo, GRACIAS. I Resumen El objetivo de este trabajo es iniciar un estudio sobre series de tiempo bas´andose en el libro de Leiva R. [1]. En el primer cap´ıtulo se introducir´an las principales definiciones y resultados te´ oricos. En los cap´ıtulos II y III se estudiar´ an y analizar´an caracter´ısticas de las diferentes representaciones de series estacionarias y no estacionarias que luego nos permitir´an identificar alg´ un modelo que pudiera haber generado una dada serie de tiempo real. En el cap´ıtulo IV se ver´an algunos m´etodos para predecir como continuar´ a el comportamiento de la serie de tiempo que se estudia, mientras que en el cap´ıtulo V se presentar´ an distintos criterios u ´tiles para determinar si el modelo elegido para la serie es adecuado o se debe proponer un nuevo modelo. Al final de cada cap´ıtulo se resuelven algunos ejercicios propuestos en la obra de Leiva [1]. Por u ´ltimo, se tomar´a una serie de tiempo real y se tratar´ a de modelar su comportamiento atendiendo a todo lo estudiado en este trabajo. Clasificaci´ on (Math. Subject Classification) 62M10, 60G10. Palabras Claves -Serie de Tiempo -Estacionaridad -Modelos ARIMA II ´Indice general Agradecimientos I Resumen II 1. Introducci´ on 1.1. Definiciones B´ asicas . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2. Estimaci´ on de las caracter´ısticas de una serie de tiempo estacionaria 1.2.1. Estimaci´ on de la media µ . . . . . . . . . . . . . . . . . . . . 1.2.2. Estimaci´ on de las autocovarianzas . . . . . . . . . . . . . . . 1.2.3. Estimaci´ on de las autocorrelaciones . . . . . . . . . . . . . . 1.2.4. Estimaci´ on de las autocorrelaciones parciales . . . . . . . . . 1.3. Ejercicios del Cap´ıtulo I . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1 3 3 3 4 4 4 2. Modelos para series de tiempo estacionarias 2.1. Representaciones de una serie de tiempo . . . . . . . . . . . . . . . . . . . . . . 2.1.1. Representaci´ on como un proceso de promedios m´oviles de orden infinito 2.1.2. Representaci´ on como un proceso autorregresivo de orden infinito . . . . 2.2. Modelos autoregresivos de orden finito . . . . . . . . . . . . . . . . . . . . . . . 2.2.1. Modelo autoregresivo de orden 1 . . . . . . . . . . . . . . . . . . . . . . 2.2.2. Modelo autoregresivo de orden p . . . . . . . . . . . . . . . . . . . . . . 2.3. Modelos de promedios m´ oviles de orden finito . . . . . . . . . . . . . . . . . . . 2.3.1. El modelo de medias m´oviles de orden 1 . . . . . . . . . . . . . . . . . . 2.4. Modelos Autorregresivos de Promedios M´oviles . . . . . . . . . . . . . . . . . . 2.4.1. Modelo ARMA(1,1). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.4.2. Modelo ARMA(p,q) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.5. Ejercicios del Cap´ıtulo II . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14 14 14 15 16 16 17 18 19 19 20 22 23 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3. Modelos para series de tiempo no estacionarias 3.1. Enfoque cl´ asico . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.1.1. Estimaci´ on de mt por el m´etodo de cuadrados m´ınimos . . . . . . . . . . . . . 3.1.2. Estimaci´ on de la tendencia mediante suavizamiento por promedios . . . . . . . 3.1.3. Estimaci´ on de la tendencia y la estacionalidad cuando la tendencia es constante durante cada per´ıodo de la componente estacional . . . . . . . . . . . . . . . . 3.1.4. Estimaci´ on de la tendencia y la estacionalidad: Caso General . . . . . . . . . . 3.2. Modelos ARIMA y SARIMA . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.2.1. Modelos ARIMA . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.2.2. Procesos SARIMA (ARIMA estacionales) . . . . . . . . . . . . . . . . . . . . . 3.3. Transformaciones . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.4. Ejercicios del Cap´ıtulo III . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . III 37 37 38 39 39 41 42 43 43 44 45 4. Predicci´ on 4.1. Predictor Lineal con Error Cuadr´atico Medio M´ınimo 4.2. Predicci´ on de un ARMA conocido todo su pasado . . 4.3. Predicci´ on de un ARIMA conocido todo su pasado . . 4.4. Ejercicios del Cap´ıtulo IV . . . . . . . . . . . . . . . . . . . . 53 53 54 56 58 5. Estimaci´ on e identificaci´ on de los modelos de series de tiempos 5.1. Identificaci´ on de procesos. Determinaci´on del orden y verificaci´on de diagn´ostico . . . 5.1.1. Criterios de Akaike: FPE, AIC, BIC . . . . . . . . . . . . . . . . . . . . . . . . 5.2. Reglas pr´ acticas para la identificaci´on de procesos . . . . . . . . . . . . . . . . . . . . 67 67 67 68 6. Aplicaci´ on 70 Bibliograf´ıa 79 IV . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Cap´ıtulo 1 Introducci´ on Muchas veces en la vida real se cuenta con observaciones ordenadas de hechos que acontecen a intervalos equiespaciados en el tiempo y que se encuentran relacionados con diferentes actividades propias del quehacer humano. Por ello, el estudio de las series de tiempo es una herramienta importante para el an´alisis de datos que no son independientes entre s´ı; por ejemplo, tasas de mortalidad y nacimiento, ´ındices de precio al consumidor, producci´on anual de soja, caudales de r´ıos, etc. El an´ alisis de series temporales se aborda usualmente desde la perspectiva que proporcionan la Probabilidad y la Estad´ıstica y permite comprender y describir los procesos generadores de los datos observados, como as´ı tambi´en realizar predicciones de valores futuros de la serie. Este trabajo apunta a realizar un primer acercamiento al estudio de las series de tiempo a partir de la teor´ıa expuesta por Leiva [1] y BrockwellDavis [2]; luego, en trabajos futuros, se espera focalizar aspectos vinculados con la generalizaci´ on a dos dimensiones de esta teor´ıa, la cual a su vez es requerida en t´opicos de an´alisis y procesamiento de im´agenes digitales. El inter´es del trabajo se centr´o en la comprensi´on de las principales definiciones y resultados te´ oricos propuestos en la bibliograf´ıa citada, poniendo especial atenci´on a la resoluci´ on de los ejercicios que se encuentran al final de cada secci´on del libro de Leiva [1]. El trabajo implic´o adem´ as una aproximaci´ on inicial al conocimiento y utilizaci´on del software estad´ıstico y de libre acceso, R, para modelar una serie de tiempo real. 1.1. Definiciones B´ asicas Definici´ on 1.1.1. Sea (Ω, A, P ) un espacio de probabilidad y sea T un conjunto no vac´ıo de ´ındices. Se llama serie de tiempo (o proceso estoc´ astico) a valores reales a una funci´ on X :T × Ω → R tal que para cada t fijo, X :t × Ω → R es una variable aleatoria Xt , que ser´ a denotada como {Xt : t ∈ T }. Debido a la dependencia que hay entre las variables aleatorias Xt , t ∈ T se hace necesario realizar algunas suposiciones sobre la forma en que se han generado los datos aleatorios, lo que motiva las siguientes definiciones: Definici´ on 1.1.2. : Una serie de tiempo {Xt : t ∈ T }, es estrictamente estacionaria si FXt1 ,Xt2 ,...,Xtn = FXt1 +k ,Xt2 +k ,...,Xtn +k para todo k, n ∈ N y para todo t1 , t2 , ..., tn , t1 + k, t2 + k, ..., tn + k en el conjunto de ´ındices T , donde FXt1 ,Xt2 ,...,Xtn denota la funci´ on de distribuci´ on conjunta de las variables Xt1 , Xt2 , ..., Xtn . Esta definici´ on nos dice que la distribuci´on conjunta de dos o m´as variables que componen la serie depende solo de la distancia entre los elementos del subconjunto de T ; sin embargo, debido a la dificultad de verificar esta condici´ on, la definici´on anterior se relaja del siguiente modo: Definici´ on 1.1.3. : Se dice que una serie de tiempo {Xt : t ∈ T } es estacionaria si se verifica que : a) E[Xt ] = µ < ∞ para todo t ∈ T 1 b) cov(Xt , Xt+k ) = E[(Xt − µ)(Xt+k − µ)] = γk < ∞ para todo t, t + k ∈ T Puede verificarse que si una serie de tiempo es estacionaria entonces es tambi´en estrictamente estacionaria, pero no viceversa (ver ejercicio 1.3.1). Enunciaremos las definiciones que nos acompa˜ naran de ahora en adelante: Definici´ on 1.1.4. : Dada una serie de tiempo estacionaria {Xt : t ∈ T } se llama funci´ on de autocovarianza te´ orica a la funci´ on γ : Z → R que a cada entero k le asigna el valor γk = cov(Xt , Xt+k ) = E[(Xt − µ)(Xt+k − µ)] Definici´ on 1.1.5. : Dada una serie de tiempo estacionaria {Xt : t ∈ T } se llama funci´ on de autocorrelaci´ on te´ orica a la funci´ on ρ : Z → R que a cada entero k le asigna el valor ρk = corr(Xt , Xt+k ) = √ cov(Xt ,Xt+k ) var(Xt )var(Xt+k ) = γk γ0 Esta u ´ltima tiene la ventaja de no estar influenciada por las unidades de medida que se usan para obtener los valores de la serie de tiempo. Algunas de las propiedades de las funciones de autocovarianza y de autocorrelaci´ on se muestran a continuaci´on: a) |γk | ≤ γ0 = var(Xt ) y |ρk | ≤ ρ0 = 1 ∀k ∈ Z b) γk = γ−k y ρk = ρ−k c) Las funciones de autocovarianza y autocorrelaci´on son semidefinidas positivas, i.e. para todo n ∈ N y para todo t1 , t2 , ..., tn ∈ T se verifica que: n X n X ai aj γ|ti −tj | ≥ 0 y i=1 j=1 n X n X ai aj ρ|ti −tj | ≥ 0 i=1 j=1 Ejemplo 1.1.1. : Dada una serie de tiempo {At : t ∈ Z} se dice que la serie es un proceso de ruido blanco si las variables aleatorias At son no correlacionadas entre si, con media µt = E[At ] y varianza 2 = var(A ) para todo t ∈ Z σA t Ejemplo 1.1.2. : Sea {At : t ∈ Z} un proceso de ruido blanco. Se dice que la serie de tiempo {Xt : t ∈ Z} es un proceso autorregresivo de orden 1 (o AR(1)) si se puede representar de la siguiente forma: Xt = φXt−1 + At para todo t ∈ Z siendo φ un n´ umero real fijo. La serie del ejemplo 1.1.2 deja en evidencia una dependencia lineal entre Xt y Xt+k . Es decir que los valores obtenidos por la varible X en el pasado influyen en el valor de X en el tiempo presente. Esta influencia se observa de forma directa entre Xt+1 y Xt+2 pues estan relacionadas de la forma Xt+2 = φXt+1 + At+2 y como a su vez Xt+1 = φXt + At+1 se ve que Xt influye de forma indirecta sobre Xt+2 . Considerando la situaci´ on general de una serie de tiempo {Xt : t ∈ T } nos gustar´ıa saber cual es la dependencia lineal directa entre dos variables Xt y Xt+k . Para ello se define el coeficiente de autocorrelaci´ on parcial de orden k como una medida de la relaci´on lineal entre observaciones separadas k per´ıodos con independencia de los valores intermedios. 2 Definici´ on 1.1.6. Se define el coeficiente de correlaci´ on parcial πk como πk = corr(Xt , Xt+k /Xt+1 , ..., Xt+k−1 ) y que es igual a 1 ρ1 ρ2 ρ1 1 ρ 1 .. .. .. . . . ρk−1 ρk−2 ρk−3 πk = 1 ρ1 ρ2 ρ1 1 ρ1 .. .. .. . . . ρk−1 ρk−2 ρk−3 . . . ρk−2 ρ1 . . . ρk−3 ρ2 .. .. .. . . . . . . ρ1 ρk . . . ρk−2 ρk−1 . . . ρk−3 ρk−2 .. .. .. . . . . . . ρ1 1 Dicha f´ormula se obtiene de realizar una regresi´on lineal de Xt+k respecto de Xt+k−1 , Xt+k−2 ,...,Xt , o sea: Xt+k = φk0 + φk1 Xt+k−1 + φk2 Xt+k−2 + ... + φkk Xt + Et+k donde {Et } es un proceso de ruido blanco y se supone que las variables Et+k no est´an correlacionadas con Xt , Xt+1 , ..., Xt+k−1 para todo t. Luego, resulta que el coeficiente φkk de Xt es el coeficiente de correlaci´on parcial te´ orico πk entre Xt+k y Xt . Definici´ on 1.1.7. Se denomina funci´ on de autocorrelaci´ on parcial te´ orica de la serie de tiempo estacionaria {Xt : t ∈ Z} a la funci´ on π :Z → R que a cada entero k le asigna el coeficiente de correlaci´ on parcial te´ orico πk entre Xt+k y Xt . 1.2. Estimaci´ on de las caracter´ısticas de una serie de tiempo estacionaria Los elementos que caracterizan a una serie de tiempo estacionaria son su media µ, su varianza σ 2 , sus autocorrelaciones ρk y sus autocorrelaciones parciales πk . Es necesario poder realizar estimaciones de estos elementos ya que en muchos casos no se puede contar con m´as de una realizaci´on de la serie de tiempo. 1.2.1. Estimaci´ on de la media µ Dada una serie de tiempo {Xt : t ∈ Z}, si se conocen solo n valores de una realizaci´on x1 , ..., xn , el estimador natural de la media com´ un es la media muestral n 1X Xn = Xt n t=1 ya que ´este es un estimador insesgado para µ y si adem´as se cumple que l´ım ρk = 0 entonces X n es k→∞ un estimador consistente en cuadrados medios para µ (ver Leiva [1]) 1.2.2. Estimaci´ on de las autocovarianzas Dadas n observaciones x1 , ..., xn de la serie de tiempo estacionaria, se puede considerar para un k fijo (0 ≤ k ≤ n − 1) los siguientes n − k pares: (x1 , x1+k ), (x2 , x2+k ), ..., (xn−k , xn ) y mirando a la 3 primera entrada de cada par como una observaci´on de una variable y a la segunda como la observaci´ on de otra variable, se puede calcular una estimaci´on de γk mediante la f´ormula: n−k 1 X (xt − x(1) )(xt+k − x(2) ) n−k t=1 donde x(1) = n−k 1 X xt n−k y n X 1 xt . n−k x(2) = t=1 t=k+1 Pn Como en general x(1) ∼ = x(2) se reemplazan ambas por xn = n1 t=1 xt y se obtiene el siguiente estimador para γk : n−k 1X γ b= (Xt − X n )(Xt+k − X n ) n t=1 Definici´ on 1.2.1. La funci´ on γ b : Z → R que a cada k ∈ Z le asigna el valor γ bk recibe el nombre de funci´ on de autocovarianza muestral. 1.2.3. Estimaci´ on de las autocorrelaciones Un estimador natural de ρk es γ bk ρbk = = γ b0 Pn−k t=1 (Xt − X n )(Xt+k − X n ) Pn 2 t=1 (Xt − X n ) k = 0, ±1, ..., ±(n − 1) para (1.2.1) Definici´ on 1.2.2. La funci´ on ρb : Z → R que a cada k ∈ Z le asigna el valor ρbk recibe el nombre de funci´ on de autocorrelaci´ on muestral, y su gr´ afico a partir de la muestra x1 , ..., xn es llamado correlograma muestral. 1.2.4. Estimaci´ on de las autocorrelaciones parciales Los estimadores naturales de πk = φkk se obtienen reemplazando en la f´ormula que expresa πk como cociente de determinantes, a los ρj por sus estimadores ρbj . Para valores muy grandes de k es conveniente usar un m´etodo recursivo para su c´alculo, comenzando por π b1 = φb11 = ρb1 y luego π bk+1 = φbk+1,k+1 = ρbk+1 − Pk b Pk b 1− bk+1−j j=1 φk,j ρ bj j=1 φk,j ρ (1.2.2) donde, φbk,j = φbk−1,j − φbk−1,k−j π bk 1.3. j = 1, ..., k Ejercicios del Cap´ıtulo I Ejercicio 1.3.1. Sean {Xt : t ∈ N ∪ {0}} variables aleatorias independientes que cumplen ( P (Xt = 1) = P (Xt = −1) = 1/2, si t es par; P (Xt = 3) = P (Xt = −3) = 1/18, si t es impar (a) ¿Es {Xt } un proceso estrictamente estacionario? 4 (1.2.3) (b) ¿Es estacionario en sentido amplio? Soluci´ on. (a) Sean n = 1 = k. Si t es par FXt (−1) = P (Xt ≤ −1) = P (Xt = −1) = 1/2 pero FXt+1 (−1) = P (Xt+1 ≤ −1) = P (Xt+1 = −3) = 1/18 por lo que FXt 6= FXt+1 y resulta que {Xt } no es estrictamente estacionario. (b) Si t es par, E[Xt ] = (−1).P (Xt = −1) + 1.P (Xt = 1) = −1/2 + 1/2 = 0 y si t es impar se puede ver del mismo modo que E[Xt ] = 0. Adem´as, si t es par V ar[Xt ] = (−1)2 .P (Xt = −1) + 12 .P (Xt = 1) = 1/2 + 1/2 = 1 (igual se ve para t impar) y si h 6= 0, cov(Xt , Xt+h ) = 0 pues Xt , Xt+h son independientes. Por lo tanto {Xt } es estacionario en sentido amplio. Ejercicio 1.3.2. Considere la serie {Xt : t ∈ Z} definida por Xt = θAt−1 + At , donde {At : t ∈ Z} es 2. un proceso de ruido blanco con varianza σA (a) ¿Qu´e tipo de proceso es {Xt }? (b) ¿Cu´ al es el valor de γk = cov(Xt , Xt+k )? (c) ¿Cu´ al es el valor de ρk = corr(Xt , Xt+k )? 2 =4 (d) Si θ = 0,7 y si σA i) ¿Cu´ al es el valor de γk ? ii) Grafique la funci´ on de autocovarianza te´ orica para k = 0, 1, 2, 3, 4, 5. 2 =4 (e) Si θ = −0,9 y si σA i) ¿Cu´ al es el valor de ρk ? ii) Grafique la funci´ on de autocorrelaci´ on te´ orica para k = 0, 1, 2, 3, 4, 5. Soluci´ on. (a) Es un proceso de medias m´oviles de orden 1 ´o MA(1). (b) γk = cov(Xt , Xt+k ) = E[Xt Xt+k ] = E[(θAt−1 + At )(θAt+k−1 + At+k )] = θ2 E[At−1 At+k−1 ] + θ(E[At−1 At+k ] + E[At At+k−1 ]) + E[At At+k ] = θ2 δk + θδk+1 + θδk−1 + δk donde δk es la autocovarianza de {At }, o sea ( 2 , si k = 0; σA δk = 0, si k 6= 0 entonces  2 2  (θ + 1)σA , 2, γk = θσA   0, 5 si k = 0; si |k| = 1; si |k| > 1 (c)   si k = 0; 1, γk 2 = θ/(θ + 1), si |k| = 1; ρk =  γ0  0, si |k| > 1 (d)   5,96, si k = 0; γk = 2,8, si |k| = 1;   0, si |k| > 1 (e)   si k = 0; 1, ρk = 0,497, si |k| = 1;   0, si |k| > 1 Ejercicio 1.3.3. Considere la serie {Xt : t ∈ Z} definida por Xt = φXt−1 + At , donde {At : t ∈ Z} 2. es un proceso de ruido blanco con varianza σA (a) ¿Qu´e tipo de proceso es {Xt }? (b) ¿Cu´ al es el valor de γk = cov(Xt , Xt+k )? (c) ¿Cu´ al es el valor de ρk = corr(Xt , Xt+k )? 2 =4 (d) Si φ = 0,7 y si σA i) ¿Cu´ al es el valor de γk ? ii) Grafique la funci´ on de autocovarianza te´ orica para k = 0, 1, 2, 3, 4, 5. 2 =4 (e) Si φ = −0,9 y si σA i) ¿Cu´ al es el valor de ρk ? ii) Grafique la funci´ on de autocorrelaci´ on te´ orica para k = 0, 1, 2, 3, 4, 5. Soluci´ on. (a) Es un proceso autorregresivo de orden 1 ´o AR(1). (b) Xt = φXt−1 + At = φ(φXt−2 + At−1 ) + At = φ2 Xt−2 + φXt−1 + At = ... = φk Xt−k + k−1 X φi At−i i=0 6 Si |φ| < 1 entonces E[(Xt − Luego, Pk−1 i=0 φi At−i )2 ] = E[Xt−k ]φk ≤ M φk −→ 0 si k −→ ∞. ∞ X Xt = φi At−i i=0 y E[Xt ] = ∞ X φi E[At−i ] = 0 i=0 Entonces, ∞ ∞ X X γk = cov(Xt , Xt+k ) = E[Xt Xt+k ] = E[( φi At−i )( φj At+k−j )] i=0 = = (c) Como ρ0 = 2 σA 1−φ2 ∞ X ∞ X φ i+j E[At−i At+k−j ] = i=0 j=0 ∞ X 2 φk σA φ2i i=0 ∞ X φ 2i+k j=0 E[A2t−i ] i=0 = = ∞ X 2 φ2i+k σA i=0 2 φk σ A 1 − φ2 entonces ρk = φk . (d) γk = 7,843(0,7)k (e) ρk = (−0,9)k Ejercicio 1.3.4. Los siguientes n´ umeros corresponden a 20 observaciones sucesivas de un proceso estoc´ astico estacionario {Xt : t ∈ N}. -0.27 0.18 -1.14 -0.86 0.47 2.09 0.63 0.48 -1.09 -0.33 0.78 -1.81 a) Grafique las observaciones en funci´ on del tiempo. 7 -0.19 0.76 1.45 -0.99 -0.94 0.45 -0.52 -0.11 b) Observando el gr´ afico de la parte (a), trate de adivinar el signo y el valor absoluto (aproximado) del coeficiente de autocorrelaci´ on de orden 1 del proceso estoc´ astico que gener´ o estos datos. c) Grafique xt versus xt+1 y observando esta nube de puntos trate nuevamente de adivinar el valor de la autocorrelaci´ on de orden 1. d) Calcule usando los datos el valor estimado ρb1 . e) Grafique xt versus xt+2 y observando esta nube de puntos trate nuevamente de adivinar el valor de la autocorrelaci´ on de orden 2. f ) Calcule usando los datos el valor estimado ρb2 . g) Calcule con los datos los valores de yt y zt siendo yt = xt − ρb1 xt+1 , zt = xt+2 − ρb1 xt+1 y grafique yt versus zt . h) Teniendo en cuenta el gr´ afico de la parte (g) adivine el valor de la autocorrelaci´ on parcial de orden 2 (π2 ) del proceso estoc´ astico que gener´ o los datos. i) Calcule con los datos una estimaci´ on de π2 . Soluci´ on. a) Serie de tiempo {Xt } b) No se puede decir nada sobre el coeficiente de autocorrelaci´on de orden 1 s´olo mirando el gr´ afico de la parte (a). Necesitar´ıamos observar como se comportan las estimaciones de la funci´ on de autocorrelaci´ on. c) xt versus xt+1 Se puede observar una tendencia lineal con pendiente negativa por lo que ρ1 < 0 8 d) De acuerdo a (1.2.1), ρb1 = −0,343 e) xt versus xt+2 No se puede apreciar ninguna tendencia por lo que ρ2 = 0. f) De acuerdo a (1.2.1), ρb2 = 0,03258 g) vt = xt − ρb1 xt+1 versus wt = xt+2 − ρb1 xt+1 h) Parece que no existe una tendencia lineal entre los puntos graficados por lo que π2 debe ser aproximadamente 0. i) De acuerdo a (1.2.2) y (1.2.3), π2 = −0,0964. Ejercicio 1.3.5. Considere la serie de tiempo {Xt : t ∈ Z} definida por Xt = 0,9At−1 + At , donde 2 = 4. {At : t ∈ Z} es un proceso de ruido blanco con varianza σA a) ¿Qu´e tipo de proceso es {Xt }?. b) ¿Cu´ al es la forma general (con el mismo orden) del tipo de proceso que identific´ o en la parte (a)?. c) ¿Cu´ al es el valor de su media? ¿Por qu´e?. d) ¿Cu´ al es el valor de γk ?. e) ¿Cu´ al es el valor de ρk ?. f ) Dados 50 datos x1 , ..., x50 generados por este proceso de ecuaci´ on Xt = 0,9At−1 + At , si se grafican los pares de puntos (xt−1 , xt ), t = 2, ..., 50, ¿Cu´ al es la ecuaci´ on de la recta te´ orica alrededor de la cual esta nube de puntos graficados tiende a concentrarse?. 9 g) Si se grafican los puntos (xt−2 , xt ), t = 3, ..., 50, ¿Cu´ al es la tendencia de esta nube de puntos?. 0,9 h) Si se grafican los puntos (yt , zt ) donde yt = xt − 1+0,9 2 xt+1 y zt = xt+2 − t = 1, ..., 48, ¿Cu´ al es la tendencia de esta nube de puntos?. Soluci´ on. 0,9 x 1+0,92 t+1 para a) {Xt } es un proceso de medias m´oviles de orden 1 (MA(1)) b) La forma general es Xt = θAt−1 + At . c) E[Xt ] = 0,9E[At−1 ] + E[At ] = 0,9 ∗ 0 + 0 = 0. d) Por lo visto en el ejercicio 1.3.2,   2   (0,9 + 1)4, si k = 0; 7,24, γk = 0,9 ∗ 4, si |k| = 1; = 3,6,     0, si |k| > 1 0, si k = 0; si |k| = 1; si |k| > 1 e)   1, ρk = 0,9/(1 + 0,92 ),   0,   si k = 0; 1, si |k| = 1; = 0,497,   si |k| > 1 0, si k = 0; si |k| = 1; si |k| > 1 f) La ecuaci´ on de la recta es z = 0,49y + 0,26. g) La nube de puntos no parece concentrarse alrededor de ninguna recta y esto es coherente con el echo de que ρ2 = 0 en un proceso MA(1). h) Los puntos siguen teniendo una tendencia lineal pero menos marcada que en los puntos anteriores y obedece a la ecuaci´ on w = −0,4v + 0,35. 10 Ejercicio 1.3.6. Considere la serie de tiempo {Xt : t ∈ Z} definida por Xt = 0,9Xt−1 + At , donde 2 = 4. {At : t ∈ Z} es un proceso de ruido blanco con varianza σA a) ¿Qu´e tipo de proceso es {Xt }?. b) ¿Cu´ al es la forma general (con el mismo orden) del tipo de proceso que identific´ o en la parte (a)?. c) ¿Cu´ al es el valor de su media? ¿Por qu´e?. d) ¿Cu´ al es el valor de γk ?. e) ¿Cu´ al es el valor de ρk ?. f ) Dados 50 datos x1 , ..., x50 generados por este proceso de ecuaci´ on Xt = 0,9Xt−1 + At , si se grafican los pares de puntos (xt−1 , xt ), t = 2, ..., 50, ¿Cu´ al es la ecuaci´ on de la recta te´ orica alrededor de la cual esta nube de puntos graficados tiende a concentrarse?. g) Si se grafican los puntos (xt−2 , xt ), t = 3, ..., 50, ¿Cu´ al es la tendencia de esta nube de puntos?. h) Si se grafican los puntos (yt , zt ) donde yt = xt − 0,9xt+1 y zt = xt+2 − 0,9xt+1 para t = 1, ..., 48, ¿Cu´ al es la tendencia de esta nube de puntos?. Soluci´ on. a) {Xt } es un proceso autorregresivo de orden 1 (AR(1)). b) La forma general es Xt = φXt−1 + At . c) De acuerdo con las cuentas hechas en el ejercicio 1.3.3, el proceso AR(1) se puede escribir como un proceso de medias m´ oviles de orden infinito y entonces E[Xt ] = 0. k 0,9 k d) Por el ejercicio 1.3.3, γk = 4 1−0,9 k = 21,05(0,9) . e) ρk = 0,9k f) La ecuaci´ on de la recta es z = 0,8y + 0,16. 11 g) La ecuaci´ on de la recta es z2 = 0,5y2 + 0,3. h) Los puntos no parecen mostrar ninguna tendencia lineal lo cual concuerda con el echo de que π2 = 0 en un proceso AR(1). Ejercicio 1.3.7. En cada uno de los siguientes casos grafique los valores dados en la tabla y especifique los modelos de serie de tiempo estacionaria {Xt } a los que podr´ıa pertenecer los pares de funciones de autocorrelaci´ on te´ orica y de autocorrelaci´ on parcial te´ orica graficadas. Establezca la relaci´ on que define en forma general a Xt en el modelo identificado y, cuando sea posible, determine el signo de cada uno de los coeficientes involucrados. a) k ρk πk 1 0.800 0.800 2 0.640 0 b) k ρk πk 1 0.497 0.497 2 0 -0.328 3 0.512 0 3 0 0.243 4 0.410 0 4 0 -0.191 5 0.328 0 6 0.262 0 5 0 0.156 6 0 -0.131 Soluci´ on. a) De acuerdo a los gr´ aficos se puede decir que los valores corresponden a un modelo AR(1) de la forma Xt = φXt−1 + At ya que en esos modelos la autocorrelaci´on te´orica es ρk = φk y la autocorrelaci´ on parcial te´ orica es πk = 0 si k 6= 1 y π1 = φ. En este caso debe ser φ > 0. 12 b) De acuerdo a los gr´ aficos se puede decir que los valores corresponden a un modelo MA(1) de la forma Xt = θAt−1 + At ya que en esos modelos la autocorrelaci´on te´orica es ρk = 0 si |k| > 1 y θ ρk = 1+θ 2 si |k| = 1. En este caso debe ser θ > 0. 13 Cap´ıtulo 2 Modelos para series de tiempo estacionarias En esta secci´ on se presentar´ an los principales modelos que pueden regir el comportamiento de una serie de tiempo estacionaria junto con las caracter´ısticas que permiten distinguir a cu´al de esos modelos corresponde una serie de datos observados. 2.1. 2.1.1. Representaciones de una serie de tiempo Representaci´ on como un proceso de promedios m´ oviles de orden infinito Un resultado importante referido a los procesos estacionarios es el probado por Wold en 1938; este dice que todo proceso estacionario {Xt : t ∈ Z} que sea no determin´ıstico puro puede expresarse como un promedio m´ ovil de orden infinito, es decir Xt = µ + ∞ X αi At−i i=0 2 y tal que con α0 = 1, {At : t ∈ Z} proceso de ruido blanco con media 0 y varianza σA P∞ 2 i=0 αi < ∞. Definici´ on 2.1.1. Se llama operador de cambio hacia atr´ as y se denota por B al operador definido por BXt = Xt−1 . Notaci´ on. Por simplicidad de usar´ a la notaci´on X˙ t = Xt − µ. P∞ Proposici´ on 1. Si {Xt } un proceso de medias m´ oviles de orden infinito Xt = µ + P i=0 αi At−i con ∞ 2 y tal que 2 α0 = 1, {At : t ∈ Z} proceso de ruido blanco con media 0 y varianza σA i=0 αi < ∞ entonces {Xt } es un proceso estacionario. Demostraci´ on. Haciendo uso del operador de cambio hacia atr´as, se puede escribir a Xt como Xt = µ + ∞ X αi B i At = µ + α(B)At i=0 donde α(B) = Es claro que P∞ i=0 αi B i y α0 = 1, o sea X˙ t = α(B)At . E[Xt ] = µ + ∞ X αi E[At−i ] = µ + 0 = µ, i=0 14 V ar(Xt ) = E[Xt2 ] − µ2 y E[Xt2 ] = E[(µ + ∞ X αi At−i )(µ + i=0 = 2µ ∞ X αi E[At−i ] + µ2 + i=0 ∞ X αj At−j )] j=0 ∞ X 2 αi αi E[At−j At−i ] = µ2 + σA i,j=0 ∞ X αi2 i=0 entonces 2 V ar(Xt ) = σA ∞ X αi2 (2.1.1) i=0 y γk = cov(X˙ t , X˙ t+k ) ∞ ∞ X X = E[( αi At−i )( αj At−j )] = = i=0 ∞ ∞ XX (2.1.2) j=0 αi αj E[At−i At+k−j ] i=0 j=0 ∞ X 2 σA αi αi+k < ∞ i=0 Como |γk | ≤ γ0 = 2 σA ∞ X αi2 < ∞, i=0 se cumplen las condiciones de la definici´on 1.1.3 por lo que resulta que {Xt } es un proceso estacionario. 2.1.2. Representaci´ on como un proceso autorregresivo de orden infinito Definici´ on 2.1.2. Se dice que una serie de tiempo {Xt : t ∈ Z} es invertible o equivalentemente es un proceso autorregresivo de orden infinito, si X˙ t = Xt − µ se puede expresar como X˙ t = At + ∞ X βi X˙ t−i i=1 con ∞ X O sea, β(B)X˙ t = At |βi | < ∞. i=1 donde β(B) = ∞ X βi B i . i=1 Proposici´ on 2. Sea {Xt } un proceso invertible como en la definici´ on 2.1.2. Si las ra´ıces de β(B) = 0 son de m´ odulo mayor que uno entonces la serie resulta ser un proceso estacionario. Demostraci´ on. Como β(x) 6= 0 cuando |x| ≤ 1, existe un ε > 0 tal que series de potencia ∞ X 1 = αj xj = α(x), |x| < 1 + ε. β(x) j=0 15 1 β(x) tiene una expansi´ on en Si (1 + 2ε ) < |x| < (1 + ε), como |αj ||x|j → 0 cuando j → ∞ resulta que αj (1 + ε/2)j → 0 cuando j → ∞ y por lo tanto, existe K ∈ (0, ∞) tal que |αj | < K(1 P + ε/2)−j ∀j = 0, 1, 2, ... En particular, P ∞ 2 un, ∞ j=0 |αj | < ∞ y α(x)β(x) ≡ 1 para todo |x| ≤ 1. Mas a´ j=0 αj < ∞. Luego, resulta que {Xt } posee una expresi´ on como un proceso de medias m´oviles de orden infinito y por la proposici´on 1, {Xt } es estacionario. An´alogamente, un proceso estacionario X˙ t = α(B)AtP ser´a invertible si las ra´ıces de α(B) = 0 son 1 de m´odulo mayor que 1 pues entonces α(B) = β(B) con ∞ i=1 |βi | < ∞. 2.2. Modelos autoregresivos de orden finito Definici´ on 2.2.1. Se dice que una serie de tiempo {Xt : t ∈ Z} es un proceso autorregresivo de orden p (se denota AR(p)), con p ∈ N fijo, si X˙ t = Xt − µ se puede expresar como X˙ t = At + p X φi X˙ t−i ∀t ∈ Z i=1 2. donde {At : t ∈ Z} es un proceso de ruido blanco con media 0 y varianza σA Al ser un caso particular de un proceso autorregresivo de orden infinito, el proceso ser´a estacionario si las ra´ıces de Φp (B) = 1 − φ1 B − φ2 B 2 − ... − φp B p = 0 son de m´odulo mayor que 1. La siguiente definici´ on se usar´ a en la demostraci´on de algunos resultados que se presentar´an a continuaci´on pero cuyas pruebas no se incluir´an en este trabajo. Definici´ on 2.2.2. Se llama ecuaci´ on lineal en diferencias con coeficientes constantes de orden k, a una ecuaci´ on de la forma yt + α1 yt−1 + ... + αk yt−k = at para t ∈ T, donde at : t ∈ T es una funci´ on real de t que no depende de y, tal que α1 , ..., αk son constantes con αk 6= 0 y donde T es un intervalo de n´ umeros enteros (o sea T ⊂ Z) que sin p´erdida de generalidad podemos suponer que es de la forma [z, ∞) o [z, z + h] o (−∞, ∞). 2.2.1. Modelo autoregresivo de orden 1 Un proceso autorregresivo de orden 1 se expresa de la forma Xt − µ = φ1 (Xt−1 − µ) + At o bien X˙ t = φ1 X˙ t−1 + At con X˙ t = Xt − µ. Entonces E[Xt ] = µ, φk1 2 γk = cov(Xt , Xt+k ) = cov(X˙ t , X˙ t+k ) = σA si |φ1 | < 1 1 − φ21 y ( φ1 , πk = 0, 16 si k = 1; si k > 1 (ver ejercicio 1.3.3). Adem´as, este proceso ser´ a estacionario si la ra´ız de Φ(B) = 1 − φ1 (B) = 0 es de m´odulo mayor que 1, o equivalentemente |φ1 | < 1. Los c´alculos de arriba se pueden realizar usando ecuaciones en diferencias: la soluci´on de la ecuaci´ on en diferencias de primer orden Xt − φXt−1 = At es la suma de una soluci´on particular y de una soluci´ on complementaria que es soluci´ on de la ecuaci´on Xt − φXt−1 = 0. Una soluci´on particular es Xt = ∞ X φi1 At−i i=0 y una soluci´on complementaria es de la forma cφt1 . Por lo que la soluci´on general es Xt = cφt1 + ∞ X φi1 At−i . i=0 Para calcular las autocorrelaciones de Xt − µ = φ1 (Xt−1 − µ) + At se multiplica a ambos miembros de la ecuaci´on por (Xt−k − µ) y se obtiene γk = φ1 γk−1 . Por lo tanto, la funci´on de autocorrelaci´ on saldr´a de resolver la ecuaci´ on ρk − φ1 ρk−1 = 0. Esta tiene soluci´on general de la forma ρk = abk donde b es soluci´on de b − φ1 = 0, o sea b = φ1 y como 1 = ρ0 = aφ01 = a entonces ρk = φk1 . 2.2.2. Modelo autoregresivo de orden p La definici´ on de que un proceso {Xt : t ∈ Z} sea autorregresivo de orden p es equivalente a que sea soluci´on de la ecuaci´ on X˙ t − φ1 X˙ t−1 − φ2 X˙ t−2 − ... − φp X˙ t−p = At donde {At : t ∈ Z} es un proceso de ruido blanco. La soluci´on general ser´ a suma de una soluci´on particular Xt∗ y de una soluci´on complementaria Xt∗∗ . Si z1 , ..., zp son las ra´ıces de φ(B) = 1 − φ1 B − φ2 B 2 − ... − φp B p = 0 entonces la soluci´on particular es Xt∗ = p p Y X 1 1 ki At = ( )At = ( )At φ(B) (1 − zi B) (1 − zi B) i=1 i=1 donde las constantes ki se eligen de tal forma que se verifique la igualdad, y si adem´as |zi | > 1 para i = 1, ..., p resulta p p ∞ ∞ X X X X j j ∗ zi B ))At = ( ki zij )At−j Xt = ( (ki i=1 j=0 i=1 j=0 Xt∗ verific´andose que es estacionaria. La soluci´on complementaria ser´a una suma de p t´erminos que depender´an del tipo de ra´ıces z1 , ..., zp que tenga su polinomio caracter´ıstico z p − φ1 z p−1 − φ2 z p−2 − ... − φp−1 z − φp . 1. Por cada ra´ız real y distinta zj hay un t´ermino de la forma cj zjt . 2. Por cada ra´ız real zi repetida mi veces hay un t´ermino de la forma (ci1 + ci2 t + ... + cimi tmi −1 )zit . 3. Por cada par de ra´ıces complejas conjugadas hay un t´ermino de la forma Rt (c1 cos(tλ)+c2 sen(tλ)). 4. Por cada par de ra´ıces complejas conjugadas repetidas h veces hay un t´ermino de la forma Rt [(c11 cos(tλ) + c12 sen(tλ)) + t(c21 cos(tλ) + c22 sen(tλ)) + ... + th−1 (ch1 cos(tλ) + ch2 sen(tλ))]. 17 (ver Leiva [1]). Para calcular la funci´ on de autocorrelaci´on de (Xt − µ) = φ1 (Xt−1 − µ) + φ2 (Xt−2 − µ) + ... + φp (Xt−p − µ) + At se multiplica a ambos miembros por Xt−k − µ y se divide todo por γ0 , obteniendose ρk = φ1 ρk−1 + φ2 ρk−2 + ... + φp ρk−p Particularizando esta u ´ltima ecuaci´ on para k = 1, ..., p se obtiene el conocido sistema de ecuaciones de Yule-Walker: ρ¯ = R.φ¯ donde φ¯t = [φ1 ...φp ], ρ¯t = [ρ1 ...ρp ] y   1 ρ1 . . . ρp−1  .. ..  .. R =  ... . . .  ρp−1 ρp−2 . . . 1 La soluci´on ser´ a una mezcla de exponenciales, debido a los t´erminos con ra´ıces reales, y sinusoidales, debidas a las ra´ıces complejas conjugadas. La funci´on de autocorrelaci´ on parcial se anula a partir del lag p+1 pues el determinante del numerador en la definici´ on 1.1.6 se hace cero por ser la u ´ltima columna combinaci´on lineal de las p primeras. 2.3. Modelos de promedios m´ oviles de orden finito Definici´ on 2.3.1. Se dice que {Xt : t ∈ Z} es un proceso de medias m´ oviles de orden q y se denota ˙ por MA(q), si Xt = Xt − µ se puede expresar como X˙ t = At + q X θi At−i ∀t ∈ Z i=1 2. donde {At : t ∈ Z} es un proceso de ruido blanco con media 0 y varianza σA Un proceso MA(q) con q fijo siempre es estacionario pues P∞ 2 i=0 αi = 1 + θ12 + ... + θq2 < ∞ Como X˙ t = θ(B)At donde θ(B) = 1 + θ1 B + ... + θq B q , la condici´on suficiente para que el proceso sea invertible es que las ra´ıces de θ(B) = 0 sean de modulo mayor que 1 (la demostraci´on es similar a la de la proposici´ on 2). Adem´as, como de vio en la demostraci´on de la proposici´on 1 (ver ecuaciones (2.1.1) y (2.1.2)) pero con αi = θi para i = 1, ..., q y αi = 0 si i > q resulta que E[Xt ] = µ, 2 V ar(Xt ) = σA (1 + θ12 + ... + θq2 ), ( 2 (θ + θ θ σA 1 k+1 + ... + θq−k θq ), k γk = 0, ( θk +θ1 θk+1 +...+θq−k θq , si k 1+θ12 +...+θq2 ρk = 0, si k si k ≤ q; si k > q ≤ q; >q Como se ve, la autocovarianza se anula para k > q y es una propiedad que caracteriza a los procesos MA(q) (ver Brockwell [2]). 18 2.3.1. El modelo de medias m´ oviles de orden 1 Un proceso de medias m´ oviles de orden 1 se expresa como Xt − µ = At + θ1 At−1 donde {At : t ∈ Z} 2. es un proceso de ruido blanco con varianza σA Como ya se vi´ o en el ejercicio 1.3.2, este proceso es estacionario con E[Xt ] = µ,  2 2  σA (1 + θ1 ), 2θ , γk = σ A 1   0, ( θ1 si 2, ρk = 1+θ1 0, si si k = 0; si k = 1; si k > 1 k = 1; k>1 En la funci´ on de autocorrelaci´ on parcial πk (definici´on 1.1.6), el determinante del numerador es θ1 k k+1 (−1) ( 1+θ2 ) pues 1 num(πk ) = 1 θ1 1+θ12 .. . 0 k+1 θ1 = (−1) 2 1 + θ1 = (−1)k+1 0 ... 0 1 .. . 0 θ1 1+θ12 θ1 1+θ12 ... .. . ... 0 .. . 0 .. . 0 1 θ1 1+θ12 θ1 1+θ12 .. . 0 θ1 1+θ12 1 .. . 0 .. . 0 0 ... 0 θ1 1+θ12 ... .. . ... 0 .. . .. . 0 θ1 1+θ12 θ1 θ1 k−1 θ1 k ( ) = (−1)k+1 ( ) 2 2 1 + θ1 1 + θ1 1 + θ12 mientras que el determinante del denominador es ∆k = usando el hecho de que ∆k = ∆k−1 − Por lo tanto, πk = θ1 1+θ12 θ1 2 ( 1+θ 2 ) ∆k−2 1 1+θ12 +...+θ12k (1+θ12 )k que se prueba f´acilmente e inducci´on en k. (−1)k+1 θ1k −(−θ1 )k −(−θ1 )k (1 − θ12 ) = = si k ≥ 1. 2(k+1) 1 + θ12 + ... + θ12k 1 + θ12 + ... + θ12k 1 − θ1 Para que este proceso sea invertible, la ra´ız de θ(B) = 1 + θ1 B = 0 debe ser de modulo mayor que 1, o equivalentemente |θ1 | < 1. 2.4. Modelos Autorregresivos de Promedios M´ oviles Definici´ on 2.4.1. Se dice que la serie de tiempo {Xt : t ∈ Z} es un proceso autorregresivo de medias m´ oviles de ´ ordenes p y q respectivamente y se lo denota por ARMA(p,q), si se verifica que X˙ t = Xt −µ se puede expresar como (1 − φ1 B − φ2 B 2 − ... − φp B p )X˙ t = (1 + θ1 B + θ2 B 2 + ... + θq B q )At , 19 donde φ1 , φ2 , ..., φp , θ1 , θ2 , ..., θq son constantes con φp 6= 0 6= θq y {At : t ∈ Z} un proceso de ruido 2. blanco con varianza σA Una condici´ on suficiente para que el proceso sea estacionario es que las ra´ıces de Φ(B) = 1 − φ1 B − 2 φ2 B − ... − φp B p = 0 sean de m´ odulo mayor que 1 (ver proposici´on 2). Mientras que una condici´ on suficiente para que el proceso sea invertible es que las ra´ıces de θ(B) = 1 + θ1 B + θ2 B 2 + ... + θq B q = 0 sean de m´odulo mayor que 1 (la demostraci´on es similar a la de la proposici´on 2). 2.4.1. Modelo ARMA(1,1). Un proceso ARMA(1,1) puede expresarse como X˙ t = φ1 X˙ t−1 + At + θAt−1 o bien como Φ(B)X˙ t = θ(B)At donde Φ(B) = 1 − φ1 B y θ(B) = 1 + θ1 B. Como se vi´o en la proposici´ on 2, este proceso es estacionario si la ra´ız de Φ(B) = 1 − φ1 B = 0 es de m´odulo mayor que 1, o equivalentemente |φ1 | < 1. En este caso se puede escribir ∞ ∞ X X θ(B) 1 X˙ t = At = (1 + θ1 B) At = (1 + θ1 B)( φi1 B i )At = (1 + (φi1 + θ1 φi−1 1 ))At Φ(B) (1 − φ1 B) i=0 i=1 que es su expresi´ on como media m´ ovil de orden infinito. Similarmente, el proceso ser´ a invertible si la ra´ız de θ(B) = 1 + θ1 B = 0 es de m´odulo mayor que 1, es decir si |θ1 | < 1. En este caso At = ∞ ∞ X X 1 Φ(B) ˙ Xt = (1−φ1 B) X˙ t = (1−φ1 B)( (−θ1 )i B i )X˙ t = (1+ (θ1i−1 (φ1 +θ1 )(−1)i ))X˙ t θ(B) (1 + θ1 B) i=1 i=0 que es su expresi´ on como un proceso autorregresivo de orden infinito. Para calcular la funci´ on de autocovarianza se puede usar su expresi´on como media m´ovil de orden infinito; entonces por (2.1.2) para k 6= 0 2 γk = σ A (αk + ∞ X αi αi+k ) i=1 2 = σA (φ1k−1 (φ1 + θ1 ) + (φ1 + θ1 )2 φk−2 1 ∞ X (φ21 )i ) i=1 2 2 k−2 = σA (φk−1 1 (φ1 + θ1 ) + (φ1 + θ1 ) φ1 2 = σA φ21 1 − φ21 (φ1 + θ1 )φ1k−1 (1 + φ1 θ1 ) 1 − φ21 y si k = 0 γ0 = 2 σA ∞ X i=0 αi2 = 2 σA (1 + ∞ X 2(i−1) (φ1 + θ1 )2 φ1 2 ) = σA (1 + i=1 2 (φ1 + θ1 )2 2 1 + 2φ1 θ1 + θ1 ) = σ ( ) A 1 − φ21 1 − φ21 Por lo tanto, γk =  2 2 ( 1+2φ1 θ1 +θ1 ), σA 2 1−φ 1 k−1 2 (φ1 +θ1 )φ1 (1+φ1 θ1 ) A 1−φ21 σ 20 si k = 0; , si k ≥ 1; y ρk = (φ1 + θ1 )φ1k−1 (1 + φ1 θ1 ) si k ≥ 1 1 + 2φ1 θ1 + θ12 Los valores de |πk | decrecer´ an de forma exponencial dependiendo del signo y magnitud de φ1 y θ1 . Los siguientes gr´ aficos representan las posibles funciones de autocorrelaci´on y autocorrelaci´on parcial de un proceso ARMA(1,1). φ1 > θ1 > 0 θ1 > φ1 > 0 θ1 < φ1 < 0 φ1 < θ1 < 0 21 φ1 < 0 < θ1 y φ1 + θ1 > 0 φ1 < 0 < θ1 y φ1 + θ1 < 0 θ1 < 0 < φ1 y φ1 + θ1 < 0 θ1 < 0 < φ1 y φ1 + θ1 > 0 2.4.2. Modelo ARMA(p,q) Como se dijo al comienzo de la secci´on 2.4, si el m´odulo de las ra´ıces de Φ(B) = 1 − φ1 B − φ2 B 2 − ... − φp B p = 0 son mayores que 1 entonces el proceso como un P∞ ser´a estacionario y se podr´ P∞ a expresar 2 proceso de medias m´ oviles de orden infinito Xt = i=0 αi At−i con α0 = 1 y i=0 αi < ∞. Con esta representaci´ on se puede calcular la funci´on de autocovarianza ya que en un proceso MA(∞) 2 γk = σA ∞ X αj αj+k para k ≥ 0. j=0 22 Lo u ´nico que hace falta es saber el valor de los αi . Para ello se puede recurrir al hecho de que ∞ X αi xi = i=0 θ(x) , Φ(x) x2 + ...)(1 − φ 2 p 2 q es decir que (1 + α1 x + α2 1 x − φ2 x − ... − φp x ) = (1 + θ1 x + θ2 x + ... + θq x ) e igualando los coeficientes de las sucesivas potencias de x se calculan los αi . Las funciones de autocorrelaci´ on y autocorrelaci´on parcial de los procesos ARMA combinan propiedades de los procesos AR y MA: en la funci´ on de autocorrelaci´on ciertos primeros valores dependen del orden de la parte MA y despu´es hay un decrecimiento dictado por la parte AR. En la funci´on de autocorrelaci´on parcial los valores iniciales dependen del orden de la parte AR y sigue un decrecimiento debido a la parte MA. Esta estructura compleja hace dif´ıcil poder identificar este modelo en la pr´actica (ver Pe˜ na [3]). El siguiente cuadro resume las caracter´ısticas de las funciones de autocorrelaci´on de un ARMA: Tabla 1 AR(p) M A(q) ARM A(p, q) 2.5. Autocorrelaci´on Muchos coeficientes no nulos que decrecen con el retardo como mezcla de exponenciales y sinoides. q primeros coeficientes no nulos y el resto cero. Decrecimiento hacia cero. Autocorrelaci´on Parcial p primeros coeficientes no nulos y el resto cero. Muchos coeficientes no nulos que decrecen con el retardo como mezcla de exponenciales y sinoides. Decrecimiento hacia cero. Ejercicios del Cap´ıtulo II Ejercicio 2.5.1. Encontrar la funci´ on de autocorrelaci´ on, la funci´ on de autocorrelaci´ on parcial y graficar ambas para k = 1, 2, 3, 4, 5. a) Xt = −0,5Xt−1 + At b) Xt = −0,2Xt−1 + 0,5Xt−2 + At c) Xt = −0,7Xt−1 + 0,1Xt−2 + At Soluci´ on. a) Como |φ1 | = 0,5 < 1 entonces ρk = (−0,5)k y ( −0,5, si k = 1; πk = 0, si k > 1; 23 b) La ecuaci´ on cuadr´ atica correspondiente es z 2 + 0,2z − 0,5 que tiene como ra´ıces a z1 = 0,614 y z2 = −0,814 entonces ρk = c1 (0,614)k + c2 (−0,814)k que con las condiciones iniciales ρ0 = 1 y ρ1 + 0,2 − 0,5ρ−1 = 0 resulta que ρk = 0,289(0,614)k + 0,71(−0,814)k y   φ 1 , si k = 1;   −0,4, si k = 1;  1−φ2 πk = φ2 , si k = 2; si k = 2; = 0,5,     0, si k > 2; 0, si k > 2; c) La ecuaci´ on cuadr´ atica correspondiente es z 2 + 0,7z − 0,1 que tiene como ra´ıces a z1 = 0,121 y z2 = −0,821 entonces ρk = c1 (0,121)k + c2 (−0,821)k que con las condiciones iniciales ρ0 = 1 y ρ1 + 0,7 − 0,1ρ−1 = 0 resulta que ρk = 0,047(0,121)k + 0,953(−0,821)k y   φ 1  −0,777, si k = 1;  1−φ2 , si k = 1;  π k = φ2 , si k = 2; si k = 2; = 0,1,     0, si k > 2; 0, si k > 2; Ejercicio 2.5.2. En cada uno de los siguientes casos grafique los valores dados en la tabla y especifique los modelos de serie de tiempo estacionaria {Xt } a los que podr´ıa pertenecer los pares de funciones de autocorrelaci´ on te´ orica y de autocorrelaci´ on parcial te´ orica graficadas. Establezca la relaci´ on que define en forma general a Xt en el modelo identificado y, cuando sea posible, determine el signo de cada uno de los coeficientes involucrados. a) k ρk πk 1 -0.833 -0.833 b) k ρk πk 1 0.714 0.714 2 0.817 0.4 2 0.657 0.300 3 -0.742 0 3 0.543 0 4 0.698 0 4 0.469 0 5 -0.645 0 5 0.397 0 6 0.602 0 6 0.339 0 24 c) k ρk πk 1 0.789 0.789 d) k ρk πk 1 -0.929 -0.929 2 0.284 -0.9 2 0.807 -0.4 3 -0.284 0 3 -0.678 0 4 -0.682 0 5 -0.767 0 6 -0.537 0 4 0.558 0 5 -0.455 0 6 0.368 0 Soluci´ on. a) De acuerdo a los gr´ aficos se puede decir que los valores corresponden a un modelo AR(2) de la forma Xt = φ1 Xt−1 + φ2 Xt−2 + At ya que en esos modelos la autocorrelaci´on parcial es  φ 1   1−φ2 , si k = 1; πk = φ2 , si k = 2;   0, si k > 2 Luego φ2 = 0,4 y φ1 = −0,5 b) De acuerdo a los gr´ aficos se puede decir que los valores corresponden a un modelo AR(2) de la forma Xt = φ1 Xt−1 + φ2 Xt−2 + At ya que en esos modelos la autocorrelaci´on parcial es  φ 1   1−φ2 , si k = 1; πk = φ2 , si k = 2;   0, si k > 2 Luego φ2 = 0,3 y φ1 = 0,43 25 c) De acuerdo a los gr´ aficos se puede decir que los valores corresponden a un modelo AR(2) de la forma Xt = φ1 Xt−1 + φ2 Xt−2 + At ya que en esos modelos la autocorrelaci´on parcial es  φ 1   1−φ2 , si k = 1; πk = φ2 , si k = 2;   0, si k > 2 Luego φ2 = −0,9 y φ1 = 1,5 d) De acuerdo a los gr´ aficos se puede decir que los valores corresponden a un modelo AR(2) de la forma Xt = φ1 Xt−1 + φ2 Xt−2 + At ya que en esos modelos la autocorrelaci´on parcial es  φ 1   1−φ2 , si k = 1; πk = φ2 , si k = 2;   0, si k > 2 Luego φ2 = −0,4 y φ1 = −1,3 Ejercicio 2.5.3. En cada uno de los siguientes casos grafique los valores dados en la tabla y especifique los modelos de serie de tiempo estacionaria {Xt } a los que podr´ıa pertenecer los pares de funciones de autocorrelaci´ on te´ orica y de autocorrelaci´ on parcial te´ orica graficadas. Establezca la relaci´ on que define en forma general a Xt en el modelo identificado y, cuando sea posible, determine el signo de cada uno de los coeficientes involucrados. a) k ρk πk 1 -0.702 -0.702 2 0.222 -0.534 3 0 -0.376 4 0 -0.180 5 0 0.027 26 6 0 0.169 b) k ρk πk 1 0.489 0.489 2 0.402 0.215 3 0 -0.357 4 0 0.052 c) k ρk πk 1 0.682 0.682 2 0.187 -0.519 3 0 0.409 4 0 -0.315 d) k ρk πk 1 -0.654 -0.654 2 0.156 -0.476 3 0 -0.364 4 0 -0.284 5 0 0.202 6 0 -0.133 5 0 0.220 5 0 -0.224 6 0 -0.124 6 0 -0.175 Soluci´ on. a) De acuerdo a los gr´ aficos se puede decir que los valores corresponden a un modelo MA(2) de la forma Xt = θ1 At−1 + θ2 At−2 + At ya que estos modelos se caracterizan por tener una funci´ on de autocorrelaci´ on que se anula para los k > 2 y es  θ1 +θ1 θ2    1+θ12 +θ22 , si k = 1; ρk = 1+θθ22+θ2 , si k = 2;  1 2  0, si k > 2 Luego θ1 < 0 y θ2 > 0 b) De acuerdo a los gr´ aficos se puede decir que los valores corresponden a un modelo MA(2) de la forma Xt = θ1 At−1 + θ2 At−2 + At ya que estos modelos se caracterizan por tener una funci´ on de autocorrelaci´ on que se anula para los k > 2 y es  θ1 +θ1 θ2   1+θ12 +θ22 , si k = 1;  ρk = 1+θθ22+θ2 , si k = 2;  1 2  0, si k > 2 Luego θ1 > 0 y θ2 > 0 27 c) De acuerdo a los gr´ aficos se puede decir que los valores corresponden a un modelo MA(2) de la forma Xt = θ1 At−1 + θ2 At−2 + At ya que estos modelos se caracterizan por tener una funci´ on de autocorrelaci´ on que se anula para los k > 2 y es  θ1 +θ1 θ2    1+θ12 +θ22 , si k = 1; ρk = 1+θθ22+θ2 , si k = 2;  2 1  0, si k > 2 Luego θ1 > 0 y θ2 > 0 d) De acuerdo a los gr´ aficos se puede decir que los valores corresponden a un modelo MA(2) de la forma Xt = θ1 At−1 + θ2 At−2 + At ya que estos modelos se caracterizan por tener una funci´ on de autocorrelaci´ on que se anula para los k > 2 y es  θ1 +θ1 θ2   1+θ12 +θ22 , si k = 1;  ρk = 1+θθ22+θ2 , si k = 2;  1 2  0, si k > 2 Luego θ1 < 0 y θ2 > 0 28 Ejercicio 2.5.4. Considere la serie de tiempo {Xt } definida por Xt −µ = θAt−1 +At , donde 0 < |θ| ≤ 1 2 . Defina el proceso {Y } como la primera y donde {At } es un proceso de ruido blanco con varianza σA t diferencia del proceso Xt , es decir Yt = Xt − Xt−1 para todo t ∈ Z. a) ¿Qu´e tipo de proceso es {Yt }? b) ¿C´ ual es el valor τk = cov(Yt , Yt+k ) en funci´ on de la funci´ on de autocovarianza γ del proceso Xt ? on de la funci´ on de autocorrelaci´ on ρ del proceso c) ¿C´ ual es el valor ρ∗k = corr(Yt , Yt+k ) en funci´ Xt ? d) ¿Es Yt un proceso invertible? Soluci´ on. a) Yt = Xt − Xt−1 = θAt−1 + At − (θAt−2 + At−1 ) = At + (θ − 1)At−1 − θAt−2 Por lo tanto, {Yt } es un proceso de medias m´oviles de orden 2 (MA(2)). b) τk = cov(Yt , Yt−k ) = E[Yt Yt−k ] = E[(Xt − Xt−1 )(Xt−k − Xt−k−1 )] = E[Xt Xt−k ] − E[Xt Xt−k−1 ] − E[Xt−1 Xt−k ] + E[Xt−1 Xt−k−1 ] = γk − γk+1 − γk−1 + γk = 2γk − γk+1 − γk−1 c) τ0 = 2γ0 − γ−1 − γ1 = 2(γ0 − γ1 ) ρ∗k = 2γk − γk−1 − γk+1 2ρk − ρk−1 − ρk+1 = 2(γ0 − γ1 ) 2(ρ0 − ρ1 ) d) Otra forma de expresar a {Yt } es Yt = (1 + (θ − 1)B − θB 2 )At = θ(B)At . Como θ(B) = (1 + (θ − 1)B − θB 2 ) = 0 tiene como ra´ıces a b1 = −1 y b2 = 1/θ pero |b1 | ≯ 1 entonces {Yt } no es invertible. 2 y considere las Ejercicio 2.5.5. Sea {At } un proceso de ruido blanco con media cero y varianza σA siguientes ecuaciones en diferencias 29 a) Xt − 1/2Xt−1 = At − 1/3At−1 . b) Xt − 1/3Xt−1 = At − 0,7At−1 + 0,1At−2 . Para cada una de estas ecuaciones considere el proceso {Xt } que es su soluci´ on. i. Halle en cada caso la funci´ on de autocovarianza del proceso {Xt }. ii. Halle en cada caso la funci´ on de autocorrelaci´ on del proceso {Xt }. iii. Calcule a funci´ on de autocorrelaci´ on parcial πk para k = 1, 2. Soluci´ on. a) i. γk =  2 2 1+θ1 +2φ1 θ1 , σA 1−φ2 ( 1,04.σ 2 , = 0,185 2A σ , si k ≥ 1; 2k−1 A si k = 0; 1 k−1 σ 2 (φ1 +θ1 )φ1 2(1+φ1 θ1 ) , A 1−φ1 ii. ρk = si k = 0; si k ≥ 1; 5 7 ∗ 2k+1 iii. π1 = ρ1 = 5/28 = 0,178 1 5/28 5/28 5/56 = 45/759 ∼ π2 = = 0,059 1 5/28 5/28 1 b) i. Como la ra´ız de la ecuaci´ on caracter´ıstica x − 1/3 = 0 es x = 1/3 entonces γk = c( 13 )k P2 k 2 donde c es tal que γk − 1/3γk−1 = σA j=k θj αj−k y αk = (1/3) 2 2,333 y Entonces para k = 0, 1 tendremos el sistema de ecuaciones γ0 − 1/3γ−1 = σA 2 2 2 γ1 − 1/3γ0 = −2σA de donde sale que γ1 = −1,375σA y γ0 = 1,875σA = c Por lo tanto, 2 1 k γk = 1,875σA ( ) 3 ii. ρk = ( 31 )k iii. π1 = ρ1 1 1/3 π2 = 1 1/3 = 1/3. 1/3 1/9 =0 1/3 1 Ejercicio 2.5.6. En cada uno de los siguientes casos grafique los valores dados en la tabla y especifique los modelos de serie de tiempo estacionaria {Xt } a los que podr´ıa pertenecer los pares de funciones de autocorrelaci´ on te´ orica y de autocorrelaci´ on parcial te´ orica graficadas. Establezca la relaci´ on que define en forma general a Xt en el modelo identificado y, cuando sea posible, determine el signo de cada uno de los coeficientes involucrados. a) k ρk πk 1 0.940 0.940 2 0.846 -0.329 3 0.762 0.129 4 0.686 -0.051 5 0.617 0.021 30 6 0.555 -0.008 b) k ρk πk 1 0.699 0.699 2 0.280 -0.408 3 0.112 0.286 4 0.045 0.218 5 0.018 0.175 6 0.007 -0.145 c) k ρk πk 1 -0.699 -0.699 2 0.280 -0.408 3 -0.112 -0.286 4 0.045 -0.218 5 -0.018 -0.175 6 0.007 -0.145 d) k ρk πk 1 -0.940 -0.940 2 0.846 -0.329 3 -0.762 -0.129 4 0.686 -0.051 5 -0.617 -0.021 6 0.555 -0.008 Soluci´ on. a) Los gr´ aficos podr´ıan corresponder a un proceso ARMA(1,1) de expresi´on Xt = φ1 Xt−1 + At + θAt−1 en donde φ1 > 0. b) Los gr´aficos podr´ıan corresponder a un proceso ARMA(1,1) de expresi´on Xt = φ1 Xt−1 + At + θAt−1 en donde φ1 > 0. c) Los gr´aficos podr´ıan corresponder a un proceso ARMA(1,1) de expresi´on Xt = φ1 Xt−1 + At + θAt−1 en donde φ1 < 0. 31 d) Los gr´aficos podr´ıan corresponder a un proceso ARMA(1,1) de expresi´on Xt = φ1 Xt−1 + At + θAt−1 en donde φ1 < 0. Ejercicio 2.5.7. Sea {Xt } un proceso ARMA(1,1) de ecuaci´ on Xt − φ1 Xt−1 = At + θ1 At−1 donde {Xt } es un proceso de ruido blanco. a) Defina el proceso {Bt } mediante la igualdad Bt = At + θ1 At−1 (2.5.1) Xt − φ1 Xt−1 = Bt . (2.5.2) exprese (2.5.1) como Suponiendo que |φ1 | < 1 y usando argumentos similares a los empleados en el caso AR(1), demuestre que {Xt } puede expresarse como Xt = ∞ X φj1 Bt−j (2.5.3) j=0 b) Usando (2.5.1) y (2.5.3) comprebe que Xt puede expresarse como Xt = At + (φ1 + θ1 ) ∞ X φj−1 1 At−j . j=1 c) Suponiendo que |θ1 | < 1 y usando (2.5.1) muestre que At puede expresarse como At = ∞ X j=0 32 (−θ1 )j Bt−j (2.5.4) d) Suponiendo que |θ1 | < 1 y usando (2.5.2) y (2.5.4) muestre que At puede expresarse como At = Xt − (φ1 + θ1 ) ∞ X (−θ)j−1 1 Xt−j . j=1 e) Considerese ahora el proceso AR(1) estacionario (causal) {Zt } soluci´ on de la ecuaci´ on Zt = φ1 Zt−1 + At . Demuestre (usando (2.5.1)y (2.5.3)) que Xt puede expresarse como Xt = Zt + θ1 Zt−1 (2.5.5) f ) Usando (2.5.5) calcule la funci´ on de autocovarianza τ del proceso {Xt } en t´erminos de la funci´ on de autocovarianza del proceso {Zt }. g) Usando el resultado de la parte (f ) calcule la funci´ on de autocovarianza γ y de autocorrelaci´ on ρ del proceso {Xt }. Soluci´ on. a) Xt = φ1 Xt−1 + Bt = φ1 (φ1 Xt−2 + Bt−1 ) + Bt = φ21 Xt−2 + φ1 Bt−1 + Bt = ... = φj1 Xt−j + j−1 X φi1 Bt−i i=0 entonces Xt − j−1 X φi1 Bt−i = φj1 Xt−j i=0 y E[Xt − j−1 X φi1 Bt−i ] = φj1 E[Xt−j ] = φj1 E[Xt ] → 0 si j → ∞ i=0 pues |φ1 | < 1. Por lo tanto, Xt = ∞ X φi1 Bt−i . i=0 b) Xt = ∞ X φi1 Bt−i i=0 = ∞ X φi1 (At−i + θ1 At−i−1 ) i=0 = ∞ X φi1 At−i i=0 = At + + ∞ X φ1i−1 θ1 At−i i=1 ∞ X φ1i−1 (θ1 + φ1 )At−i i=1 33 c) At = Bt − θ1 At−1 = Bt − θ1 (Bt−1 − θ1 At−2 ) = ... = j−1 X (−θ1 )j At−j + (−θ1 )j At−j i=0 Luego, E[At − j−1 X (−θ1 )j At−j ] = (−θ1 )j E[At−j ] = (−θ1 )j E[Xt ] → 0 si j → ∞ i=0 pues |θ1 | < 1. Por lo tanto, ∞ X At = (−θ1 )j Bt−j . j=0 d) ∞ X At = (−θ1 )j Bt−j = j=0 ∞ X (−θ1 )j (Xt−i − φ1 Xt−i−1 ) j=0 = Xt + ∞ X (−θ1 )j−1 (−θ1 − φ1 )Xt−j j=1 = Xt − (φ1 + θ1 ) ∞ X (−θ)j−1 1 Xt−j j=1 e) Como Zt es soluci´ on de Zt = φ1 Zt−1 + At entonces Zt = Zt + θ1 Zt−1 = At + ∞ X P∞ i i=0 φ1 At−i . φi1 At−i + θ1 i=1 = At + (φ1 + θ1 ) ∞ X Y adem´as, φi−1 1 At−i i=1 ∞ X φ1i−1 At−i i=1 = Xt φk 2 1 si |φ1 | < 1. f) Como Zt = φ1 Zt−1 + At es un AR(1), entonces la autocorrelaci´on es τk = σA 1−φ21 Luego, γk = cov(Xt , Xt−k ) = E[(Zt + θ1 Zt−1 )(Zt−k + θ1 Zt−k−1 )] = E[Zt Zt−k ] + θ1 (E[Zt Zt−k−1 ] + E[Zt−1 Zt−k ]) + θ12 E[Zt−1 Zt−k−1 ] = τk + θ1 (τk+1 + τk−1 ) + θ12 τk g) 2 γk = σ A 2 = σA k+1 k−1 φk1 φk1 2 φ1 2 φ1 2 2 + θ (σ + σ ) + θ σ 1 A A 1 A 1 − φ21 1 − φ21 1 − φ21 1 − φ21 φk−1 1 (φ1 + θ1 (φ21 + 1) + θ12 φ1 ) 1 − φ21 34 Y ρk = φk1 Este ejercicio nos muestra una forma simple de calcular la funci´on de autocovarianza de un proceso ARMA(1,1) a partir de conocer la expresi´on de la funci´on de autocovarianza de un proceso AR(1). Ejercicio 2.5.8. Sea {Xt } un proceso estacionario con funci´ on de autocovarianza γ. Defina la funci´ on G mediante la igualdad ∞ X G(z) = γk z k (2.5.6) k=−∞ suponiendo que la serie converge para todo z en un anillo de la forma {z : 1/r < |z| < r} para alg´ un r > 1. La funcion G definida por (2.5.6) se la llama funci´ on generadora de autocovarianzas. a) ¿Cu´ al es la funci´ on generadora de autocovarianzas correspondiente a un proceso de ruido blanco? b) Sea {Xt } el proceso de ecuaci´ on ∞ X Xt = αj At−j (2.5.7) j=−∞ donde {At } es un proceso de ruido blanco. Suponiendo que existe r > 1 donde ∞ X |αj |xj j=−∞ converge para todo x ∈ {z : 1/r < |z| < r}, demuestre que la funci´ on G puede expresarse como 2 G(z) = σA α(z)α(z −1 ) donde ∞ X α(x) = αk xk . (2.5.8) (2.5.9) k=−∞ c) Considere ahora el caso en que {Xt } es un proceso ARMA(p,q) de ecuaci´ on Φ(B)Xt = θ(B)At (2.5.10) donde Φ(z) es distinto de 0 para |z| = 1. Expresando apropiadamente (2.5.10) y usando (2.5.7), (2.5.8) y (2.5.9), probar que la funci´ on generadora de autocovarianzas de este proceso esta dada por −1 2 θ(z)θ(z ) G(z) = σA para 1/r < |z| < r. (2.5.11) Φ(z)Φ(z −1 ) d) Use (2.5.6) y (2.5.11) para hallar la expresi´ on de la funci´ on de autocovarianza de un proceso MA(2). Soluci´ on. a) Como ( 2, σA γk = 0, si k = 0; si k = 6 0 entonces 2 G(z) = γ0 z 0 = σA 35 b) 2 σA α(z)α(z −1 ) = 2 σA 2 = σA = ∞ X k=−∞ ∞ X αk z ∞ X k αi z −i i=−∞ ∞ X αk αk+i z i k=−∞ i=−∞ ∞ X γi z i i=−∞ = G(z) c) Si las ra´ıces de Φ(B) son de m´ odulo mayor que 1, entonces Xt = 2 2 G(z) = σA α(z)α(z −1 ) = σA θ(B) Φ(B) At = α(B). Luego, θ(z)θ(z −1 ) Φ(z)Φ(z −1 ) d) 2 G(z) = σA θ(z)θ(z −1 ) 2 = σA (1 + θ1 z + θ2 z 2 )(1 + θ1 z −1 + θ2 z −2 ) 2 = σA [(1 + θ12 + θ22 ) + θ2 z −2 + (θ1 + θ1 θ2 )z −1 + (θ1 + θ1 θ2 )z + θ2 z 2 ] ∞ X γk z k = k=−∞ Entonces, igualando los coeficientes de cada potencia se obtiene que: 1. si |k| > 2, γk = 0 2. γ−2 = θ2 3. γ−1 = θ1 + θ1 θ2 4. γ0 = 1 + θ12 + θ22 5. γ1 = θ1 + θ1 θ2 6. γ2 = θ2 36 Cap´ıtulo 3 Modelos para series de tiempo no estacionarias No todas las series de tiempo pueden ser modeladas mediante procesos ARMA, sin embargo, en la pr´actica la mayor´ıa de las series reales no estacionarias pueden ser transformadas en series que responden al patr´ on de alg´ un proceso ARMA(p,q). En este cap´ıtulo abordaremos el estudio de m´etodos que producen transformaciones sobre una serie de tiempo para convertirla en otra susceptible de ser modelada mediante un proceso ARMA. 3.1. Enfoque cl´ asico Una serie de tiempo {Zt } puede ser expresada como Zt = mt + st + Xt (3.1.1) donde * mt es una funci´ on determin´ıstica de variaci´on suave conocida como la tendencia. * st es una funci´ on determin´ıstica con per´ıodo conocido denominada estacionalidad. Se la entiende como una pauta regular de comportamiento per´ıodico de la serie. * Xt es la componente aleatoria denominada ruido aleatorio o componente estacionaria. Ejemplo 3.1.1. Un ejemplo de una serie de tiempo real es el ´ındice diferencial normalizado de vegetaci´ on (tambi´en llamado NDVI) de la Regi´ on Amaz´ onica cuyo gr´ afico se observa abajo: 37 A continuaci´ on se mostrar´ an algunos m´etodos llamados “cl´asicos” para estimar y luego extraer las componentes determin´ısticas mt y st esperando que la componente aleatoria residual Xt pueda ser modelada por un proceso ARMA (ver Leiva [1]). 3.1.1. Estimaci´ on de mt por el m´ etodo de cuadrados m´ınimos Asumamos que en la ecuaci´ on (3.1.1) la componente de estacionalidad st no est´a presente; entonces (3.1.1) queda: Zt = mt + Xt con E[Xt ] = 0. Se trata primero de estimar la tendencia suponiendo que mt es un polinomio en t con grado fijo de la forma mt = a0 + a1 t + ... + ak tk . Luego sePaplica el m´etodo de m´ınimos cuadrados para aproximar a0 , ..., ak , o sea, elegirlos tal que minimizen nt=1 (Zt − mt )2 . Se tendr´a entonces m bt = b a0 + b a1 t + ... + b ak tk donde b a0 , b a1 , ..., b ak son las estimaciones obtenidas por el m´etodo de m´ınimos cuadrados (ver Leiva [1]). Por u ´ltimo, a los valores residuales x bt = zt − m b t se los estima apropiadamente y con dicha estimaci´ on se puede mejorar la de mt . Ejemplo 3.1.2. Se aplic´ o a la serie del ejemplo 3.1.1 el m´etodo de m´ınimos cuadrados suponiendo que mt es un polinomio de grado 1 de la forma mt = a0 + a1 t. De ´esta estimaci´ on de la tendencia result´ o que a0 = −3,809 y a1 = 0,033 y su gr´ afico se muestra a continuaci´ on: Luego se estim´ o el residuo resultando una serie como se observa en el siguiente gr´ afico: 38 3.1.2. Estimaci´ on de la tendencia mediante suavizamiento por promedios En este caso tambi´en se supone que la componente de estacionalidad st en el ejemplo 3.1.1 no est´a presente y que mt es aproximadamente lineal en intervalos de tiempo de longitud L = 2h + 1 con h real y fijo. Una estimaci´ on de mt es m b t = Pt = h X 1 Zt+i para t ∈ [h + 1, n − h] 2h + 1 i=−h ya que h h h X X X 1 1 1 (mt+i + Xt+i ) = mt+i + Xt+i ∼ = mt + 0 = mt 2h + 1 2h + 1 2h + 1 Pt = i=−h i=−h i=−h pues como mt ∼ = At + b y E[Xt ] = 0, entonces h h h X X X 1 1 1 mt+i ∼ A(t + i) + b = (At + b) + i∼ = = mt . 2h + 1 2h + 1 2h + 1 i=−h i=−h i=−h Para los valores extremos de mt se pueden usar los promedios ponderados unilaterales con α ∈ (0, 1): m bt = n−t X α(1 − α)j Zt+j si t = 1, ..., h y m bt = j=0 t−1 X α(1 − α)j Zt−j si t = n − h + 1, ..., n j=0 Esta ponderaci´ on se realiza asumiendo que los valores extremos de la serie tendr´an una menor influencia en su comportamiento. Estas u ´ltimas estimaciones no son altamente sensibles a los valores de α. Emp´ıricamente se ha encontrado que con α ∈ (0,1; 0,3) las estimaciones son razonables [2]. 3.1.3. Estimaci´ on de la tendencia y la estacionalidad cuando la tendencia es constante durante cada per´ıodo de la componente estacional Consideremos que Zt se puede descomponer como Zt = mt + st + Xt (notar que ahora est´ an presente tanto la tendencia como la componente estacional) donde E[Xt ] = 0 y st = st+d ∀t ∈ Z con Pd s = 0; siendo d un n´ umero natural fijo que responde al tiempo en que la serie vuelve a repetir j=1 j P su comportamiento, o sea, el per´ıodo de la serie. La condici´on dj=1 sj = 0 no es restrictiva pues si Pd j=1 sj = s 6= 0 entonces d X d X sj = s + (sj − s) j=1 j=1 d X s =s+ (sj − ) d j=1 =s+ d X s∗j j=1 P con s∗j = sj − ds . Luego, Zt = mt + ds + s∗t + Xt = m∗t + s∗t + Xt con dj=1 s∗j = 0. Para eliminar la tendencia y la estacionalidad simult´aneamente, supondremos que los cambios en la tendencia son peque˜ nos durante cada uno de los per´ıodos de longuitud d de la estacionalidad. 39 P as Suponiendo que en el k-´esimo per´ıodo de longitud d se verifica que dj=1 sj+(k−1)d = 0 y que adem´ m(k−1)d+1 ∼ ¯ k (debido a que la tendencia es muy peque˜ na en cada = m(k−1)d+2 ∼ = ... ∼ = m(k−1)d+d−1 ∼ =m uno de los per´ıodos), entonces un estimador de mt con t ∈ [(k − 1)d + 1, kd] es m b t = Pk pues d Pk = 1X Z(k−1)d+j d j=1 = d d d j=1 j=1 j=1 1X 1X 1X m(k−1)d+j + s(k−1)d+j + x(k−1)d+j d d d ∼ ¯k +0+0 =m =m ¯k donde r es el n´ umero de per´ıodos de longitud d. Ahora, una buena estimaci´ on para sj es r−1 sbj = 1X (Zkd+j − m b kd+j ) para j = 1, ..., d r k=0 pues r−1 sbj = 1X (Zkd+j − m b kd+j ) r k=0 = = = r−1 1X r 1 r k=0 r−1 X (mkd+j + skd+j + Xkd+j − m ¯ k+1 ) k=0 r−1 1X r (Zkd+j − Pk+1 ) (mkd+j − m ¯ k+1 ) + k=0 1 ∼ =0+0+ r r−1 X r−1 r−1 k=0 k=0 1X 1X Xkd+j + skd+j r r skd+j k=0 r−1 = 1X sj r k=0 = sj ya que mkd+j ∼ ¯ k+1 , E[Xt ] = 0 y sj = sj+d = ... = sj+(r−1)d . O sea, la estimaci´on de st est´ a dada =m por sbt = sbj para t = kd + j ∀k = 0, 1, ..., r − 1. Finalmente, la estimaci´ on de la componente estacionaria viene dada por: b t = Zt − m X b t − sbt = Zt − Pk − sbj si t = kd + j, para k = 0, ..., r − 1 y j = 1, ..., d. Ejemplo 3.1.3. Se aplic´ o a los primeros 100 valores de la serie del ejemplo 3.1.1 el m´etodo estudiado en esta secci´ on. Suponiendo que el per´ıodo es d = 12, el gr´ afico de la serie junto con la tendencia resulta ser: 40 Luego se estimaron la componente estacional junto con el residuo lo que result´ o en: 3.1.4. Estimaci´ on de la tendencia y la estacionalidad: Caso General Este m´etodo es m´ as general que el anterior ya que no se supone que la tendencia sea aproximadamente constante en cada per´ıodo. Primero se hace una estimaci´ on inicial de la tendencia mediante un suavizamiento con un promedio, de acuerdo con el siguiente procedimiento: Si d = 2h + 1 entonces h X 1 Pt = Zt+i para h + 1 ≤ t ≤ n − h 2h + 1 i=−h 41 y si d = 2h entonces Pt = h−1 X 1 1 1 [ Zt−h + ( Zt+i ) + Zt+h ] para h + 1 ≤ t ≤ n − h 2h 2 2 i=−h+1 La influencia de las componentes estacionales en cada promedio es igual a cero pues la longitud de los intervalos de tiempo donde se analizan los promedios es igual a la longitud del per´ıodo de la componente estacional; mientras que por ser E[Xt ] = 0, el promedio de las componentes estacionarias es aproximadamente igual a cero. Luego, Pt es una buena estimaci´on de mt para t = h + 1, ..., n − h. Partiendo de esta estimaci´ on inicial de mt , se consideran las desviaciones δt = zt − Pt para h + 1 ≤ t ≤ n − h y se calculan promedios de estas desviaciones qj : j = 1, ..., d. As´ı para el caso en que n = rd, si d = 2h + 1,  1 Pr−1  δj+kd , si j = 1, ..., h;  r−1 Pr−1k=1 1 qj = r k=0 δj+kd , si j = h + 1;   1 Pr−2 k=0 δj+kd , si j = h + 2, ..., d r−1 y si d = 2h, ( qj = 1 r−1 1 r−1 Pr−1 k=1 δj+kd , si j = 1, ..., h; Pr−2 k=0 δj+kd , si j = h + 1, ..., d Notemos que en cada promedio s´ olo est´an involucradas las desviaciones referidas a una misma componente estacional. Las componentes estacionales sj se estiman mediante d sbj = qj − 1X qi . d i=1 Se reestima la tendencia mediante un suavizado con promedios de los datos desestacionalizados zt o por ajuste mediante cuadrados m´ınimos. Finalmente, se calculan las estimaciones de la componente b t = Zt − m estacionaria Xt mediante X b t − sbt . 3.2. Modelos ARIMA y SARIMA El m´etodo que se ver´ a en esta secci´on se basa en la utilizaci´on de diferencias de los datos originales para la eliminaci´ on de la tendencia y la estacionalidad y poder luego modelar la componente estacionaria resultante. Definici´ on 3.2.1. Sea {Zt : t ∈ Z} una serie de tiempo, se define el operador diferencia ∆ por ∆Zt = (1 − B)Zt = Zt − Zt−1 Se puede notar que si aplicamos este operador k veces a un polinomio de grado k se obtiene una constante (ver ejercicio 3.4.6). Definici´ on 3.2.2. Sea {Zt : t ∈ Z} una serie de tiempo y sea d ∈ N, se define el operador diferencia de tama˜ no d, ∆d por ∆d Zt = (1 − B d )Zt = Zt − Zt−d Se puede verificar que si se aplica este operador a un proceso Zt = mt + st + Xt donde st es la componente estacional de per´ıodo d se tendr´a ∆d Zt = m∗t + Xt∗ donde m∗t = mt − mt−d y Xt∗ = Xt − Xt−d . 42 3.2.1. Modelos ARIMA Los procesos ARIMA generalizan a los procesos ARMA. Estos procesos ARIMA sirven para modelar series de tiempo que luego de diferenciarlas un n´ umero finito de veces, se reducen a un proceso ARMA. Definici´ on 3.2.3. Si k ∈ N ∪ {0}, se dice que {Zt : t ∈ Z} es un proceso Autorregresivo Integrado de Promedios M´ oviles de ´ ordenes p,k,q y se lo indica por ARIMA(p,k,q), si la serie de tiempo {Xt } que se obtiene de diferenciar k veces a {Zt }, es decir Xt = ∆k Zt , es un proceso ARMA(p,q). Equivalentemente, se puede decir que {Zt } es un proceso ARIMA(p,k,q) si satisface la ecuaci´ on en diferencias Φ(B)(1 − B)k Zt = θ(B)At con Φ(B) = 1 − φ1 B − ... − φp B p y θ(B) = 1 + θ1 B + ... + θq B q donde {At } es un proceso de ruido blanco y donde se supone que los polinomios Φ(x) y θ(x) no tiene ra´ıces comunes y que Φ(x) 6= 0 si x = 1. En especial, es interesante que la serie {Xt } que resulta de la definici´on 3.2.3 defina un proceso estacionario e invertible, por lo que se requerir´a que las ra´ıces de Φ(x) y θ(x) sean todas de m´ odulo mayor que 1. Notemos que el proceso ARIMA(p,k,q) se reduce a un: a) Proceso ARMA(p,q) cuando k = 0 b) Proceso MA(q) cuando p = 0 = k c) Proceso AR(p) cuando q = 0 = k d) Proceso ARI(p,k) cuando q = 0 e) Proceso IMA(k,q) cuando p = 0 Estos procesos ARIMA se caraterizan por tener una funci´on de autocorrelaci´on te´orica con decaimiento positivo muy lento. Sin embargo, esta condici´on es suficiente pero no necesaria ya que por ejemplo los procesos ARMA con una ra´ız del polinomio de su parte autorregresiva cercana a 1 tambi´en presentan un comportamiento similar. Un problema que tienen estos procesos ARIMA(p,k,q) es que modelan solamente la no estacionariedad que tienen una ra´ız igual a 1 en su parte autorregresiva, es decir, ψ(B)Zt = θ(B)At donde ψ(B) = φ(B)(1 − B)k con φ(B) = 1 − φ1 B − ... − φp B p tal que Φ(x) 6= 0 si −1 < x ≤ 1. 3.2.2. Procesos SARIMA (ARIMA estacionales) En los modelos ARIMA se incorporaba en la descomposici´on cl´asica (3.1.1) la funci´on per´ıodica st con el fin de modelar la presencia de estaciones en una serie de tiempo. Esta componente estacional es r´ıgida en el sentido que su contribuci´on se repite en cada uno de los ciclos. Los modelos SARIMA permiten cierta aleatoriedad en el comportamiento estacional de un ciclo al siguiente, lo que parece m´as razonable en la pr´ actica. Definici´ on 3.2.4. Si k, h ∈ N ∪ {0}, se dice que {Zt : t ∈ Z} es un proceso ARIMA Estacional de per´ıodo d y ´ ordenes (p, k, q) y (r, h, s)d y se lo denota por SARIMA (p, k, q), (r, h, s)d si la serie diferenciada Xt = ∆k ∆hd Zt = (1 − B)k (1 − B d )h Zt , satisface la ecuaci´ on Φ(B)ψ(B d )Xt = θ(B)η(B d )At donde {At } es un proceso de ruido blanco y donde Φ(x) = 1−φ1 x−...−φp xp , ψ(x) = 1−ψ1 x−...−ψr xr , θ(x) = 1 + θ1 x + ... + θq xq y η(x) = 1 + η1 x + ... + ηs xs . 43 Para identificar estos procesos primero se determinar´an los valores k y h tal que los datos diferenciados Xt = (1 − B)k (1 − B d )h Zt sean aparentemente estacionarios. Despu´es, se seleccionar´ an r y s tales que los valores de autocorrelaci´ on ρbid con i ∈ N de la serie {Xt } pueda ser modelado por un proceso ARMA(r,s). Por u ´ltimo, se elegir´an p y q tal que los valores ρbj con i = 1, ..., d − (q + 1) sean compatibles con las de un proceso ARMA(p,q). 3.3. Transformaciones Con frecuencia se verifica que var(Zt ) = f (µt ) donde f es una funci´ on no negativa. A veces, este problema se puede resolver aplicando una transformaci´on T a la serie {Zt } para obtener una serie {Yt } donde Yt = T (Zt ) que va a tener varianza constante a lo largo del tiempo. Una forma de encontrar T es usar el desarrollo de Taylor de orden 1 alrededor de µt como sigue: T (Zt ) ∼ = T (µt ) + T 0 (µt )(Zt − µt ). Como supongo que T (µt ) es no aleatoria, resulta que var(T (Zt )) ∼ = T 0 (µt )2 f (µt ). Y como quiero var(T (µt )) = k constante entonces √ k T (µt ) = p f (µt ) 0 o sea, √ k p dµt f (µt ) Z T (µt ) = (3.3.1) Box y Cox en 1964 introdujeron para estabilizar la varianza de una serie de tiempo {Zt } a la familia T (Zt ) = Ztλ − 1 λ donde λ es llamado el par´ ametro de la transformaci´on. Notar que para λ = 0, Ztλ − 1 = log(Zt ) λ→0 λ l´ım p Ejemplo 3.3.1. Si var(Zt ) = mµt , o sea, la desviaci´ on est´ andar de Zt var´ıa proporcionalmente al 2 2 nivel µt , entonces f (µt ) = m µt y √ √ Z √ Z 1 k k p T (µt ) = dµt = log(µt ) (3.3.2) dµt = k 2 2 mµt m m µt Por lo que la transformaci´ on que estabiliza la varianza en este caso es la transformaci´ on logar´ıtmica. 44 3.4. Ejercicios del Cap´ıtulo III Ejercicio 3.4.1. Considere el siguiente conjunto de datos que corresponden al consumo bimestral de gas natural de una familia durante el per´ıodo de enero de 1983 a diciembre de 1991. ˜ ANO 1983 1984 1985 1986 1987 1988 1989 1990 1991 1 96 90 68 111 103 133 139 114 126 2 200 223 240 223 309 341 246 135 127 BIMESTRES 3 4 5 332 162 118 333 232 86 367 238 135 460 321 128 486 379 136 595 342 132 378 378 252 258 394 272 241 454 251 6 65 68 106 124 101 77 165 167 153 a) Grafique los datos en funci´ on de t. b) Calcule y grafique las funciones muestrales de autocorrelaci´ on y autocorrelaci´ on parcial de estos datos. c) Suponga que los datos han sido generados por un proceso que responde a la forma Zt = mt + st + Xt , donde mt es la tendencia, st es la estacionalidad con per´ıodo d y donde {Xt } es un proceso estacionario con media 0. Bas´ andose en las partes (a) y (b), determine el valor del per´ıodo d. d) Suponga que la tendencia es constante dentro cada per´ıodo y calcule y grafique las estimaciones mt , st y xt para t = 1, ..., 54. e) Utilice las estimaciones de st halladas en la parte (d) para calcular los datos desestacionalizados, es decir calcule yt = zt − sbt para t = 1, ..., 54. Con estos valores, que responden aproximadamente a la forma yt = zt − sbt ∼ = mt + Xt , calcule otra estimaci´ on de la tendencia ajust´ andoles por m´ınimos cuadrados una recta de ecuaci´ on r(t) = a0 + a1 t. Dibuje en un mismo gr´ afico la recta obtenida, la estacionalidad estimada y los residuos ybt = zt − sbt − mt . f ) Calcule y grafique las funciones muestrales de autocorrelaci´ on y autocorrelaci´ on parcial de los residuos obtenidos en la parte (e). g) Determine un modelo apropiado para Xt en base a lo observado en la parte (f ) y estime sus coeficientes. Soluci´ on. a) La serie de tiempo es 45 b) Las funciones de autocorrelaci´ on y de autocorrelaci´on parcial son respectivamente: c) De acuerdo a los gr´ aficos de las partes (a) y (b) se puede decir que el per´ıodo es d=6. d) Las estimaciones de mt , st y xt son respectivamente: e) Con la nueva estimaci´ on de mt = r(t) = 173,4 + 1,65t, los gr´aficos de la tendencia, la estacionalidad y el residuo queda: 46 f) Las funciones de autocorrelaci´ on y de autocorrelaci´on parcial son respectivamente: g) De acuerdo a lo observado en la parte (f) y atendiendo al criterio de parsimonia podemos decir que el residuo xt puede ser modelado por un proceso AR(2) y usando la funci´on de autocorrelaci´ on parcial te´ orica de estos procesos podemos estimar los coeficientes como φ1 = 0,485 y φ2 = −0,409. Ejercicio 3.4.2. Considere el conjunto de datos zt con t = 1, ..., 54 del ejercicio 3.4.1. que corresponden al consumo bimestral de gas natural de una familia durante el per´ıodo de enero de 1983 a diciembre de 1991. a) Para corregir la no estacionariedad que se observa en la varianza use la transformaci´ on logar´ıtmica. Es decir calcule yt = log(zt ) para t = 1, ..., 54. b) Calcule y grafique las funciones muestrales de autocorrelaci´ on y autocorrelaci´ on parcial. c) En base a lo observado en la parte (b), diferencie a los datos transformados con un retardo (lag) apropiado y luego r´esteles la media de los datos diferenciados. d) Grafique los datos finales obtenidos en la parte (c) y calcule y grafique las funciones muestrales de autocorrelaci´ on y autocorrelaci´ on parcial. e) Bas´ andose en lo observado en la parte (d), determine posibles modelos para los datos graficados en (d) y sus correspondientes modelos para los datos originales. 47 Soluci´ on. a) Los siguientes son los valores de zt del ejercicio anterior a los cuales le aplicamos logaritmo. ˜ ANO 1983 1984 1985 1986 1987 1988 1989 1990 1991 1 4.564 4.499 4.219 2.397 4.634 4.890 4.934 4.736 4.836 2 5.298 5.407 5.480 5.407 5.733 5.831 5.505 4.905 4.844 BIMESTRES 3 4 5.805 5.087 5.808 5.446 5.905 5.472 6.131 5.771 6.186 5.937 6.388 5.834 5.934 5.934 5.552 5.976 5.484 6.118 5 4.770 4.454 4.905 4.852 4.912 4.882 5.529 5.605 5.525 6 4.174 4.219 4.663 4.820 4.615 4.343 5.105 5.117 5.030 b) Las funciones de autocorrelaci´ on y autocorrelaci´on parcial resultan ser: c) Con un lag = 6 la media de los valores diferenciados resulta ser 0.044558. Por lo que los valores de los datos diferenciados menos su media son: ˜ ANO 1983 1984 1985 1986 1987 1988 1989 1990 1991 1 2 -0.109 -0.325 -1.866 2.192 0.211 -0.001 -0.243 0.056 0.064 0.029 -0.118 0.282 0.054 -0.371 -0.645 -0.106 BIMESTRES 3 4 -0.042 0.053 0.181 0.011 0.158 -0.499 -0.426 -0.113 0.315 -0.019 0.255 0.126 -0.147 0.056 -0.003 0.097 5 6 -0.361 0.406 -0.098 0.0161 -0.075 0.602 0.032 -0.125 0.001 0.399 0.112 -0.250 -0.316 0.718 -0.033 -0.132 d) Los datos y sus funciones de autocorrelaci´on y autocorrelaci´on parcial resultan ser ahora: 48 e) Como las autocorrelaciones de los datos diferenciados se comportan como ruido blanco, entonces el modelo correspondiente a los datos de consumo es ((1 − B 6 )log(Zt ) − 0,045) = At con {At } ruido blanco. Ejercicio 3.4.3. En cada uno de los siguientes casos los valores dados en la tabla son los primeros seis valores de las funciones te´ oricas de autocorrelaci´ on y autocorrelaci´ on parcial de la primera diferencia Xt = Zt − Zt−1 de un proceso {Zt }. Especifique los modelos de serie de tiempo estacionaria {Xt } a los que podr´ıan pertenecer los pares de funciones te´ oricas de autocorrelaci´ on y autocorrelaci´ on parcial dadas y establezca luego la relaci´ on que define en forma general a {Zt } en el modelo que le corresponde de acuerdo al modelo especificado para Xt . Cuando sea posible determine el signo y/o valores de los coeficientes involucrados. a) k ρk πk 1 0 0 2 0 0 3 0 0 4 0 0 5 0 0 b) k ρk πk 1 0.213 0.213 2 -0.284 -0.3455 c) k ρk πk 1 0.933 0.933 2 0.807 -0.500 d) k ρk πk 1 0.727 0.727 2 0.655 0.267 6 0 0 3 0 0.184 3 0.663 0 3 0.589 0.105 4 0 -0.187 4 0.524 0 4 0.530 0.042 5 0 0.141 5 0.403 0 5 0.477 0.017 6 0 -0.128 6 0.302 0 6 0.429 0.007 49 Soluci´ on. a) Dado que no hay correlaci´on ni correlaci´on parcial entre los datos, decidimos que {Xt } es un proceso de ruido blanco. Por lo que (1 − B)Zt = At . b) Dado que la correlaci´ on de {Xt } se anula para k > 2 podemos decir que es un proceso MA(2) con θ2 < 0. Por lo que (1 − B)Zt = At + θ1 At−1 + θ2 At−2 . c) Dado que la correlaci´ on parcial de {Xt } se anula para k > 2 podemos decir que es un proceso AR(2) con φ2 = −0,5 < 0 y φ1 = 1,4 > 0 Por lo que (1 − B)Zt = At + φ1 (1 − B)Zt−1 + φ2 (1 − B)Zt−2 , o sea (1 − (1 + φ1 )B + (φ1 − φ2 )B 2 + φ2 B 3 )Zt = At . d) Dado que ambas funciones de correlaci´on y correlaci´on parcial son positivas y decresen a cero, diremos que se trata de un proceso ARMA(1,1), o sea (1 − φ1 B)(1 − B)Zt = (1 + θ1 B)At . Ejercicio 3.4.4. Sea {Xt } un proceso ARIMA(p,2,q) que satisface la ecuaci´ on en diferencias Φ(B)(1 − B)2 Xt = θ(B)At , (3.4.1) donde {At } es un proceso de ruido blanco. Demuestre que cualquiera sean las variables aleatorias V y W, el proceso {Yt } definido por Yt = V + tW + Xt , tambi´en satisface la ecuaci´ on (3.4.1). Soluci´ on. Φ(B)(1 − B)2 Yt = Φ(B)(1 − B)2 (V + tW + Xt ) = Φ(B)(1 − B)2 V + Φ(B)(1 − B)2 tW + Φ(B)(1 − B)2 Xt = Φ(B)(V − 2V + V ) + W Φ(B)(t − 2(t − 1) + (t − 2)) + Φ(B)At = Φ(B)At . Por lo tanto, {Yt } satisface la ecuaci´ on (3.4.1). Ejercicio 3.4.5. Sean {Yt } y {Zt } dos procesos. Demuestre que ∆k (aYt + cZt ) = a∆k Yt + c∆k Zt ∀n ∈ Z Soluci´ on. Por inducci´ on, Si k = 1, ∆(aYt + cZt ) = aYt + cZt − aYt−1 + cZt−1 = a(Yt − Yt−1 ) + c(Zt − Zt−1 ) = a∆Yt + c∆Zt . Supongamos cierto para k − 1, veamos para k ∆k (aYt + cZt ) = ∆k−1 (∆(aYt + cZt )) = ∆k−1 (a∆Yt + c∆Zt ) = a∆k−1 (∆Yt ) + c∆k−1 (∆Zt ) = a∆k Yt + c∆k Zt . Ejercicio 3.4.6. Si Pk es un polinomio en t de grado k dado por Pk (t) = a0 + a1 t + ... + ak tk con ak 6= 0, demostrar que 50 (3.4.2) a) ∆i Pk (t) es un polinomio de grado ≤ k − i, es decir ∆i Pk (t) = b0 + b1 t + ... + bk−i tk−i , para todo i = 0, 1, ..., k. b) ∆k Pk (t) = k!ak . Soluci´ on. a) Se puede probar f´ acilmente usando la f´ormula del binomio sobre (t − 1)j que ∆tj = tj − (t − 1)j es un polinomio de grado j − 1 con coeficiente principal j. Luego, usando la linealidad de ∆ y lo anterior se puede ver que si Pk (t) = a0 + a1 t + ... + ak tk es un polinomio de grado k, entonces ∆Pk (t) es un polinomio de grado k − 1 con coeficiente principal kak . Por otro lado, se ve tambi´en usando la f´ormula del binomio que si i ≤ j, ∆i tj es un polinomio j! de grado j − i con coeficiente principal (j−i)! . Combinando estos resultados y usando la linealidad de ∆ resulta que ∆i Pk (t) es un polinomio k! . de grado menor o igual que k − i con coeficiente principal ak (k−i)! k! b) De la parte (a) se deduce inmediatamente que ∆k Pk (t) = ak (k−k)! = ak k!. Ejercicio 3.4.7. Sea ∆d el operador diferencia de tama˜ no d definido por ∆d = 1 − B d . a) Demuestre que ∆d es un operador lineal. Concretamente, siendo {Yt } y {Zt } dos procesos demuestre que ∆d (aYt + cZt ) = a∆d Yt + c∆d Zt . b) Demuestre que ∆kd es un operador lineal para todo k ∈ N. c) Pruebe que ∆d tk es un polinomio de grado k − 1 con coeficiente principal igual a kd. d) Demuestre, usando el resultado de las partes (a) y (c) anteriores, que si Pk (t) es un polinomio de grado k como el de la ecuaci´ on (3.4.2) entonces ∆id Pk (t) es un polinomio de grado ≤ k − i, k! i es decir ∆d Pk (t) = b0 + b1 t + ... + bk−i tk−i , con coeficiente principal igual a ak di (k−i)! , para todo i = 0, 1, ..., k. En particular, ∆kd Pk (t) = k!dk ak . Soluci´ on. a) ∆d (aYt + cZt ) = aYt + cZt − aYt−d + cZt−d = a(Yt − Yt−d ) + c(Zt − Zt−d ) = a∆d Yt + c∆d Zt . b) Se prueba del mismo modo que se prob´o en el ejercicio 3.4.5 usando inducci´on en k y la parte (a) para el caso k = 1.  P k k−i ti que es un polinomio de grado k − 1 con coeficiente c) ∆d tk = tk − (t − d)k = − k−1 i=0 i (−d) k! principal (k−1)! d = kd. d) Se puede ver usando la f´ ormula del binomio que si i ≤ k, ∆id tk es un polinomio de grado k − i k! con coeficiente principal di (k−i)! . Si Pk (t) = a0 + a1 t + ... + ak tk , por la linealidad de ∆id y k! lo anterior, ∆id Pk (t) es un polinomio de grado k − i con coeficiente principal ak di (k−i)! . En particular, ∆kd Pk (t) = ak dk k!. Ejercicio 3.4.8. Demuestre que si {Xt } es un proceso estacionario con E[Xt ] = µ y autocovarianza γk entonces ∗ )= a) Xt∗ = ∆Xt tal que t ∈ Z es un proceso estacionario con media cero y autocovarianza cov(Xt∗ , Xt+k 2γk − γk−1 − γk+1 . 51 ∗ )= b) Xt∗ = ∆d Xt tal que t ∈ Z es un proceso estacionario con media cero y autocovarianza cov(Xt∗ , Xt+k 2γk − γ|k−d| − γk+d . Soluci´ on. a) Xt∗ = Xt − Xt−1 E[Xt∗ ] = E[Xt ] − E[Xt−1 ] = µ − µ = 0 ∗ cov(Xt∗ , Xt+k ) = E[Xt Xt+k − Xt Xt+k−1 − Xt−1 Xt+k + Xt−1 Xt−1+k ] = E[Xt Xt+k ] − E[Xt Xt+k−1 ] − E[Xt−1 Xt+k ] + E[Xt−1 Xt−1+k ] = 2γk − γk−1 − γk+1 b) Xt∗ = Xt − Xt−d E[Xt∗ ] = E[Xt ] − E[Xt−d ] = µ − µ = 0 ∗ cov(Xt∗ , Xt+k ) = E[Xt Xt+k − Xt Xt+k−d − Xt−d Xt+k + Xt−d Xt−d+k ] = E[Xt Xt+k ] − E[Xt Xt+k−d ] − E[Xt−d Xt+k ] + E[Xt−d Xt−d+k ] = 2γk − γ|k−d| − γk+d Ejercicio 3.4.9. Sean {Vt } y {Wt } dos procesos estacionarios con medias µV y µW y autocovarianzas τk y δk respectivamente. Si Zt = Vt + tWt para todo t ∈ Z demostrar que a) {Zt } no es estacionario pues no cumple ninguna de las dos propiedades de la definici´ on de estacionariedad. b) Yt = ∆Zt no es estacionario a menos que Wt = W para todo t. Soluci´ on. a) E[Zt ] = E[Vt ] + tE[Wt ] = µV + tµW (depende de t). cov(Zt , Zt+k ) = τk +(t+k)E[Vt Wt+k ]+tE[Wt Vt+k ]+(t2 +tk)δk −µ2V −(k+2t)µV µW −(t2 +kt)µ2W que depende de t. b) E[Yt ] = E[Zt ] − E[Zt−1 ] = µW cov(Zt , Zt+k ) = E[(Zt −Zt−1 )(Zt+k −Zt+k−1 )]−µ2W que depende de t pues cov(Zt , Zt+k ) depende de t. 52 Cap´ıtulo 4 Predicci´ on Uno de los motivos m´ as importantes del estudio de las series de tiempo es el hecho de poder predecir los valores futuros, conocidos sus datos pasados. Esto es, se quiere predecir el valor Xn+m conocidos los valores X1 , X2 , ..., Xn y a este valor predicho se lo designar´a por X1,n (m). El criterio que usaremos para comparar predictores del valor Xn+m es mediante la evaluaci´on del Error Cuadr´ atico Medio de cada predictor X1,n (m) dado por ECM (X1,n (m)) = E[(Xn+m − X1,n (m))2 ] y eligiendo como el mejor a aqu´el que minimiza este valor. Si P est´a definida como la clase de todos los predictores de Xn+m donde se realiza la b´ usqueda y si adem´as se agrega la condici´on de que estos predictores sean funciones de X1 , X2 , ..., Xn , entonces el mejor predictor de Xn+m , que se denota b1,n (m), es X b1,n (m) = E[Xn+m /Xt : 1 ≤ t ≤ n]. X (4.0.1) Esto se puede entender como que el mejor estimador para el tiempo n + m es el promedio de los Xn+m correspondientes a todas las posibles realizaciones de Xt cuyos valores hasta el tiempo n coinciden con las observaciones x1 , ..., xn . Sin embargo, este resultado no es u ´til ya que requiere conocer la estructura probabil´ıstica de {Xt }. Por ello se restringir´a la clase P de predictores admitidos de Xn+m a los que son una combinaci´ on lineal de X1 , X2 , ..., Xn y de esta familia el que tenga menor Error Cuadr´ atico e1,n (m). En general, este predictor Medio se lo denominar´ a mejor predictor lineal, que se denotar´a por X no coincidir´a con el dado en (4.0.1), pero igual hay casos donde se da la coincidencia como se ver´ a en la secci´on (4.2). 4.1. Predictor Lineal con Error Cuadr´ atico Medio M´ınimo Sea {Xt } un proceso estacionario con funci´on de autocovarianza γ y media cero. El mejor predictor lineal de Xn+m est´ a dado por 0 (m) e1,n (m) = φ(m) Xn + ... + φ(m) X , n,n X1 = X n φ n,1 con X 0n = (Xn , Xn−1 , ..., X1 ), (m) (m) 0 φm = (φn,1 , ..., φ(m) n,n ) (m) y donde φn,1 , ..., φn,n se obtienen de resolver Γn φm = γ n,m donde γ n,m = (γm , γm+1 , ..., γm+n−1 )0 53 (4.1.1) y  γ0 γ1 .. .   Γn =   γ1 γ0 .. . . . . γn−1 . . . γn−2 .. .. . . . . . γ0 γn−1 γn−2      Cuando Γn es singular, el sistema puede tener infinitas soluciones pero todas ellas dan lugar al mismo e1,n (m). En el caso donde Γn es no singular, la soluci´on del sistema es predictor X φ(m) = Γ−1 n γ n,m (4.1.2) y en este caso e1,n (m)) = E[(Xn+m − X e1,n (m))2 ] ECM (X = E[(Xn+m − X 0n φ(m) )2 ] 0 0 = E[Xn+m Xn+m ] + φ(m) E[X n X 0n ]φ(m) − 2E[Xn+m X 0n ]φ(m) 0 = γ0 + φ(m) Γn φ(m) − 2γ 0n,m Γ−1 n γ n,m = γ0 − γ 0n,m Γ−1 n γ n,m Desde el punto de vista pr´ actico es poco u ´til cuando se posee una gran cantidad de datos, por eso existen algunos algoritmos que permiten superar estos problemas. Uno de ellos es el algoritmo de Durbin-Levinson que establece un procedimiento recursivo para calcular el vector de coeficientes φ. e1,n (m) como combinaci´on lineal de las innovaciones Otro algoritmo se basa en la representaci´on de X e Xn−j+m − X1,n−j (m) para j = m, ..., n + m − 1, es decir mediante la determinaci´on de los coeficientes θn−1+m,j con j = m, ..., n + m − 1 tales que se verifica e1,n (m) = X n+m−1 X e1,n−j (m)) θn−1+m,j (Xn−j+m − X j=m y tiene como ventaja que sirve tambi´en para procesos no estacionarios cuya funci´on de autocovarianza γij no depende u ´nicamente de la distancia |i − j| entre los tiempos. 4.2. Predicci´ on de un ARMA conocido todo su pasado El mejor predictor lineal de Xn+m conocido los valores de {Xt } con −∞ < t ≤ n se lo denota en (m). Los resultados te´ por X oricos que se obtengan de esta secci´on servir´an para estimaciones que se hagan cuando se conocen solo los valores de X1 , X2 , ..., Xn . Consideremos una serie {Xt } que corresponde a un proceso ARMA(p,q) estacionario e invertible dado por Xt − φ1 Xt−1 − ... − φp Xt−p = At + θ1 At−1 + ... + θq At−q , donde {At } es un proceso de ruido blanco formado por variables aleatorias independientes entre s´ı. Como es un proceso estacionario (causal) se puede representar como un proceso MA(∞) de la forma Xt = ∞ X αj At−j ∀t ∈ Z (4.2.1) j=0 con α0 = 1 y P∞ 2 j=0 αj < ∞. En particular, Xn+m = ∞ X j=0 54 αj An+m−j (4.2.2) Entonces, en (m) = X ∞ X ∀t ∈ Z βj An−j (4.2.3) j=0 en (m) es funci´ pues X on lineal de {Xt } con −∞ < t ≤ n y que a la vez por (4.2.1) son funciones de {At } con t ≤ n. Para calcular los βi se usan (4.2.2) y (4.2.3) para llegar a la expresi´on en (m) = Xn+m − X m−1 X αi An+m−i + i=0 ∞ X (αm+j − βj )An−j j=0 y como las variables del ruido blanco est´an no correlacionadas entre s´ı resulta que 2 en (m))2 ] = σA E[(Xn+m − X m−1 X 2 αi2 + σA i=0 ∞ X (αm+j − βj )2 (4.2.4) j=0 que ser´a m´ınimo cuando la segunda sumatoria se anule, es decir cuando βj = αm+j para j = 0, 1, .... Por ello la f´ormula (4.2.3) se reduce a en (m) = X ∞ X αm+j An−j (4.2.5) j=0 y (4.2.4) queda 2 en (m)) = σA ECM (X m−1 X αi2 (4.2.6) i=0 que claramente aumenta a medida que sePdesea predecir un Xn+m m´as alejado de Xn pero igualmente 2 este crecimiento va a estar acotado por ∞ i=0 αi < ∞. en (m) dado por (4.2.5) coincide con el predictor o´ptimo X bn (m). Probemos que X Dado que ( 0, si j < m; E[An+m−j /Xt : −∞ < t ≤ n] = An+m−j , si m ≤ j pues los valores de Xt : −∞ < t ≤ n dependen solo de At : t ≤ n, entonces bn (m) = E[Xn+m /Xt : −∞ < t ≤ n] = X ∞ X αj E[An+m−j /Xt : −∞ < t ≤ n] j=0 = ∞ X αj An+m−j = j=m ∞ X αm+i An−i i=0 en (m) =X Si a Xt : −∞ < t ≤ n se le agrega con el devenir del tiempo la observaci´on Xn+1 puede volver Pse ∞ e a calcular el predictor lineal de Xn+m . Por la f´ormula (4.2.5) para m = 1, Xn (1) = j=0 α1+j An−j y por (4.2.2) para m = 1, Xn+1 = An+1 + ∞ X en (1) αi An+1−i = An+1 + X i=1 55 por lo que en (1) = An+1 Xn+1 − X (4.2.7) Luego, en+1 (m − 1) = X ∞ X αm−1+i An+1−i i=0 = αm−1 An+1 + ∞ X αm+j An−j j=0 en (1)) + X en (m) = αm−1 (Xn+1 − X en+1 (m − 1) modifica al predictor X en (m) ya que que puede interpretarse como que el predictor X incorpora informaci´ on y este cambio tiene que ver con el error cometido al estimar Xn+1 a partir de Xt : t ≤ n. 4.3. Predicci´ on de un ARIMA conocido todo su pasado Sea {Zt } un proceso ARIMA(p,k,q) donde φ(B)(1 − B)k Zt = θ(B)At . Supongamos que el proceso de ruido blanco {At } est´ a compuesto por variables aleatorias independientes entre s´ı y los polinomios p φ(x) = 1−φ1 x−...−φp x y θ(x) = 1+θ1 x+...+θq xq no tienen ra´ıces comunes, siendo las ra´ıces de θ(B) en m´odulo mayor que 1, y adem´ as φ(1) 6= 0. Supongamos tambi´en que conocemos Zt : −∞ < t ≤ n y trataremos de encontrar el predictor ´ optimo Zbn (m) de Zn+m . Por ser un proceso invertible, admite una expresi´on como un proceso AR(∞) de la forma β(B)Zt = At , donde (4.3.1) ∞ ∞ X X φ(B)(1 − B)k i β(B) = = βi B con |βi | < ∞ θ(B) i=0 que es lo mismo que Zt = At − ∞ X i=0 βi Zt−i (4.3.2) i=1 En particular, para t = n + m Zn+m = An+m − ∞ X βi Zn+m−i (4.3.3) i=1 Multiplicamos ambos miembros de (4.3.3) por el operador τ (B) = 1 + τ1 B + ... + τm−1 B m−1 y elegimos los τ1 , ..., τm−1 de tal forma que se hagan cero los coeficentes correspondientes a Zn+m−1 , ..., Zn+1 . O 56 sea, despejando Zn+m y agrupando los distintos Zt se tiene Zn+m = τ (B)An+m − m−1 X τi Zn+m−i − i=1 = τ (B)An+m − = τ (B)An+m − m−1 X βi (1 + i=1 j=1 = τ (B)An+m − ∞ X m−1 X i=1 j=1 i=j+1 m−1 X ∞ X i−1 X (τi + βi )Zn+m−i − (τi + βi )Zn+m−i − i−1 X i=2 j=1 (τi + βi + m−1 X ( i X βi−j τj Zn+m−i − τj B j )Zn+m−i j=1 βi Zn+m−i ∞ X βi Zn+m−i i=m ∞ X i−1 X (βi + i=m i=1 j=0 ∞ X m−1 X i=m βi−j τj Zn+m−i − βi−j τj )Zn+m−i − τj βi−j )Zn+m−i − βi (1 + i=m i=2 j=1 m−1 X ∞ X τj B j )Zn+m−i − m−1 X i=1 = τ (B)An+m − m−1 X βi−j τj )Zn+m−i − Zn+m−1 (τ1 + β1 ) j=1 ∞ m−1 X X ( τj βi−j )Zn+m−i i=m j=0 donde τ0 = β0 = 1 que es lo mismo que Zn+m = τ (B)An+m − m−1 X (m) βi Zn+m−i − i=1 donde (m) βi ∞ X (m) βi Zn+m−i (4.3.4) i=m (P i si 1 ≤ i ≤ m − 1; j=0 τj βi−j , = Pm−1 j=0 τj βi−j , si i ≥ m y como se quiere que los coeficientes de Zn+m−1 , ..., Zn+1 sean cero, i X τj βi−j = 0 para i = 1, ..., m − 1 j=0 Es decir, se eligen τ1 , ..., τm−1 del operador τ (B) = 1 + τ1 B + ... + τm−1 B m−1 tal que resuelvan      1 0 0 ... 0 τ1 β1  β1     1 0 ... 0     τ2   β2   β2     β1 1 . . . 0   τ3  =  β3     ..  .   .  .. .. . . ..   . . .   ..   ..  . . βm−1 βm−2 βm−3 βm−4 . . . 1 τm−1 Entonces, (4.3.4) se reduce a Zn+m = An+m + m−1 X τj An+m−j − j=1 ∞ X (m) βi Zn+m−i (4.3.5) i=m y como E[As /Zt : t ≤ n] = 0 si s > n entonces bm = E[Zn+m /Zt : t ≤ n] = − Z ∞ X (m) βi Zn+m−i (4.3.6) i=m y tiene ECM de predicci´ on: bm ) = E[(Zn+m − Zbm )2 ] = (1 + τ 2 + ... + τ 2 )σ 2 ECM (Z 1 m−1 A 57 (4.3.7) 4.4. Ejercicios del Cap´ıtulo IV Ejercicio 4.4.1. Sea {Xt } un ARMA(p,q) estacionario (causal) de ecuaci´ on φ(B)Xt = θ(B)At , donde 2 y sea X = α(B)A la expresi´ {At } es un proceso de ruido blanco con varianza σA on de este proceso t t P∞ P∞ 2 θ(B) i como un M A(∞), es decir α(B) = φ(B) = i=0 αi B con α0 = 1 y i=0 αi < ∞. Considere el en (m) de Xn+m basado en Xt : −∞ < t ≤ n dado por predictor lineal ´ optimo X (P ∞ i=0 αm+i An−i , si m > 0; en (m) = X Xn+m , si − (n − 1) < m ≤ 0 a) Demuestre que (P q−m en (m − p) = en (m − 1) − ... − φp X en (m) − φ1 X X i=0 θm+i An−i , 0, si q ≥ m; si q < m (use la expresi´ on v´ alida para k ∈ N, φ(B)αk = θk donde se supone que αj = 0 si j < 0 y θk = 0 si k > q). en−i−1 (1) b) En la expresi´ on obtenida en la parte (a), reemplace An−i por la innovaci´ on Xn−i − X y compruebe que la f´ ormula resultante es la misma que se obtiene de la ecuaci´ on φ(B)Xn+m = en (k) si k > 0, An−k por 0 si k > 0 y An+k por la innovaci´ θ(B)An+m reemplazando Xn+k por X on e Xn+k − Xn+k−1 (1) si k > 0. en (k) para m = 1 y m = 2, en el c) Use la f´ ormula obtenida en la parte (a) para hallar la expresi´ on X caso en que {Xt } es un AR(p) estacionario (causal) de ecuaci´ on φ(B)Xt = At siendo X1 , ..., Xn las observaciones con n > p. Soluci´ on. a) Sea m > p, en (m) − φ1 X en (m − 1) − ... − φp X en (m − p) = X ∞ X (αm+i − i=0 Como αk − p X φi αk−i i=1 entonces la primera ecuaci´ on queda (P q−m i=0 0, p X φj αm+i−j )An−i j=1 ( 0, si k ≥ q + 1 = θk , si 0 ≤ k ≤ q θm+i An−i , si q ≥ m si q < m b) en (m) − φ1 X en (m − 1) − ... − φp X en (m − p) φ(B)Xn+m = X = θ(B)An+m = An+m + q X θj An+m−j j=0 q−m X = θm+j An−j j=−m que es igual a (P q−m i=0 0, en−i−1 (1)), si q ≥ m θm+i (Xn−i − X si q < m como se quer´ıa. 58 en (m) − φ1 X en (m − 1) − ... − φp X en (m − p) = 0 c) En un AR(p), X en (1 − p) = φ1 Xn + ... + φp Xn−p+1 en (0) + ... + φp X en (1) = φ1 X Si m = 1, X Si m = 2, en (2) = φ1 X en (1) + ... + φp X en (2 − p) X en (1) + ... + φp Xn−p+2 = φ1 X = φ1 (φ1 Xn + ... + φp Xn−p+1 ) + ... + φp Xn−p+2 = (φ21 + φ2 )Xn + (φ1 φ2 + φ3 )Xn−1 + ... + (φ1 φp + φp−1 )Xn−p+1 + φp Xn−p+2 Ejercicio 4.4.2. Sea {Xt } un proceso de promedios m´ oviles de orden q que satisface la ecuaci´ on Xt = 2 θ(B)At , donde {At } es un proceso de ruido blanco con varianza σA y donde θ(B) = 1+θ1 B +...+θq B q con θq 6= 0. a) Demostrar que (P q−m en (m) = X i=0 si m ≤ q si m > q θm+i An−i , 0, b) Mostrar que si q ≥ m, en (m) = X q−m X en−i−1 (1)), θm+i (Xn−i − X i=0 en (m) lo que constituye un procedimiento iterativo que permite calcular el mejor predictor lineal X de Xn+m a partir de Xt : t ≤ n. Soluci´ on. a) Este es un caso particular (p P = 0) de lo que se vio para un proceso ARMA donde el en (m) = ∞ αm+i An−i con αj = θj si 0 ≤ j ≤ q y αj = 0 si estimador ´ optimo resulta ser X i=0 j > q. Entonces, si m + i ≤ q con i = 0, ..., q − m, o sea m ≤ q queda en (m) = X q−m X αm+i An−i = i=0 q−m X αm+i An−i i=0 ; y si m + i > q con i = 0, ..., q − m, es decir q < m resulta que en (m) = 0 X como se quer´ıa. en (1) = Pq−1 θ1+i An−i y Xt = Pq θi At−i entonces b) Como X i=0 i=0 Xn−i = q X θj An−i−j = An−i + j=0 An−i q X j=1 en−i−1 (1) = An−i + X en−i−1 (1) = Xn−i − X en−i−1 (1)). en (m) = Pq−m θm+i (Xn−i − X Por lo tanto, X i=0 59 θj An−i−j Ejercicio 4.4.3. Sea {Xt } un proceso autorregresivo de orden p estacionario (causal) de ecuaci´ on Φ(B)Xt = At , (4.4.1) 2 y donde Φ(B) = 1 + φ B + ... + φ B p con donde {At } es un proceso de ruido blanco con varianza σA 1 p φp 6= 0. a) Halle la expresi´ on de los coeficientes αi : i ∈ N que corresponden a la representaci´ on de este proceso como un M A(∞), es decir a la representaci´ on donde α(B) = 1 Φ(B) Xt = α(B)At P P∞ 2 i = ∞ i=0 αi < ∞. i=0 αi B con α0 = 1 y (4.4.2) b) Suponga que por un error de especificaci´ on se ha considerado que el proceso {Xt } admite la representaci´ on θ(B)Xt = Et , (4.4.3) donde θ(B) = 1 − θ1 B − ... − θq B q con θq 6= 0 y q 6= p. En consecuencia, {Et } no es un proceso de ruido blanco, pero Xt admite a´ un la representaci´ on Xt = δ(B)Et , (4.4.4) P∞ 2 P 1 i donde δ(B) = θ(B) = ∞ i=0 δi < ∞. Mostrar que en este caso {Et } es i=0 δi B con δ0 = 1 y en realidad un proceso ARMA(p,q); para ello halle la expresi´ on Et en funci´ on de θ(B), Φ(B) y {At }. c) Suponga que en el modelo equivocado (4.4.3) se calcula el mejor predictor lineal de Xn+m . Si designamos a este predictor con la notaci´ on Xn∗ (m) muestre que se verifica que Xn∗ (m) = θ1 Xn∗ (m − 1) + ... + θq Xn∗ (m − q). (4.4.5) d) Use (4.4.3) y (4.4.5) para mostrar que el error de predicci´ on Yn∗ (m) = Xn+m − Xn∗ (m) satisface la ecuaci´ on en diferencias Yn∗ (m) − θ1 Yn∗ (m − 1) − ... − θq Yn∗ (m − q) = En+m . (4.4.6) Yn∗ (m) = En+m + δ1 En+m−1 + ... + δm−1 En+1 (4.4.7) e) Compruebe que es soluci´ on de la ecuaci´ on en diferencias (4.4.6) hallada en la parte (d). f ) Use (4.4.4) y (4.4.7) para mostrar que el predictor lineal Xn∗ (m) est´ a dado por Xn∗ (m) = ∞ X δm+i En−i , i=0 y en consecuencia se verifica que Yn∗ (m) = Y ∗ n − 1(m + 1) − δm En ∗ (1)). (recordar que En se puede expresar como En = Xn − Yn−1 60 en (m) el mejor predictor lineal de Xn+m usando el modelo correcto (4.4.1) y sea g) Sea ahora X en (m) el error de predicci´ Yn (m) = Xn+m − X on usando este modelo. Compruebe que el error de predicci´ on cometido usando el modelo equivocado (4.4.3) se puede expresar como en (m)) + (X en (m) − X ∗ (m)), Yn∗ (m) = (Xn+m − X n siendo los dos sumandos entre par´entesis no correlacionados entre s´ı y, en consecuencia, en (m) − X ∗ (m)). var(Yn∗ (m)) = var(Yn (m)) + var(X n Soluci´ on. a) Se buscan los αi tal que α(B)φ(B) = (1 + α1 B + α2 B 2 + α3 B 2 + ...)(1 − φ1 B − ... − φp B p ) = 1. Distribuyendo e igualando potencias se tiene: α0 φ0 = 1 ∗ 1 = 1 α1 φ0 − α0 φ1 = 0 entonces α1 = φ1 α0 −φ2 α0 − φ1 α1 + α2 φ0 = 0 entonces α2 = α0 φ2 + α1 φ1 −φ3 α0 − φ2 α1 − φ1 α2 + α3 φ0 = 0 entonces α3 = α0 φ3 + α1 φ2 + α2 φ1 ... Pk Pk−1 α φ = siguiendo as´ ı se llega a que si 1 ≤ k ≤ p, α = i k−i k j=1 φj αk−j y si k > p, i=0 P αk = pj=1 φj αk−j . 1 b) Et = θ(B)Xt = θ(B) φ(B) At por lo tanto {Et } satisface la ecuaci´on φ(B)Et = θ(B)At por lo que es un proceso ARMA(p,q). P c) Con el modelo equivocado Xn∗ (s) =P ∞ j=0 δs+j En−j ∀s ∈ N donde los δi resultan de resolver la q ecuaci´on δ(B)θ(B) = 1 o sea δk = j=1 θj δk−j si k > q. Luego, Xn∗ (m) − q X θi Xn∗ (m ∞ X − i) = i=1 j=0 ∞ X = δm+j En−j − i=1 (δm+j − j=0 pues m + j ≥ m > q y entonces δm+j = q X q X θi ∞ X δm+j−i En−j j=0 θi δm−i+j )En−j = 0 i=1 Pq i=1 θi δm−i+j . d) Yn∗ (m) = Xn+m − Xn∗ (m) entonces Yn∗ (m) − q X θi Yn∗ (m − i) = Xn+m − Xn∗ (m) i=1 − q X θi (Xn+m−i − Xn∗ (m − i)) i=1 = Xn+m − q X θi Xn+m−i − i=1 Xn∗ (m) + q X θi Xn∗ (m − i) i=1 q = (1 − θ1 B − ... − θq B )Xn+m = En+m e) Yn∗ (m) − q X i=1 θi Yn∗ (m − i) = En+m + m−1 X δj En+m−j − j=1 q X j=1 61 θi (En+m−i + m−i−1 X j=1 δj En+m−i−j ) = θ(B)En+m + = θ(B)En+m + ∞ X j=1 ∞ X δj En+m−j − Xn∗ (m) − q X ∞ X θi δj En+m−j−i + i=1 j=1 δj (En+m−j − j=1 q X q X θi Xn∗ (m − i) i=1 θi B i En+m−j ) i=1 ∞ X = θ(B)En+m + θ(B)( δj En+m−j − En+m ) j=0 = θ(B)δ(B)En+m = En+m f) Xn∗ (m) = Xn+m − Yn∗ (m) = = ∞ X δj En+m−j − j=0 ∞ X m−1 X δj En+m−j = j=m j=0 ∞ X δj En+m−j δi+m En−i i=0 Por lo tanto, ∗ Yn∗ (m) − Yn−1 (m + 1) = Xn+m − ∞ X δi+m En−i − (Xn+m − i=0 =− ∞ X ∞ X δi+m+1 En−i−1 ) i=0 δi+m En−i + i=0 ∞ X δj+m En−j = −δm En j=1 g) en (m)) + (X en (m) − X ∗ (m)) Yn∗ (m) = Xn+m − Xn∗ (m) = (Xn+m − X n Por un lado, en (m) = Xn+m − X ∞ X αj An+m−j − = αm+j An−i j=0 j=0 ∞ X ∞ X αi+m An−i − i=−m ∞ X j=0 αm+j An−i = −1 X i=−m = α0 An+m + α1 An+m−1 + ... + αm−1 An+1 62 αi+m An−i por otro lado, en (m) − X ∗ (m) = X n = = ∞ X j=0 ∞ X j=0 ∞ X αm+j An−j − αm+j An−j − αm+j An−j − j=0 ∞ X i=0 ∞ X i=0 ∞ X δi+m En−i δi+m δi+m i=0 ∞ X j=0 ∞ X βj An−i−j βk−i An−k k=i en (m) involucra a An+1 , ..., An+m−1 , An+m resulta que involucra a los At : t ≤ n y como Xn+m − X que los dos t´erminos que est´ an entre par´entesis no estan correlacionados. Por lo tanto, en (m))+var(X en (m)−Xn∗ (m)) = var(Yn (m))+var(X en (m)−Xn∗ (m)). var(Yn∗ (m)) = var(Xn+m −X Ejercicio 4.4.4. Sea {Zt } un proceso autorregresivo integrado de promedios m´ oviles de ´ ordenes (p,k,q) invertible de ecuaci´ on Φ(B)(1 − B)k Zt = θ(B)At , 2 y sea β(B)Z = A la expresi´ donde {At } es un proceso de ruido blanco con varianza σA on de este t t proceso como un AR(∞), es decir que ∞ X Φ(B)(1 − B)k = βi B i , β(B) = θ(B) i=0 P∞ bn (m) de Zn+m basado en los datos optimo Z con β0 = 1 y i=0 |βi | < ∞. Considere el predictor ´ Zt : t ≤ n dado por ∞ X (m) en (m) = − X βj Zn+m−j si m > 0, j=m donde (m) βj = m−1 X τi βj−i para j ≥ m, i=0 donde τ0 = 1 y donde τ1 , ..., τm−1 se obtienen del sistema   1 0 0 ... 0 τ1  β1   τ2 1 0 . . . 0    β2  β1 1 ... 0     τ3  ..  .. .. . . .. .    . . .  .. . . βm−2 βm−3 βm−4 . . . 1 τm−1         =     β1 β2 β3 .. .        βm−1 bn (m), que mide el error de predicci´ a) Muestre que la variable Yn (m) = Zn+m − Z on, est´ a dado por Yn (m) = m−1 X τi An+m−i si m > 0. i=0 b) Calcule la varianza de la variable Yn (m). 63 c) Compruebe que para m fijo mayor que 1 se verifica que {Yn (m) : n ∈ N} es un proceso de promedios m´ oviles de orden m − 1 y que cuando m = 1 resulta ser un proceso de ruido blanco. d) Demuestre que para m ≥ 1 se verifica que ( 0, cov(Yk (m), Yh (m)) = Pm−1−|k−h| 2 σA τj τj+|k−h| , j=0 si m ≤ |k − h| si m > |k − h| e) Demuestre que para m, s ≥ 1 se verifica que min{m,s}−1 cov(Yn (m), Yn (s)) = X 2 σA τj τj+|m−s| . j=0 f ) Compruebe que se verifica que 2 var(Yn (m + k)) = var(Yn (m)) + σA k−1 X 2 τn+i . i=0 g) Deduzca del resultado de la parte (f ) que var(Yn (m)) como funci´ on de m es mon´ otona creciente y que, en consecuencia, la longitud L(m) del intervalo de predicci´ on del (1 − α) × 100 % de confianza para Zn+m es tambi´en una funci´ on mon´ otona creciente de m. bn (m), Yn (m)) = 0 y que por lo tanto se verifica que var(Zn+m ) = var(Yn (m))+ h) Compruebe que cov(Z bn (m)). var(Z i) Compruebe que Yn−1 (m + 1) = Yn (m) + τm An y deduzca de esto que bn (m) = Z bn−1 (m + 1) + τm (Zn − Z bn−1 (1)). Z (4.4.8) j) Notar que el resultado (4.4.8) obtenido en la parte (i) anterior, constituye un procedimiento que bn−1 (m + 1) es el predictor de Xn+m a partir de los sirve para actualizar predictores, ya que Z b datos hasta el tiempo n − 1 y Zn (m) es el predictor de Zn+m cuando se cuenta con el dato adicional Zn . Analice como act´ ua este procedimiento cuando el nuevo valor observado Zn es b menor (mayor) que lo que Zn−1 (1) predijo. (m) j=m βj Zn+m−j P∞ bn (m) = Zn+m + Soluci´ on. a) Yn (m) = Zn+m − X vio en las p´ aginas 49 y 50. = Pm−1 i=0 τi An+m−i como se b) m−1 X 2 var(Yn (m)) = E[Yn (m) ] = E[( τi An+m−i )2 ] i=0 m−1 X = E[ = τi2 A2n+m−i + i=0 m−1 X τi2 E[A2n+m−i ] i=0 X = i 6= j τi τj An+m−i An+m−j ] m−1 X 2 τi2 σA i=0 P c) Si m > 1, por la parte (a), Yn (m) = m−1 i=0 τi An+m−i que claramente es un proceso M A(m − 1). Si m = 1, Yn (1) = An+1 o sea es un proceso de ruido blanco. 64 Pm−1 P d) Se sabe que cov(Yk (m), Yh (m)) = E[Yk (m)Yh (m)] = m−1 j=0 τi τj E[Ak+m−i Ah+m−j ] i=0 Si m ≤ |k − h|, |k + m − i − h − m + j| = |(k − h) + (j − i)| ≥ |k − h| − |j − i| > m + (1 − m) = 1 entonces E[Ak+m−i Ah+m−j ] = 0. Por otra parte, supongamos sin p´erdida de generalidad que k < h. Si |k + m − i − h − m + j| = |(k − h) + (j − i)| = 0 entonces 0 ≤ i = j + k − h ≤ m − 1, es decir 0 ≤ i ≤ m − 1 − (h − k) y 2 como 0 ≤ m − 1 − (h − k) ≤ m − 1 (o sea |k − h| < m). Por lo tanto, E[Ak+m−i Ah+m−j ] = σA se quer´ıa. Pm−1 P e) Como se sabe cov(Yn (m), Yn (s)) = m−1 j=0 τi τj E[An+m−i An+s−j ]. i=0 Supongamos primero que m > s entonces |n + m − i − (n + s − j)| = |m − s + j − i| = 0 y luego m − s ≤ m − 1 por lo que 0 ≤ j ≤ s − 1 y se tiene cov(Yn (m), Yn (s)) = Ps−1 0 ≤ i = j + 2 j=0 τj τj+(m−s) σA . Supongamos ahora que s ≥ m entonces |n + m − i − (n + s − j)| = |m − s + j − i| = 0 y luego 0 ≤ j = i − m + s ≤Pm − 1 por lo que m − s ≤ 0 ≤ j ≤ m − 1 + m − s ≤ m − 1 y se tiene 2 cov(Yn (m), Yn (s)) = m−1 i=0 τi τi−(m−s) σA . f) var(Yn (m + k)) = cov(Yn (m + k), Yn (m + k)) = m+k−1 X 2 τj2 σA j=0 = m−1 X 2 τj2 σA + j=0 m+k−1 X 2 τj2 σA j=m = var(Yn (m)) + 2 σA k−1 X 2 τi+m i=0 g) Si m < l, existe un k ∈ N tal que l = m + k. Entonces, 2 var(Yn (l)) = var(Yn (m + k)) = var(Yn (m)) + σA k−1 X 2 τi+m ≥ var(Yn (m)). i=0 Por lo tanto, es mon´ otona creciente. La parte de este ´ıtem que corresponde a longitud del intervalo de predicci´on no se desarrollar´ a ya que es necesario conocer sobre ruido blanco gaussiano que no se desarrollo en este trabajo. Pm−1 h) Por (a), Yn (m) = an involucrados An+1 , ..., An+m . Mientras que i=0 τi An+m−i , o sea est´ P∞ P∞ Pm−1 (m) b Zn (m) = − j=m βj Zn+m−j = − j=m i=0 τi βj−i Zn+m−j que involucra a los At : t ≤ n bn (m) no est´an correlacionados y cov(Zbn (m), Yn (m)) = 0. pues β(B)Zt = At por lo que Yn (m) y Z Adem´as, bn (m)) var(Zn+m ) = var((Zn+m − Zbn (m)) + Zbn (m)) = var(Yn (m) + Z bn (m)). = var(Yn (m)) + var(Z i) Yn (m) − Yn−1 (m + 1) = m−1 X τi An+m−i − i=0 m X i=0 65 τi An+m−i = −τm An . bn−1 (m + 1) y por el resultado Adem´as, como Yn (m) = Zn+m − Zbn (m), Yn−1 (m + 1) = Zn+m − Z anterior resulta que bn (m) = Zbn−1 (m + 1) + τm An Z que debido a que Zbn−1 (1) = − ∞ X (1) βj Zn−j j=1 =− ∞ X βj Zn−j = −( j=1 ∞ X βj Zn−j − β0 Zn ) j=0 = −(An − Zn ) = Zn − An la f´ormula queda bn−1 (1)) Zbn (m) = Zbn−1 (m + 1) + τm (Zn − Z bn (m) = Z bn−1 (m + 1) + τm (Zn − Z bn−1 (1)) nos dice que el predictor Zbn−1 (m + 1) de j) La relaci´ on Z Zn+m se modifica al agregar el conocimiento del dato Xn y esta modificaci´on depender´ a de τm y del error al predecir Zn pero no se puede analizar el cambio de acuerdo a si Zn es mayor (o menor) que Zbn−1 (1) pues este tambi´en depender´a de como sea τm . 66 Cap´ıtulo 5 Estimaci´ on e identificaci´ on de los modelos de series de tiempos 5.1. Identificaci´ on de procesos. Determinaci´ on del orden y verificaci´ on de diagn´ ostico Las t´ecnicas de identificaci´ on permiten decidir cual es el proceso que mejor se ajusta a los datos observados. Entre estas t´ecnicas est´ a el estudio del comportamiento de las funciones de autocorrelaci´ on muestral y autocorrelaci´ on parcial muestral. Una manera de verificar que el modelo diagnosticado es el correcto es analizando el residuo bt = Xt − X bt A bt es la estimaci´ bt se comporta como un proceso de ruido donde X on hecha con el modelo elegido. Si A blanco es porque el modelo elegido es el adecuado. Por otro lado, si con los datos X1 , ..., Xn se ha bt = Xt − X bt son estimado un proceso ARMA(p,q) de ecuaci´on Φp (B)Xt = θq (B)Et y los residuos E tales que se comportan como un proceso ARMA(r,s), entonces el modelo adecuado para los datos originales es un ARMA(p+r,q+s). Un test basado en la funci´ on de autocorrelaci´on muestral ρbAb(k) de los residuos que sirven para verificar la hip´ otesis nula H0 : ρA (1) = ρA (2) = ... = ρA (m) = 0 es el llamado Test de Portmanteau de falta de ajuste, que usa la estad´ıstica P = n(n + 2) m X k=1 1 ρb2 (k) n − k Ab que tiene una distribuci´ on aproximadamente Chi-cuadrado con m − e grados de libertad, donde e es el n´ umero de par´ ametros estimados en el modelo. 5.1.1. Criterios de Akaike: FPE, AIC, BIC Akaike propuso en 1969 el m´etodo del Error de Predicci´on Final (FPE) que sirve para determinar el orden del modelo AR que mejor ajusta a los datos. Consiste en calcular: F P E(k) = n+k 2 σ b n−k A (5.1.1) 2 es la estimaci´ donde n es el n´ umero de observaciones y σ bA on m´axima veros´ımil de la varianza residual 2 b b que est´a dada por σ bA = γ b0 − φ1 γ b1 − ... − φp γ bp y se elige el orden p como aquel valor de k que minimiza F P E(k). Otro m´etodo propuesto por Akaike en 1971 y 1974, es el Criterio de Informaci´on de Akaike y se 67 lo indica por AIC. Si el n´ umero de par´ametros estimados independientemente es e entonces el AIC correspondiente est´ a definido por b + 2e AIC(e) = −2log(verosimilitud maximada) + 2e = −2L (5.1.2) En particular, para modelos ARMA con ruido blanco gaussiano, m 2 b = − m log(b σA )− L 2 2 (5.1.3) para m suficientemente grande, donde m es el n´ umero efectivo de datos considerados. En 1976, Shibata demostr´ o que el AIC tiende a sobreestimar el orden de la parte autorregresiva y en 1978 y 1979, Akaike desarroll´ o una extens´on bayesiana del AIC a la que denomin´o Criterio de Informaci´on Bayesiana (que se indica BIC) y demostr´o mediante simulaci´on que este m´etodo es menos proclive a sobreestimar el valor de p. Este est´a definido por 2 BIC(e) = nlog(b σA ) − (n − e)log( σ b2 − σ b2 n−e ) + elog(n) + elog( X 2 A ) n eb σA (5.1.4) 2 es la donde e es el n´ umero de par´ ametros estimados, n el n´ umero total de variables observadas, σ bX 2 es la estimaci´ varianza muestral de la serie y σ bA on m´aximo veros´ımil de la varianza residual. El m´etodo consiste en elegir el modelo que minimiza BIC(e). Hurvich y Tsai propusieron en 1989 una versi´on corregida del criterio de Akaike original (se lo de2e.n nomin´o AICC) que toma el sumando 2e de AIC y lo reemplaza por n−e−1 que penaliza a´ un m´ as la incorporaci´on de nuevos par´ ametros al modelo, eliminando de esa forma la sobreestimaci´on del criterio de Akaike. 5.2. Reglas pr´ acticas para la identificaci´ on de procesos En la pr´actica, los m´etodos de identificaci´on y estimaci´on propuestos en la secci´on anterior puede convertirse en una tarea larga y complicada. A continuaci´on se esbozar´an los pasos a seguir para simplificar la b´ usqueda: 1. Es importante entender que es lo miden los datos, cual es el significado dentro de la disciplina en la que se encuadran y conocer dentro de lo posible los m´etodos con que ellos se obtuvieron. Esto sirve para conocer ciertos aspectos sobre el comportamiento de los datos y tambi´en es u ´til como elemento de juicio entre modelos diferentes igualmente aceptables seg´ un los criterios de identificaci´ on. 2. El an´alisis de la representaci´ on gr´afica de los datos nos permitir´a observar si se satisfacen las condiciones b´ asicas de estacionariedad: media y varianza constantes en el tiempo. Tambi´en se puede descubrir la presencia de periodicidades propias de las componentes estacionales. Si existe una tendencia debe observarse si la dispersi´on en los datos var´ıa con la tendencia. En ese caso, deben transformarse los datos para estabilizar la varianza. 3. Una vez estabilizada la varianza se debe observar si hay tendencia y/o estacionalidad, en cuyo caso se aplicar´ a algunos de los m´etodos discutidos en la secci´on 3: m´etodo cl´asico o m´etodo de la diferenciaci´ on. En este u ´ltimo caso se aplicar´a el operador diferencia de tama˜ no d hasta que se elimine la estacionalidad y luego se usar´a el operador diferencia tantas veces como sea necesario para deshacerse de la tendencia. 68 4. Cuando los datos no presenten desviaciones evidentes de estacionariedad, se intentar´a ajustar por un modelo ARMA. Las funciones de autocorrelaci´on muestral y autocorrelaci´on parcial muestral permitir´ an seleccionar de forma descriptiva los posibles modelos que ajustan observando los coeficientes de correlaci´ on que son significativamente distintos de cero. Por lo general se considera que un valor es significativo si en m´odulo supera el valor √2n . Esto proviene de suponer que si las variables et no estuvieran correlacionadas entre s´ı entonces los coeficientes de correlaci´ on tendr´ıan una distribuci´ on muestral normal con media cero y una cierta desviaci´on est´ andar, esto es, que en el 95 % de los casos los coeficientes caer´an entre la media (= 0) m´as menos 1.96 de desviaci´ on est´ andar (que se suele redondear a 2), y solamente el 5 % de los coeficientes de correlaci´ on se apartar´ an de una media cero en ±1,96 de desviaci´on est´andar. Es decir, el coeficiente de correlaci´ on es estad´ısticamente significativo si con un 95 % de confianza la media de las correlaciones es cero. Una vez elegidos, se estiman los par´ametros de cada uno de ellos. Es conveniente establecer intervalos de confianza para los coeficientes, pues cuando el intervalo de alg´ un coeficiente incluye el cero se puede considerar como modelo alternativo el que corresponde a forzar en el modelo original que dicho coeficiente sea cero. 5. La decisi´ on sobre el modelo que mejor ajusta a los datos se realiza usando los criterios de identificaci´ on de procesos vistos en la secci´on anterior. El criterio corregido de Akaike y el criterio de informaci´ on bayesiana son los que mejor desempe˜ no tienen. Sin embargo, cuando hay m´as de un modelo adecuado seg´ un estos criterios se recurre a tomar consideraciones como el an´alisis de los residuos, simplicidad del modelo y facilidad en su interpretaci´on desde el punto de vista de la disciplina en la que se est´a estudiando la serie de tiempo en consideraci´on. 69 Cap´ıtulo 6 Aplicaci´ on Las im´agenes satelitales proporcionan diversas herramientas para discriminar entre masas vegetales en una regi´on fitogeogr´ afica, siendo el ´ındice diferencial normalizado de vegetaci´on, NDVI, una de ellas. De este indicador es posible obtener numerosas variables, tales como contenido de agua en las hojas, productividad neta de la vegetaci´ on, contenido de clorofila en la hoja, evapotranspiraci´on potencial, etc. El NDVI var´ıa entre -1 y +1, aunque es frecuente el tratamiento de valores reparametrizados del ´ındice; cuanto m´ as pr´ oximo a 1, mayor resulta el vigor vegetal (Chuvieco [6]). Diferentes autores han estudiado series temporales de este ´ındice. Por ejemplo, Barbosa et al [5] emplearon una serie mensual correspondiente a 20 a˜ nos del NDVI en la regi´on noreste de Brasil para estudiar el comportamiento temporal del ´ındice en la regi´ on. En este cap´ıtulo nos propusimos encontrar un modelo ARIMA para explicar el comportamiento del NDVI en la selva Misionero-Brasile˜ na. Esta zona corresponde a la regi´on fitogeogr´afica Paranaense (Cabrera [4]), que comprende la provincia de Misiones y gran parte de la provincia de Corrientes. El estudio se bas´ o en una serie de tiempo mensual de 17 a˜ nos de NDVI de la zona (per´ıodo 1982-1998). Los datos fueron provistos por CONAE (Comisi´on Nacional de Actividades espaciales) y fueron , producidas por el sensor AVHRR (Advanced Very High Resolution Radiometer) que se encuentra a bordo de sat´elites meteorol´ ogicos de la serie NOAA (National Oceanic and Atmospheric Administration). Este sensor genera im´ agenes en cinco canales o bandas, cuyas combinaciones posibilitan la determinaci´on de diferentes indicadores de propiedades biof´ısicas, entre ellos el NDVI que se calcula con datos de los canales 1 y 2 del sensor, seg´ un: N DV I = Ch2 − Ch1 Ch1 + Ch2 donde Ch1 y Ch2 son los valores provistos por las bandas 1 y 2 del AVHRR, respectivamente. El Gr´afico 6.1 muestra los datos de la serie de NDVI objeto de an´alisis y el 6.2, el gr´afico de las funciones de autocorrelaci´ on simple y parcial de la serie. 70 Gr´ afico 6.1: Serie de tiempo de NDVI Gr´ afico 6.2: Funci´ on de autocorrelaci´ on simple de la serie (izq.). Funci´ on de autocorrelaci´ on parcial de la serie (der.) De la simple observaci´ on de ρk y πk (Gr´afico 6.2) se dificulta la propuesta de un modelo ARIMA que ajuste razonablemente los datos de la serie, ya que ambas funciones muestran correlaciones significativas (fuera de las bandas de confianza) para rezagos altos. Por otro lado, en la Gr´afico 6.1 se puede apreciar a partir del dato 115 un aumento en la tendencia que se detiene aproximadamente en el dato correspondiente a k = 150. Ello nos alent´o a analizar por separado dos sub-series de la serie original; por un lado los datos entre 1 y 115, y por otro, los datos entre 116 y 228. Denotaremos a estas series como {z1t } y {z2t } respectivamente. A continuaci´ on presentamos el gr´ afico de {z1t }: Gr´ afico 6.3: Serie de tiempo {z1t } de NDVI en la regi´ on Paranaense. 71 Restamos a {z1t } su media para conseguir una serie de media cero y estabilizamos su varianza aplicando la transformaci´ on logaritmo (ver secci´on 3.3). Resulta entonces que las funciones de autocorrelaci´ on simple y parcial muestrales de la serie corregida corresponden a las exhibidas en el Gr´afico 6.4. A los efectos de no complicar en exceso la notaci´on, de ahora en m´as {z1t } designar´a la versi´on corregida de {z1t }. Gr´ afico 6.4: Funciones de autocorrelaci´ on simple y parcial de {z1t } A trav´es del criterio de Akaike Corregido (AICC) y usando el m´etodo de Yule-Walker, el programa R propuso un modelo AR(1) para la serie. Pero atendiendo a las funciones de autocorrelaci´ on se propone finalmente un modelo ARMA(1,1): Z1t = φZ1t−1 + θA1t−1 + A1t donde {A1t } es un proceso de ruido blanco y φ y θ son par´ametros a estimar. Estimamos φ y θ por el m´etodo de M´ınimos Cuadrados Condicionales: φb = 0, 3206 θb = 0,1198 El gr´afico 6.5, muestra el patr´ on de los residuos estandarizados correspondientes al modelo ajustado a la serie. r Gr´ afico 6.5 Residuos estandarizados para el ARMA(1,1) ajustado a la serie {z1t }. Si el modelo fuese adecuado, la serie generada a partir de la diferencia entre la serie observada y la estimada deber´ıa comportarse, una vez estandarizada, como un proceso de ruido blanco con media 0 y 72 varianza 1. Esta serie recibe el nombre de “serie de residuos estandarizados”, o simplemente “residuos estandarizados”. Un an´ alisis descriptivo de los residuos estandarizados se puede apreciar a trav´es de los siguientes gr´ aficos: Gr´ afico 6.6: Cuantiles de los residuos vs. cuantiles de una normal est´ andar. Gr´ afico 6.7: histograma de los residuos estandarizados. Gr´ afico 6.8: Funciones de autocorrelaci´ on simple y parcial muestrales de los residuos estandarizados. Se observa que los residuos estandarizados no est´an correlacionados, aunque el supuesto de normalidad parece no satisfacerse (es evidente la asimetr´ıa en el gr´afico 6.7). Teniendo en cuenta que {z1t } constituye el pasado “ m´ as remoto”de la serie original, concentraremos nuestra atenci´on en el an´ alisis de la serie {z2t }, la historia m´ as reciente de la serie de NDVI objeto de estudio (valores de la serie 73 de t = 116 a t = 228). El siguiente gr´ afico muestra el perfil de estos datos corregidos (se estabiliz´ o su varianza con la ayuda de la funci´ on logaritmo y se rest´o su media). Igual que procedimos con {z1t }, de ahora en adelante {z2t } denotar´ a a la serie {z2t } estabilizada o corregida: Gr´ afico 6.8: Serie {z2t } corregida. Gr´ afico 6.9: Funciones de autocorrelaci´ on simple y parcial de la serie {z2t } corregida. Analizando descriptivamente las autocorrelaciones, una opci´on ser´ıa quedarnos con un modelo autorregresivo y proponer un modelo AR(15). Cuando le sugerimos a R este modelo, en base al criterio de Akaike Corregido (AICC) y usando el m´etodo de Yule-Walker, el programa R propone un modelo AR(2). Sin embargo la funci´ on de autocorrelaci´on parcial de la serie no muestra evidencias de que la autocorrelaci´ on parcial de orden 2 sea estad´ısticamente significativa, por lo que no se justificar´ıa incluir un t´ermino autorregresivo de orden 2 en el modelo. Se pens´o entonces en un modelo ARMA(1,4), atendiendo a las autocorrelaciones simples y parciales m´as notables fuera de las bandas en el gr´ afico 6.9: Z2t = φZ2t−1 + θ1 A2t−1 + θ2 A2t−2 + θ3 A2t−3 + θ4 A2t−4 + A2t donde {A2t } es un proceso de ruido blanco y φ y θ1 , θ2 , θ3 , θ4 son par´ametros a estimar. Cabe mencionar en este punto que un modelo constituye una representaci´on de la realidad. Pueden formularse diferentes modelos y el hecho de quedarnos con uno, no implica que no exista otro que pueda superarlo sobre la base de ciertos criterios o requisitos que se desee satisfacer, como por ejemplo la parsimonia. 74 Las estimaciones de los par´ ametros involucrados en el modelo ARMA(1,4) propuesto se muestran a continuaci´on: φb = 0,6715, θb1 = 0,2416, θb2 = 0,0413, θb3 = 0,0649, θb4 = 0,0495 El gr´afico 6.10 muestra los residuos estandarizados una vez que se ajust´o a la serie z2t un modelo ARMA(1,4). Gr´ afico 6.10: Residuos estandarizados para el ARMA(1,4) ajustado a la serie {z2t }. El an´alisis descriptivo de los residuos muestra que el modelo es aceptable. Este estudio puede apreciarse a trav´es de los gr´ aficos 6.11, 6.12 y 6.13. Gr´ afico 6.11: Cuantiles de los residuos vs. cuantiles de una normal est´ andar. 75 Gr´ afico 6.12: histograma de los residuos estandarizados. Gr´ afico 6.13: Funciones de autocorrelaci´ on simple y parcial muestrales de los residuos estandarizados. Debido a que la serie original de NDVI (Gr´afico 6.1) presenta un cambio en la tendencia entre t = 115 y t = 150 aproximadamente, y teniendo presente que lo importante ser´a poder predecir valores futuros de la serie a partir de su pasado m´as reciente, consideramos s´olo los u ´ltimos 7 a˜ nos de la serie. A esta nueva serie la llamamos {zt } y sus valores corregidos y funciones de autocorrelaci´ on simple y parcial se muestran en los siguientes gr´aficos: Gr´ afico 6.14: Datos correspondientes a la serie {zt }. Gr´ afico 6.15: Funciones de autocorrelaci´ on simple y parcial de la serie {zt } corregida. 76 A partir de las correlaciones proponemos un modelo AR(15) que comparando a trav´es del criterio AICC y el m´etodo de Yule-Walker con modelos autorregresivos de menor orden, es desplazado por un modelo AR(1): Zt = φZt−1 + At donde {At } es un proceso de ruido blanco y φ un par´ametro a estimar. Se estim´o el valor de φ por el m´etodo de Yule-Walker, obteni´endose φb = 0,4243 Los gr´aficos 6.16 a 6.19, muestran los residuos estandarizados una vez que se ajust´o la serie {zt } un modelo AR(1) y un an´ alisis descriptivo de los mismos. Gr´ afico 6.16: Residuos estandarizados para el modelo AR(1) ajustado a la serie {zt }. Gr´ afico 6.17: Cuantiles de los residuos vs. cuantiles de una normal est´ andar. 77 Gr´ afico 6.18: histograma de los residuos estandarizados. Gr´ afico 6.19: Funciones de autocorrelaci´ on simple y parcial de los residuos para el ajuste AR(1) obtenido para la serie {zt } corregida. El an´alisis previo evidencia que el modelo propuesto es adecuado ya que se cumplen sobre los residuos los supuestos asumidos, con lo cual ser´ıa factible su utilizaci´on para realizar predicciones a corto y mediano plazo del NDVI de la regi´on Paranaense. 78 Bibliograf´ıa [1] Leiva, R. 1995. Introducci´ on al An´ alisis de las series de tiempo. Facultad de Ciencias Econ´ omicas de la Universidad Nacional de Cuyo. ISBN 9-503-90057-3 e-ISBN 978-9-503-90057-4. [2] Brockwell, P., Davis, R. 2002. Introduction to Time Series and Forecasting. Springer Texts in Statistics. ISBN 0-387-95351-5. Springer New York Dordrecht Heidelberg London. [3] Pe˜ na S´anchez de Rivera D. 1993. Estad´ıstica, modelos y m´etodos: Fundamentos 2. Modelos lineales y series temporales. Alianza Universal Textos. ISBN 84-206-8110-5 e-ISBN 978-8-420-689937. Alianza Editorial. ´ (1976) Regiones fitogeogr´aficas argentinas. En Kugler WF (Ed.) Enciclopedia ar[4] Cabrera, AL gentina de agricultura y jardiner´ıa. Tomo 2. 2a edici´on. Acme. Buenos Aires. Argentina. Fasc´ıculo 1. pp. 1-85. [5] Barbosa, HA, Huetea, AR, Baethgenb, WE (2006) A 20-year study of NDVI variability over the Northeast Region of Brazil. Journal of Arid Environments 67 288-307. [6] Chuvieco Salinero E (1996) Teledetecci´on Ambiental Capitulo 2, 57-59. Editorial Ariel SA. 79