Note introduttive all econometria

Note introduttive all’econometria Domenico Suppa∗ 25 maggio 2016

Il modello di regressione lineare Date n osservazioni campionare relative alla variabile Y e a k − 1 variabili X i (con k < n): Y = [Yi ]   X = 1,[X i, j ]

per i = 1,..., n per i = 1,..., n e per j = 2,..., k

nel modello di regressione lineare si ipotizza che tra Y e X sussista una relazione lineare alterata solo dalla presenza di una componente erratica ε. Le condizioni imposte, per ipotesi, sono le seguenti:1 Y = Xβ + ε

ipotesi di linearità della relazione tra Y e X

(I1)

E[ε|X] = 0

ipotesi di esogeneità dei regressori X

(I2)

ˆ si Ci si pone l’obiettivo di stimare il vettore dei parametri incogniti β. Tale stima, β, ˆ può ottenere minimizzando la somma dei quadrati dei residui2 u = Y − Xβ: min u0 u βˆ

(1)

Si tratta, quindi, di determinare βˆ di modo che la lunghezza del vettore u sia la più piccola possibile. Ciò equivale a calcolare βˆ imponendo la condizione che u sia ∗

Università di Napoli "Federico II". Gli elementi della prima colonna della matrice X sono posti tutti pari a 1, in tal modo il primo elemento del vettore β, indicato con β0 , è il termine costante (l’intercetta) della funzione lineare che determina Y. Riguardo alla condizione (I2) essa risulta generalmente troppo restrittiva e quindi, quando è possibile, si richiede solo l’incorrelazione tra l’errore e i regressori della medesima osservazione: E[ε|X = xi ] = 0. 2 È chiaro che u 6= ε: l’errore ε resterà sempre incognito mentre i residui calcolati u possono consentire di inferire le caratteristiche di ε solo se sono soddisfatte le ipotesi (I1), (I2) e quelle di seguito specificate nel testo. L’ipotesi che le k colonne di X siano linearmente indipendenti garantisce che la matrice X0 X sia invertibile (unicità della stima) e che essa sia definita positiva (in tal modo la somma dei quadrati dei residui è minima). Per valutare se una matrice simmetrica definita positiva è quasi-singolare (cioè molto vicina ad una matrice non invertibile) si utilizza l’indice di condizionamento che è pari al rapporto tra il massimo e il minimo autovalore della matrice (se questo rapporto è molto distante da 1, allora la matrice è quasi-singolare e i risultati dell’inversione dipendono drasticamente dalla precisione di macchina del calcolatore). V. appendice A3.

1

1

ortogonale al sottospazio di Rn generato dalle colonne di X. Tale sottospazio avrà dimensione k, se le k colonne di X sono linearmente indipendenti (questa è un’ipotesi ulteriore che può essere verificata facilmente), e pertanto il vettore u giacerà nel sottospazio di Rn di dimensione n− k (Rn−k è il complemento ortogonale di Rk rispetto a Rn ). Imposta la condizione di ortogonalità X0 u = 0, determiniamo βˆ risolvendo il seguente sistema di equazioni: Y = Xβˆ + u

(2)

−1 βˆ = (X0 X) X0 Y

(3)

dal quale otteniamo:

−1

Lo stimatore dei minimi quadrati (X0 X) X0 Y è costruito sulla validità della ipotesi (I1). Se è valida anche l’ipotesi (I2), tale stimatore è corretto. Infatti, data la (I1), possiamo riscrivere la (3) nel modo seguente: −1 −1 βˆ = (X0 X) X0 (Xβ + ε) = β + (X0 X) X0 ε

(4)

e, grazie alla (I2), verificare che lo stimatore dei minimi quadrati è non distorto:3 −1 0

ˆ E[β|X] = E[(X0 X)

X (Xβ + ε)|X] −1 0

= β + E[(X0 X)

X ε|X]

−1 0

= β + (X0 X)

X E[ε|X]

= β +0 = β

(5)

La distorsione da variabile omessa L’importante risultato di correttezza delle stime ottenute con il metodo dei minimi quadrati (espresso nella (5)) richiede (oltre alle ipotesi (I1) e (I2)) che tra i regressori siano incluse tutte le variabili esplicative dalle quali dipende probabilisticamente Y; in altre parole non devono esservi variabili omesse (si tratta di un requisito che talvolta viene esplicitato come una ulteriore ipotesi).4 Le conseguenze dell’omissione di una variabile esplicativa, importante nella determinazione di Y, possono essere illustrate nel modo seguente. Poniamo di dividere i regressori in due gruppi X1 e X2 , in modo da avere il seguente processo che genera Y: Y = X1 β 1 + X2 β 2 + ε

(6)

ma, supponiamo di omettere i regressori X2 nella specificazione del modello che vogliamo stimare, allora la stima di β1 sarà:5 −1 βˆ1∗ = (X01 X1 ) X01 Y

(7) −1

Per effettuare i passaggi è necessario tener presente il fatto che (X0 X) X0 è funzione deterministica di −1 −1 X, per cui E[(X0 X) X0 |X] = (X0 X) X0 . 4 Come vedremo, aggiungere questa ulteriore ipotesi non è strettamente necessario e in pratica risulterebbe ridondante in quanto è già contenuta nella (I2). 5 Premoltiplichiamo la specificazione scorretta Y = X1 βˆ1∗ + u per X0 , assumiamo la condizione di ortogonalità X0 u = 0 (che non è corretta), ed esplicitiamo rispetto a βˆ1∗ . 3

2

ma, dato il modello vero, sostituendo la (6) nella (7) e prendendo il valore medio, si ha: −1 0 X1 X2 β 2

E[βˆ1∗ |X1 ] = β1 + (X01 X1 )

(8)

Quindi la stima βˆ1∗ è distorta (sistematicamente diversa da β1 ) se non si verifica almeno una delle seguenti condizioni: −1

• (X01 X1 ) X01 X2 = 0 (i regressori X2 sono incorrelati con - cioè ortogonali ai regressori X1 ); • β2 = 0 (i regressori X2 sono incorrelati con - cioè ortogonali a - Y). L’ultima condizione (β2 = 0) ci dice che l’omissione di regressori (X2 ) irrilevanti nel processo di generazione di Y non altera i risultati della regressione di Y su X1 (l’intuizione è ovvia). La prima condizione, invece, indica che quando i regressori sono tra loro ortogonali (cioè incorrelati), abbiamo la possibilità di condurre separatamente le regressioni di Y su X1 e di Y su X2 avendo la certezza di ottenere gli stessi parametri che avremmo ottenuto dalla regressione condotta su tutti i regressori in un’unica soluzione. Più in generale, però, la distorsione causata dall’omissione di una variabile −1 rilevante (quantificata dal termine (X01 X1 ) X01 X2 β2 della (8)), che non sia ortogonale ai regressori utilizzati, conduce a stime distorte. Inoltre, la distorsione da variabile omessa non diminuisce all’aumentare della numerosità campionaria (le stime quindi risultano distorte e anche inconsistenti). In effetti, si dovrebbe riflettere sul fatto che questa omissione invalida l’ipotesi (I2): la variabile omessa, correlata con i regressori utilizzati, viene inglobata nei residui della regressione che, a questo punto, non possono più essere incorrelati con i regressori utilizzati (viene meno l’ipotesi (I2)). Ora risultà più chiaro il senso della seguente affermazione: sotto le ipotesi (I1) e (I2) sul modello di regressione lineare, • gli stimatori dei minimi quadrati producono stime non distorte dei parametri β.

Sfericità degli errori e teorema di Gauss-Markov Riguardo agli errori ε è stata già introdotta l’ipotesi (I2) di ortogonalità rispetto ai regressori X del modello di regressione lineare. Si è visto che le stime dei parametri βˆ sono non distorte se è soddisfatta questa ipotesi (insieme alla (I1)). Per indagare la variabilità di βˆ è necessario formulare una ulteriore ipotesi inerente la distribuzione degli errori.6 In condizioni ottimali gli errori dovrebbero essere omoschedastici e non autocorrelati: Var[ε|X] = E[εε0 |X] = σ2 I

ipotesi di sfericità degli errori

(I3)

Dove σ2 è la varianza costante di ε e I è la matrice identità n × n. Se i regressori X non sono stocastici (perché, ad esempio, scaturiscono da un disegno sperimentale), allora (tenuto conto della (4)) le varianze e le covarianze del vettore βˆ possono essere calcolate nel modo seguente: ˆ VarX [β|X] = E[(βˆ − β)(βˆ − β)0 |X] 6

Al riguardo, si veda ad esempio Greene (2012, cap. 4).

3

−1 0

= E[(X0 X)

−1

X εε0 X(X0 X)

−1 0

= (X0 X)

−1

X E[εε0 |X]X(X0 X)

−1 0

= (X0 X)

|X]

−1

X (σ2 I)X(X0 X) −1

= σ2 (X0 X)

(9)

Se i regressori X sono stocastici (cioè sono delle variabili casuali), allora (v. appendice A1):     ˆ = EX VarX [β|X] ˆ ˆ Var[β] + VarX E[β|X] −1

= EX [σ2 (X0 X)

−1

] = σ2 EX [(X0 X)

]

(10)

−1 ˆ Tale risultato tenendo presente che per la (9) VarX [β|X] = σ2 (X0 X) e per  si ottiene  ˆ la (5) VarX E[β|X] = VarX [β] = 0.

Dalle ipotesi (I1), (I2) e (I3) si ricava il teorema di Gauss-Markov: • gli stimatori dei minimi quadrati dei parametri β del modello di regressione lineare hanno varianza minima tra tutti gli stimatori lineari e non distorti di β (sono BLUE, best linear unbiased estimator). Per provare il teorema consideriamo un qualunque altro generico stimatore di β −1 che sia, per ipotesi, lineare e non distorto ma diverso da (X0 X) X0 Y. La stima di β −1 sarà data da βˆC = CY (dove C è una matrice di dimensione k × n al pari di (X0 X) X0 ). Il requisito della non distorsione consiste nella seguente uguaglianza: E[βˆC |X] = E [CY|X] = E [(CXβ + Cε)|X] = β

(11) −1

che è soddisfatta solo se CX = I. In generale avremo C = D + (X0 X) X0 (dove D è −1 definita implicitamente) e quindi CX = DX + (X0 X) X0 X = DX + I = I solo se DX = 0. Assumendo che sia soddisfatta questa condizione, calcoliamo la varianza di βˆC :7 VarX [βˆC |X] = E[(βˆC − β)(βˆC − β)0 |X] = E[Cε(Cε)0 ] = E[Cεε0 C0 ] = CE[εε0 ]C0 = σ2 CC0 ” — 2 0 −1 0 0 −1 0 0 = σ (D + (X X) X )(D + (X X) X ) −1

= σ2 DD0 + σ2 (X0 X)

ˆ = σ2 DD0 + VarX [β|X]

(12)

l’ultimo passaggio riprende il risultato della (9). La (12) stabilisce che VarX [βˆC |X] ˆ poiché DD0 e (X0 X)−1 sono entrambe definite raggiunge il minimo solo quando βˆC = β, positive. È bene ribadire che il teorema di Gauss-Markov può essere applicato solo se il modello di regressione lineare è correttamente specificato (ipotesi (I1) e (I2)) e se gli errori sono omoschedastici e non autocorrelati (ipotesi (I3)).8 In condizioni alternative, 7 8

−1 Sostituendo C nella (4) al posto di (X0 X) X0 , abbiamo (βˆC − β) = Cε. Un ipotesi molto più forte della (I3) potrebbe richiedere che gli errori siano variabili casuali indipendenti e identicamente distribuite (iid).

4

più generali, altri stimatori (non lineari e magari distorti) possono essere più efficienti9 degli stimatori dei minimi quadrati nell’inferenza sul legame funzionale esistente tra Y e X.

La varianza dell’errore Siccome σ2 è (e resta) ignota, la definizione della varianza degli stimatori dei minimi quadrati, fornita dalla (9) o dalla (10), non è direttamente utilizzabile nei calcoli. Pertanto, è necessario stimare la varianza degli errori per determinare la varianza degli stimatori dei minimi quadrati. La soluzione è basata sulla definizione dei residui. Dalla (2) e dalla (3), abbiamo: u = Y − Xβˆ ” — ” — −1 −1 = I − X(X0 X) X0 Y = I − X(X0 X) X0 [Xβ + ε] ” — −1 = I − X(X0 X) X0 ε = MX ε

(13)

−1

dove con MX è stata indicata la matrice [I − X(X0 X) X0 ].10 Notiamo che, per l’ipotesi (I2), dalla (13) segue che E[u|X] = E[MX ε|X] = 0. Sfruttando le proprietà della matrice MX , si può dimostrare che: • la varianza campionaria dei residui: Su2 =

u0 u n−k

(14)

è uno stimatore corretto della varianza degli errori. Omettendo per semplicità il pedice deponente assegnato alla matrice MX : E[u0 u] = E[(Mε)0 (Mε)] = E[ε0 M0 Mε] = E[ε0 Mε] = E[t r(ε0 Mε)] = E[t r(Mεε0 )] = t r(E[Mεε0 ]) = t r(ME[εε0 ]) = t r(Mσ2 I) = t r(M)σ2 = (n − k)σ2

(15)

Questo risultato proviene dal fatto che ε0 Mε è uno scalare (la cui traccia t r(ε0 Mε) coincide con il suo valore) e dalla definizione di M (che è una funzione deterministica di X).11 La matrice delle covarianze di u è data da: Var[u] = E[uu0 ] = E[(Mε)(Mε)0 ] = E[Mεε0 M0 ] = M(σ2 I)M0 = σ2 M

9

(16)

Uno stimatore è più efficiente rispetto a un altro se ha una varianza più bassa. Una generica matrice MX (denominata residual maker) trasforma un qualunque vettore (o variabile) nei rispettivi residui ortogonali allo spazio generato dalle componenti (colonne) della matrice X. È semplice dimostrare che MX è simmetrica (M = M0 ), idempotente (M2 = MM = M) e semidefinita positiva (i suoi autovalori sono non negativi). La matrice HX = I − MX (hat matrix) che proietta un qualunque vettore sullo spazio generato dalle colonne di X, gode delle stesse proprietà di MX . Gli autovalori di entrambe le matrici assumono valore 0 o 1 e si può facilmente verificare che il prodotto HX MX ha come risultato una matrice nulla. Dato che H è simmetrica e idempotente, abbiamo: H2 x = λ2 x = Hx = λx e quindi λ2 − λ = 0, da cui segue che λ può assumere solo i valori 0 e 1. 11 Si tenga presente anche che, quando le dimensioni delle matrici€permettonoŠla commutazione € Šdel

−1 −1 prodotto, t r(AB) = t r(BA). Pertanto: t r ε0 Mε = t r Mεε0 e t r X(X0 X) X0 = t r X0 X(X0 X) . 10

5

La (16) mostra che anche se gli errori ε sono tra loro incorrelati e hanno varianza costante (rispettando l’ipotesi (I2)), in generale, i residui u non hanno varianza costante e presentano un certo grado di autocorrelazione.12

Ipotesi di normalità degli errori Per concludere l’elenco delle proprietà finite di cui può godere la stima dei minimi quadrati dei parametri del modello di regressione lineare aggiungiamo l’ipotesi di normalità dei residui. La motivazione è dovuta al fatto che per condurre dei test statistici sui parametri stimati è necessario conoscere la distribuzione dei rispettivi stimatori. Se gli errori del modello di regressione lineare si distribuiscono in modo normale e indipendente (fermo restando le altre ipotesi (I1), (I2) e (I3)), applicando la proprietà riproduttiva della v. c. normale alla (4), risulta che anche i (gli stimatori dei) parametri sono distribuiti in modo normale (essendo questi ultimi il risultato di unaŠcombinazione € 2 0 −1 ˆ lineare degli errori come mostrato dalla (4)): β|X ∼ N β,σ (X X) . L’ipotesi di normalità degli errori può essere formalizzata nel modo seguente:
ε ∼ N 0,σ2 I distribuzione normale degli errori (I4) −1

Sotto questa ipotesi, date la (4) e la (9), indicando con (X0 X)ii il termine diagonale −1 della matrice (X0 X) corrispondente al parametro βi , possiamo affermare che: βˆi − βi ∼ N (0,1) = z q σ (X0 X)−1 ii

per i = 1,..., k

(17)

la v. c. definita dalla (17) si distribuisce come la normale standardizzata, indicata con z. Applicando la proprietà riproduttiva della variabile casuale normale alla (13), si può dedurre che, se per l’errore ε vale l’ipotesi di normalità (I4), anche i residui hanno
2 una distribuzione normale u ∼ N 0,σ MX . Inoltre, si può dimostrare, che la variabile casuale

u0 u σ2

si distribuisce come una v. c. Chi-quadrato con n − k gradi di libertà:13 2

 ε  S (n − k) u0 u  ε 0 u 2 = MX = ∼ χ(n−k) 2 2 σ σ σ σ

(18)

Rapportando la (17) alla radice quadrata della (18) (dopo aver diviso quest’ultima per n − k), deduciamo che tale rapporto si distribuisce come una variabile casuale t di 12

Nei successivi paragrafi emergeranno ulteriori considerazioni sulla relazione tra errori e residui in base all’esame delle matrici MX e HX . In particolare, esaminando i termini diagonali della matrice MX si può mostrare che la variabilità dei residui è sempre minore di quella degli errori. Ciò è evidente anche intuitivamente, in quanto i residui sono calcolati in base ad un criterio di distanza minima (l’ortogonalità alle colonne della matrice X) mentre per gli errori l’ipotesi di ortogonalità è generalmente soddisfatta Pn solo con una certa approssimazione. 13 2 Ricordiamo che χ(n) = i=1 zi2 , dove zi sono variabili casuali normali standardizzate e tra loro indipendenti. Il risultato presentato nella (18) deriva dal fatto che se x, ad m componenti, si distribuisce come una variabile casuale normale multivariata con vettore delle medie pari a µ e matrice delle covarianze pari a S, allora la forma quadratica (x − µ)0 S−1 (x − µ) si distribuisce come una variabile casuale chi-quadrato con m gradi di libertà. (Piccolo, 1998, pp. 423 e 879), Greene (2012), v. appendice A2.

6

Student con n − k gradi di libertà: βˆ −β qi i

σ (X0 X)−1 βˆi − βi z ii ∼È = t (n−k) =r q 2 Su2 (n−k) χ(n−k) Su (X0 X)−1 ii 2 (n−k)σ

(19)

n−k

Grazie alla (19), se è valida la (I4), si può codurre il test t di Student sul confronto tra 2 βˆi e βi . Nella (19) il numeratore z deve essere indipendente dal denominatore χ(n−k) , cioè lo stimatore del parametro deve essere indipendente dallo stimatore dell’errore standard degli errori.14

Cenni alle proprietà asintotiche dei minimi quadrati Circa la consistenza dello stimatore dei minimi quadrati di β, osserviamo che la sua varianza, ottenuta dalla (9), può essere riscritta nel modo seguente:  ‹ σ2 1 0 −1 ˆ XX (20) Var[β] = n n Il secondo membro della (20) è costituito dal prodotto di due termini: il primo, σ2 /n, tende a zero al crescere di n, il secondo (grazie all’applicazione della legge dei grandi numeri, se la matrice X0 X non è singolare e i regressori presentano alcuni requisiti di regolarità) tende (in probabilità) ad un valore limite costante.15 Pertanto, siccome la varianza dello stimatore dei minimi quadrati tende a zero al crescere di n, tale stimatore è consistente oltre che corretto.16 Nel caso gli errori non siano normalmente distribuiti, al fine di ottenere stimatori dei minimi quadrati la cui distribuzione sia asintoticamente normale, è necessario disporre di una elevata numerosità campionaria e formulare l’ipotesi di indipendenza seriale:17 (xi ,εi ) per

i = 1,..., n

sono serialmente indipendenti

(I5)

(questa ipotesi rafforza l’incorrelazione richiesta dall’ipotesi di sfericità dei residui formulata dalla (I3)). In un certo senso, il venir meno dell’indipendenza, a causa degli errori non più distribuiti normalmente, richiede un assunto esplicito sull’indipendenza degli errori. Permane, invece, immutata l’ipotesi di varianza costante degli errori. Sotto queste condizioni è possibile applicare il teorema limite centrale.18 14

V. Davidson e MacKinnon (2003, cap. 4). V. appendice A2. V. Greene (2012, p. 103 e ss.). 16 Uno stimatore θˆn è consistente per il parametro θ se, al crescere della numerosità campionaria n, tende a 1 la probabilità che θˆn assuma il valore θ . 17 Indichiamo con xi il vettore riga trasposto della matrice dei regressori relativo all’osservazione i-esima. 18 Le due formulazioni del teorema limite centrale più utilizzate in econometria sono: quella di LindebergLevy e, in particolare, quella di Lindeberg-Feller; riguardo alla prima, la formulazione (più semplice) è la seguente (vedi appendice A9 a p. 25): Se x 1 ,..., x n è un capione casuale (cioè a componenti indipendenti) estratto da una popolazione la cui Pn p 2 i=1 x i funzione di distribuzione ha media µ e varianza σ finite e x n = n , allora, per n → ∞, n(x n − µ) 2 si distribuisce come una v. c. normale N (0,σ ) . Il secondo teorema limite centrale, quello di Lindeberg-Feller, invece, pur conservando l’ipotesi di indipendenza, ammette che le v. c. x i presentino diverse varianze (e anche diverse medie). È proprio questo il teorema limite necessario per ottenere la (21), dato che le varianze dei termini xi εi non possono essere ritenute costanti anche se è costante la varianza di εi al variare di i. Greene (2012). 15

7

Tenendo presente la (4), possiamo scrivere:19 p

X0 X n(βˆ − β) = n 

‹−1 

X0 ε p n

‹

a



∼ N 0, σ2



X0 X n

‹−1 

(21)

Il membro di sinistra si distribuisce asintoticamente come una v. c. normale € 0 Š−1con media zero e varianza finita. Si è già detto che, al crescere di n il termine XnX tende in probabilità ad una matrice con valori finiti. La stima della varianza σ2 può essere basata ancora sulla varianza campionaria dei residui Su2 definita dalla (14).20 Svolgiamo in modo più formale e dettagliato questo discorso. Nel seguito indichiap

mo con “plim” o → la convergenza in probabilità.21 Nell’approccio alla regressione lineare su dati di elevata numerosità campionaria si assume che:22  0 ‹ XX plim =Q con Q definita positiva e finita (I6) n Date le ipotesi (I1), (I2), (I3), (I5), (I6), abbandonando l’ipotesi (I4), dalla (4) si deduce che:23 • i minimi quadrati stimano in modo consistente i parametri del modello di regressione lineare. Dato che:  0 ‹−1 0 Xε XX 0 −1 0 ˆ per la (4) (22) β = β + (X X) X ε = β + n n la dimostrazione della precedente proposizione è basata sulle seguneti convergenze in probabilità:24 

X0 X n

‹−1

p

→ Q−1

per la (I6)

n

p X0 ε 1 X = xi εi → E[xi εi ] = Exi [xi E[εi |xi ]] = 0 n n i=1 p

βˆ → β + Q−1 0 = β

per la LDGN e la (I2) per le precedenti.

(23)

Inoltre, si dimostra che: • gli stimatori dei minimi quadrati si distribuiscono in modo asintoticamente normale. 19

a

Il simbolo ∼ indica che il termine di sinistra “si distribuisce asintoticamente come” il termine di destra. 20 Per le dimostrazioni di questi enunciati si veda Greene (2012, pp. 103 e ss. e app. D). 21 Una successione di variabili casuali X n converge alla variabile casuale X se e solo se, per qualsiasi δ > 0: lim P r ob(|X n − X |) < δ = 1 n→∞

in altri termini, al crescere di n, tende ad essere nulla la probabilità che X n assuma valori che si discostano da X . 22 Nella prassi si richiede che siano rispettate le condizioni (minimali) di Grenander sul comportamento regolare della matrice dei regressori. Greene (2012, p. 105). 23 Naturalmente, l’ipotesi (I4) potrebbe essere ancora valida senza in alcun modo compromettere il discorso svolto in questa sede. 24 LDGN è l’acronimo di legge debole dei grandi numeri.

8

Consideriamo il termine

0 X pε n

=

Pn i=1

xi εi p n

della (21), abbiamo E[xi εi ] = 0 e:

Var[xi εi ] = E[xi εi εi 0 xi 0 ] n  0  0 0 X   xi εi xi εi p = E xi E[εi εi 0 |xi ]xi 0 = σ2 E[xi xi 0 ] → σ2 Q = p p n n i=1 Pn l’ultimo passaggio è giustificato all’identità X0 X = i=1 xi xi 0 ; infatti, dividendo ambo i membri di quest’ultima per n e considerando l’ipotesi (I6), segue: Pn 0 X0 X p i=1 xi xi E[xi xi 0 ] = = →Q n n Pertanto, applicando il teorema limite centrale di Lindeberg-Feller, il termine 2

0 X pε
n

della

(21) tende a distribuirsi asintoticamente come una v. c. normale N 0, σ Q . Gra25 zie a questi risultati, e alle (23), applicando il teorema di Slutsky, deduciamo che € Š −1 −1 2 −1 0 la distribuzione asintotica del primo membro della (21) è: N Q 0, Q σ QQ . p ˆ Quindi, come riportato
nella (21), la v. c. n(β − β) tende a distribuirsi come una v. c. 2 −1 normale N 0, σ Q . Infine: • la varianza dei residui Su2 è uno stimatore consistente della varianza degli errori σ2 . Dalla (14), ricordando la definizione della matrice M: ” — 0 0 −1 0 0 0 ε I − X(X X) X ε uu ε Mε Su2 = = = n−k n−k n−k  0  0 ‹ 0 ‹−1  0 ‹ εX XX Xε n εε (24) = − n−k n n n n n Il fattore costante n−k tende a 1 al crescere di n, il secondo termine della differenza nelle parentesi quadre, grazie alla (I6) e alle convergenze in probabilità già calcolate 0 nella (23), tende ad annullarsi. Il termine εnε converge, sotto condizioni minimali, a σ2 per n tendente all’infinito. Quindi: plim Su2 = σ2 e, pertanto, la varianza asintotica −1 di βˆ è uguale a Su2 (X0 X) . In conclusione, grazie alla teoria asintotica che permette di determinare la (21), possiamo applicare il test sulla media di una v. c. normale standard z quando i dati sono numerosi e per essi non è possibile sostenere l’ipotesi di normalità degli errori.26 Per quanto si è detto in questo paragrafo:

p Su

n(βˆi − βi ) = r
X0 X −1 n ii

βˆi − βi ∼z q −1 0 Su (X X)ii

p

25

(25)

p

Se g è una funzione continua che non dipende da n, allora g(x n ) → g(x) se x n → x. 26 D’altra parte è anche ben noto che la variabile casuale t (n) tende velocemente alla v. c. normale standard all’aumentare di n. Quindi il test sui parametri può essere basato sulla v. c. normale standard ogni volta che vi è una considerevole numerosità campionaria: prossima o superiore a 100.

9

A0. Il teorema di Frisch-Waugh-Lovell e i valori anomali Il teorema di Ragnar Frisch, Frederick V. Waugh27 e Michael C. Lovell28 afferma che ogni regressione lineare può essere condotta in più stadi. La procedura può essere descritta nel modo seguente: i regressori vengono divisi in due gruppi, si calcola la regressione della variabile dipendente rispetto al primo gruppo, poi si calcolano le regressioni del secondo gruppo di variabili esplicative rispetto al primo gruppo di regressori, si effettua infine la regressione dei residui della prima regressione rispetto ai residui delle altre regressioni. I parametri ottenuti nell’ultimo stadio di questo procedimento sono uguali a quelli che si ottengono da un unica regressione condotta su tutte le variabili esplicative (in un’unica soluzione).29 Il modello di partenza è il seguente:30 Y = X1 β 1 + X2 β 2 + u (A0-1) Sotto la condizione di ortogonalità dei residui rispetto ai regressori (implicita nel metodo dei minimi quadrati), si possono determinare i parametri β1 , ottenendo: −1 0 X1 Y −

β1 = (X01 X1 )

−1 0 X1 X2 β 2

(X01 X1 )

(A0-2)

Sostituendo β1 nella (A0-1) e riorganizzando i termini: −1 0 X1 Y −

Y = X1 (X01 X1 )

−1 0 X1 ]Y =

[I − X1 (X01 X1 )

−1 0 X1 ],

Posto MX 1 = [I − X1 (X01 X1 )

−1 0 X1 X2 β 2 +

X1 (X01 X1 )

X2 β 2 + u

−1 0 X1 ]X2 β2 +

[I − X1 (X01 X1 )

u

(A0-3) (A0-4)

la (A0-4) può essere riscritta in sintesi:

M X 1 Y = M X 1 X2 β 2 + u

(A0-5)

Teorema FWL: I parametri β2 e i residui u ottenuti da quest’ultima equazione di regressione sono identici a quelli prodotti dal modello espresso dall’equazione (A0-1). La figura 1 illustra l’interpretazione grafica del teorema in uno spazio a tre dimensioni. Per giungere alla dimostrazione del teorema basta notare che la retta passante per r2 è parallela alla perpendicolare a X 1 passante per X 2 e che le parallele passanti per a1 , a2 e X 1 individuano dei segmenti corrispondenti e proporzionali (in base al teorema di Talete) sulle rette passanti per r1 e r2 . Quindi, l’univocità di β2 è garantita dalla corrispondenza tra le proiezioni di Y lungo le rette per r1 e r2 . Naturalmente, i residui sono individuati, in ogni caso, dalla perpendicolare rispetto al piano generato da X 1 e X 2 passante per Y .31 Un’applicazione del teorema FWL è costituita dall’analisi dell’influenza dei cosiddetti valori anomali sulla regressione. Si tratta di punti (o osservazioni) lontani dalla funzione di regressione lineare. Tale distanza è valutata in un duplice senso: verticale (gli outliers) e orizzontale (i punti di leverage). La combinazione di queste due distanze 27

Frisch e Waugh (1956). Lovell (1963). 29 Anche i residui sono gli stessi. 30 Rispetto alla simbologia utilizzata in precedenza, con u indichiamo ancora i residui mentre β rappresenta il vettore dei parametri stimati. Il simbolo Y identifica ancora la cosiddetta variabile dipendente e X sono sempre le variabili esplicative o regressori (compresa la costante). 31 Vedi Davidson e MacKinnon (2003). 28

10

Y

r2 X2

r1 a2

a1

0

X1 figura 1: Il teorema FWL in uno spazio a tre dimensioni determina l’influenza del punto preso in considerazione sulla pendenza della funzione di regressione lineare. Indicando con ξ i residui di una generica regressione: Y = Xβ + ξ

(A0-6)

˜ che assume valore 1 in corriintroduciamo una variabile (vettore colonna) dummy e spondenza dell’osservazione per la quale si vuole valutare l’influenza e ha tutte le altre componenti uguali a zero. Questa operazione comporta una variazione dei parametri ˜ e l’equazione (A0-6) diventa: (da β a β) ˜α + u Y = Xβ˜ + e

(A0-7)

˜ premoltiplichiamo la precedente equazione per HX Per misurare la differenza (β − β), ˆ = HX Y = Xβ, che HX X = X0 HX = X e che i residui u (ricordando che per definizione Y ˜); abbiamo i seguenti ovvi passaggi: sono ortogonali sia a X, sia a e ˜ Xβ = Xβ˜ + αHX e

(A0-8)

˜ = αHX e ˜ X(β − β)

(A0-9)

˜ = αX HX e ˜ X X(β − β)

(A0-10)

−1 ˜ β − β˜ = α(X0 X) X0 e

(A0-11)

0

0

˜0 u = 0, possiamo sviluppare Grazie al teorema FWL, considerato che MX X = 0 e che e l’equazione (A0-7) nel modo seguente: ˜+u MX Y = αMX e 11

(A0-12)

α=

˜0 MX Y = α˜ ˜+e ˜0 u e e0 M X e

(A0-13)

˜0 MX Y = α˜ ˜ e e0 M X e

(A0-14)

˜0 M X Y ˜0 M X Y e e ˜ u = = 0 0 ˜ MX e ˜ [I − HX ]˜ ˜ e e e 1 − ˜h

(A0-15)

˜ e ˜h indicano rispettivamente la componente del vettore dei Nell’ultima espressione, u residui u e l’elemento diagonale della matrice HX corrispondenti all’osservazione della quale si vuole misurare l’influenza. Determinato il valore del parametro α, possiamo sostituirlo nell’equazione (A0-11), ottenendo: β − β˜ =

˜ u 1 − ˜h

−1 0

(X0 X)

˜ Xe

(A0-16)

La (A0-16) mostra che l’effetto dell’osservazione presa in esame è tanto più grande ˜ e quanto più elevato è ˜h. Quest’ultimo valore, infatti, appartiene quanto maggiore è u ˜ (un vettore all’intervallo [0,1], in quanto è pari al quadrato della proiezione di e di norma unitaria) sul sottospazio generato dalle colonne di X. In altri termini, per l’idempotenza di HX : ˜h = e ˜0 H X e ˜= e ˜0 H X H X e ˜ = kHX e ˜ k2

(A0-17)

In realtà ˜h non può essere minore di 1/n (dove n è la numerosità delle osservazioni). Se la matrice dei regressori fosse costituita dalla sola costante con tutte le componenti uguali a 1, indicando con H1 la matrice delle proiezioni, avremo: 1 ˜ k2 = kH1 e n

(A0-18)

ma se oltre alla costante vi sono altri regressori, avremo H1 HX = HX H1 = H1 , perché il vettore costante appartiene a X, e quindi: 1 ˜k2 = kH1 HX e ˜k2 ≤ kHX e ˜ k2 = kH1 e n

(A0-19)

˜ su X non può essere in modulo minore della proiezione di e ˜ sul La proiezione di e vettore con tutte le componenti pari a 1 (dato che questo vettore è una colonna della matrice X; infatti, la proiezione ortogonale di un qualunque punto P su una retta r non può essere maggiore della proiezione dello stesso punto su un sottospazio che contiene la retta, in quanto la proiezione ortogonale del punto P sulla retta r non è altro che la proiezione ortogonale sulla retta del punto del sottospazio, contenente la retta, che è proiezione del punto P sul sottospazio al quale appartiene la retta; in un triangolo rettangolo ogni cateto non può avere lunghezza maggiore dell’ipotenusa). Inoltre, il valore medio dei termini diagonali di HX è k/n, poiché vi sono n termini diagonali e la traccia di HX è pari al numero k dei regressori. I valori ˜h tendono ad essere sempre più elevati via via che le componenti dei regressori si allontanano dalle rispettive medie ˜0 X, abbiamo: (generando l’effetto leverage). Premoltiplicando la (A0-16) per e   ˜h ˜ u −1 0 0 0 0 ˜ = ˜ X(β − β) ˜ X(X X) X e ˜= ˜ e e u (A0-20) 1 − ˜h 1 − ˜h ˜ tale valore è l’effetto sulla stima dell’elemento di Y individuato dalla componente di e ˜ dall’equazione (A0-9), premoltiplicando ambo i membri ˜ = Xβ, diversa da zero. Posto Y 12

per se stessi, dividendo per kSu2 e tenendo conto della (A0-15), si ottiene la distanza di Cook:32 ˜h˜ ˆ − Y) ˜ 0 (Y ˆ − Y) ˜ ˆ − Y)k ˜ 2 (Y k(Y α2 0 u2 ˜ ˜ = = e H e = X 2 kSu2 kSu2 kSu2 kSu2 (1 − ˜h)

(A0-21)

Nella figura 2 è riportata una simulazione del modello: y = β0 + β1 x + u

1 2

3 4 5

0

2

4

6

8

10

0

6

2

4

1

6

8

10

6

8

10

∆β1

Cook’s distance

2 3 4 5

∆β0

6

3 6 21

45 0

2

4

6

8

10

0

2

4

figura 2: Distanza di Cook Il primo grafico rappresenta le 100 osservazioni (con x in ascissa e y in ordinata) e la retta di regressione. Il secondo grafico, a destra, riporta la relazione tra x (in ascissa) Pn

32

u2

ˆ 0 (Y−Y) ˆ (Y−Y)

i=1 i Con Su2 si è indicata la varianza dei residui della regressione: Su2 = n−k = . V.: Weisberg n−k ˜ u p ˜s = (2005) e anche Belsey et al. (1980). Ogni valore u è detto residuo studentizzato, con questa

notazione la (A0-21) può scriversi come

˜h˜ u2s k(1−˜h)

=

Su 1−˜h ˜2s ˜h u k 1−˜h = outlier × leverage = influence.

13

e ˜h (in ordinata). La linea orizzontale continua indica il valore medio h = k/n = 2/100, la linea tratteggiata corrisponde al doppio di tale media. Il grafico illustra l’effetto leverage esercitato da ogni x: valori più distanti dalla media x (corrispondente alla linea verticale) causano un valore più elevato di ˜h (i primi 6 valori sono numerati). Nel grafico in basso a sinistra sono riportati in ordinata i valori che risultano dalla formula (A0-21) per ogni osservazione (x è riportato in ascissa come negli altri grafici e la linea orizzontale tratteggiata corrisponde al valore 4/n). Dal grafico si può notare come i valori che influiscono di più sulla stima di y (con maggiore ordinata) non siano ˜ (l’eventuale condizione in generale tutti quelli con un più alto leverage. Il valore di u ˜ di outlier), infatti, modifica l’effetto sintetizzato da h (leverage). L’ultimo grafico, in basso a destra, riporta l’equazione (A0-16) per entrambi i parametri del modello. Si può notare come gli effetti siano differenziati per i diversi parametri.

A1. Medie e varianze condizionate La variabile casuale doppia (X,Y) ha la seguente distribuzione di probabilità: Definizioni di probabilità marginali e probabilità congiunta

Probabilità marginali e Probabilità congiunta

pX i = pX i ,Y1 + pX i ,Y2

Y X

X1 X2

Y1

Y2

pX 1 ,Y1 pX 2 ,Y1 pY1

pX 1 ,Y2 pX 2 ,Y2 pY2

pY j = pX 1 ,Y j + pX 2 ,Y j

pX 1 pX 2 1

pX i ,Y j = P(X i ∩ Y j )

Le distribuzioni di probabilità marginali (pX 1 ,pX 2 ) e (pY1 ,pY2 ), rispettivamente nell’ultima riga e nell’ultima colonna della tabella, sono date dalla somma delle probabilità congiunte (pX i ,Y j ). Consideriamo ora le variabili casuali condizionate, del tipo (X|Yi ), e le loro distribuzioni di probabilità pX|Y j = P(X|Y j ). Data la definizione di probabilità congiunta: pX i ,Y j = P(X i ∩ Y j ) = pX i P(Y j |X i ) = pX j pYi |X j = pY j P(X i |Y j ) = pY j pX j |Yi

(A1-1)

avremo le seguenti probabilità condizionate: pX i |Y j =

pX i ,Y j

pY j |X i =

pY j

pX i ,Y j pX i

e il riepilogo dettagliato riportato nelle seguenti tabelle: Probabilità condizionate di X Y Y1 Y2 X 1 pX 1 |Y1 pX 1 |Y2 X X 2 pX 2 |Y1 pX 2 |Y2 1 1

Probabilità condizionate di Y Y Y1 Y2 X 1 pY1 |X 1 pY2 |X 1 1 X X 2 pY1 |X 2 pY2 |X 2 1

14

(A1-2)

Utilizzando le probabilità condizionate si possono calcolare i valori medi condizionati di X e Y: E[X|Y1 ] = E[X|Y2 ] = E[Y|X 1 ] = E[Y|X 2 ] =

pX 1 ,Y1 pY1 pX 1 ,Y2 pY2 pX 1 ,Y1 pX 1 pX 2 ,Y1 pX 2

X1 + X1 + Y1 + Y1 +

pX 2 ,Y1 pY1 pX 2 ,Y2 pY2

pX 1 ,Y2 pX 1 pX 2 ,Y2 pX 2

X 2 = pX 1 |Y1 X 1 + pX 2 |Y1 X 2 X 2 = pX 1 |Y2 X 1 + pX 2 |Y2 X 2 Y2 = pY1 |X 1 Y1 + pY2 |X 1 Y2 Y2 = pY1 |X 2 Y1 + pY2 |X 2 Y2

Calcolando le medie delle medie condizionate, si dimostra che vale la legge delle aspettative iterate (LIE, law of iterated expectations). Con riferimento alle precedenti due tabelle, tenuto conto delle identità (A1-2), (posto E[X|Y] = (E[X|Y1 ], E[X|Y2 ]) e E[Y|X] = (E[Y|X 1 ], E[Y|X 2 ])), calcoliamo: EY [E[X|Y]] = pY1 E[X|Y1 ] + pY2 E[X|Y2 ]

= pY1 pX 1 |Y1 X 1 + pX 2 |Y1 X 2 + pY2 pX 1 |Y2 X 1 + pX 2 |Y2 X 2 = pX 1 ,Y1 X 1 + pX 2 ,Y1 X 2 + pX 1 ,Y2 X 1 + pX 2 ,Y2 X 2 = pX 1 X 1 + pX 2 X 2 = E[X]

(A1-3)

EX [E[Y|X]] = pX 1 E[Y|X 1 ] + pX 2 E[Y|X 2 ]

= pX 1 pY1 |X 1 Y1 + pY2 |X 1 Y2 + pX 2 pY1 |X 2 Y1 + pY2 |X 2 Y2 = pX 1 ,Y1 Y1 + pX 2 ,Y1 Y2 + pX 1 ,Y2 Y1 + pX 2 ,Y2 Y2 = pY1 Y1 + pY2 Y2 = E[Y]

(A1-4)

In definitiva, la media delle medie condizionate EY [E[X|Y]] è uguale alla media non condizionata E[X]. Ogni realizzazione della variabile casuale media condizionata E[X|Yi ] è funzione di Yi (la media condizionata varia al variare di Yi , per = 1,..., n e viene detta funzione di regressione). La variabile casuale media condizionata E[X|Y], quindi, eredita la propria distribuzione da quella della variabile casuale condizionante Y e lo stesso vale anche per la media di ogni funzione deterministica di X.  condizionata  2 Pertanto, avremo, ad esempio, EY E[X |Y] = E[X2 ]. Il risultato della (A1-3) può essere utilizzato per scomporre la varianza della variabile X nella somma della varianza (condizionata) della funzione di regressione e della varianza (condizionata) dei residui di X rispetto alla funzione di regressione. Premesso che la varianza può essere sempre espressa come differenza tra la media dei quadrati e il quadrato della media, per una generica variabile casuale Z:       Var[Z] = E (Z − E[Z])2 = E Z 2 + (E[Z])2 − 2Z E[Z] = E Z 2 − (E[Z])2 (A1-5) calcoliamo il valore atteso della varianza condizionata di X rispetto a Y:   EY [VarY [X|Y]] = EY E[X2 |Y] − (E[X|Y])2     = EY E[X2 |Y] − EY (E[X|Y])2 15

  Sommando e sottraendo (E[X])2 (ricordando la (A1-3) e che EY E[X2 |Y] = E[X2 ]):   EY [VarY [X|Y]] = E[X2 ] − EY (E[X|Y])2 − (E[X])2 + (E[X])2  
= E[X2 ] − (E[X])2 − EY (E[X|Y])2 − (EY [E[X|Y]])2 applicando la (A1-5) possiamo scrivere: EY [VarY [X|Y]] = Var[X] − VarY [E[X|Y]] e infine otteniamo: Var[X] = VarY [E[X|Y]] + EY [VarY [X|Y]]

(A1-6)

Quindi, la varianza di X, Var[X], è data dalla somma della varianza condizionata della funzione di regressione VarY [E[X|Y]] e della media condizionata delle varianze di X rispetto alla funzione di regressione EY [VarY [X|Y]]. Questo risultato resta così dimostrato in termini generali, ma, per fissare le idee, proviamo ad applicarlo all’esemplificazione riportata nelle tabelle precedenti, utilizzate nel calcolo dei valori attesi condizionati espressi dalla (A1-3) e dalla (A1-4). Siamo condotti alla seguente esemplificazione dei termini della (A1-6):   VarY [E[X|Y]] = EY (E[X|Y])2 − (EY [E[X|Y]])2 = = pY1 (E[X|Y1 ])2 + pY2 (E[X|Y2 ])2 − (E[X])2

(A1-7)

EY [VarY [X|Y]] =pY1 VarY1 [X|Y1 ] + pY2 VarY2 [X|Y2 ] =
=pY1 pX 1 |Y1 (X 1 − E[X|Y1 ])2 + pX 2 |Y1 (X 2 − E[X|Y1 ])2 +
pY2 pX 1 |Y2 (X 1 − E[X|Y2 ])2 + pX 2 |Y2 (X 2 − E[X|Y2 ])2 = =pX 1 ,Y1 (X 1 − E[X|Y1 ])2 + pX 2 ,Y1 (X 2 − E[X|Y1 ])2 + pX 1 ,Y2 (X 1 − E[X|Y2 ])2 + pX 2 ,Y2 (X 2 − E[X|Y2 ])2 =
=pX 1 ,Y1 X 1 2 + (E[X|Y1 ])2 − 2X 1 E[X|Y1 ] +
pX 2 ,Y1 X 2 2 + (E[X|Y1 ])2 − 2X 2 E[X|Y1 ] +
pX 1 ,Y2 X 1 2 + (E[X|Y2 ])2 − 2X 1 E[X|Y2 ] +
pX 2 ,Y2 X 2 2 + (E[X|Y2 ])2 − 2X 2 E[X|Y2 ] =   =E[X2 ] + pY1 (E[X|Y1 ])2 + pY2 (E[X|Y2 ])2 − 2EY (E[X|Y])2

(A1-8)

Sommando la (A1-7) e la (A1-8): VarY [E[X|Y]] + EY [VarY [X|Y]] = E[X2 ] − (E[X])2 = Var[X].

A2. Indipendenza tra trasformazioni di vettori di v. c. normali Sia z un vettore ad n componenti di variabili casuali normali standardizzate e incorrelate (quindi indipendenti) e sia q(z) una forma quadratica con argomento z: q(z) = z0 Az

z ∼ N(0,I)

16

(A2-1)

dove A è una matrice simmetrica n × n. Dal teorema spettrale33 sappiamo che A è diagonalizzabile mediante una similitudine ortogonale. Sia C la matrice ortogonale (C0 C = I) che realizza la diagonalizzazione di A:34 C0 AC = L

⇐⇒

A = CLC0

(A2-2)

dove L è la matrice diagonale degli autovalori {λi } di A.35 Sostituendo A con CLC0 nella definizione di q(z): q(z) = z0 CLC0 z

(A2-3)

Definito il vettore y = C0 z, abbiamo: E[y] = C0 E[z] = 0

Var[y] = C0 Var[z]C = C0 IC = I

y ∼ N(0,I)

(A2-4)

Quindi, la forma quadratica q(z) può essere rappresentata nel modo seguente: 0

q(z) = y Ly =

n X

λi yi2

(A2-5)

i=1

Se la matrice A, associata a q(z), è anche idempotente (oltre che simmetrica), i suoi autovalori possono assumere solo i valori 1 o 0; in questo caso, grazie alla (A2-5), possiamo stabilire che:36 0

0

q(z) = z Az = y Ly =

n−k X

2 y 2j ∼ χ(n−k)

(A2-6)

j=1

dove j indicizza gli n − k autovalori non nulli di A che compongono ordinatamente la diagonale di L. In particolare, se la matrice associata alla forma quadratica qA(z) = z0 Az ha rango pari a n (essendo A simmetrica ed idempotente: A = I), avremo: 2 q(z) = z0 Az ∼ χ(n)

(A2-7)

Siano qA(z) = z0 Az e qB (z) = z0 Bz due forme quadratiche idempotenti, esse sono indipendenti solo se AB = 0.37 Infatti, posto yA = Az e yB = Bz, risulta qA(z) = yA0 yA e qB (z) = yB 0 yB (per la simmetria e l’idempotenza di A e B). Queste due funzioni, rispettivamente di yA e yB , sono indipendenti solo se sono indipendenti yA e yB . In altri termini, dato che E[yA] = 0 e E[yB ] = 0, la matrice delle covarianze deve soddisfare la condizione E[yAyB 0 ] = 0. Ma: E[yAyB 0 ] = E[Azz0 B0 ] = AE[zz0 ]B0 = AIB0 = AB = 0. Questo risultato permette di stabilire l’indipendenza tra numeratore e denominatore della definizione della v. c. F di Fisher-Snedecor e della v. c. t di Student.38 33

Lang (1984, p. 245). Le colonne della matrice C costituiscono una base ortonormale dello spazio Rn . 35 Per costruzione le matrici A e L hanno: · gli stessi autovalori e, quindi, · la stessa traccia (pari alla somma degli autovalori), · lo stesso determinante (pari al prodotto degli autovalori), · lo stesso rango (pari al numero degli autovalori non nulli). 36 V. nota 10 e Piccolo (1998, p. 508). 37 Questa affermazione è nota come teorema di Craig e consente di dimostrare che media e varianza di un campione casuale estratto dal una v. c. normale sono tra loro indipendenti (e viceversa, se media e varianza campionarie sono indipendenti allora il campione è estratto da una popolazione che si distribuisce in modo normale) Piccolo (1998, p. 508). 38 Riguardo alla v. c. t di Student, la forma lineare Bz e la forma quadratica z0 Az sono indipendenti solo se AB = 0. La dimostrazione si ottine allo stesso modo di quella fornita per due forme quadratiche. 34

17

A3. Indice di condizionamento Dal punto di visto operativo, il calcolo dei parametri βˆ richiede l’inversione della matrice X0 X, come mostrato dalla (3). Il punto è che la “precisione” dell’operazione di inversione di una matrice, su qualunque calcolatore, è segnata dai limiti posti dalla precisione di macchina. Pertanto, quando una matrice è quasi-singolare (cioè non molto diversa da una matrice singolare), l’approssimazione per mezzo dei numeri di macchina dei risultati nelle varie fasi della procedura di calcolo dell’inversa può generare un output molto diverso da quella che è la determinazione "esatta".39 Questo tipo di errore, che non è contemplato dalla struttura probabilistica del modello di regressione lineare, determina la differnza tra la determinazione esatta di βˆ e quella ottenuta dai calcoli βˆ∗ . In termini percentuali, possiamo esprimere tale differenza come: e=

ˆ ||βˆ∗ − β|| ˆ ||β||

(A3-1)

Considerando le condizioni più semplici possibili, per misurare e, posto A = X0 X e c = X0 Y, introduciamo una perturbazione nella variabile Y (diciamo δY, cosi che δc = X0 δY) per la quale, mentre si intende risolvere il sistema Aβˆ = c, in realtà il calcolatore risolve il sistema Aβˆ∗ = c + δc. Essendo βˆ = A−1 c e βˆ∗ = A−1 c + A−1 δc, avremo: ˆ ||βˆ∗ − β|| ||A−1 δc|| = e= (A3-2) ˆ ||A−1 c|| ||β|| Il massimo valore del numeratore ||A−1 δc|| può essere determinato considerando che la matrice A, di dimensione k × k, è simmetrica e definita positiva e che, pertanto, può essere diagonalizzata (come nel paragrafo precedente) mediante una base ortonormale di vettori vi che generano Rk . Ogni vettore s ∈ Rk può essere rappresentato come combinazione lineare dei vettori vi mediante opportuni ed univoci coefficienti t i , Pk ottenendo: s = i=1 t i vi . Visto che ogni vi è un autovettore di A (per cui Avi = λi vi ), che vi 0 vi = 1 e vi 0 v j6=i = 0, abbiamo: ‚ ||As||2 = (As)0 As =

k X

Œ0 ‚ t i Avi

k X

i=1

‚ =

k X

Œ‚ t i λi vi 0

i=1

Essendo ||s||2 =

Pk

i=1 t i 2

2

Πt i Avi

i=1 k X

Œ t i λi vi =

k X

i=1

t i 2 λi 2

(A3-3)

i=1

, indicando con λmax il massimo autovalore di A, abbiamo:

||As|| =

k X

t i 2 λi 2 ≤ λmax 2 ||s||2 ⇔ ||As|| ≤ λmax ||s||

(A3-4)

i=1

Siccome la matrice A−1 è anch’essa simmetrica definita positiva e ha come spettro i reciproci degli autovalori di A, analogamente alla (A3-4): ||A−1 s|| ≤ 39

1 λmin

||s||

(A3-5)

l’indice Variance inflation factor può essere utilizzato per misurare l’effetto della collinearità sulla stima dei parametri. Greene (2012).

18

Quindi, per il numeratore della (A3-2) abbiamo ||A−1 δc|| ≤

||δc|| λmin

e, sostituendo:

ˆ ||βˆ∗ − β|| 1 ||δc|| ≤ e= ˆ ˆ λmin ||β|| ||β|| Il minimo valore del denominatore della (A3-2) sarà dato da quindi, sostituendolo nella (A3-6), abbiamo: e=

(A3-6) ||c|| λmax

ˆ e, ≤ ||A−1 c|| = ||β||

ˆ λmax ||δc|| ||βˆ∗ − β|| ≤ ˆ λmin ||c|| ||β||

(A3-7)

λ

L’indice di condizionamento λmax esprime il rapporto massimo possibile tra l’errore min relativo iniziale sui dati (rappresentati mediante numeri di macchina) e l’errore relativo finale commesso sui risultati (anche questi espressi in numeri di macchina). Il valore ideale minimo dell’indice di condizionamento è 1, al suo aumentare i risultati dell’elaborazione sono sempre meno accurati e se tende all’infinito l’elaborazione è impossibile (la matrice A è singolare). Il reciproco dell’indice di condizionamento, quindi, rappresenta, in termini relativi, la distanza della matrice A dalla matrice singolare ad essa più vicina, una distanza che va da 0 (se A è singolare) a 1 (nessun errore λ nell’elaborazione). L’indice λmax può essere interpretato come l’elasticità massima (o min massima reattività) dell’errore nei risultati rispetto agli errori in fase di immissione (e in qualunque altra fase dell’elaborazione) dei dati. Questi risultati possono essere generalizzati in vari modi, ad esempio considerando perturbazioni anche nella matrice A, oltre che nel vettore c, oppure analizzando il caso delle matrici non simmetriche, ottenendo, ovviamente, per il calcolo dell’indice, formule diverse dalla precedente.40

A4. Intervallo di previsione Date le ipotesi del teorema di Gauss-Markov, Yˆi = xi βˆ è la stima BLUE di EX [Yi |xi ] = xi β, l’ampiezza dell’intervallo nel quale ricade l’errore di previsione della funzione di regressione può essere stimata utilizzando la varianza campionaria del residuo ui = Yˆi − Yi = xi (βˆ − β) + εi : ”

0 — VarX [ui |X ] = E xi (βˆ − β) + εi xi (βˆ − β) + εi € Š −1 = σ2 1 + xi (X0 X) xi 0 (A4-1) −1

Quindi Su2 = Su2 (1 + xi (X0 X) i

xi 0 ).

Se accettiamo anche l’ipotesi di normalità degli errori, abbiamo che Pertanto l’intervallo di previsione di Yi sarà dato da: q Yˆi ± t (n−k,1−α/2) Sui = Yˆi ± t (n−k,1−α/2) Su (1 + xi (X0 X)−1 xi 0 )

Yi −Yˆi Sui

∼ t (n−k) .

(A4-2)

dove α ∈ [0,1] è un opportuno livello di fiducia assegnato. Nell’approccio della teoria asintotica la v. c. di riferimento non sarà più la t ma la z. 40

Nell’esposizione si è seguito Quarteroni e Saleri (2006, pp. 134 e ss), per una trattazione generale dell’indice di condizionamento Murli (2007, pp. 221 e ss.).

19

Se si vuole determinare l’ampiezza dell’intervallo di confidenza per l’errore che si commette prevedendo Yi con il suo valore medio E[Yi |xi ], allora Su2 deve essere i

−1 Su2 (xi (X0 X) xi 0 )

stimato con (poiché il valore atteso di εi è zero). Questo l’intervallo di confidenza sarà molto meno ampio rispetto a quello determinato con la (A4-2).

A5. Trasformazione di Box e Cox Talvolta è possibile applicare una opportuna trasformazione alla variabile y per fare in modo di ricondurla alla “normalità” quando questa condizione non è soddisfatta. La trasformazione di Box e Cox (1964) può essere impiegata a questo fine, essa ha il vantaggio di richiedere la stima di un unico parametro aggiuntivo λ. Poniamo che la variabile yλ ( y), trasformata della y, abbia una distribuzione normale e indichiamo con f ( yλ ) la sua funzione di densità di probabilità (la pdf, cosi che F ( yλ ) è la sua primitiva, la cdf), allora avremo che (considerando la funzione composta F ( yλ ( y))) dy la funzione di densità di y sarà data da f ( yλ ) d yλ . In pratica la trasformazione in y λ −1

questione è yλ = λ ma se λ = 0 allora si pone yλ = ln( y) (che è il limite per λ → 0 della trasformazione originaria). Con questo accorgimento la trasformazione è definita su tutto l’asse reale rispetto a λ ma solo per valori non negativi di y. Si tratta, quindi, di rendere massima la verosimiglianza rispetto al fatto che il campione sia stato estratto da una popolazione (rappresentata da yλ ) distribuita in modo normale. La funzione di verosimiglianza è la seguente: L (λ) =

n Y

f ( yi ) =

i=1

=

n Y

n Y

f ( yλ i )

i=1

f

( yλ i ) yiλ−1

d yλi d yi

= (2πσλ2 )

−n 2

e

−

Pn

i=1

( yλi −µλ )2 2σ2 λ

n Y

yiλ−1

(A5-1)

i=1

i=1

e passando al logaritmo: n n 2 X n log 2π n log σλ X ( yλi − µλ )2 l o g L (λ) = − − − + (λ − 1) log yi 2 2 2σλ2 i=1 i=1

(A5-2)

raggruppando i termini costati in C, bisogna risolvere il problema: n n X X 1 2 max log L (λ) = C − log ( yλi − µλ ) + (λ − 1) log yi λ 2 i=1 i=1

(A5-3)

per determinare il valore di λ.

A6. Test di Wald Cosideriamo due modelli econometrici: un modello ristretto con m regressori è un modello completo con k > m regressori. Nei due modelli la variabile dipendente Y è la stessa e i primi m regressori sono gli stessi. Date le caratteristiche geometriche dei minimi quadrati, il modello ristretto produce dei residui che hanno una maggiore varianza rispetto al modello completo. Il test di Wald si pone l’obiettivo di verificare il potere esplicativo dei k − m regressori aggiunti nel modello completo rispetto a quello 20

ristretto. Il test è quindi basato sul confronto tra le varianze dei residui dei due modelli. Se i k − m regressori aggiunti nel modello completo hanno un buon potere esplicativo nei confronti della variabile Y , allora la varianza dei residui u del modello completo deve essere significativamente inferiore alla varianza dei residui u∗ prodotti dal modello ristretto. Il grafico della figura 3 (nel quale z = Y e y è il regressore aggiunto) illustra il significato geometrico di queste proposizioni: z

u

u*

y Hz

x

figura 3: Test di Wald In formule:41 u∗ = Y − Xβ∗ = Y − Xβ∗ + Xβ − Xβ = Y − Xβ + X(β − β∗ ) = u + X(β − β∗ ) u∗ 0 u∗ = u0 u + (β − β∗ )0 X0 X(β − β∗ ) ≥ u0 u perché il doppio prodotto u0 X(β − β∗ ) è nullo e u∗ 0 u∗ − u0 u = (β − β∗ )0 X0 X(β − β∗ ) = ||X(β − β∗ )||2 ≥ 0

(A6-1)

Sotto l’ipotesi H0 che le varianze di u∗ e u siano uguali, considerato che lo stimatore al primo membro della (A6-1) (diviso per σ2 ) si distribuisce come una v. c. χ 2 che ha (n − m) − (n − k) = k − m gradi di libertà42 e che lo stimatore u0 u (diviso per σ2 ) si distribuisce come una variabile casuale χ 2 con n − k gradi di libertà, trattandosi di due stimatori indipendenti, una volta divisi entrambi per i rispettivi gradi di libertà, il loro rapporto si distribuisce come una v. c. F di Fisher-Snedecor con k − m e n − k gradi di libertà:43 u∗ 0 u∗ −u0 u k−m u0 u n−k

∼ F(k−m,n−k)

41

(A6-2)

Ancora una volta, per non appesantire la notazione, indichiamo con β i parametri stimati e non le v. c. da cui sono generati. 42 Per il teorema di Fisher e Cochran, v. Piccolo (1998, p. 508), 43 Per la dimostrazione di queste affermazioni v. Piccolo (1998, pp. 507 e ss.). Il test di Wald può essere costruito con il metodo del rapporto di verosimiglianza, v. Piccolo (1998, pp. 639 e ss.).

21

Utilizziamo le seguenti uguaglianze per definire una misura del fitting della regressione: 0 ˆ ˆ (Y−Y) (Y−Y) 44 ˆ = Xβ e l’indice di determinazione, indicato con R2 := . Ricordando che Y Y = Xβ:

(Y−Y)0 (Y−Y)

45

Y − Y = (X − X)β + u ˆ − Y)0 (Y ˆ − Y) + u0 u (Y − Y)0 (Y − Y) = β 0 (X − X)0 (X − X)β + u0 u = (Y ˆ − Y) ˆ − Y)0 (Y (Y (Y − Y)0 (Y − Y) u0 u (Y − Y)0 (Y − Y)

+

u0 u (Y − Y)0 (Y − Y)

=1

= 1 − R2

(A6-3)

Dividendo numeratore e denominatore della (A6-2) per (Y − Y)0 (Y − Y), possiamo verificare la seguente uguaglianza: u0∗ u∗ −u0 u k−m u0 u n−k

=

(R2 − R2∗ )/(k − m) (1 − R2 )/(n − k)

∼ F(k−m,n−k)

(A6-4)

dove R2∗ è l’R-quadro (l’indice di determinazione) del modello ristretto. Se consideriamo il modello ristretto con la sola costante (l’intercetta verticale), per qualunque modello completo con k > 1 regressori, la (A6-4) diventa:46 0 ˆ ˆ (Y−Y) (Y−Y) k−1 u0 u n−k

=

R2 /(k − 1) ∼ F(k−1,n−k) (1 − R2 )/(n − k)

(A6-5)

Infatti, nel modello restretto (con la sola intercetta verticale) l’indice R2∗ è pari a zero. In questo caso il test di Wald consente di esprimere un giudizio statistico complessivo su tutti i parametri del modello di regressione (escluso β0 ). Naturalmente, se k = 2, questo test produce gli stessi risultati del test t condotto sotto l’ipotesi H0 che β1 sia 2 nullo. In generale F(1,n−2) = t (n−2) , ma si può anche verificare che, per il modello di regressione semplice (con soli due parametri), la statistica al primo membro della (A6-5) è pari al quadrato della statistica che compare al primo membro della (19). 44

Questa definizione dell’R2 è frutto della scomposizione della devianza totale della Y nella somma della ˆ e della devianza dei residui u e non rispetta il criterio generale di parsimonia devianza della stima Y che dovrebbe sovrintendere alla specificazione del modello di regressione. Una misura corretta, per tener conto della perdita di gradi di libertà all’aumentare del numero dei regressori, è data dall’indice S2

di determinazione corretto R2c = 1 − S 2u . Altre misure di fitting proposte per selezionare la migliore Y

specificazione del modello sono l’Akaike Information Criterion (AIC), il Bayesian Information Criterion (BIC) e il Criterio di Schwarz. In ogni caso questi indici cercano di coniugare l’obiettivo di ridurre la varianza dei residui con la necessità di rendere semplice il percorso esplicativo che mira ad interpretare la variabile dipendente, riducendo al minimo il numero di regressori (si tratta, in ogni caso, di proporre una soluzione accettabile al trade-off che sussite tra questi due obiettivi). 45 I simboli Y e X indicano rispettivamente la media campionaria di Y e il vettore riga delle medie campionarie delle colonne della matrice X. 46 Siccome in questo caso u0∗ u∗ = (Y − Y)0 (Y − Y), dal secondo passaggio della (A6-3) segue che: ˆ − Y)0 (Y ˆ − Y) u0∗ u∗ − u0 u = (Y − Y)0 (Y − Y) − u0 u = (Y

22

Una ulteriore applicazione del test di Wald è il test di Chow per il confronto di due identici modelli definiti su set di dati alternativi e tra loro indipendenti. L’ipotesi H0 è che i due modelli presentano gli stessi coefficienti contro l’ipotesi alternativa che i coefficienti sono diversi. Il test di Chow può essere impiegato per verificare la stabilità strutturale del modello (che consegue all’accettazione di H0), oppure per individuare un break strutturale verificatosi nel tempo o per quantificare l’impatto di un trattamento (ad esempio un intervento di policy) sulla variabile dipendente (in questi due casi l’ipotesi statistica è confermata se si rifiuta H0). Il modello vincolato stima i parametri sulla totalità delle osservazioni (producendo, così, una maggiore devianza dei residui u∗ ). Il modello non vincolato stima gli stessi parametri su dei sottoinsiemi di dati (tipicamente ripartiti in due sottoinsiemi). La devianza dei residui del modello non vincolato u0 u = u1 0 u1 + u2 0 u2 è minore di quella del modello vincolato per costruzione ed è ottenuta come somma dei residui stimati nei diversi sottoinsiemi nei quali è stato diviso il campione. Supponendo che i dati costituiti da n osservazioni siano stati ripartiti in due sottoinsiemi di numerosità n1 e n2 , con riferimento alla (A6-2), i gradi di libertà del numeratore saranno dati da (n − k) − [(n1 − k) + (n2 − k)] = k e quelli del denominatore da (n1 − k) + (n2 − k) = n1 + n2 − 2k = n − 2k: u∗ 0 u∗ −u0 u k u0 u n−2k

∼ F(k,n−2k)

(A6-6)

A7. Stima numerica di massima verosimiglianza Sia `(β) il logaritmo della funzione di verosimiglianza L (x;β), la massimizzazione ∂ `(β)

di quest’ultima si effettua risolvendo l’equazione di verosimiglianza ∂ β = 0. Molto spesso questa equazione non ha una forma risolvibile rispetto a β, allora è necessario utilizzare dei metodi numerici applicati allo sviluppo in serie di Taylor dell’equazione di verosimiglianza (nel punto β0 ): ˆ ∂ `(β0 ) ∂ 2 `(β0 ) ∂ `(β) ' + (βˆ − β0 ) '0 ∂β ∂β ∂ β2  2 −1 ∂ `(β ) ∂ `(β0 ) 0 βˆ ' β0 − 2 ∂β ∂β

(A7-1)

ˆ Quindi, è possibile utilizzare un processo iterativo per trovare la soluzione β, basta sostituire βˆ con βn+1 e β0 con βn ad ogni iterazione fino a quando |βn+1 − βn | non risulta minore ad un valore piccolo fissato a piacere (eventualmente, il doppio della precisione di macchina). Tale procedimento è detto metodo di Newton, la sua f (x) formulazione generale è: x n+1 = x n − f 0 (x) . La convergenza a βˆ non è sempre garantita, tutto dipende dalle caratteristiche della funzione di cui bisogna determinare gli zeri. Non è detto che le derivate prime della funzione f (x) siano esprimibili in forma analitica, in tal caso è necessario una loro approssimazione (ottenuta calcolando la pendenza della retta che interseca la f (x) in due punti vicini). In ogni caso il metodo può fallire per vari motivi.

23

A8. Varianza dei parametri e log-verosimiglianza R∞ Sia f (x;β) una funzione di densità di probabilità, per definizione −∞ f (x;β)d x = 1 è costante rispetto a β. Quindi, se vi sono le condizioni per applicare il teorema di Libeniz, differenziando rispetto a β: ∂ ∂β

∞

∞

∂f f f (x;β)d x = dx = −∞ ∂ β f −∞

Z

Z

∞

∂ log f f dx = −∞ ∂ β

Z

Z

`0 f = E[`0 ] = 0 (A8-1)

dove ` è la funzione di log-verosimiglianza e `0 è la sua derivata rispetto a β (tale derivata è detta anche score). La somma delle varianze di ogni componente di `0 (una per ogni osservazione campionaria) è definita informazione di Fisher ed è indicata con (nI) := Var[`0 ]. Il valore I indica di quanto varia la pendenza del piano tangente alla funzione di log-verosimiglianza nell’intorno di β. Quanto maggiore è questa varianza, tanto più le sezioni superiori di ` saranno “vicine” alla perpendicolare passante per β. Consideriamo la derivata seconda `00 di ` rispetto a β, per essa si ha: `00 =

f 00 f − f 0 f 0 f 00 ∂ `0 ∂ f0 = = − (`0 )2 = ∂β ∂β f f2 f

Prendendo il valore atteso di ambo i membri della (A8-2): Z 00 Z  0 2  0 2  0 2 f ∂2 00 E[` ] = f − E (` ) = f − E (` ) = −E (` ) = −Var[`0 ] 2 f ∂β

(A8-2)

(A8-3)

  Quindi: Var[`0 ] = E (`0 )2 = −E[`00 ] (la varianza di `0 è pari all’opposto della matrice hessiana della funzione di log-verosimiglianza). Se nella (A7-1) si pone β0 uguale al valore vero di β, allora: ∂ 2 `(β) βˆ − β = − ∂ β2 

−1

∂ `(β) `0 = − 00 ∂β `

e, moltiplicando ambo i membri della (A8-4) per naria): p

p

(A8-4)

n (dove n è la numerosità campio-

 ‹−1 `0 n `0 `00 ˆ n(β − β) = − p 00 = − p n` n n

(A8-5)

Essendo `0 pari alla sommma di tutte le sue componenti (log f (X i ;β)), per il teorema 0 p 0 limite centrale, la distribuzione del termine p` n = n `n tende ad essere quella di v. c. normale con media zero è varianza nI: N (0,nI). Per la legge dei grandi numeri il 00 termine `n (pari alla media delle varianze delle singole componenti di `0 ) si approssima p a nI all’aumentare di n. Pertanto, n(βˆ − β) tende ad essere distribuito come una v. c. normale N (0,(nI)−2 (nI)) = N (0,(nI)−1 ). In definitiva, all’aumentare di n, la distribuzione di βˆ risulta approssimarsi a quella di v. c. normale N (β,(nI)−1 ) dove (nI)−1 è l’inversa della matrice delle varianze-covarianze dello score (o matrice inversa dell’informazione di Fisher).47 I termini diagonali della matrice (nI)−1 approssimano ˆ La varianza di ogni stimatore non distorto le varianze delle stime dei parametri β. 47

La matrice (nI) è più spesso indicata con I, inglobando il fattore n.

24

ottenuto con il metodo della massima verosimiglianza ha un limite inferiore stabilito ˆ ≥ 1 = 1 0 . Infatti, per uno stimatore dalla disuguaglianza di Cramér-Rao: Var[β] (nI) Var[` ] non distorto: Z Z Z ∂ ˆ 0 f = E[β` ˆ 0] = 1 βˆ f = β ⇒ βˆ f = β` (A8-6) ∂β ricordando che cov[x, y] = E[x y] − E[x]E[ y] e che E(`0 ) = 0, abbiamo: ˆ 0 ] = cov[β` ˆ 0] = 1 E[β` e, per la disuguaglianza di Schwarz:48  0 ˆ 0 ] 2 ≤ Var[β]Var[` ˆ 1 = cov[β` ] che dimostra la disuguaglianza di Cramér-Rao. I precedenti risultati forniscono fondamento teorico a gran parte della statistica inferenziale e, in particolare, consentono di capire la logica di uno dei test statistici più generali, il test del rapporto di verosimiglianza. Tale test è basato sul confronto di due massimi della funzione di verosimiglianza, il primo relativo ad un modello ristretto (nel quale, ad esempio, alcuni parametri sono posti pari a zero) e il secondo relativo ad un modello senza restrizioni (che, ad esempio, include tutti i parametri). Se la funzione di verosimiglianza del modello ristretto raggiunge un valore massimo significativamente maggiore del massimo della funzione di verosimiglianza del modello senza restrizioni, allora non è possibile rifiutare l’ipotesi H0 che il modello ristretto sia più verosimile del modello senza restrizioni, viceversa è necessario accettare l’ipotesi alternativa (favorevole al modello con restrizioni) e, quindi, rifiutare H0. Si tratta di L prendere una decisione in base al valore assunto dal seguente rapporto: Q = L0 dove la A massima verosimiglianza al numeratore è funzione dei parametri del modello ristretto, mentre quella al denominatore è funzione dei parametri del modello senza restrizioni. Espandendo in serie di Taylor il logaritmo del numeratore di Q (nell’intorno di βA), otteniamo: 1 `0 = `A + (β0 − βA)`0 + (β0 − βA)2 `00 2 0 considerato che (data H0) E[` ] = 0 in βA e che, quindi, in tale punto del dominio di `A, risulta `0 ' 0, possiamo scrivere: −2logQ = −2(`0 − `A) = −(β0 − βA)2 `00 Dai paragrafi precedenti sappiamo che −`00 ' −E[`00 ] è l’informazione di Fisher (nI), pertanto −2(`0 −`A) si distribuisce come (β0 − βA)2 (nI). Ma, per quanto si è dimostrato nel paragrafo precedente, (β0 − βA) si distribuisce come un v. c. normale N (0,(nI)−1 ) e 1

ciò implica che (β0 − βA)(nI) 2 si distribuisce come una v. c. normale standard N (0,1). Quindi, abbiamo che (β0 − βA)2 (nI) si distribuisce come una variabile casuale χ g2 , dove g è la differenza tra il numero di parametri del denominatore e il numero di parametri del numeratore. Si è cosi determinato che la statistica −2logQ si distribuisce come una v. c. χ g2 . Tale risultato è noto come teorema di Wilks (Piccolo, 1998, p. 647). 48

R R R 2 Per dimostrare la disuguaglianza di Cauchy-Schwarz, ( g f ) ≤ ( g 2 )( f 2 ), basta considerare la funzione di t definita come che l’integrale del suo quadrato è sempre non negativo R R R g + t f e verificare R (g + t f )2 = ( g 2 ) + 2t( g f ) + t 2 ( f 2 ) ≥ 0. Ciò può verificarsi solo se il discriminante di questa R R R 2 equazione di secondo grado in t è minore di zero: ( g f ) −( g 2 )( f 2 ) ≤ 0 (cvd). Tale discriminante è uguale a zero solo se g e f sono linearmente dipendenti, per cui esiste un t tale che la funzione R R R 2 g + t f = 0. In questo caso ( g f ) = ( g 2 )( f 2 ).

25

A9. Teorema limite centrale La formulazione più semplice del teorema limite centrale riguarda la convergenza della somma standardizzata di n v. c. X i ∼ iid verso la distribuzione normale standardizzata. Indicando con σ2 la varianza e con µ la media delle v. c. X i : Pn  Pn Pn X − E i i i Xi i X i − nµ q Pn  = p 2 nσ Var i Xi P p  1 Pn n 1 n n i Xi −µ n i Xi −µ ∼ N (0,1) (A9-1) = = σ p σ n Quindi, si ha che all’aumentare di n: p  1 Pn n n i Xi −µ p ∼ N (0,1) ⇔ n(X n − µ) ∼ N (0,σ2 ) σ dove X n è la media campionaria delle n v. c. X i .

26

(A9-2)

B0. Test statistici premessa Nel test delle ipotesi statistiche α è la probabilità di commettere l’errore del Io tipo: rifiutare l’ipotesi H0 quando essa è vera. Fissato α, rifiuto H0, al livello α, se il valore stimato a partire dalla realizzazione campionaria ricade nella Regione Critica (la regione dello spazio parametrico che, dato il campione e il valore prescelto per α, implica il rifiuto dell’ipotesi nulla H0). Nel seguito n è la numerosità del campione, e: n

n 1 X S = (x i − x)2 n − 1 i=1

1X x= xi n i=1

2

(1)

sono rispettivamente la media e la varianza campionarie. In generale faremo riferimento a popolazioni distribuite in modo normale e ipotizzeremo che i campioni siano casuali (rappresentativi), cioè costituiti da realizzazioni di variabili casuali indipendenti ed identicamente distribuite (distribuite come la popolazione).49

B1. Test sulla media 1. Se la varianza σ2 della popolazione dalla quale è stato estratto il campione è nota, allora la seguente statistica (calcolata sotto l’ipotesi H0: la media della popolazione è uguale a µ0 ) è la realizzazione di uno stimatore che si distribuisce come una variabile casuale normale standardizzata: z ∼ N (0,1), zc =

x − µ0 p σ/ n

(2)

l’ipotesi H0, può essere rifiutata al livello di significatività α, se |zc | > zα/2 . Dove, zα/2 è il valore teorico di z che delimita la Regione Critica (corrispondente ai valori di z più distanti dalla media che possono verificarsi con probabilità α).50 2. Se la varianza della popolazione è incognita, lo stimatore della media, che si distribuisce come una variabile casuale t di Student con n − 1 gradi di libertà51 , da luogo al calcolo della seguente statistica: tc =

x − µ0 p S/ n

(3)

l’ipotesi H0, può essere rifiutata al livello di significatività α, se |t c | > t n−1,α/2 . Dove, t n−1,α/2 è il valore teorico di t n−1 che delimita la Regione Critica (corrispondente ai valori di t n−1 più distanti dalla media che possono verificarsi con probabilità α). L’ipotesi potrebbe essere formulata in modo unidirezionale e, in questo caso, il confronto dovrebbe essere effettuato con il valore t n−1,α . 49

V. Piccolo D., 1998, Statistica, Il Mulino, Bologna. V. figura 4. 51 In quanto al punto precedente si è sostituita al denominatore della statistica σ con S. 50

27

Il test sulla differenza tra i valori medi di due campioni casuali indipendenti (di numerosità n ed m) generati da variabili casuali normali è basato sull’ipotesi che le varianze delle due popolazioni, X e Y , siano uguali anche se incognite. L’ipotesi H0 afferma che la differenza tra le medie delle due popolazioni è d0 contro l’ipotesi alternativa H1 che tale differenza sia maggiore di (minore o diversa da) d0 . La statistica da calcolare è: v t nm(n + m − 2) x − y − d0 (4) tc = q n+m (n − 1)S 2 + (m − 1)S 2 x

y

essa dovrà essere confrontata con il valore teorico della variabile casuale t di Student in corrispondenza ad n + m − 2 gradi di libertà ed opportuni livelli di α/2 o α (a secondo che l’ipotesi sia bidirezionale, |t c | > t n+m−2,α/2 , o unidirezionale t c > t n+m−2,α ). In sostanza questo test verifica se le due popolazioni sono somiglianti. Per campioni grandi la precedente statistica può essere semplificata nel modo seguente: zc = q

x − y − d0

(5)

S 2x /n + S 2y /m

da confrontare con i valori della variabile casuale normale standardizzata zα o zα/2 secondo i casi. Se le varianze sono diverse (problema di Behrens-Fisher), un test asintotico (test di Welch), può essere basato sulla seguente statistica: tc = q

x − y − d0

(6)

S 2x /n + S 2y /m

e i gradi di libertà della variabile casuale t di Student saranno dati dall’espressione: S 2x /n − S 2y /m

(7)

S 2y S 2x + n2 (n−1) m2 (m−1)

Il test sulla differenza tra valori medi per dati appaiati è riferito a coppie di variabili casuali. L’ipotesi nulla è specificata come H0 : µ x − µ y = v0 contro l’ipotesi alternativa H1 : µ x − µ y > v0 . La statistica del test è: tc =

d − v0 p Sd / n

(8)

dove: n

di = x i − yi

d=

1X (x i − yi ) n i=1

Sd2 =

n 1 X (di − d)2 n − 1 i=1

(9)

Fissato α, l’ipotesi H0 viene rifiutata, al livello α, se t c risulta maggiore al valore teorico t n−1,α . Anche in questo caso l’ipotesi potrebbe essere bidirezionale e il valore teorico della v. c. t di Student da prendere in considerazione sarà t n−1,α/2 . Un caso particolare di test sulla media è quello del test asintotico sulla proporzione. Nel caso dell’ipotesi semplice unidirezionale, H0 è formulata come p = p0 28

contro H1 : p > p0 . La distribuzione del carattere in esame è di tipo bernoulliano: ogni elemento assume l’attributo 1 con probabilità p e valore 0 con probabilità (1 − p). La media della popolazione è, quindi, pari a p e la varianza della popolazione è p(1 − p). La statistica da calcolare per sottoporre a test l’ipotesi H0 è: zc = p

ˆp − p0

(10)

p0 (1 − p0 )/n

dove ˆp è la media campionaria (la somma dei valori 0 e 1 divisa per n). L’ipotesi H0 verrà rifiutata se zc > zα . Quest’ultimo test, in definitiva, sottopone a verifica un’ipotesi sull’andamento delle frequenze di un attributo. In questa stessa tipologia rientra anche il test sulla differenza tra le frequenze relative di un attributo detto anche test sul confronto tra proporzioni relative a variabili casuali indipendenti. L’ipotesi H0 sostiene che la differenza tra le due proporzioni è pari a d0 e l’ipotesi alternativa H1 è espressa da pˆx − pˆy > d0 . La statistica da calcolare è: zc = p

pˆx − pˆy − d0

(11)

ˆp(1 − ˆp)(1/n + 1/m)

dove: ˆp =

n pˆx + m pˆy

(12)

n+m

Il valore calcolato dovrà essere posto a confronto con la variabile casuale normale standardizzata zα o zα /2 a secondo che l’ipotesi sia rispettivamente unidirezionale (tenuto conto del verso) o bidirezionale.

B2. Test sulla varianza Se non è nota la media della popolazione dalla quale è stato estratto il campione, l’ipotesi H0 che può essere sottoposta a test è la seguente: H0 : σ2 = σ02 contro H1 : σ2 > σ02 . Calcolata la seguente statistica: χc2 =

n−1 σ02

S2

(13)

2 essa deve risultare maggiore del valore teorico χn−1,α corrispondente al livello di significatività α e ad n−1 gradi di libertà. Nel caso sia nota la media bisogna considerare n gradi di libertà. Se l’ipotesi è formulata come: H0 : σ2 = σ02 contro H1 : σ2 6= σ02 allora, il rifiuto di H0, al livello di significatività α, sarà possibile solo se la statistica χc2 assume valori esterni all’intervallo che ha come estremi, rispettivamente, inferiore e superiore i valori teorici: 2 χn−1,1−α/2

2 χn−1,α/2

(14)

Il test sul rapporto tra le varianze di variabili casuali normali e indipendenti sottopone a verifica l’ipotesi H0 : σ2x = σ2y contro l’ipotesi alternativa H1 : σ2x > σ2y . La statistica da calcolare è: Pn 1 2 S 2x n−1 i=1 (x i − x) Fc = 2 = 1 Pm (15) 2 Sy ( y − y) i m−1 i=1 29

e i valori critici teorici della v. c. Fn−1,m−1 di Fisher-Snedecor per n − 1 e m − 1 gradi di libertà sono quelli corrispondenti al livello di significatività α o 1 − α (secondo la direzione), nel caso unidirezionale, oppure α/2 e 1 − α/2 nel caso bidirezionale. È utile tener presente la seguente relazione: Fn−1,m−1,1−α =

1 Fn−1,m−1,α

(16)

Questo test permette di rifiutare l’ipotesi nulla di uguaglianza delle varianze solo quando il numeratore della statistica e molto diverso dal denominatore.

B3. Test non parametrici Il test sulla bontà di adattamento formula l’ipotesi H0 che le frequenze osservate siano uguali alla frequenze teoriche di un preciso modello probabilistico. La statistica necessaria per effettuare questo confronto tra frequenze osservate e frequenze teoriche è la seguente: χc2

k X (ni − npi )2 = npi i=1

(17)

dove ni sono le frequenze osservate per le k classi la cui unione costituisce il supporto delle probabilità pi , quindi i valori npi sono le frequenze teoriche calcolate in base al modello probabilistico scelto dall’ipotesi H0. Fissato il livello α, la regione critica per il 2 rifiuto di H0 è data da χc2 > χk−1,α . Il test sulla indipendenza formula l’ipotesi H0 di indipendenza di due variabili casuali, X e Y . La statistica del test è costruita in base allo scostamento tra frequenze ˆ i j calcolate sotto l’ipotesi di indipendenza: osservate ni j e frequenze teoriche n χc2

=

k X h X ˆ i j )2 (ni j − n i=1 j=1

ˆi j n

(18)

La regione critica per H0 è χc2 > χ(k−1)(h−1),α .

B4. ANOVA L’ANOVA (analisi della varianza) analizza le medie in più di due gruppi mediante il confronto delle devianze. Questa metodologia è basata sulla decomposizione della devianza totale, D, di un fenomeno multivariato nel quale un determinato effetto y può dipendere da k cause (o trattamenti): nj nj k X k X k X X X 2 2 D= ( yi j − y) = ( yi j − y j ) + n j ( y j − y)2 = DW + DB j=1 i=1

j=1 i=1

(19)

j=1

dove y è la media globale, y j è a media del gruppo j soggetto alla causa j, D è la devianza globale, DW è la devianza all’interno dei gruppi (within) e DB è la devianza tra

30

i gruppi (between). Dividendo le devianze per i rispettivi gradi di libertà, si ottengono la varianza within e la varianza between, rispettivamente: 2 SW =

DW n−k

SB2 =

DB k−1

Pn j

( yi j − y j )2

(20)

Considerato che: k X n−k = (n j − 1)

S 2j

=

j=1

i=1

nj −1

(21)

risulta: Pk 2 SW

2 j=1 (n j − 1)S j

= Pk

j=1 (n j − 1)

(22)

per cui la varianza within è una media ponderata delle varianze all’interno di ogni gruppo. Il test statistico F di Fisher-Snedecor può essere utilizzato per stabilire se la varianza between è significativamente maggiore della varianza within. L’ipotesi nulla è che le due varianze siano uguali, il suo rigetto significa che almeno una causa ha agito, modificando la variabilità dell’effetto medio tra i gruppi oltre il livello medio di variabilità registrato all’interno dei gruppi. I gradi di libertà del numeratore della statistica sono k − 1 e quelli del denominatore son n − k, pertanto il valore calcolato dalla statistica: Fc =

SB2 2 SW

(23)

verrà confrontato con i valore teorico, dipendente da α, della v.c. Fk−1,n−k,α . Perché sia valida questa procedura è necessario verificare che le componenti erratiche individuali siano distribuite normalmente.

31

B5. Esempio di individuazione della Regione Critica Il grafico seguente mostra un esempio di individuazione della regione critica per un test unidirezionale sulla media nel caso della v. c. z.

0.3 0.2

RC(α)=zc > zα

0.1

Densita della v.c. normale standardizzata

0.4

Individuazione della Regione Critica per H0

0.0

α −3

−2

−1

0

1

zα 2

3

z

figura 4: La regione critica per H0, RC(α), va da zα a ∞.

B6. Simulazione della regressione lineare in R Per studiare la regressione può essere conveniente simulare il DGP (data generating process) e poi inferire dai dati le caratteristiche del DGP che sono, in realtà, già note. Ciò permette di verificare i risultati e di constatarne la coerenza rispetto alle prescrizioni teoriche del modello di regressione lineare. Poniamo che le osservazioni ( yi , x i ) siano generate da una relazione lineare con errori che si distribuiscono come una v.c. normale standardizzata ε ∼ N (0,1): y i = β0 + β1 x i + ε i dove β0 = 5 e β1 = 1.2 [v. R1 e figura 5]. Uno dei test utilizzati per verificare la corretta specificazione funzionale del modello di regressione è il test di Ramsey (reset test).52 Se il modello è stato specificato correttamente non si può rifiutare l’ipotesi H0: «la forma funzionale è ben specificata» dato il valore del p-value (che deve risultare maggiore di 0.05 per accettare H0). Inoltre, in base alle ipotesi sul modello di regressione, i residui dovrebbero distribuirsi in modo "normale" (vedi grafico in alto a destra della figura 5). L’analisi 52

Per i test citati è necessario caricare alcune librerie - ad. es. lmtest - come specificato negli esempi.

32

15 5

10

y

2

4

6

8

10

x

figura 5: Scatter plot (x i , yi ) grafica consente di avere un’idea della distribuzione degli errori e l’ipotesi di normalità può essere verificata applicando il test di Shapiro-Wilk sulla normalità dei residui ˆei −e (preventivamente standardizzati con la trasformazione σ ), in tal caso l’ipotesi H0: «i e residui sono distribuiti come una v.c. normale» può essere accettata solo se il p-value è sufficientemente elevato. Con obiettivi analoghi, dal grafico della figura 7 (Normal Q-Q Plot), si può verificare se gli errori standardizzati si distribuisco come la v. c. normale standardizzata: i punti dovrebbero disporsi lungo la retta a 45o , se è verificata l’ipotesi di normalità. Riguardo alla corretta specificazione del modello di regressione lineare, è utile esaminare il grafico dei residui rispetto ai valori stimati della variabile dipendente: i punti dovrebbero disporsi orizzontalmente rispetto alle ascisse e simmetricamente rispetto alla linea orizzontale con intercetta zero (v. grafico in basso a sinistra nella figura 6). Lo stesso tipo di grafico può essere tracciato rispetto ad ogni variabile esplicativa: se la disposizione dei residui non è casuale ciò indica che il modello "vero" non è lineare rispetto a quel particolare regressore. Per testare l’eteroschedasticità dei residui si può utilizzare il test di Breusch-Pagan per l’ipotesi H0: «la varianza degli errori e costante» (omoschedasticità) che verrà accettata se il p-value è elevato. Soprattutto per le serie storiche viene utilizzato il test Durbin-Watson per verificare che non vi sia autocorrelazione nei residui sotto l’ipotesi H0: «i residui sono incorrelati» che può essere accettata se il valore del p-value è abbastanza alto. Nello stesso contesto seriale è possibile utilizzare la funzione di autocorrelazione totale (acf) e la funzione di autocorrelazione parziale (pacf) per inferire sul processo ARMA che potrebbe aver generato i residui. In tutti i software statistici l’output del modello di regressione prevede il calcolo 33

0.2 0.0

5

0.1

10

y

Density

15

0.3

0.4

Istogramma degli errori

2

4

6

8

10

−4

−2

x

0

2

4

2

4

e

4 2 −2

0

et

0 −2

Residui

2

4

Residual vs fitted

6

8

10

12

14

16

−2

Fitted

0 et−1

figura 6: Grafici della regressione e dei residui Var(e)

dell’indice di determinazione R2 = 1 − Var( y) e della sua versione corretta (non distorta): (N −1)Var(e)

R2c = 1− (N −k−1)Var( y) (dove k è il numero delle variabili esplicative). l’indice R2 (come la sua versione corretta) indica quanta parte della variabilità complessiva (TSS: Total Sum of Squares) della variabile dipendete y è spiegata dalla regressione ed è, quindi, basato sulla decomposizione della varianza: la variabilità totale è uguale per definizione alla somma della variabilità della regressione (ESS: Explained Sum of Squares della stima ˆy ) e della variabilità degli errori (RSS: Residual Sum of Squares). Un test globale su tutti i parametri del modello di regressione può essere condotto con il test F di Fisher che confronta la variabilità della stima ˆy con la variabilità dei residui ˆe. L’ipotesi nulla è H0:«le variabili esplicative non spiegano la variabile dipendente» (tutti i coefficienti del modello non sono significativamente diversi da zero). Per rifiutare tale ipotesi (ed attribuire qualche validità al modello) è necessario che la Fc calcolata sotto H0 sia sufficientemente maggiore della F t eor ica , ottenuta in corrispondenza di k e n − k − 1 gradi di libertà (e cioè che il p-value sia sufficientemente piccolo). Il test sui singoli coefficienti βˆj formula l’ipotesi H0: «il coefficiente βˆj è, in media, pari a zero» e confronta la statistica test ˆt j =

βˆj ˆβj σ

con la t t eor ica che si distribuisce come la

variabile casuale t di Student ed è ben approssimata dalla normale standardizzata al crescere della numerosità delle osservazioni. La t calcolat a deve essere significativamente maggiore della t t eor ica e, quindi, il p-value deve essere abbastanza piccolo per poter rifiutare H0. Utili informazioni al riguardo possono essere dedotte anche dagli intervalli di confidenza che scaturiscono da questo test (ad es. la posizione dello zero).

34

1 0 −2

−1

Sample Quantiles

2

3

Normal Q−Q Plot

−2

−1

0

1

2

Theoretical Quantiles

figura 7: Adattamento dei residui alla Normale Codice in R In sequenza, cercando di rispettare l’ordine degli argomenti trattati nel testo: R1

#Il DGP dell’esempio n<-500 # numerosit` a del campione e<-rnorm(n,0,1) # genera una v.c. normale x<-runif(n,1,10) # genera la variabile esplicativa x y=5+1.2*x+e # genera la variabile dipendente y plot(x,y) # il grafico cor.test(x,y,method="pearson") # test sulla correlazione R2

# La stima del modello library(lmtest) # carica la libreria lmtest r<-lm(y ~ x) # stima del modello lineare summary(r) # mostra i risultati della stima confint(r) # intervalli di confidenza dei coefficienti names(r) # mostra le informazioni contenute in r e<-residuals(r) # salva i residui stimati in e resettest(r) # reset test sulla regressione t.test(e) # test sulla media dei residui shapiro.test( 35

(e-mean(e))/sd(e)) # test di normalit` a dei residui R3: analisi grafica

# Regressione e distribuzione errori par(mfrow=c(2,2)) # divide il grafico plot(x,y) # nube dei punti abline(r, col="red", lwd=2) # linea di regressione hist(e,freq=F, ylim=c(0,.45), main="Istogramma",cex.main=.92) # istogramma rug(e,lwd=.25,side=1) # evidenziare la concentrazione lines(density(e), lwd=2, col="blue") # densit` a della distribuzione z<-seq(min(e), max(e),.1) # genera il supporto lines(z,dnorm(z, mean=mean(e), sd=sd(e)), col="red", lty=2) # curva della normale standardizzata plot(fitted(r), e, ylab="Residui", xlab="Fitted", main= "Residual vs fitted") # grafico errori verso y stimato abline(h=0, col="red") # linea sulla media (zero) dei residui plot(e[-n], e[-1], xlab=expression(e[t-1]), ylab= expression(e[t])) # grafico della correlazione dei residui abline(0,1,col="red") # bisettrice par(mfrow=c(1,1)) # ripristina la finestra del grafico qqnorm((e-mean(e))/sd(e)) abline(0,1, col="red") # Normal Q-Q plot plot(x, e, ylab="Residui", xlab="x", main="Residual vs x"); abline(h=0, col="red") # grafico errori rispetto ad x pairs(data.frame(x, fitted(r), r$residuals)) # grafico pairs R4: Test sulla regressione e stima robusta

# Altri test sulla regressione vcov(r) # matrice var-covarianze dei coefficienti bptest(r) # Breusch-Pagan test (omoschedasticit` a residui) dwtest(r) # test di Durbin-Watson (autocorrelazione residui) library(lmtest) # carica la libreria lmtest coeftest(r) # test sui coefficienti coeftest(r, df=Inf, vcov=vcovHAC(r)) # test robusto sui coefficienti library(sandwich) # carica la libreria sandwich print(vcovHC(r)) # varcov dei parametri con eteroschedasticit` a print(vcovHAC(r)) # con eteroschedasticit` a ed autocorrelazione library(robust) # libreria per la stima robusta lmRob(y ~ x) # regressione robusta summary(lmRob(y ~ x)) # riepilogo dei risultati della stima

36

Ulteriori comandi in R: Per salvare il grafico corrente nel file "scatterxy" in formato eps:

par(family="Times") # da porre all’inzio dello script dev.copy(postscript, file="scatterxy.eps", height=6, width=6, horizontal=F, onefile=F, fonts="Times") dev.off()

Riferimenti bibliografici Belsey, D. A., Kuh, E. e Welsch, R. E. (1980), Regression diagnostics, Wiley. Box, G. E. P. e Cox, D. R. (1964), ‘An analysis of transformations’, Journal of the Royal Statistical Society Vol. 26(2), 211–252. Davidson, R. e MacKinnon, J. G. (2003), Econometric Theory and Methods, Oxford University Press, New York. Frisch, R. e Waugh, F. V. (1956), ‘Partial time regressions as compared with individual trends’, Econometrica Vol. 1(4), 387–401. Greene, W. H. (2012), Econometric Analysis - 7/E., Prentice Hall. Lang, S. (1984), Algebra lineare - 7/E., Boringhieri, Torino. Lovell, M. (1963), ‘Seasonal adjustment of economic time series’, Journal of the American Statistical Association Vol. 58, 993–1010. Murli, A. (2007), Matematica numerica: Metodi, Algoritmi e Software. Parte prima, Liguori Editore, Napoli. Piccolo, D. (1998), Statistica, Il Mulino, Bologna. Quarteroni, A. e Saleri, F. (2006), Introduzione al calcolo scientifico, Springer Verlag. Weisberg, S. (2005), Applied Linear Regression - ed. III, Wiley.

37

Indice Il modello di regressione lineare

1

La distorsione da variabile omessa

2

Sfericità degli errori e teorema di Gauss-Markov

3

La varianza degli errori

5

L’ipotesi di normalità degli errori

6

Cenni alle proprietà asintotiche dei minimi quadrati

7

A0. Il teorema di Frisch-Waugh-Lovell e i valori anomali

10

A1. Medie e varianze condizionate

14

A2. Indipendenza tra trasformazioni di vettori di v. c. normali

16

A3. Indice di condizionamento

18

A4. Intervallo di previsione

19

A5. Trasformazione di Box e Cox

20

A6. Test di Wald

20

A7. Stima numerica di massima verosimiglianza

23

A8. Varianza dei parametri e log-verosimiglianza

24

A9. Teorema limite centrale

26

B0. Test statistici premessa

27

B1. Test sulla media

27

B2. Test sulla varianza

29

B3. Test non parametrici

30

B4. ANOVA

30

B5. Esempio di individuazione della Regione Critica

32

B6. Simulazione della regressione lineare in R

32

Bibliografia

37

38