STATISTICA MATEMATICA II · UNIVERSITA’ DEGLI STUDI DI MILANO STATISTICA MATEMATICA II Serie...

47
UNIVERSITA’ DEGLI STUDI DI MILANO STATISTICA MATEMATICA II Serie storiche Prof. Vincenzo CAPASSO Dipartimento di Matematica F. Enriques Via C. Saldini 50, 20133 Milano

Transcript of STATISTICA MATEMATICA II · UNIVERSITA’ DEGLI STUDI DI MILANO STATISTICA MATEMATICA II Serie...

UNIVERSITA’ DEGLI STUDI DI MILANO

STATISTICA MATEMATICA II

Serie storiche

Prof. Vincenzo CAPASSO

Dipartimento di Matematica F. Enriques

Via C. Saldini 50, 20133 Milano

Indice

1 Processi stocastici e serie temporali 3

2 Stazionarieta 4

3 Previsione 123.1 Predittori lineari . . . . . . . . . . . . . . . . . . . . . . . . . 133.2 Rappresentazione di serie storiche . . . . . . . . . . . . . . . . 16

3.2.1 Modelli MA . . . . . . . . . . . . . . . . . . . . . . . . 183.2.2 Modelli AR . . . . . . . . . . . . . . . . . . . . . . . . 19

3.3 Previsione nei modelli AR . . . . . . . . . . . . . . . . . . . . 243.3.1 AR(1) . . . . . . . . . . . . . . . . . . . . . . . . . . . 243.3.2 AR(p) . . . . . . . . . . . . . . . . . . . . . . . . . . . 26

3.4 Previsione nei modelli MA . . . . . . . . . . . . . . . . . . . . 293.5 Previsione nei modelli ARMA . . . . . . . . . . . . . . . . . . 323.6 La varianza . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34

4 Stima 364.1 Stima della media . . . . . . . . . . . . . . . . . . . . . . . . . 364.2 Stima della funzione di autocovarianza . . . . . . . . . . . . . 394.3 Stima della funzione di autocorrelazione . . . . . . . . . . . . 404.4 Stima dei parametri . . . . . . . . . . . . . . . . . . . . . . . . 41

4.4.1 AR(1) . . . . . . . . . . . . . . . . . . . . . . . . . . . 414.4.2 AR(p) . . . . . . . . . . . . . . . . . . . . . . . . . . . 414.4.3 MA(q) . . . . . . . . . . . . . . . . . . . . . . . . . . . 45

2

1 Processi stocastici e serie temporali

Una serie temporale e una collezione di osservazioni Xt : t = 1, 2, . . ., dovet indica l’istante nel quale la variabile Xt viene misurata. Per poter dare unadefinizione matematicamente rigorosa, e necessario introdurre la teoria deiprocessi stocastici.Sia (Ω,F , P ) uno spazio di probabilita. Una serie temporale si puo definirecome una realizzazione di un processo stocastico. Ricordiamo che

Definizione 1.1. Un processo stocastico e un insieme di variabili aleatorie(Xt)t∈T definite su (Ω,F , P ) a valori in (R,B(R)), dove T e un insieme diindici.

Se T e un intervallo di R, allora il processo stocastico si dice continuo esi indica con X(t) : t ∈ T; se T ⊆ N (Z), allora il processo stocastico sidice discreto e si indica con Xt : t ∈ T. (Per semplicita descriveremo soloprocessi stocastici a tempo discreto, anche se le definizioni e i risultati cheverranno presentati valgono anche nel caso continuo).Fissato t ∈ T , Xt(·) e una variabile aleatoria chiamata stato del processoall’istante t; fissato ω ∈ Ω, l’applicazione t ∈ T → Xt(ω) e una traiettoria (orealizzazione) del processo, cioe una serie temporale (si denota con X·(ω)).Dal punto di vista matematico, conoscere probabilisticamente un proces-so stocastico significa conoscere per ogni insieme di istanti (t1, . . . , tn) ladistribuzione congiunta di (Xt1 , . . . , Xtn).

Definizione 1.2. Sia T l’insieme di tutti i vettori t = (t1, . . . , tn)′ ∈ T n :t1 < . . . < tn, n = 1, 2, . . .. Si dicono funzioni di distribuzione (finito-dimensionali) del processo stocastico Xt, t ∈ T le funzioni Ft(·), t ∈ T definite per t = (t1, . . . , tn)′ da

Ft(x) = P(Xt1 ≤ x1, . . . , Xtn ≤ xn), x = (x1, . . . , xn)′ ∈ Rn.

A questo punto, la domanda naturale che ci si pone e la seguente: conoscen-do le funzioni di distribuzione Ft, ∀ t ∈ T , conosciamo probabilisticamenteil processo Xtt∈R?Si dimostra che, se specifichiamo in modo compatibile le distribuzioni di(Xt1 , . . . , Xtn) ∀ (t1, . . . , tn) ∈ T , allora il processo stocastico e completa-mente caratterizzato da questa famiglia.

3

Definizione 1.3. Siano u,v ∈ T con u = (u1, . . . , uk) ⊂ v = (v1, . . . , vn) esia πu,v : Rdimv → Rdimu la proiezione da Rdimv in Rdimu, definita da

πu,v(Fv)(x1, . . . , xk) = P ((Xv1 , . . . , Xvn) ∈ Rn : Xui≤ xi, i = 1, . . . , k).

Si dice che la famiglia Ftt∈T verifica la condizione di compatibilita se

Fu = πu,v(Fv) ∀ u,v ∈ T con u ⊂ v.

Dal momento che R e uno spazio polacco (cioe metrico, completo e sepa-rabile), possiamo applicare il

Teorema 1.1 (di Kolmogorov-Bochner). Si supponga di conoscere perogni t ∈ T la distribuzione Ft definita su Rn, dove n e il numero di elementidi t; si supponga che Ftt∈T verifichi la condizione di compatibilita. Alloraesiste un’unica distribuzione di probabilita FT : RT → [0, 1] t.c. ∀t ∈ TFt = πt,T (FT ).

Da questo teorema segue che le distribuzioni finito-dimensionali determi-nano univocamente la distribuzione dell’intero processo.

2 Stazionarieta

Una classe molto importante di processi stocastici e quella dei processi sta-zionari.

Definizione 2.1. Una serie temporale (Xt)t∈T si dice strettamente staziona-ria se ∀h ∈ R, ∀t1, . . . , tn ∈ T , ∀n ∈ N, ∀x = (x1, . . . , xn) ∈ Rn

FXt1 ,...,Xtn(x1, . . . , xn) = FXt1+h,...,Xtn+h

(x1, . . . , xn).

In altre parole, slittare gli istanti temporali di una quantita h, non alterala distribuzione congiunta.Vediamo ora le principali proprieta di una serie strettamente stazionaria:

1. prendendo n = 1 nella definizione precedente, e immediato osservareche la stretta stazionarieta implica che la distribuzione di Xt sia lastessa per ogni t. Infatti:∀t ∈ T , ∀h ∈ R, ∀x ∈ R

FXt(x) = FXt+h(x),

cioe le variabili Xt sono identicamente distribuite.

4

2. In particolare dalla proprieta precedente segue chese Xt ∈ L1(P) ∀t ∈ T , allora E[Xt] = E[Xt′ ] = µ ∀t, t′ ∈ T ,se Xt ∈ L2(P) ∀t ∈ T , allora Var(Xt) = σ2 ∀t ∈ T .Dunque, se un processo e strettamente stazionario, le funzioni valoreatteso e varianza sono costanti, cioe non dipendono dal valore di t.

3. Per n = 2, si ricava che la distribuzione congiunta di Xt1 e Xt2 dipendesolo da (t2 − t1):

FXt1 ,Xt2(x1, x2) = FXt1+h,Xt2+h

(x1, x2) =

= FX0,Xt2−t1(x1, x2) =

= Ft2−t1(x1, x2)

per h = −t1.

Esempio 2.1. Caso i.i.d.Consideriamo un processo stocastico di variabili aleatorie i.i.d.

Xt1 , . . . , Xtn , . . .

Allora tale processo e strettamente stazionario. Infatti ∀h ∈ R, ∀n ∈ N,∀x = (x1, . . . , xn) ∈ Rn, ∀t1, . . . , tn ∈ T :

FXt1 ,...,Xtn(x1, . . . , xn) = FXt1

(x1) . . . FXtn(xn) =

=n∏

j=1

FXt1(xj) =

= FXt1+h,...,Xtn+h(x1, . . . , xn).

Esempio 2.2. Processo gaussianoIl processo stocastico reale (Xt)t∈R+ si dice processo gaussiano se ∀n ∈ N∗ e∀ (t1, t2, . . . , tn) ∈ Rn

+ il vettore aleatorio n-dimensionale (Xt1 , Xt2 , . . . , Xtn)ha distribuzione gaussiana multivariata, cioe la sua densita di probabilita e:

fXt1 ,Xt2 ,...,Xtn(x) =

1

(2π)n/2√

det Kexp

(−1

2(x− b)T K−1(x− b)

)

dove bi = E[Xti ] per i=1,. . . ,n e K = (Kij) = (Cov(Xti , Xtj)) per i,j=1,. . . ,n.Si osservi che questo processo e univocamente determinato assegnando le

5

funzioni

b : t ∈ T → E[Xt] ∈ RK : (t1, t2) ∈ T × T → Cov(Xt1 , Xt2) ∈ R,

dove K e una funzione simmetrica e semidefinita positiva.Per la dimostrazione di questo enunciato utilizzeremo due proprieta delledistribuzioni gaussiane:

A. e sempre possibile costruire un vettore aleatorio gaussiano corrispon-dente ad una data matrice (n× n) simmetrica semidefinita positiva Ke un punto b∈ Rn;

B. se X e un vettore gaussiano n-dimensionale e Z = AX dove A e unamatrice m × n, allora Z e un vettore gaussiano m-dimensionale. Inparticolare, se X1, . . . , Xn e un vettore gaussiano, allora lo e ancheX1, . . . , Xm, m ≤ n.

Torniamo alla dimostrazione dell’enunciato. Dati t1, . . . , tn ∈ T , t1 < . . . < tn,sia (Xt1 , . . . , Xtn) un vettore aleatorio gaussiano con matrice di covarianza[K(tr, ts)], r, s = 1, . . . , n e valore atteso (b1, . . . , bn) (tale processo esiste perla proprieta A). Allora, per la B, anche (Xti1

, . . . , Xtij) e un vettore gaus-

siano con matrice di covarianza [K(tir , tis)], r, s = 1, . . . , j e valore atteso(bi1 , . . . , bij). Questo si accorda con la definizione iniziale della distribuzionecongiunta di (Xti1

, . . . , Xtij). Ma allora l’ipotesi di compatibilita richiesta

dal teorema di Kolmogorov-Bochner e soddisfatta, e il risultato segue imme-diatamente.

Si osservi infine che se b = cost = µ e se K(t1, t2) = K(t1 − t2), allora ilprocesso e fortemente stazionario.

L’esempio del processo gaussiano ci suggerisce di introdurre la seguente

Definizione 2.2. Un processo (Xt)t∈T in L2 si dice debolmente stazionariose

i. E[Xt] = µ, µ costante indipendente da t,

ii. E[Xt1Xt2 ] = E[Xt1+hXt2+h], ∀ t1, t2 ∈ T, ∀ h ∈ R.

6

Si osservi che dalle i. e ii. segue che:

E[Xt] = E[X0] ∀ t ∈ T

E[X2t ] = E[X2

0 ] ∀ t ∈ T

E[Xt1Xt2 ] = E[X0Xt2−t1 ] ∀ t1, t2 ∈ T.

Quindi ∀ t1, t2 ∈ T

Cov(Xt1 , Xt2) = E[Xt1Xt2 ]− E[Xt1 ]E[Xt2 ] =

= E[X0Xt2−t1 ]− µ2 =

= Cov(X0, Xt2−t1).

Si osservi che quest’ultima quantita dipende soltanto da t2 − t1. Inoltre

V ar(Xt) = E[X2t ]− (E[Xt])

2 =

= E[X20 ]− (E[X0])

2 =

= V ar(X0).

Osservazione 2.1. Per un processo in L2 la stazionarieta stretta implicaquella debole. Il viceversa non e vero in generale.

Osservazione 2.2. Un processo gaussiano che sia debolmente stazionario,cioe per il quale

b = µ = cost ∀ t ∈ T

K(t1, t2) = K(t1 − t2), ∀ t1, t2 ∈ T × T,

e anche fortemente stazionario.

La stazionarieta debole richiede che i primi due momenti siano invariantirispetto al tempo; piu in generale si richiede questa proprieta per i primi kmomenti:

Definizione 2.3. Una serie temporale (Xt)t∈T si dice debolmente stazionariadi ordine k se, per ogni t1, . . . , tn ∈ T e per ogni h ∈ R, tutti i primi k mo-menti di Xt1 , Xt2 , . . . , Xtn esistono e sono invarianti rispetto a cambiamentidi tempo, ovvero:

E[Xk1t1 Xk2

t2 . . . Xkntn ] = E[Xk1

t1+hXk2t2+h . . . Xkn

tn+h],

∀ h ∈ R e ∀ k1, k2, . . . , kn ∈ N tali che k1 + k2 + . . . + kn = k.

7

Osservazione 2.3. All’interno di questa nuova definizione, una serie debol-mente stazionaria corrisponde ad una serie debolmente stazionaria di ordine2.

Osservazione 2.4. Di solito i processi reali non sono stazionari ma possonoavere dei ’trend’ rispetto al tempo.L’analisi delle serie temporali tende ad identificare questi trend con metodiper certi versi analoghi a quelli della teoria generale della regressione. Inseguito si sottrae il trend e si analizza la serie stazionaria ottenuta.

Sia ora (Xt)t∈T un processo in L2 debolmente stazionario. Le nozioniintrodotte precedentemente ci portano in modo naturale all’introduzione didue funzioni di fondamentale importanza, le funzioni di autocovarianza e diautocorrelazione.Si consideri la quantita

Cov[Xt, Xt+h] = E[(Xt − µ)(Xt+h − µ)].

Dalla discussione precedente sappiamo che questa quantita dipende solo dalvalore di h e non dipende da t. Possiamo dunque scrivere:

γ(h) = E[(Xt − µ)(Xt+h − µ)]. (1)

Questa funzione prende il nome di funzione di autocovarianza, mentre laquantita h viene abitualmente chiamata lag. Si osservi che nel caso continuola funzione γ e definita per tutti i valori di h, mentre nel caso discreto edefinita solo per h = 0,±1,±2, . . ..Per ogni h, γ(h) misura la covarianza tra coppie di variabili del processoseparate da un intervallo di lunghezza h.Il grafico della funzione γ prende il nome di covariogramma.La funzione di autocorrelazione, invece, e definita da

ρ(h) =Cov[Xt, Xt+h]√

V ar(Xt)V ar(Xt+h). (2)

Poiche γ(0) = V ar(Xt), ∀ t, possiamo scrivere

ρ(h) =γ(h)

γ(0). (3)

8

La rappresentazione grafica della funzione di autocorrelazione prende il nomedi correlogramma.Per ogni h, ρ(h) rappresenta il coefficiente di correlazione tra due variabilidel processo (Xt) separate da un intervallo di lunghezza h. Al crescere dih, ci aspettiamo che la correlazione tra Xt e Xt+h diminuisca. Infatti, se he grande in generale il processo all’istante t + h avra ’dimenticato’ il valoreassunto all’istante t. Di conseguenza, ci aspettiamo che sia γ(h) che ρ(h)tendano a zero per |h| → +∞. Il modo in cui ρ(h) decresce a zero puoessere interpretato come una misura della ’memoria’ del processo. Si ricordicomunque che stiamo parlando di casi ’tipici’; ci sono infatti processi per iquali ρ(h) non decresce a zero per |h| → +∞, ma si comporta in manierasinusoidale.

Esempio 2.3. Rumore biancoUn processo εtt∈Z si dice rumore bianco o white noise, se εt e una succes-sione di variabili aleatorie non correlate, identicamente distribuite con medianulla e varianza σ2

ε :

E[εt] = 0

Var[εt] = σ2ε > 0

Cov(εt, εt′) = 0, t 6= t′.

Questo processo e debolmente stazionario; infatti, la sua funzione di auto-covarianza e

γ(h) =

0 h 6= 0σ2

ε h = 0.

Segue che la sua funzione di autocorrelazione e

ρ(h) =

0 h 6= 01 h = 0

.

Se si assume che εt abbia distribuzione normale N(0, σ2ε), allora il processo e

fortemente stazionario e prende il nome di rumore bianco gaussiano.

Esempio 2.4. Consideriamo un processo stocastico del tipo

Xt = β1 + β2t,

doveβ = (β1, β2)

T ∼ N(β,Σ)

9

e

β = (β1, β2)T , Σ =

(σ11 σ12

σ21 σ22

).

Questo processo potrebbe essere visto come il modello stimato di una regres-sione lineare semplice a rumore gaussiano.Calcoliamone la media:

E[Xt] = E[β1] + E[β2]t =

= β1 + β2t 6= cost.

La covarianza risulta:

Cov[Xt, Xt+h] = E[(Xt − E(Xt))(Xt+h − E(Xt+h))] =

= E[(β1 − β1) + (β2 − β2)t][(β1 − β1) + (β2 − β2)(t + h)] =

= σ11 + σ12(t + h) + σ12t + σ22t(t + h).

E chiaro che la covarianza non dipende solo dall’incremento h, ma anche dagliistanti t e t + h, che diventano rilevanti. Questo processo evidentemente none stazionario.

Vediamo ora le principali proprieta della funzione di autocovarianza γ.

Proposizione 2.1. La funzione di autocovarianza di una serie storica realedebolmente stazionaria e una funzione pari di h.

Dimostrazione. Supponiamo per semplicita E[Xt] = 0 ∀t. Allora

γ(h) = E[XtXt+h] = E[Xt−hXt] = γ(−h).

Osservazione 2.5. Ci possiamo limitare a studiare la dipendenza di γ da|h|, ossia dai valori non negativi di h.

Proposizione 2.2. La funzione di autocovarianza γ di una serie temporaledebolmente stazionaria (Xt) e semidefinita positiva.

Dimostrazione. Supponiamo per semplicita E[Xt] = 0 ∀t. Siano

tj ∈ R, j = 1, . . . , n

aj ∈ R, j = 1, . . . , n

γ(tj − tk) = Cov(Xtj , Xtk).

10

Allora avremo

0 ≤ V ar

(n∑

j=1

ajXtj

)=

= E

[n∑

j=1

n∑

k=1

ajakXtjXtk

]=

=n∑

j=1

n∑

k=1

ajakE[XtjXtk ] =

=n∑

j=1

n∑

k=1

ajakγ(tj − tk).

Prendendo n = 2 nella dimostrazione precedente possiamo ricavare un’al-tra interessante proprieta:

0 ≤ a21γ(0) + a2

2γ(0) + 2a1a2γ(t1 − t2),

da cui segue1

2(a2

1 + a22) ≥ −a1a2

γ(t1 − t2)

γ(0).

Ponendo h = t1 − t2 e scegliendo successivamente

−a1 = a2 = 1

a1 = a2 = 1

si ottiene la nota relazione

|γ(h)| ≤ γ(0).

Come conseguenza diretta delle proprieta della funzione di autocovarianzaderivano le proprieta della funzione di autocorrelazione:

1. ρ(h) = ρ(−h) ∀h ∈ R;

2. ρ(0) = 1;

3. |ρ(h)| ≤ 1 ∀h ∈ R;

4. ρ(h) e una funzione semidefinita positiva (questo risultato segue daquello analogo per la funzione di autocovarianza; basta dividere tuttoper γ(0) nella dimostrazione).

11

3 Previsione

Uno dei principali problemi dell’analisi delle serie storiche e quello della pre-visione di un valore futuro del processo a partire da un insieme di osservazionidelle sue realizzazioni passate.Consideriamo un processo in L2 stazionario. Supponiamo di aver osservato isuoi valori fino all’istante n:

. . . , Xn−2, Xn−1, Xn

e di voler prevedere i valori futuri del processo agli istanti n + 1,. . . , n + k(k ≥ 1).Chiaramente, per fare cio, vorremmo utilizzare le informazioni contenutenelle osservazioni che abbiamo a disposizione. Indichiamo la previsione conXn+k; vorremmo dunque calcolare Xn+k come una qualche funzione delleXt, t ≤ n, ovvero:

Xn+k = f(. . . , Xn−2, Xn−1, Xn).

La questione che si pone e: come scegliere la funzione f? Per risponderedobbiamo innanzittutto scegliere un ’criterio’ per misurare l’accuratezza diXn+k come previsione di Xn+k. Il criterio piu semplice e anche quello piulargamente utilizzato consiste nel minimizzare l’errore medio quadratico:

E[(Xn+k − Xn+k)2]. (4)

Il problema diventa allora trovare una funzione f che minimizzi la quantita(4). Sappiamo che in questo contesto la miglior scelta per f e il valore attesocondizionato di Xn+k dati Xt, t ≤ n:

Xn+k = E[Xn+k|Xt, t ≤ n]

o meglioXn+k = E[Xn+k|Fn], (5)

dove Fn = σ(Xt, t ≤ n) e la σ-algebra generata dai dati.Si osservi che in L2 Xn+k e la proiezione ortogonale di Xn+k ∈ L2(Ω,Fn+k, P )sul sottospazio L2(Ω,Fn, P ). In queste condizioni

(Xn+k − Xn+k) ⊥ L2(Ω,Fn, P )

12

eE[(Xn+k − Xn+k)

2] = minZ∈L2(Ω,Fn,P )

E[(Xn+k − Z)2].

Una buona proprieta di questo stimatore e la non distorsione. Infatti, per leproprieta del valore atteso, abbiamo:

E[Xn+k] = E[E[Xn+k|Fn]] = E[Xn+k].

A questo punto sappiamo che lo stimatore esiste e che e unico q.c., ma comelo costruiamo? Lo stimatore (5) e un importante risultato teorico, ma hauna scarsa utilita pratica. Per poter effettivamente calcolare la previsione,dobbiamo considerare un modello, cioe dobbiamo ridurci a casi particolari,per esempio ai modelli lineari.

3.1 Predittori lineari

Siano date le osservazioni (in numero finito)

X1, . . . , Xn.

Supponiamo che la previsione Xn+k sia una funzione lineare dei dati:

Xn+k =n∑

i=1

β(k)i Xi = XT β(k) ≡ η,

dove

X =

X1...

Xn

β(k) =

β(k)1...

β(k)n

.

Per calcolare Xn+k, dobbiamo solo determinare i coefficienti β(k)1 , . . . , β

(k)n che

minimizzano l’errore quadratico medio (4), ovvero:

β(k)

= arg minβ(k)∈Rn

E[(Xn+k −XT β(k))2].

In questo modo abbiamo ottenuto una riduzione dello spazio L2(Ω,Fn, P ) aduno spazio finito dimensionale.

Introduciamo ora la definizione di matrice inversa generalizzata di Moore-Penrose, che ci servira nel prossimo teorema:

13

Definizione 3.1. L’inversa generalizzata di Moore-Penrose di una matriceA (m× n) di ordine p e una matrice A+ (n×m) tale che

AA+A = A

A+AA+ = A+.

Si dimostra cheA+ = (AT A)−1AT .

Se A e quadrata e det A 6= 0 allora A+ = A−1.

Teorema 3.1. Se gli Xj, j = 1, . . . , n, derivano da un processo debolmentestazionario a media nulla e funzione di covarianza nota γ, allora il predittoreMQ lineare dell’osservazione Xn+s, s ≥ 1 e

Xn+s = XT bs,

dovebs = V +

nnVns,

essendo Vnn = E[XXT ], V Tns = [γ(n + s − 1), γ(n + s − 2), . . . , γ(s)] e V +

nn

l’inversa generalizzata di Moore-Penrose di Vnn.

Dimostrazione. Cerchiamo il

minβ(s)∈Rn

E[(Xn+s −XT β(s))2] (6)

con

η = XT β(s) =n∑

i=1

β(s)i Xi. (7)

Derivando rispetto ai β(s)i (i = 1, . . . , n) si trova che il problema (6) e

equivalente alla ricerca di η del tipo (7) tale che

E[(Xn+s − η)Xi] = 0 i = 1, . . . , n,

ovvero

E

[(Xn+s −

n∑j=1

β(s)j Xj

)Xi

]= 0 i = 1, . . . , n

Da cui si ricava

E[Xn+sXi] =n∑

j=1

β(s)j E[XjXi],

14

che equivale a

γ(n + s− i) =n∑

j=1

β(s)j γ(i− j), i = 1, . . . , n.

Dunque, i parametri (β(s)j )1≤j≤n si possono scegliere tra le soluzioni del siste-

ma lineare:

γ(n + s− 1) = β(s)1 γ(0) + β

(s)2 γ(1) + . . . + β(s)

n γ(n− 1)

γ(n + s− 2) = β(s)1 γ(1) + . . . + β(s)

n γ(n− 2)...

...

γ(s) = β(s)1 γ(n− 1) + . . . + β(s)

n γ(0).

Possiamo riscrivere il sistema nella forma

Vns = Vnnβ(s)

(8)

da cui

β(s)

= V +nnVns.

Calcoliamo infine l’errore quadratico medio, che sappiamo essere minimo:

σ2t = E[(Xt − Xt)

2] =

= E[(Xt − Xt)(Xt − Xt)] =

= E[(Xt − Xt)Xt]− E[(Xt − Xt)Xt] =

= E[X2t ]− E[XtXt] =

= γ(0)−n∑

j=1

βjE[XjXt] =

= γ(0)−n∑

j=1

βjγ(t− j).

Osservazione 3.1. Il sistema (8) sara ben posto se si dispone di informazioniper γ per tutte le coppie i e j, anche per i valori futuri. Servira dunque unmodello parametrico per γ.

Nel seguito ci limiteremo ad analizzare la famiglia di processi AR e MAper i quali il problema della previsione si riconduce ad un numero finito dioperazioni sui dati disponibili.

15

3.2 Rappresentazione di serie storiche

Vediamo innanzittutto un teorema fondamentale di rappresentazione di unprocesso stocastico debolmente stazionario.

Teorema 3.2 (Wold). Supposto che Xtt∈Z sia un processo stocastico amedia nulla debolmente stazionario, esso puo essere sempre rappresentatonella forma

Xt =∞∑

j=0

ψjεt−j, (9)

dove εtt∈Z e una successione di rumori bianchi, cioe di variabili aleatorienon correlate, a media nulla e con varianza comune σ2

ε e ψjj∈N e unasuccessione di numeri reali tale che

∞∑j=0

|ψj|2 < +∞. (10)

Osservazione 3.2. La rappresentazione (9) va intesa in media quadratica,ovvero:

E

(Xt −

n∑j=0

ψjεt−j

)2 n−→ 0.

Questa rappresentazione e detta rappresentazione lineare generale della serie.

Osservazione 3.3. La condizione (10) garantisce la stazionarieta del pro-cesso (9). Infatti:

E[Xt] = limn↑+∞

E

[n∑

j=0

ψjεt−j

]=

= limn↑+∞

n∑j=0

ψjE[εt−j] = 0,

Var(Xt) = limn↑+∞

n∑j=0

ψ2j E[ε2

t−j] =

= σ2ε

∞∑j=0

ψ2j

16

e

γ(h) = E[XtXt+h] =

= limn↑+∞

limm↑+∞

n∑j=0

m∑i=0

ψjψiE[εt−jεt+h−i] =

= σ2ε

∞∑i=0

ψiψi+h,

poiche E[εrεs] = δrsσ2ε .

Di conseguenza la funzione di autocorrelazione risulta:

ρ(h) =

∑∞i=0 ψiψi+h∑∞

i=0 ψ2i

.

Chiaramente le funzioni di autocovarianza e autocorrelazione dipendono solodallo shift h. Comunque, dal momento che coinvolgono somme infinite, perla stazionarieta dobbiamo mostrare che γ(h) e finita per ogni h:

|γ(h)| = |E[XtXt+h]| ≤ [Var(Xt)Var(Xt+h)]1/2 = σ2

ε

∞∑j=0

ψ2j .

Introduciamo ora l’operatore backshift :

BjXt = Xt−j

e chiamiamo ψ(B) l’operatore formale:

ψ(B) =∞∑

j=0

ψjBj.

Possiamo quindi riscrivere la (9) nella seguente forma:

Xt = ψ(B)εt.

Un’altra utile rappresentazione per il processo Xt e nota come rappresen-tazione autoregressiva (AR):

Xt = π1Xt−1 + π2Xt−2 + . . . + εt,

17

dove εt e un rumore bianco. L’equazione precedente si puo riscrivere nellaseguente forma

π(B)Xt = εt (11)

dove

π(B) = 1−∞∑

j=1

πjBj

con

1 +∞∑

j=1

|πj| < +∞.

Secondo la definizione di Box & Jenkins un processo si dice invertibile se sipuo scrivere nella forma (11).Veniamo ora ai modelli finiti.

3.2.1 Modelli MA

Un modello moving average di ordine q e un modello del tipo (9) ridotto alfinito.

Definizione 3.2. Sia εt un rumore bianco avente media nulla e varianzaσ2

ε . Il processo Xt e un processo a media mobile di ordine q (MA(q)) se

Xt = εt + ψ1εt−1 + . . . + ψqεt−q, (12)

dove ψi sono costanti reali.

Vediamo subito che

E[Xt] = 0

Var(Xt) = σ2ε(1 + ψ2

1 + . . . + ψ2q ).

Calcoliamo ora la funzione di autocovarianza:

γ(k) = Cov(XtXt+k) =

= Cov(εt + . . . + ψqεt−q, εt+k + . . . + ψqεt+k−q) =

=

0 k > q

σ2ε

∑q−ki=0 ψiψi+k k = 0, 1, . . . , q

γ(−k) k < 0

,

18

dove abbiamo posto ψ0 = 1. Osserviamo che γ non dipende da t, dunqueil processo e debolmente stazionario. Se εt ∼ N(0, σ2

ε), allora Yt e unprocesso gaussiano (e somma di normali indipendenti) e quindi e strettamentestazionario.La funzione di autocorrelazione risulta

ρ(k) =

1 k = 0Pq−ki=0 ψiψi+kPq

i=0 ψ2i

k = 1, . . . , q

0 k > qρ(−k) k < 0

Esempio 3.1. Si consideri la serie definita da

Xt =4∑

i=0

εt−i.

Allora

ρh =

5−|h|

5|h| ≤ 4

0 altrimenti.

Proposizione 3.1. La serie storica Xt definita da (12) e invertibile se leradici dell’equazione caratteristica

mq + ψ1mq−1 + ψ2m

q−2 + . . . + ψq = 0

soddisfano la condizione |m| < 1.

Dimostrazione. Si veda [3].

3.2.2 Modelli AR

Definizione 3.3. Sia εt un rumore bianco avente media nulla e varianzaσ2

ε . Il processo Xt e un processo autoregressivo di ordine p (AR(p)) se

Xt + π1Xt−1 + . . . + πpXt−p = εt, (13)

dove πi sono costanti reali.

Proposizione 3.2. La serie storica Xt definita da (13) e stazionaria se leradici dell’equazione caratteristica

mp +

p∑j=1

πjmp−j = 0

soddisfano la condizione |m| < 1.

19

Dimostrazione. Si veda [3].

Esempio 3.2. AR(1)Per studiare i modelli AR cominciamo con lo studio di un modello autore-gressivo del primo ordine debolmente stazionario. Esso si puo esprimere nellaforma

Xt = ρXt−1 + εt (14)

oppure equivalentemente(1− ρB)Xt = εt.

Si osservi che, per la proposizione precedente, la condizione che ci assicurala debole stazionarieta della serie e |ρ| < 1.Sostituendo quindi nella (14) la quantita

Xt−1 = ρXt−2 + εt−1

otteniamo

Xt = ρ(ρXt−2 + εt−1) + εt =

= ρ2Xt−2 + ρεt−1 + εt

e per induzione

Xt = ρkXt−k +k−1∑j=0

ρjεt−j.

Ora calcoliamo

E

(Xt −

k−1∑j=0

ρjεt−j

)2 = E[(ρkXt−k)

2] =

= ρ2kE[X2t−k] =

= ρ2kcostk−→ 0 per |ρ| < 1.

(E[X2t−k] e costante per la debole stazionarieta). Dunque possiamo esprimere

Xt nella forma

Xt =+∞∑j=0

ρjεt−j, (15)

dove l’uguaglianza va intesa in media quadratica.Si noti come cio che abbiamo ottenuto sia coerente con il teorema di Wold e

20

in particolare con la condizione di stazionarieta che in esso compare.Dalla rappresentazione (15) segue che

E[Xt] = limk↑+∞

E

[k−1∑j=0

ρjεt−j

]=

= limk↑+∞

k−1∑j=0

ρjE[εt−j] = 0.

Inoltre

E[X2t ] = lim

k↑+∞

k−1∑j=0

E[(ρjεt−j)2]

e poiche gli εt sono non correlati, otteniamo

E[X2t ] = lim

k↑+∞

k−1∑j=0

ρ2jE[ε2t−j] =

= σ2ε lim

k↑+∞

k−1∑j=0

ρ2j =

=σ2

ε

1− ρ2.

Per quanto riguarda la covarianza si osservi che

Xt+h = ρhXt +h−1∑j=0

ρjεt+h−j,

per cui

γ(h) = E[XtXt+h] = ρhE[X2t ] + E

[Xt

h−1∑j=0

ρjεt+h−j

]=

= ρh σ2ε

1− ρ2+ lim

l→+∞

l−1∑m=0

h−1∑j=0

ρm+jE[εt−mεt+h−j] =

= ρh σ2ε

1− ρ2,

per h = 0, 1, 2, . . .. Come ci aspettavamo γ(h) dipende solo da h. Inoltreossevando l’espressione finale per la γ vediamo che

21

- per ρ > 0, γ(h) e monotona decrescente;

- per ρ < 0, γ(h) e altalenante con valore assoluto decrescente.

Dunque se siamo in presenza di dati la cui funzione di autocovarianza edecrescente possiamo modellizzarli utilizzando un modello AR(1) con ρ > 0.Viceversa se la funzione γ e altalenante potremmo utilizzare un modelloAR(1) con ρ < 0.

Esercizio. A partire da un processo AR(1) possiamo generare una mar-tingala. Si consideri il processo

Xk = θXk−1 + εk, k = 1, 2, . . .

con θ ∈ R, X0 ≡ 0, εk ∼ N(0, 1) indipendenti. Costruiamo la funzione diverosimiglianza (applicando la regola P (A ∩B) = P (A|B)P (B)):

Lt(θ) =t∏

k=1

pθX(Xk|Xk−1)p

θX(X0),

dove pθX(X0) = 1 e pθ

X(Xk|Xk−1) e la distribuzione condizionata di Xk datoXk−1. Poiche Xk e Xk−1 sono normali, anche la distribuzione condizionatasara normale (si veda [4] p.85) con

E[Xk|Xk−1] = E[θXk−1 + εk|Xk−1] = θXk−1

Var[Xk|Xk−1] = E[(Xk − θXk−1)2|Xk−1] = E[ε2

k] = 1,

ovvero

pθX(Xk|Xk−1) =

1√2π

exp[−1

2(Xk − θXk−1)

2].

Quindi

Lt(θ) =

(1

)t/2

exp

[−1

2

t∑

k=1

(Xk − θXk−1)2

].

Consideriamo ora il logaritmo della funzione di verosimiglianza:

lt(θ) = ln Lt(θ) = cost− 1

2

t∑

k=1

(Xk − θXk−1)2 ,

22

e calcoliamone la derivata:

lt(θ) =∂

∂θln Lt(θ) =

t∑

k=1

(Xk − θXk−1) Xk−1.

Si verifica che lt(θ)t∈N e una martingala a media nulla rispetto alla fil-trazione Ftt∈N generata dagli esperimenti. Infatti:

E[lt(θ)− lt+1(θ)|Ft] = E[Xt(Xt+1 − θXt)|Ft] =

= E[XtXt+1|Ft]− θE[X2t |Ft] =

= XtE[Xt+1|Ft]− θX2t =

= θX2t − θX2

t = 0,

e

E[lt(θ)] =t∑

k=1

E[Xk−1(Xk − θXk−1)] =

=t∑

k=1

(E[Xk−1Xk]− θE[X2

k−1])

=

=t∑

k=1

(γ(1)− θγ(0)) =

=t∑

k=1

1

1− θ2− θ

1− θ2

)= 0.

Esempio 3.3. AR(2)Consideriamo un modello autoregressivo del secondo ordine debolmente sta-zionario:

Xt + α1Xt−1 + α2Xt−2 = εt t ∈ Z,

con α2 6= 0. Per la debole stazionarieta richiediamo che le radici di

m2 + α1m + α2 = 0 (16)

siano in valore assoluto minori di 1.Dette m1 ed m2 le radici dell’equazione (16), si dimostra che (cfr. [3] pp.

23

54-58) per radici distinte (reali o complesse) la funzione di autocovarianza edata da

ρ(h) =mh+1

1 (1−m22) + mh+1

2 (1−m21)

(m1 −m2)(1 + m1m2)

per h = 0, 1, 2, . . .. Se le radici sono complesse (e coniugate), la funzione ρ sipuo anche esprimere nella forma:

ρ(h) =rh sin(hθ + δ)

sin δ,

dove

r = α1/22 = (m1m2)

1/2

cos θ =m1 + m2

2(m1m2)1/2= − α1

2α1/22

tan δ =1 + α2

1− α2

tan θ.

Infine per radici reali e coincidenti si ha che:

ρ(h) =

[1 + h

(1−m2

1 + m2

)]mh, h ∈ N.

3.3 Previsione nei modelli AR

3.3.1 AR(1)

Consideriamo un modello AR(1) stazionario:

Xt = ρXt−1 + εt,

dove εt sono variabili aleatorie indipendenti tali che E[εt] = 0 ∀ t, V ar[εt] =σ2

ε ∀ t. Per la stazionarieta del processo richiediamo che |ρ| < 1.Supponiamo ora di aver identificato il modello (cioe di conoscere ρ e σ2

ε) e diavere a disposizione i dati fino all’istante n:

X1, . . . , Xn.

Ci proponiamo di effettuare la previsione all’istante n + 1 e, piu in generale,agli istanti n+s con s > 1. Sappiamo che la miglior previsione, cioe quella che

24

minimizza l’errore quadratico medio, e il valore atteso condizionato. Dunqueavremo:

Xn+1 = E[Xn+1|Fn]

dove Fn = σ(X1, . . . , Xn). Dal momento che conosciamo il modello per Xn+1

avremo:Xn+1 = E[ρXn + εn+1|Fn] = ρXn

(Xn e Fn misurabile e εn+1 e indipendente da Fn).Calcoliamo ora lo scarto quadratico, che sappiamo essere minimo:

E[(Xn+1 − Xn+1)2] = E[ε2

n+1] = σ2ε .

Esso coincide dunque con la varianza del rumore.Per effettuare la previsione all’istante generico n + s ricordiamo che, tramitesostituzioni successive, il processo si puo riscrivere nella seguente forma:

Xn+s = ρsXn +s∑

j=1

ρs−jεn+j, s ∈ N∗.

La previsione risulta dunque essere

Xn+s = E[Xn+s|Fn] = ρsXn.

Si osservi come, essendo |ρ| < 1, la sensibilita della previsione all’ultimo datodiminuisca all’aumentare di s. Calcoliamo ora l’errore quadratico medio:

E[(Xn+s − Xn+s)2] = E

(s∑

j=1

ρs−jεn+j

)2 =

= E

[s∑

j=1

ρ2(s−j)ε2n+j

]=

=s∑

j=1

ρ2(s−j)E[ε2n+j] =

= σ2ε

s∑j=1

ρ2(s−j) =

= σ2ε

1− ρ2s

1− ρ2→ σ2

ε

1− ρ2per s →∞.

Si osservi che la quantita E[(Xn+s − Xn+s)2] e crescente rispetto ad s e per

s tendente all’infinito converge alla varianza del processo Xt.

25

3.3.2 AR(p)

Consideriamo ora un processo AR(p) stazionario:

Xt =

p∑j=1

θjXt−j + εt, (17)

dove εt sono variabili aleatorie indipendenti aventi media nulla e varianzaσ2

ε . Introducendo il polinomio:

ϑ(B) = 1−p∑

j=1

θjBj,

dove B e l’operatore di shift, possiamo riscrivere il modello (17) nella seguenteforma:

ϑ(B)Xt = εt. (18)

Affinche il processo sia stazionario richiediamo che le soluzioni dell’equazionepolinomiale associata a (17), ovvero di:

ϑ(m−1) = mp −p∑

j=1

θjmp−j = 0,

soddisfino la condizione |m| < 1.Analogamente al caso precedente, supponiamo di aver identificato il model-lo e di avere a disposizione i dati fino all’istante n. Effettuiamo quindi laprevisione all’istante n + 1:

Xn+1 = E

[p∑

j=1

θjXn+1−j + εn+1|Fn

]=

p∑j=1

θjXn+1−j.

L’errore quadratico medio risulta:

E[(Xn+1 − Xn+1)2] = E[ε2

n+1] = σ2ε .

26

E possibile effettuare la previsione all’istante generico n + s di un processoAR(p):

Xn+s = E[Xn+s|Fn] =

= E

[p∑

j=1

θjXn+s−j + εn+s|Fn

]=

=

∑pj=1 θjXn+s−j s = p + 1, p + 2, . . .∑s−1j=1 θjXn+s−j +

∑pj=s θjXn+s−j s = 2, 3, . . . , p

La previsione viene dunque effettuata ricorsivamente partendo da Xn+1, ot-tenuto precedentemente.Per poter calcolare l’errore quadratico medio commesso e necessario scrivereil processo al passo n+s in funzione dei dati noti. Si procede ricorsivamente.Riscriviamo la (17) nella seguente forma:

Xn+s = εn+s + θ1Xn+s−1 + . . . + θpXn+s−p; (19)

all’istante precedente avremo

Xn+s−1 = εn+s−1 + θ1Xn+s−2 + . . . + θpXn+s−1−p. (20)

Sostituendo la (20) nella (19) otteniamo:

Xn+s = εn+s + θ1(εn+s−1 + θ1Xn+s−2 + . . . + θpXn+s−1−p) + . . . +

+θpXn+s−p =

= εn+s + θ1εn+s−1 + (θ21 + θ2)Xn+s−2 + . . . + θ1θpXn+s−1−p.

Ripetendo questo procedimento per s− 1 volte arriviamo all’espressione:

Xn+s = εn+s + δ∗1εn+s−1 + . . . + δ∗s−1εn+1 + α∗s−1,1Xn + α∗s−1,2Xn−1 + . . . +

+α∗s−1,pXn−p+1. (21)

Sostituendo poiXn = εn + θ1Xn−1 + . . . + θpXn−p

nell’espressione precedente si ottiene

Xn+s = εn+s + δ∗1εn+s−1 + . . . + δ∗s−1εn+1 + α∗s−1,1εn + (α∗s−1,2 + θ1α∗s−1,1)×

× Xn−1 + . . . + (α∗s−1,p + θp−1α∗s−1,1)Xn−p+1 + α∗s−1,1θpXn−p.

27

Da quest’ultima equazione si ricava la relazione ricorsiva per i coefficienti:

δ∗s = α∗s−1,1

α∗s,j = α∗s−1,j+1 − θjα∗s−1,1 j = 1, . . . , p− 1

α∗s,p = θpα∗s−1,1.

Questo modo di procedere puo essere espresso in termini piu formali. Siconsideri il modello AR(p) scritto in forma (18), allora formalmente

Xt = ϑ(B)−1εt,

dove

ϑ(B)−1 =∞∑

j=0

δjBj;

in altre parole i δj sono i coefficienti dell’espressione polinomiale

(1−

p∑j=1

θjzj

)−1

=∞∑

j=0

δjzj.

Si verifica che δ∗j = δj ∀ j, cioe i δ∗j sono i coefficienti dei primi termini delmodello AR espresso in forma MA.A partire dalla (21) possiamo ricavare una nuova espressione per la previsioneall’istante n + s:

Xn+s = E[Xn+s|Fn] = α∗s−1,1Xn + α∗s−1,2Xn−1 + . . . + α∗s−1,pXn−p+1.

Il valore atteso dell’errore quadratico risulta:

E[(Xn+s − Xn+s)2] = E[(εn+s + δ∗1εn+s−1 + . . . + δ∗s−1εn+1)

2] =

=s−1∑i=0

δ∗2i E[ε2n+s−i] =

= σ2ε

s−1∑i=0

δ∗2i ,

dove δ∗0 = 1.Si osservi che al crescere di s l’errore quadratico aumenta e per s →∞ tendealla varianza del processo Xt.Per la costruzione dell’intervallo di fiducia, si veda la sezione 3.6.

28

3.4 Previsione nei modelli MA

Consideriamo ora un modello a media mobile di ordine q:

Yt =

q∑j=1

βjεt−j + εt (22)

e imponiamo le seguenti ipotesi:

- le radici dell’equazione caratteristica

mq +

q∑j=1

βjmq−j = 0

siano minori di 1 in valore assoluto (condizione di invertibilita);

- i coefficienti βj, j = 1, . . . , q, siano noti;

- i rumori εt siano variabili aleatorie indipendenti a media nulla e varian-za σ2

ε .

Supponiamo, inoltre, di avere a disposizione i dati fino all’istante n:

Y1, . . . , Yn.

Se conoscessimo i rumori εt, t = n, n− 1, . . . , n− q + 1, per il processo (22),allora, utilizzando argomenti analoghi ai precedenti, la miglior previsionelineare sarebbe:

Yn+s = E[Yn+s|ε1, . . . , εn] =

∑qj=s βjεn+s−j, 1 ≤ s ≤ q,

0, s > q.

In realta gli εt non li conosciamo; quindi, dobbiamo costruire dei predittoribasati su degli stimatori per gli εt, dove questi stimatori saranno funzionidelle osservazioni Yt. Si utilizza una procedura ricorsiva. Partiamo dallastima di ε1. Formalmente avremmo:

Y1 =

q∑j=1

βjε1−j + ε1.

29

Da quest’ultima equazione vorremmo ricavare ε1, ma non abbiamo a dis-posizione ε−q+1, . . . , ε0. Si risolve allora questo problema scegliendo comestimatore per questi rumori la loro media, ovvero:

ε−j+1 = 0, j = 1, . . . , q.

In questo modo ricaviamo:ε1 = Y1.

In modo analogo consideriamo

Y2 =

q∑j=1

βj ε2−j + ε2,

ovveroY2 = β1ε1 + ε2,

da cuiε2 = Y2 − β1ε1.

Sulla base di questo ragionamento iterato otteniamo:

ε3 = Y3 − β1ε1 − β2ε2

...

εn = Yn −q∑

j=1

βj εn−j.

Indichiamo ora con εn(Y−m, . . . , Yn) lo stimatore di εn calcolato sulla basedelle osservazioni Y−m, . . . , Yn. Si puo dimostrare che per un processo a mediamobile invertibile, gli stimatori εn convergono ad εn in media quadratica perm → +∞.

Esempio 3.4. Vediamo il calcolo esplicitamente per un modello MA(1) deltipo

Yn = ϑεn−1 + εn.

Per tale modello risultaεn = Yn − ϑεn−1

e sostituendo all’indietro gli ε otteniamo

εn = Yn − ϑ(Yn−1 − ϑεn−2)

= Yn − ϑYn−1 + ϑ2(Yn−2 − ϑεn−3).

30

Facendo il limite per m → +∞ otteniamo:

limm→+∞

E[(εn(Y−m, . . . , Yn)− εn)2

]=

= E

(n∑

i=−∞(−1)n−iϑn−iYi − εn

)2 =

= E

(n∑

i=−∞(−1)n−iϑn−i(εi + ϑεi−1)− εn)

)2 =

= E

(n∑

i=−∞(−1)n−iϑn−iεi +

n∑i=−∞

(−1)n−iϑn−i+1εi−1 − εn

)2 =

= E

(n∑

i=−∞(−1)n−iϑn−iεi −

n−1∑j=−∞

(−1)n−jϑn−jεj − εn

)2 =

= E[(εn − εn)2

]= 0.

Tornando al caso generale, calcoliamo la previsione al passo s per ilmodello (22):

Yn+s =

∑qj=s βj εn+s−j, 1 ≤ s ≤ q,

0, s > q.

Si osservi che l’errore di previsione per Yn+1 calcolato in media quadraticatende a σ2

ε quando il numero di osservazioni m diventa grande. Infatti:

limm→+∞

E

[(Yn+1 − Yn+1(Y−m, . . . , Yn)

)2]

=

= limm→+∞

E

(Yn+1 −

q∑j=1

βj εn+1−j

)2 =

= E(ε2n+1) + lim

m→+∞

q∑j=1

βjE[(εn+1−j − εn+1−j)

2] = σ2ε ,

dove nella penultima uguaglianza e stata utilizzata l’indipendenza degli εi.

31

3.5 Previsione nei modelli ARMA

I metodi visti per i modelli AR e MA si estendono facilmente ai processiautoregressivi a media mobile. Si consideri un processo ARMA(p,q):

Yt =

p∑j=1

θjYt−j +

q∑i=1

βiεt−i + εt, (23)

che soddisfi le ipotesi di stazionarieta ed invertibilita, ovvero le radici deipolinomi

mp −p∑

j=1

θjmp−j = 0

e

mq +

q∑i=1

βimq−j = 0

siano minori di 1 in valore assoluto. Inoltre supponiamo gli εt variabili aleato-rie indipendenti (0, σ2

ε).A partire dai dati Y1, . . . , Yn, la previsione per il modello (23) e data da

Yn+s =

∑pj=1 θjYn−j+1 +

∑qi=1 βiεn−i+1, s = 1,∑s−1

j=1 θjYn−j+s +∑p

j=s θjYn−j+s +∑q

i=s βiεn−i+s, s = 2, 3, . . . ,

(24)dove

εt =

0, t = p, p− 1, . . . ,Yt −

∑pj=1 θjYt−j −

∑qi=1 βiεt−i, t = p + 1, p + 2, . . . , n,

e si intende che θj = 0 per j > p e

q∑i=s

βiεn−i+s = 0, s = q + 1, q + 2, . . . ,

p∑j=s

θjYn−j+s = 0, s = p + 1, p + 2, . . . .

Esempio 3.5. Nella Tabella 1 sono riportate osservazioni Y1, . . . , Y7 del se-guente modello ARMA(1,1):

Yt = 0.9Yt−1 + 0.7εt−1 + εt,

32

dove gli εt sono variabili aleatorie indipendenti normali (0, 1).Calcoliamo ora le previsioni per questo processo utilizzando il metodo delTeorema 3.1. Le prime undici autocorrelazioni per questa serie storica sono,approssimativamente, 1.00, 0.95, 0.85, 0.77, 0.69, 0.62, 0.56, 0.50, 0.45, 0.41e 0.37.Per costruire il predittore lineare per Y8 ricaviamo i sette coefficienti necessaririsolvendo il seguente sistema di equazioni:

1.00 0.95 0.85 . . . 0.560.95 1.00 0.95 . . . 0.620.85 0.95 1.00 . . . 0.69

......

......

0.56 0.62 0.69 . . . 1.00

b11

b12

b13...

b17

=

0.500.560.62

...0.95

Questo sistema ci fornisce i sette coefficienti per la previsione all’istante n + 1.Il vettore dei coefficienti bT

1 = (0.06,−0.18, 0.32,−0.50, 0.75,−1.10, 1.59)viene poi applicato alle osservazioni (Y1, . . . , Y7). La previsione all’istanten + 2 si ottiene sostituendo il secondo membro del sistema precedente con(0.45, 0.50, 0.56, 0.62, 0.77, 0.85)T . Le previsioni ottenute utilizzando questicoefficienti sono riportate nella Tabella 1.

t Osservazioni Previsioni

1 -1.58 -2 -2.86 -3 -3.67 -4 -2.33 -5 0.42 -6 2.87 -7 3.43 -8 - 3.049 - 2.7310 - 2.46

Tabella 1: Osservazioni e previsioni per la serie storica Yt = 0.9Yt−1 +0.7εt−1 +εt

Vediamo ora come costruire la previsione con il metodo appena introdot-

33

to. Imponiamo ε0 = ε1 = 0 e calcoliamo:

ε2 = Y2 − 0.9Y1 − 0.7ε1 = Y2 − 0.9Y1 = −1.438

ε3 = Y3 − 0.9Y2 − 0.7ε2 = −0.089

ε4 = Y4 − 0.9Y3 − 0.7ε3 = 1.035

ε5 = Y5 − 0.9Y4 − 0.7ε4 = 1.793

ε6 = Y6 − 0.9Y5 − 0.7ε5 = 1.237

ε7 = Y7 − 0.9Y6 − 0.7ε6 = −0.003.

Quindi, la previsione per Y8 risulta:

Y8 = E[Y8|Y1, . . . , Y7] = 0.9Y7 + 0.7ε7 = 3.085

e

Y9 = 0.9Y8 = 2.776

Y10 = 0.9Y9 = 2.498.

In questo caso la previsione ottenuta con il metodo degli ε e abbastanzavicina a quella ottenuta con il metodo del Teorema 3.1, anche se la dimensionedel campione e piccola. Se il coefficiente di εt−1 fosse stato piu vicino ad uno,la differenza sarebbe stata maggiore.

3.6 La varianza

Come noto, possiamo sempre scrivere un modello stazionario e invertibileARMA (che contiene come casi particolari i processi AR e MA), come

Yt =+∞∑j=0

vjεt−j, (25)

dove v0 = 1,∑+∞

j=0 |vj| < +∞ e gli εt sono variabili aleatorie (0, σ2ε) non

correlate.Quando la serie e scritta in questa forma, si ottengono facilmente le varianzedegli errori di previsione.

Teorema 3.3. (2.9.4 Fuller)Dato un modello ARMA invertibile e stazionario nella forma (25) e il predit-tore

Yn+s =n+s−1∑

j=s

vj εn−j+s,

34

dove εt, t = 1, . . . , n, sono definiti come nella (24), allora

limn→+∞

E

[Yn+s − Yn+s −

s−1∑j=0

vjεn−j+s

]2 = 0.

Utilizzando questo teorema possiamo costruire la matrice di covarianzaper i predittori di Yn+1, Yn+2, . . . , Yn+s. Siano YT

n,s = [Yn+1, Yn+2, . . . , Yn+s] e

YTn,s = [Yn+1, Yn+2, . . . , Yn+s]. Allora

limn→+∞

E[(Yn,s − Yn,s )( Yn,s − Yn,s)T ] =

=

1 v1 . . . vs−1

v1

∑1j=0 v2

j . . .∑1

j=0 vjvs+j−2

......

...

vs−1

∑1j=0 vjvs+j−2 . . .

∑s−1j=0 v2

j

σ2

ε .

Si osservi come gli elementi diagonali della matrice precedente ci forniscanodelle informazioni sulla varianza per n → +∞. Per n finito essi non costi-tuiscono la vera varianza. In ogni caso ci danno un’idea di come la varianzacresca quando ci si allontana dai dati. Se il processo Yt e stazionario, alloraper s → +∞ la varianza della previsione converge alla varianza del processostesso.Possiamo adottare gli elementi diagonali della matrice precedente come sti-matori della varianza. Si noti che tale scelta comporta due problemi: in-nanzittutto si tratta di quantita asintotiche e in secondo luogo, in generale,non conosciamo σ2

ε . La varianza del rumore dovra quindi essere stimata;nella prossima sezione introdurremo uno stimatore per σ2

ε .

Costruiamo ora un intervallo di fiducia per la previsione. Ricordiamo chei predittori che abbiamo considerato sono non distorti, ovvero

E[Yn+s − Yn+s] = 0.

Dunque l’errore quadratico medio dei predittori coincide con la varianza.Al fine di poter utilizzare queste informazioni per costruire un intervallo diconfidenza per la previsione, e necessario conoscere la distribuzione dei rumoriεt. Se la serie e normale, un intervallo di confidenza a livello 1−α e dato da

Yn+s ±Nα/2

[Var(Yn+s − Yn+s)

]1/2

35

ovvero

P

Yn+s ∈ Yn+s ±Nα/2

[s−1∑j=0

v2j

]1/2

σε

= 1− α,

dove Nα/2 e il valore che soddisfa

P (N > Nα/2) = α/2.

Si osservi come la banda si allarghi man mano che ci si allontana dai datidisponibili. Questo vuol dire che al crescere di s, la previsione e sempre menoaccurata.

4 Stima

4.1 Stima della media

Quando avevamo come dati n realizzazioni i.i.d. Y1, . . . , Yn di una variabilealeatoria Y ∈ L2, utilizzavamo come stimatore del suo valor medio:

Y =1

n

n∑j=1

Yj.

Per studiare questo stimatore utilizzavamo le proprieta di indipendenza e diidentica distribuzione dei dati:

E(Y ) =1

n

n∑j=1

E(Yj)id=

1

nnE(Y ) = E(Y )

Var(Y )i=

1

n2

n∑j=1

Var(Yj)id=

Var(Y )

n.

La situazione che si presenta con le serie storiche e molto diversa. Infattiin questo contesto i dati si rivelano una sola volta e non si hanno replichedella stessa quantita. In particolare, le osservazioni sono, in generale, cor-relate tra loro; inoltre solo quando il processo e strettamente stazionarioesse presentano la stessa distribuzione di probabilita. Dunque non potremopiu comportarci come nel caso i.i.d., ma dovremo introdurre qualche altroconcetto.

36

Sia Yt un processo debolmente stazionario, cioe un processo in L2 taleche:

E(Yt) = µ ∀ t

Var(Yt) = σ2 ∀ t

Cov(Yt, Yt+h) = γ(h) = σ2ρ(h) ∀ t, h

Supponiamo di avere a disposizione le osservazioni Y1, . . . , Yn. Anche se nonabbiamo posto alcuna ipotesi sulla distribuzione di queste variabili, sappiamoche esse hanno tutte lo stesso valore atteso µ. Questo ci suggerisce cheil solito stimatore, ovvero la media campionaria, possa ancora fornirci unastima ragionevole per µ. Verificheremo questo fatto attraverso lo studio delleproprieta dello stimatore. Proponiamo dunque il seguente stimatore per lamedia:

Y (n) =1

n

n∑t=1

Yt.

Abbiamo

E(Y (n)) =1

n

n∑t=1

E(Yt) =1

nnµ = µ,

dove abbiamo utilizzato l’ipotesi di stazionarieta debole. Quest’ultima ugua-glianza ci dice che lo stimatore Y (n), ottenuto facendo la media storica deidati, e non distorto.Calcoliamo ora la varianza dello stimatore:

Var(Y (n)) =1

n2

n∑t=1

n∑s=1

Cov(Yt, Ys) =

=σ2

n2

n∑t=1

n∑s=1

ρ(t− s) (26)

A questo punto effettuiamo un cambio di variabile da s, t a s e k = (t− s).(La doppia sommatoria (26) si estende su un reticolo quadrato di punti.Anziche sommare prima sulle righe e poi addizionare tra loro i risultati ot-tenuti, sommiamo sulle diagonali e poi addizioniamo le somme diagonali).La sommatoria su k va da −(n−1) a n−1; per k > 0, s va da 1 fino a n−k,mentre per k < 0, s va da −k a n. Si osservi che l’addendo nella (26) e una

37

funzione solo di k e quindi sommando su s otteniamo l’espressione (N − |k|)sia per k positivo che per k negativo. Abbiamo quindi

Var(Y (n)) =σ2

n2

(n−1)∑

k=−(n−1)

(N − |k|)ρ(k) =

=σ2

n

(n−1)∑

k=−(n−1)

(1− |k|

n

)ρ(k). (27)

Si osservi come, contrariamente al caso i.i.d., la sommatoria non ci forniscan volte la varianza del processo.Per poter affermare che Y (n) e uno stimatore consistente per µ dobbiamodimostrare che la quantita (27) tende a zero per n → +∞. Si dimostra chevale la seguente

Proposizione 4.1. Se ρ(k) → 0 per k → +∞, allora

limn→+∞

1

n

(n−1)∑

k=−(n−1)

ρ(k) = 0.

Dimostrazione. Si veda [6].

La proposizione precedente implica

Var(Y (n)) → 0 per n → +∞.

Ma allora

limn→+∞

1

n

n∑t=1

Yt = µ

in media quadratica, ovvero

Y (n) L2−→ µ per n → +∞.

E poiche la convergenza in L2 implica quella in probabilita si ha che Y (n) euno stimatore debolmente consistente per µ.Tra Y (n) e µ sussiste un’interessante relazione. Infatti per ogni t, µ e ilvalore atteso di Yt e rappresenta la media spaziale. Y (n) invece e la media delprocesso Yt rispetto agli istanti temporali t = 1, . . . , N e prende il nome di

38

media temporale. Utilizzando Y (n) per stimare µ abbiamo in effetti sostituitouna media spaziale con la corrispondente media temporale. In generale, se lostimatore e consistente in media quadratica allora la media temporale su unperiodo di n punti converge in L2, per n → +∞, alla corrispondente mediaspaziale. I processi che soddisfano questa proprieta vengono detti ergodici.

4.2 Stima della funzione di autocovarianza

Procedendo in modo analogo a quanto visto per la media, stimiamo la fun-zione di autocovarianza facendo la media temporale tra tutte le coppie divariabili osservate che si trovano alla stessa distanza. Vengono proposti duestimatori:

γ(k) =1

n

n−|k|∑t=1

(Yt − Y (n))(Yt+|k| − Y (n))

ˆγ(k) =1

n− |k|n−|k|∑t=1

(Yt − Y (n))(Yt+|k| − Y (n))

per k = 0,±1,±2, . . . ,±(n− 1). Si dimostra che (cfr. Wei):

E(γ(k)) ' γ(k)− |k|n

γ(k)−(

n− |k|n

)Var(Y (n))

E(ˆγ(k)) ' γ(k)− Var(Y (n)).

E chiaro che entrambi questi stimatori sono distorti. Comunque, se ρ(k) → 0per k → +∞ e quindi il processo e ergodico nella media, e limn→+∞ Var(Y (n)) =0, allora entrambi gli stimatori γ(k) e ˆγ(k) sono asintoticamente non distorti.Vediamo brevemente i vantaggi e gli svantaggi di questi due stimatori:

- la distorsione su γ(k) e maggiore di quella su ˆγ(k), specialmente agrandi k;

- la varianza su γ(k) e spesso minore che su ˆγ(k). Infatti si dimostra cheVar(γ(k)) = O(1/n) per ogni k, mentre Var(ˆγ(k)) = O(1/(n − |k|)).Di conseguenza quando |k| si avvicina a n − 1, Var(ˆγ(k)) aumentasostanzialmente;

- γ(k) condivide con γ(k) la proprieta di essere semidefinita positiva,mentre ˆγ(k) non e necessariamente semidefinita positiva.

39

Per queste ragioni si preferisce utilizzare γ(k) come stimatore della funzionedi autocovarianza.Infine, ci chiediamo quando il processo e ergodico nell’autocovarianza. Aquesto proposito vale la seguente

Proposizione 4.2. Se∑+∞

i=−∞ |γi| < +∞ allora

limn→+∞

Var(γ(k)) = 0.

Dalla proposizione segue che, sotto opportune ipotesi, γ(k) converge inL2 a γ(k), ovvero:

limn→+∞

γ(k) = γ(k),

in media quadratica e questo significa che il processo e ergodico nella funzionedi autocovarianza. Ovviamente dalla convergenza in L2 segue la convergenzain probabilita, cioe lo stimatore γ(k) e debolmente consistente.

4.3 Stima della funzione di autocorrelazione

Poiche la funzione di autocorrelazione ρ(k) e legata alla funzione di autoco-varianza γ(k) dalla relazione

ρ(k) =γ(k)

γ(0),

uno stimatore naturale per ρ(k) e

ρ(k) =γ(k)

γ(0)=

∑n−|k|t=1 (Yt − Y (n))(Yt+|k| − Y (n))∑n

t=1(Yt − Y (n))2,

dove k = 0,±1,±2, . . . ,±(n−1) e Y (n) = 1n

∑nt=1 Yt. Anche in questo caso lo

stimatore ρ(k) e distorto. Se Xt e un processo gaussiano stazionario, alloraper grandi valori di n si ha che

ρ(k) ∼ AN(ρ(k), Var(ρ(k))),

dove

Var(ρ(k)) ' 1

n

+∞∑i=−∞

(ρ2(i)+ρ(i+k)ρ(i−k)−4ρ(k)ρ(i)ρ(i−k)+2ρ2(k)ρ2(i)).

Se ρ(k) = 0 per k > m, allora la varianza si riduce a:

Var(ρ(k)) ' 1

n(1 + 2ρ2(1) + 2ρ2(2) + . . . + 2ρ2(m)).

40

4.4 Stima dei parametri

4.4.1 AR(1)

Si consideri la serie storica stazionaria definita da

Yt − µ = ρ(Yt−1 − µ) + εt,

dove εt sono variabili aleatorie normali (0, σ2) indipendenti e |ρ| < 1. Ab-biamo gia discusso il problema della stima di µ. Sappiamo che nei modelliAR(1) ρ(1) = ρ; quindi, per stimare ρ possiamo utilizzare la stimatore delcoefficiente di autocorrelazione del primo ordine:

ρ = ρ(1) =γ(1)

γ(0)=

∑n−1t=1 (Yt − Y (n))(Yt+1 − Y (n))∑n

t=1(Yt − Y (n))2,

dove appunto

µ = Y (n) =1

n

n∑t=1

Yt.

Osservazione 4.1. Si noti l’analogia con la stima dei parametri nei modellidi regressione lineare.

4.4.2 AR(p)

In questa sezione studiamo gli stimatori per il processo autoregressivo diordine p. Consideriamo la serie storica

Yt +

p∑i=1

αiYt−i = θ0 + εt, (28)

dove le radici del polinomio

mp +

p∑i=1

αimp−i = 0

sono minori di uno in valore assoluto e εt sono variabili aleatorie (0, σ2) noncorrelate.Possiamo usare anche qui tecniche di regressione lineare per il modello (28);per cui e conveniente riscriverlo nella seguente forma:

Yt = θ0 + θ1Yt−1 + θ2Yt−2 + . . . + θpYt−p + εt, (29)

41

dove θi = −αi, i = 1, . . . , p.Cominciamo dal caso p = 2:

Yt = θ0 + θ1Yt−1 + θ2Yt−2 + εt.

Moltiplichiamo l’equazione precedente per XTt = (1, Yt−1, Yt−2)

T e sommiamoper avere

n∑t=3

XTt Yt =

n∑t=3

XTt Xtθ +

n∑t=3

XTt εt, (30)

dove θ = (θ0, θ1, θ2)T .

Calcoliamo ora il valore atteso per il terzo termine:

E

[n∑

t=3

XTt εt

]= E

[n∑

t=3

εt,

n∑t=3

Yt−1εt,

n∑t=3

Yt−2εt

]T

= 0,

poiche εt e Yt−i (i > 0) sono quantita non correlate. Per cui nella stimapossiamo limitarci a considerare l’equazione:

n∑t=3

XTt Yt =

n∑t=3

XTt Xtθ.

Lo stimatore per θ risulta allora essere:

θ =

[n∑

t=3

XTt Xt

]−1 [n∑

t=3

XTt Yt

].

Si osservi che tutto quanto e stato fatto vale anche per p > 2. Per cui sedefiniamo

Xt = (1, Yt−1, . . . , Yt−p)

An =1

n− p

n∑t=p+1

XTt Xt

vn =1

n− p

n∑t=p+1

XTt Yt,

lo stimatore e dato da:θ = A−1

n vn. (31)

42

Non abbiamo informazioni sulla varianza di questo stimatore. Si di-mostra, pero, che tale stimatore e anche stimatore ai minimi quadrati, ovveroe quello che minimizza la quantita:

Q(θ0, θ1, . . . , θp) =n∑

t=p+1

ε2t =

=n∑

t=p+1

(Yt − θ0 − θ1Yt−1 − . . .− θpYt−p)2.

Infatti, derivando rispetto ai parametri θ0, θ1, . . . , θp otteniamo le equazioni:

n∑t=p+1

(Yt − θ0 − θ1Yt−1 − . . .− θpYt−p) = 0

en∑

t=p+1

(Yt − θ0 − θ1Yt−1 − . . .− θpYt−p)Yt−j = 0, j = 1, . . . , p.

Si vede chiaramente che questo sistema coincide con (31).Ricordiamo che lo stimatore ai minimi quadrati ha una buona proprieta. In-fatti quando i rumori εt sono non correlati e hanno tutti la stessa varianza,questo metodo fornisce stimatori ottimali nella classe degli stimatori linearinon distorti. Si dimostra che lo stimatore MQ e efficiente all’interno di que-sta classe, ovvero ha varianza minima.Le relazioni (31) includono le covarianze della variabile Y , che noi non conos-ciamo. Nei paragrafi precedenti abbiamo pero imparato a stimarle. Viste leproprieta asintotiche della stima della funzione di autocovarianza possiamodire che per n À 1 (31) e equivalente a

θ∗

= (A∗n)−1v∗n, (32)

dove

A∗n =

γ(0) γ(1) . . . γ(p− 1)γ(1) γ(0) . . . γ(p− 2)

......

...γ(p− 1) γ(p− 2) . . . γ(0)

e

v∗n =

γ(1)...

γ(p)

.

43

Lo stimatore ai minimi quadrati per σ2 e dato da

σ2 =1

n− 2p− 1

n∑t=p+1

ε2t ,

doveε2

t = Yt −Xtθ.

Le proprieta asintotiche di θ e σ2 sono date nel seguente

Teorema 4.1. (8.2.1 Fuller)Sia Yt il processo definito da

Yt = θ0 +

p∑j=1

θjYt−j + εt, t = p + 1, p + 2, . . . .

Si supponga che (Y1, Y2, . . . , Yp) sia fissato, oppure che (Y1, Y2, . . . , Yp) siaindipendente da εt per t > p e E[|Yi|2+ν ] < ∞ per i = 1, 2, . . . , p e qualcheν > 0. Le radici del polinomio

mp −p∑

j=1

θjmp−j = 0

siano minori di uno in valore assoluto. Si supponga inoltre che la successioneεt∞t=1 sia una successione di variabili aleatorie (0, σ2) con

E[(εt, ε2t )|At−1] = (0, σ2) q.c.

eE[|εt|2+ν |At−1] < L < ∞ q.c.

per ogni t e qualche ν > 0, dove At−1 = σ(εt : j ≤ t− 1). Se (Y1, Y2, . . . , Yp)e aleatorio, la σ-algebra e generata da εt : j ≤ t − 1 e da (Y1, Y2, . . . , Yp).Allora

(a)√

n(θ − θ)L−→ N(0, A−1σ2),

dove

A = limn→+∞

n−1

n∑t=p+1

E(XTt Xt),

44

(b) σ2 P−→ σ2.

Osservazione 4.2. Anche se abbiamo ottenuto dei buoni risultati asintotici,il comportamento degli stimatori ottenuti con poche osservazioni puo differireconsiderevolmente da quanto visto nella teoria asintotica. Se le radici delpolinomio caratteristico sono vicine a zero, l’avvicinamento alla distribuzionenormale e piuttosto rapido. Per esempio, se ρ = 0 in un modello AR(1), l’ap-prossimazione normale e molto buona per n > 30. D’altra parte, se le radicisono vicine a uno in valore assoluto, per avere una buona approssimazionenormale occorre che n sia molto grande.

Rimane aperto il problema della determinazione dell’ordine p. Si puo ot-tenere, per esempio, utilizzando l’F-test per verificare l’ipotesi p = p1 control’ipotesi p = p2 < p1.

4.4.3 MA(q)

Purtroppo, la stima dei parametri nei processi a media mobile e piu compli-cata rispetto al caso AR.Dai risultati del paragrafo 3.2.1 sappiamo che c’e una relazione tra la fun-zione di autocorrelazione di una serie storica a media mobile e i parametridel modello che la descrive.Vediamo questo metodo nel dettaglio per il seguente modello MA(1):

Xt = εt + βεt−1.

In questo caso abbiamo

ρ(1) =β

1 + β2,

mentre tutte le autocovarianze di ordine superiore sono nulle. Utilizzandoil calcolo elementare, e facile vedere che ρ(1) raggiunge un massimo di 0.5per β = 1 e un minimo di −0.5 per β = −1. Dunque, ρ(1) per un modelloMA(1) cadra sempre nell’intervallo [−0.5, 0.5].L’idea per stimare β e quella di partire dallo stimatore di ρ(1):

ρ(1) =

∑n−1t=1 (Xt − X(n))(Xt+1 − X(n))∑n

t=1(Xt − X(n))2.

45

Per lo stimatore di β allora troviamo:

β =

1−[1−4(ρ(1))2]1/2

2ρ(1)0 ≤ |ρ(1)| ≤ 0.5

−1 ρ(1) < −0.51 ρ(1) > 0.50 ρ(1) = 0.

Ovviamente, se ρ(1) giace significativamente al di fuori dell’intervallo (−0.5, 0.5),allora un modello a media mobile del primo ordine e sospetto.

Si dimostra che β e uno stimatore consistente per β, ma non e efficienteper β 6= 0. Stimatori piu efficienti si possono ottenere con i metodi ai minimiquadrati o di massima verosimiglianza (si veda [3]).

46

Riferimenti bibliografici

[1] Ash, R.B. e Gardner, M.F. (1975) Topics in Stochastic Processes,Academic Press, New York.

[2] Brockwell, P. and Davis R. (1991), Time Series: Theory and Methods,Springer-Verlag, New York.

[3] Fuller, W.A. (1996) Introduction to Statistical Time Series, John Wiley& Sons, New York.

[4] Priestley, M.B. (1981) Spectral Analysis and Time Series, Volume 1,Academic Press, London.

[5] Priestley, M.B. (1981) Spectral Analysis and Time Series, Volume 2,Academic Press, London.

[6] Wei, W.W. (1990) Time Series Analysis, Addison-Wesley PublishingCompany.

47