Wygodny jest zapis macierzowy zależności (1): Przy zapisie wektora w postaci kolumny oraz oznaczeniu transpozycji przez z’ wektora z zależność (1) można zapisać w postaci

ANALIZA ZALEŻNOŚCI WIELU ZMIENNYCH.

REGRESJA LINIOWA WIELOKROTNA.

I. Model liniowy regresji wielokrotnej.

^{Yi 01xi12xi2...p1xi,p1i}, (1)

n

i 1,..., , ^{1,2,...,n} są niezależnymi zmiennymi

losowymi o takim samym rozkładzie z wartością średnią 0 i wariancją ^2. ^i są błędami losowymi.

Założenia:

(i) Obserwujemy wartości zmiennych ^{Y ,...,1 Yn} (zmiennych objaśnianych).

(ii) ^{xi1,...,xi,p1,i1,...,n}, są znane ( zmienne objaśniające )

(iii) ^{0,1,...,p1,2} są nieznanymi parametrami modelu

(iv) ^{i ~ N(0,)}, ^{i 1,...,n} (losowe błędy ).

Cel eksperymentu – wnioskowanie na temat parametrów modelu

Wygodny jest zapis macierzowy zależności (1):

Przy zapisie wektora w postaci kolumny oraz

oznaczeniu transpozycji przez z’ wektora z zależność (1) można zapisać w postaci

Y = X^{ }, (2) gdzie Y = ^{(Y1,Y2,...,Yn)'} jest wektorem zmiennych

objaśnianych, ^{'(0,1,...,p1)} jest wektorem nieznanych współczynników, a ^{'(1,2,...,n)} wektorem błędów

losowych. Ponadto X jest macierzą wymiaru ^{n p} postaci, zawierającą zmienne objaśniające:

X =



























































1 , 1

1 , 2 21

1 , 1 11

1 1 1

p n n

p p

x x

x x

x x

.

Równanie (2) z przyjętymi założeniami nazywamy liniowym modelem regresji wielokrotnej.

Uwaga. Szczególnymi przypadkami modelu (2) są:

(a) model regresji jednokrotnej (liniowej), gdy ^p2 Y =































 





































n n

x x









1

1 0

1 , 1 , 1

1 1

.

Wyraz wolny^0 można traktować jako współczynnik odpowiadający dodatkowej zmiennej objaśniającej

(b) prosta próba losowa^:



























Yn

Y₁

=



























1 1

 ₊



























n

₁

,

gdzie ^{E(Yi) }, ^{i Yi} są niezależnymi „błędami”,

. ,...,

1 n

i 

Własności wektora losowego Y = ^(Y1,...,Yn)'.

 ) (Y_i

E E(₀₁x_i1..._p1x_i,p1)E(_i), skąd

Yi

 = ^{0 1xi1...p1xi,p1} = ^'xⁱ, gdzie ^{'(0,1,...,p1)}, x’ⁱ= ^{(1,xi1,...,xi,p1)}.

Niech ^{Y (Y1,...,Yn)'}. Wówczas

Y 

 X^.

Var(^Yi) = Var(^{0 1xi1...p1xi,p1i)} =

Var(^i)2

Cov(^Yi,Yj) = 0 dla ^{i j}, gdyż ^{Y ,...,1 Yn} są niezależne.

Stąd, definiując macierz kowariancji wektora losowego Y :

_Y =





















































) , ( )

, ( )

, (

) , ( )

, ( )

, (

) , ( )

, ( )

, (

2 1

,

2 2

2 1

2

1 2

1 1

1

n n n

n

n n

Y Y Cov Y

Y Cov Y

Y Cov

Y Y Cov Y

Y Cov Y

Y Cov

Y Y Cov Y

Y Cov Y

Y Cov

otrzymujemy

 _Y = ^2I,

gdzie I jest macierzą jednostkową wymiaru ^nn, tzn.

mającą na przekątnej 1, a poza przekątną 0.

II. Metoda najmniejszych kwadratów.

Niech b = ^{(b0,b1,...,bp1)'} będzie ustalonym wektorem, a y

= ^(y1,...,yn)' realizacją wektora zmiennych objaśnianych Y

= ^(Y1,...,Yn)'.

Niech Q (b) będzie kwadratem odległości wektora y od wektora Xb.

Wówczas

Q (b) = ^_^{n      }

i yi b bxi bp xi p 1

1 2 , 1 1

1

0 ... )

(

( =

= (y – Xb)’(y – Xb).

Definicja. Wartością estymatora wektora współczynników ^ wyznaczonym metodą

najmniejszych kwadratów (MNK) nazywamy wektor b minimalizujący funkcję Q(^).

Funkcja Q(^) osiąga minimum w punkcie b, w którym zerują się pochodne cząstkowe :

( b_i Q



 b) = 0, ^{i 1,...,n.} (3)

Q(^) jest funkcją kwadratową, stąd (3) jest układem równań liniowych, który w postaci macierzowej przyjmuje postać:

X’Xb = X’y. (4)

Załóżmy, że macierz X’X jest odwracalna ( kolumny są liniowo niezależne ). Wtedy rozwiązaniem równania (4) jest wektor

b = (X’X)^1X’y. (5) Zastępując w (5) y przez Y otrzymujemy estymator MNK wektora współczynników regresji wielokrotnej ^ postaci:

^b= (X’X)^1X’Y. (6)

Własności estymatora MNK

Stwierdzenie. Niech U będzie r – wymiarowym

wektorem losowym o wartości średniej ^U i macierzy kowariancji ^U oraz niech A będzie macierzą rozmiaru

r

s . Wówczas dla s – wymiarowego wektora losowego V = AU mamy

V = A^U oraz ^{ V} A^U A’.

D. ^{Vi } _j^_^r₁^aijUj dla ^{i 1,...,s.} Stąd, obliczając wartość średnią obu stron mamy

 

 r j ij U

V_i a _j

1 

 , czyli ^V = A^U. Analogicznie, otrzymujemy

  

  r j

r

k ij hk j k

h

iV a a E U U

V E

1 1 ( )

)

( oraz

  

  r j

r

k ij hk j k

h

i E V a a E U E U

V E

1 1 ( ) ( )

) ( )

( .

Zatem

Cov(^Vi,Vh) = ^{E(ViVh)E(Vi)E(Vh)} =

= _j^{ r}_{ 1 1k}^{r aijahkCov(Uj,Uk)}.

Stąd ^{ V} A^U A’. c.k.d.

Twierdzenie. Estymator ^bjest nieobciążonym estymatorem ^, tzn. ^{b } X^ oraz

 _b^ _²(X’X)^1,

D. Wiemy, że ^b= (X’X)^1X’Y, ^{Y } X^ . Podstawiając w poprzednim twierdzeniu A = (X’X)^1X’ otrzymujemy ^{b  } .

Wykorzystując wzór na macierz kowariancji wektora, własność macierzy: (AB)’ = B’A’, oraz ^{ Y} = ^2I mamy

 _b^ (X’X)^1X’(^2I) ((X’X)^1X’)’ =

(X’X)^1X’(^2I) X((X’X)^1)’ = ^2 (X’X)^1,

gdyż macierz (X’X)^1 jest symetryczna. c.k.d.

W szczególności

Var^{(bi) b2i 2}(X’X)^ii1, ^{i 0,1,...,p1.}

Np. w przypadku regresji jednokrotnej ( p=2) mamy:

Var^{(b1) b21 2} (X’X)^221 = ^_^{n }

i xi x

1

2

2 / ( )

 .

Wartość przewidywana dla i-tej obserwacji:

i  Y

1 , 1 2

2 1 1

0b x_i b x_i ...b_p x_{i p}

b   

= x’^{i b}. Wektor wartości przewidywanych:

)' ,..., ,

(Y₁ Y₂ Y_n

Y     = X^b= X(X’X)^1X’Y = HY, gdzie H = X(X’X)^1X’.

Uwaga. Macierz H jest symetryczna ( H = H’ ) oraz H²y = Hy dla każdego wektora y.

Wartości resztowe (rezydua).

e = ^{(e1,e2,...,en)'} Y - ^Y = (I - H)Y = wektor rezyduów

Stwierdzenie.

(i) E(e)^{(0,0,...,0)'}, (ii) ^{ e 2}( I – H ).

D. (i) ^E(Y) = E(X^b) = X E(^b) = X^ = ^E(Y) E(e) = E(Y -^Y) = ^E(Y) - ^E(Y) = ^(0,0,...,0)'.

(ii) ^{ e} (I - H)^2I(I – H)’ = ^2(I – 2H + H²) = = ^2(I – H), gdyż H² = H.

Niech

SSE = ^_^{n }

i ei 1

2 e’e.

Można pokazać, że

E(e’e) = ^{(n p)2}.

Stąd błąd średniokwadratowy (zdefiniowany podobnie jak dla regresji jednokrotnej)

2 

S 

 

  n

i i

p e n

p

n 1

2 1

1 e’e = _n¹ _p^SSE

jest nieobciążonym estymatorem wariancji ^2.

Liczbę ^{n p} nazywamy liczbą stopni swobody sumy kwadratów błędów = liczba niezależnych obserwacji n pomniejszona o liczbę więzów nakładanych na ^{Yi,i 1,...,n} , równą p.

Stąd, wobec ^{ b 2}(X’X)^1 oraz ^E(S2)2, otrzymujemy błędy standardowe estymatorów ^bi współczynników ^i jako pierwiastki z

(^{SEbi)2 S2}(X’X)^ii1, ^{i0,1,..., p1.}

Określimy współczynnik determinacji wielokrotnej.

Ocena „dobroci” dopasowania modelu regresji wielokrotnej.

 

  n

i yi y SST

1

)2

( = całkowita suma kwadratów

( Total Sum of Sqaures ) ( miara zmienności samych ^{y1,..., yn)}.

 

  n

i yi y SSR

1

)2

(  = regresyjna ( modelowa ) suma

kwadratów ( Regression ( Model ) Sum of Squares

( miara zmienności ^{y1,..., yn)}. Można pokazać:

 

     



 

n

i i

n

i i

n

i i

i y y y y y

y

1

2 1

2 1

2 ( ) ( )

)

(   .

^SST = ^SSE + ^SSR

R² = _SST^{SSR  1}_SST^SSE =

 

 



 n

i i

n

i i

y y

y y

1

2 1

2

) (

) ( 

= współczynnik determinacji wielokrotnej = zmienność wyjaśniona przez model/zmienność całkowita Im mniejsze ^SSE tym model bardziej adekwatny.

Współczynnik determinacji jest miarą stopnia

dopasowania modelu do obserwacji ( ocenia jakość tego dopasowania ).

Testy dla wektora współczynników ^ . (A)

0 ...

: _{1 2 1}

0    _p 

H    ,

1:

H co najmniej jeden ze współczynników ^{1,...,p1} jest różny od 0.

Niech:

SSE = ^_^{n }

i Yi Yi 1

)2

( 

, SSR = ^_^{n }

i Yi Y

1

)2

(

. Jeśli ^H0 jest prawdziwa, to

(a) ^{2SSR~ 2p1}, ^{2SSE~n2p} oraz

zmienne losowe SSR i SSE są niezależne.

(b) Statystyka

^{F }_^22_SSE^SSR_/(^/(_n^p_^1_p⁾₎ = _SSE^SSR_/(^/(_n^p_^1_p⁾₎

ma rozkład F Snedecora z ^p1 i ^{n p} stopniami swobody.

Zbiór krytyczny testu hipotezy ^H0 przeciw ^H1 na poziomie istotności ^1 ma postać:

} :

{f f F_obl f_{1 ,p 1,n p}

C    _{  }

(B) Niech ^{i {0,1,...,p1}} - ustalone.

0

0: _i 

H  , ^{H1:i 0} Wiemy, że

bi i i i

SE T b



 

 ~^tnp .

W szczególności, jeśli ^H0 jest prawdziwa, to

p n b

i i t

SE T b

i

 ~ 



 .

Stąd zbiór krytyczny ma postać:

} :

{t t t_{1 /2,n p}

C   _{ } .

Prognoza wartości ^Y na podstawie x⁰ Obserwowane ^{Y ,...,1 Yn}:

i p i p i

i

i x x x

Y ₀_{1 1}_{2 2}... _{1 , 1} , Nieobserwowane

0 1 , 0 1 02

2 01 1 0 0

0 Y(x )  x  x ..._p x _p 

Y ,

gdzie ^{1,2,...,n,0} są niezależnymi zmiennymi losowymi o rozkładach ^N(0,).

W notacji wektorowej

Y(x⁰) = x⁰’^ + ^0 gdzie x⁰ = (^{x01,x02,...,x0,p1)'},

Zadanie:

(a) ocena ( estymacja ) wartości średniej

) (x₀

Y = ^E[Y(x⁰ )] zmiennej objaśnianej w sytuacji, gdy wektorem zmiennych objaśniających jest x⁰

(b) przewidywanie ( prognoza ) wartości Y(x⁰).

(a) Estymacja ^Y(x0):

) (x₀

Y = E(x⁰’^ + ^0) = E(x⁰’^ ) + E(^0) = x⁰’^. Niech ^Y(x⁰) = x⁰’^b - estymator ^Y(x0).

x E

Y_{( )} 

 0

 (x⁰’^b) = x⁰’E(^b) = x⁰’^ = ^Y(x0).

Zatem ^Y(x⁰) jest nieobciążonym estymatorem ^Y(x0).

2 ) (x₀ Y

 = x⁰’^b x⁰ = ^2 x⁰’ (X’X)^1 x⁰ Stąd błąd standardowy estymatora ^Y(x⁰)

1 0 0'

)

( ₀ S x (X'X) x SE_Y^ _x  ^ ,

co pozwala otrzymać granice przedziału ufności dla ^Y(x0) na poziomie ufności ^1 jako realizacje zmiennych

) , (

2 / 1

0) ₀

(x t _{n p}SE_{Y x}

Y  ^

  .

(b) Prognoza Y(x⁰) = x⁰’^ + ^0 przy pomocy ^Y(x⁰).

Podobnie jak dla regresji jednokrotnej obliczamy

2

) ( ) (x₀ Y x₀ Y 

 = ^2(1 + x⁰’ (X’X)^1 x⁰)

Stąd błąd standardowy estymatora

1 0 0'

) ( )

( _{0 0} S 1 x (X'X) x SE_Y^ _{x Y x}   ^ ,

co pozwala otrzymać granice przedziału ufności Y(x⁰) dla na poziomie ufności ^1 jako realizacje zmiennych

) ( ) , (

2 / 1

0) _{0 0}

(x t _{n p}SE_{Y x Y x} Y  _{ } ^ _



Diagnostyka modelu regresji

(a) Wykres rezyduów pozwala wykryć odstępstwa od modelu, podobnie jak w przypadku regresji

jednokrotnej, takie jak: nieliniowość równania regresji, skorelowanie i niejednakowa wariancja błędów, rozkład błędów różny od normalnego.

(b) Identyfikacja obserwacji odstających – realizacji zmiennych, które nie spełniają zależności (1):

i p i p i

i

i x x x

Y ₀_{1 1}_{2 2}... _{1 , 1} .

Możliwe powody: błędny zapis danych lub zależność (1) prawdziwa tylko w pewnym zakresie zmiennych

objaśniających.

Wiemy: ^{ e 2}( I – H ).

Stąd błąd standardowy i – go rezyduum

i

e S h

SE _i  1 , gdzie ^hi = Hⁱⁱ= i – ty element diagonalny macierzy H,

Studentyzowana wartość resztowa:

ei i

i e SE

r  /

niweluje różną zmienność rozkładów rezyduów.

Wykres {(

i , r

), i  1 ,..., n }

wartości, które prawdopodobnie odpowiadają niektórym obserwacjom odstającym, za wyjątkiem tych dla

których wartość ^{Yi Yi} jest mała.

Identyfikację obserwacji odstających poprawimy rozpatrując modyfikację rezyduów:

) (i i i

i Y Y

d    ,

gdzie ^Yi(i) jest wartością przewidywaną zmiennej objaśnianej dla x = xⁱ w modelu regresji, w którym usunęliśmy obserwację ^Yi, tzn. skonstruowanym dla danych:

Jⁱ =

{(

di = rezyduum modyfikowane

di i

i d SE

t  / = studentyzowane rezyduum modyfikowane

Można pokazać, że

2 / 1

) 2

1 (

1 













 

i i i

i SSE h e

p e n

t ~ ^{tn p1}.

Duża wartość ^ti wskazuje, że obserwacja i – ta jest odstająca

(a) Testujemy n hipotez:

0i :

H Obserwacja i – ta nie jest odstająca przeciw

1i :

H Obserwacja i – ta jest odstająca.

(b) ^H0: żadna obserwacja nie jest odstająca przeciw

1:

H są obserwacje odstające

Przyjmujemy ^H1, jeśli przyjmiemy co najmniej jedną hipotezę ^H1i. Wówczas poziom istotności takiego testu ustalamy z zależności (przy założeniu, że ^H0 jest

prawdziwa):

P(_i^_ⁿ₁ {^H0i nie odrzucone}) = 1 - ^P(_i^_ⁿ₁ {^H0i odrzucone})

 

 n i

HP i 1 ({ 0

1

stąd ^P(H0 odrzucone ) ^{1n '} = ograniczenie na poziom istotności testu z (b), zatem ^{ }_n^' powinno być poziomem istotności indywidualnych testów w (a).

Rzeczywisty poziom takiego testu jest znaczne niższy niż

' ( ze względu na grube oszacowanie ), zatem test

znajduje mniej obserwacji odstających niż test dokładnie na poziomie istotności ^'.

Identyfikacja obserwacji wpływowych.

Obserwacja wpływowa, to taka, której usunięcie ze zbioru danych powoduje duża zmianę wektora

estymatorów MNK. Podejrzane są o to:

(i) obserwacje odstające

(ii) obserwacje, dla których wektor zmiennych objaśniających różni się znacznie od wektora średnich ^{x (1,x1,...,xp1)'}. Miarą odstępstwa xⁱ od ^x jest i – ty wyraz diagonalny macierzy H : ^hi, ponieważ wiadomo, że

  n

i hi p

1 oraz dla każdego i ^{1hi }_n¹, zatem można przyjąć, że typowa wartość ^hi nie przekracza znacznie wartości ^{p /n}. W praktyce przyjmujemy, że obserwacja (x^i,Yi), dla której

n

h_i 2p może być potencjalnie obserwacją wpływową.

Wówczas usuwamy ją ze zbioru danych i sprawdzamy na ile zmienił się wektor estymatorów MNK.

(iii) Odległość Cooke’a definiujemy

2 2

2 2

1

) 2 (

) 1 ( )

(

i i i

n

j j j i

i h

h pS

e pS

Y Y

D  

 

 ^





,

gdzie ^Yj(i) jest wartością przewidywaną dla j – tej obserwacji ^Yj na podstawie danych z usuniętą i – tą obserwacją.

Wartość ^Di odpowiada wpływowi, jaki na prognozę znanych wartości zmiennej objaśnianej ma usunięcie ze zbioru danych i – tej obserwacji.

Duża wartość ^Di wskazuje, że obserwacja i – ta jest wpływowa.

Współliniowość występuje, gdy niektóre zmienne są liniowo zależne, np.

i i i

i x x

Y 5 ₁2 ₂ oraz ^{xi2 2xi1}.

Wówczas – nie ma jednoznacznego modelu, można zredukować liczbę zmiennych objaśniających.

Wykrywamy współliniowość lub zależność bliską współliniowości następująco:

(i) ^r(x¹,x²⁾ jest bliski 1.

(ii) Wartość współczynnika determinacji

wielokrotnej ^Ri2 obliczonego dla hipotetycznego modelu, w którym xⁱ jest zmienną objaśnianą a pozostałe x ^j, ^{j i}, są zmiennymi objaśniającymi, jest bliska 1.

Równoważnie, wartość tzw. współczynnika podbicia (ang. – variance inflation factor ):

1 2) 1

(  ^

 _i

i R

VIF jest duża.

Wybór zmiennych objaśniających w liniowym modelu regresji

Cel – selekcja zmiennych objaśniających aby otrzymać model najprostszy.

Metody selekcji sekwencyjnej: