12. Regrese Jedním z hlavních úkolů matematické statistiky je hledání a studium závislostí mezi dvěma či více proměnnými. Závisle proměnná se zpravidla označuje Y a nezávisle proměnné X1 , … , Xi ,i ¥ 1. Závislosti mezi Y a skupinou proměnných X mohou být funkční nebo statistické. V případě přímé funkční závislosti je náhodná část ve vyjádření této závislosti nulová, proto se jí nebudeme samozřejmě zabývat. V běžných situacích však nastává situace , kdy Y má povahu náhodné veličiny . Nezávisle proměnné X1 , … , Xi mohou být nenáhodnými ( pevnými ) veličinami nebo náhodnými veličinami. Statistickou závislostí mezi náhodnou veličinou Y a veličinami X1 , … , Xi rozumíme předpis, který každé uspořádané i – tici x1 , … , xi přiřazuje právě jedno podmíněné rozdělení náhodné veličiny Y . Náhodná veličina Y je tedy statisticky závislá na veličinách X1 , … , Xi , jestliže změny v hodnotách x1 , … , xi mají za následek změnu podmíněného rozdělení náhodné veličiny Y.
12.1 Teoretické základy V praxi není většinou zákon rozdělení ( typ náhodné veličiny Y ) znám. Máme většinou k dispozici pouze data ve formě uspořádaných i+1 – tic , kde prvních i složek odpovídá hodnotám x a poslední složka odpovídá hodnotě y. Na základě experimentálních údajů se provede výběr funkce , která má co nejlépe vystihovat rozložení skutečných údajů , tuto funkci nazýváme teoretickou regresní funkcí a je uváděna ve tvaru y = f ( x1 ,… , xi ; γ 0 ,… , γ p ) (12.1),
kde hodnoty γ 0 ,… , γ p jsou nazývány teoretickými regresními parametry ( koeficienty ) . Tyto parametry mají povahu konstant, které nejsou známy . V modelovaných vztazích je budeme nahrazovat jejich bodovými odhady a0 , … , ap tak , aby funkce Y = f ( x1 ,… , xi ; a0 ,… , a p ) (12.2) co nejlépe aproximovala naměřené hodnoty Y .Funkce (12.2) se nazývá empirickou regresní funkcí a veličiny a0 , … , ap empirickými regresními parametry. Empirické regresní parametry jsou náhodné veličiny. Jednou z metod jak je určit je použít metodu nejmenších čtverců. V dalším textu bude dále závisle proměnná Y ( vysvětlovaná proměnná ) náhodnou veličinou , naopak nezávisle proměnné ( vysvětlující proměnné ) X1 , … , Xi jsou nenáhodnými veličinami. Regresní závislostí bude pro naše potřeby jednostranná závislost dané náhodné veličiny ( např. Y ) na jiné veličině – ne nutně náhodné veličině ( např. X1 , … , Xi ). Pro další účely třídění a studování regresních závislostí se rozlišují následující situace: a. Náhodná veličina závisí na jedné další veličině – párová regrese. Například : y = γ 0 .x +
γ1 x +γ3
b. Náhodná veličina –závisí na dvou a více veličinách – mnohonásobná regrese γ .x .Například : y = γ 1.x1 + 2 2 + γ 0 γ 3 + x3 Z pohledu tvaru funkce (12.1) můžeme dále dělit jednotlivé případy regresních závislostí na: 1. Lineární regresi – uvažovaná funkce f je lineární vzhledem k parametrům γ 0 ,… , γ p . Například : y = γ 1.x1 + γ 2 .x2 + γ 0
2. Nelineární regresi – uvažovaná funkce f je nelineární funkcí vzhledem γ .x + γ .x k parametrům γ 0 ,… , γ p . Například : y = 1 1 2 2 . γ 3 .x3 + γ 4 .x4 Regresní analýzou nazveme tedy tu část matematické statistiky , která se zabývá studiem a konstrukcí regrese. Nejdříve určíme podmínky regresního modelu. Definice 12.1 Určujeme regresi Y na i nenáhodných veličinách x1 , … , xi . Nechť je dále dáno n ( n > p+1 ) uspořádaných i+1 – tic ( x1k ,… xik , yk ) , kde k = 1 , … , n . Tyto údaje byly
postupně získány empiricky. Na základě těchto údajů byla zvolena teoretická regresní funkce s neznámými parametry γ 0 ,… , γ p . Dále předpokládejme , že pro náhodné veličiny Yk , k=1,…,n platí:
1. Yk = f ( x1 ,… , xi ; γ 0 ,… , γ p ) + ek
(12.2)
Veličiny ek představují náhodné chyby měření a jsou zřejmě také náhodnými veličinami. Tyto náhodné veličiny musí splňovat následující podmínky: 2. E( ek ) = 0 (12.3) 2 3. Var( ek ) = s (12.4) 4. cov( ek , el ) = 0 (12.5) Dále platí 5. E(Yk) = f ( x1k ,… , xik ; γ 0 ,… , γ p ) 6. VAR(Yk) = s2 . K řešení výše popsaného modelu používáme nejčastěji metodu nejmenších čtverců. Součástí řešení je nalezení odhadů parametrů γ 0 ,… , γ p . Tyto odhady budeme dále označovat symboly c1 , … , cp . Metoda spočívá v nalezení minima čtverců hodnot ek Obecněji budeme studovat tuto metodu aplikovanou na případ lineární regrese. Hledáme tedy minimum funkce n
n
k =1
k =1
(
)
S r (γ 0 ,… , γ p ) = ∑ ek2 = ∑ Yk − f ( x1 ,… , xi ; γ 0 ,… , γ p ) , 2
(12.6)
požadujeme tedy , aby součet čtverců odchylek naměřených hodnot Yk od teoretických hodnot f ( x1 ,… , xi ; γ 0 ,… , γ p ) byl minimální. Najdeme tedy nejdříve podezřelé body ( stacionární body S ) : ∂S (γ 0 ,… , γ p ) = 0, k = 0,… , p (12.7) ∂γ k Poté ověříme pomocí klasických metod, zda v nalezených bodech se nachází minimum. Soustava rovnic (12.6) se nazývá soustavou normálních rovnic . Řešení této soustavy je jednoznačné , jestliže alespoň p + 1 z uspořádaných k – tic [ x1k ,..., xik ] , k=1,…,n je navzájem různých. Řešením soustavy (12.6) budou potom hodnoty a1 , … , ap , které budou bodovým odhadem parametrů γ 0 ,… , γ p . Odtud tedy získáme bodové odhady teoretické regresní funkce (12.1) ve tvaru Yˆk = f ( x1k ,...xik ; a0 ,...a p ) , Rozdíly eˆk = Yk − Yˆk se nazývají rezidui a veličina
k = 1,..., n (12.8)
n
n
k =1
k =1
S r = ∑ eˆk2 = ∑ (Yk − Yˆk ) 2
(12.9)
se nazývá reziduálním součtem čtverců. Studiem řešení soustavy normálních rovnic a reziduálního součtu čtverců se budeme zabývat dále v následujících částech této kapitoly a v kapitole následující.
12.2 Model lineární regrese V této části se budeme zabývat nejjednodušším případem lineární regrese – případem , kdy odhadujeme jen dva parametry a a b. Takovýto speciální vztah označujeme jako regresní přímku a je vyjádřen Y = a + b.x +e (12.10) K odvození bodových odhadů skutečných parametrů a a b použijeme klasickou metodu nejmenších čtverců. Bodové odhady označíme pro naše potřeby a a b .Tedy podle (12.6) je n
n
k =1
k =1
S ( a, b ) = ∑ ek2 = ∑ (Yk − a − b.xk )
2
(12.11)
Tento výraz budeme derivovat podle proměnných a , b . Výsledná soustava normálních rovnic je potom následující n ∂S ( a, b ) = −2.∑ (Yk − a − b.xk ) = 0 (12.12) ∂a k =1 n ∂S ( a, b ) (12.13) = −2.∑ xk . (Yk − a − b.xk ) = 0 ∂b k =1 Přepíšeme – li tuto soustavu klasickým způsobem máme: n
n
n.a + b.∑ xk = ∑ Yk , k =1
k =1
n
n
k =1
k =1
n
a.∑ xk + b.∑ x = ∑ xk Yk 2 k
,
(12.14)
k =1
Z této soustavy vyplývá , že řešení existuje , je – li determinant soustavy různý od 2
n nuly. Tedy n.∑ x − ∑ xk ≠ 0 , ale tento vztah nastane, jestliže hodnoty xk budou k =1 k =1 navzájem různé ( stačí dokonce , aby aspoň dva byly navzájem různé ). Soustava tedy bude mít právě jedno řešení. Nejprve upravíme první rovnici v (12.14) tím , že ji vydělíme hodnotou n , získáme tak a + b.x = Y (12.15) Budeme – li znát jednu z hodnot a , b pomocí tohoto vztahu dopočteme druhou . Určíme hodnotu b , dosazením (12.15) do druhého vztahu v (12.14). n
2 k
n
n
(Y − b.x ) .n.x + b.∑ x = ∑ x Y k =1
2 k
k =1
k k
⇒
n n b. ∑ xk2 − n.x 2 = ∑ xk Yk − n.x .Y ⇒ k =1 k =1 n
b=
∑x Y k =1 n
k k
− n.x .Y
2 2 ∑ xk − n.x k =1
(12.16)
Ze vztahů (12.15) a (12.16) určíme samozřejmě hodnoty a , b. Příklad 12.2 V 15 za sebou následujících dnech byla zaznamenána hodnota relativní vlhkosti v % a hodnota orosení v mm3. Naměřené údaje jsou uvedeny v tabulce dále Relativní vlhkost Hodnota orosení
46 53 29 61 36 39 47 49 52 38 55 32 57 54 44 12 15
7 17 10 11 11 12 14
9 16
8 18 14 12
Provedeme – li výpočet získáme následující hodnoty a = -2,51 a hodnota b = 0,32. Výsledek můžeme také zobrazit graficky 18
16
14
Hodnota orosení
12
teoretická hodnota
10
8
6 25
30
35
40
45
50
55
60
65
Z výše uvedeného grafu je patrno , že žádná ze skutečných naměřených hodnot Y na přímce neleží , přesto však je vyjádření vztahu mezi x a Y pomocí lineárního vztahu dobré. Podle předpokladů lineárního modelu předpokládáme , že střední hodnota chyb ek je rovna nule, jejich rozptyl je konstantní a roven s2 . Dále předpokládáme, že hodnoty Yk ( tzv. vysvětlované proměnné ) jsou typu normální rozdělení. Odtud tedy vyplývá , že E (Yk ) = α + β .xk a dále VAR (Yk ) = σ 2 . Celkově tedy Yk ∼ N (α + β .xk , σ 2 ) . Budeme dále vyšetřovat vlastnosti bodových odhadů parametrů a a b. Tedy n n n ∑ xk Yk − n.x .Y ∑ ( xk − x ) E (Yk ) ∑ ( xk − x ) . (α + β .xk ) = k =1 E (b) = E k =1n = k =1 n = n 2 2 2 2 2 2 ∑ xk − n.x ∑ xk − n.x ∑ xk − n.x k =1 k =1 k =1
n 2 2 ∑ xk − n.x k =1 = β . k =1 =β. = k =1 n n 2 2 2 2 ∑ xk − n.x ∑ xk − n.x k =1 k =1 Protože je E(b) = b , je podle definice 11.12 z kapitoly 11. tento odhad nevychýlený. Zjistíme dále hodnotu rozptylu tohoto odhadu n n n 2 2 x Y n . x . Y ( x x ) . VAR Y σ . ( xk − x ) 2 − − ( k) ∑ k k ∑ k ∑ k =1 = k =1 VAR(b) = VAR k =1n = = 2 2 n n 2 2 2 2 2 2 ∑ xk − n.x ∑ xk − n.x ∑ xk − n.x 1 = k k =1 k =1 n
n
α .∑ ( xk − x ) + β . ∑ xk . ( xk − x )
σ2
=
(12.17)
n 2 2 ∑ xk − n.x k =1 Při výpočtu jsme využili následující rovnost n
∑ (x k =1
k
n
− x ) 2 = ∑ xk2 − n.x 2
(12.18)
k =1
Podobně nyní zjistíme střední hodnotu a rozptyl pro odhad a. n ∑ Yk n α + β .x k E (a ) = E k =1 − b.x = ∑ − x .E (b) = α + β .x − x .β = α , n n k =1 tedy i bodový odhad parametru a je nevychýlený. Rozptyl hodnoty a je roven n
∑VAR(Yk )
n
+ ( x 2 ) .VAR ( b ) =
σ2
σ 2 .x 2
∑x
2 k
k =1 = σ 2. n n 2 2 n. ∑ xk2 − n.x 2 ∑ xk − n.x k =1 k =1 Závěrem můžeme tedy konstatovat , že náhodné veličiny a a b jsou typu normální rozdělení ( jsou lineární kombinací nezávislých náhodných veličin typu normálního rozdělení ) a můžeme navíc určit jejich parametry tedy n 2 σ . xk2 ∑ 2 σ k =1 a dále b ∼ N β , . a ∼ N α , n n 2 2 2 2 n. ∑ xk − n.x ∑ xk − n.x k =1 k =1 Tyto dva vztahy budou klíčové pro další určení intervalových odhadů skutečných parametrů a a b.
VAR (a ) =
k =1
n2
n
+
Věta 12.3
Nechť Sr je reziduální součet čtverců definovaný vztahem (12.6) . Potom Důkaz: Proveden například v [1] .
Sr
σ
2
∼ χ n2−2 .
Poznámka 12.4 Výše uvedený důkaz budeme motivovat následující úvahou: Náhodné veličiny Yk nezávislé , jsou dále normálního typu. Pokud je znormujeme máme Yk − E (Yk ) Yk − α − β .xk = (12.19) σ VAR (Yk )
jsou
Proto tedy 2
Yk − α − β .xk 2 ∼ χn ∑ σ k =1 Protože ale nahradíme skutečné hodnoty a a b jejich nevychýlenými odhady a a b , sníží se počet stupňů volnosti o 2 – tedy o počet parametrů, které jsme odhadovali. n
Pro další práci označme n
n
k =1
k =1
S xY = ∑ ( xk − x ) . (Yk − Y ) = ∑ xk Yk − n.x .Y n
(12.20)
n
S xx = ∑ ( xk − x ) = ∑ xk2 − n.x 2
(12.21)
SYY = ∑ (Yk − Y ) = ∑ Yk2 − n.Y 2
(12.22)
2
k =1 n
2
k =1
k =1 n
k =1
Jestliže použijeme tuto notaci můžeme koeficienty a a b zapsat také takto S b = xY ; a = Y − b.x . S xx Dále odvodíme následující výpočetní identitu :
S xx .SYY − S xY 2 Sr = S xx
(12.23)
n
S r = ∑ (Yk − a − b.xk ) k =1
2
2
S S =∑ Yk − Y + xY .x − xY .xk = S xx S xx k =1 n
S xx .(Yk − Y ) − S xY .( xk − x ) S xx 2 .SYY − 2.S xx .S xY .S xY + S xY 2 .S xx = ∑ = = 2 S xx S k =1 xx S .S − S xY 2 = xx YY S xx 2
n
Ze všech výše uvedených rovností vyplývají následující tvrzení:
1.
Sr
σ
2
∼ χ n2−2
(12.24)
(b − β ) σ
2.
S xx Sr 2 σ .( n − 2 )
=
( n − 2 ).S xx . Sr
( b − β ) ∼ tn − 2
(12.25)
3.
n.( n − 2 ) .S xx n
∑ x .S k =1
2 k
.( a − α ) ∼ tn − 2
(12.26)
r
Toto jsou základní vztahy , pomocí nichž je možno konstruovat jak intervalové odhady neznámých parametrů a a b , ale i neznámé hodnoty s. Dále jsou na nich založeny klasické statistické hypotézy o těchto parametrech. Často jsme postaveni před otázku nalezení odhadu ( prognózy ) hodnot v nějakém bodě x0 , na základě znalosti vstupních dat {xk , Yk }. Podle způsobu naší konstrukce je zřejmé, že přirozeným bodovým odhadem je v tomto případě hodnota a + b.x0 ( je opět nevychýleným odhadem skutečné neznámé hodnoty – dokažte! ) . Podobně jako v předchozích krocích nás bude především zajímat typ náhodné veličiny . n
a + b.x0 =
∑ Yk k =1
n
n
−
∑ Y .( x k =1
k
k
− x ) . ( x − x0 ) S xx
n 1 ( x − x ) . ( x − x0 ) = ∑ Yk . − k (12.27) S xx k =1 n
Z předchozího vztahu (12.25) vyplývá ( proč ? ), že náhodná veličina a + b.x0 je typu normální rozdělení . Je tedy nutno zjisti jeho střední hodnotu a rozptyl. A. E ( a + b.x0 ) = E (a ) + x0 .E (b) = α + β .x0 B. n
VAR(a + b.x0 ) = ∑ k =1
2 2 2 n 1 ( xk − x ) . ( x − x0 ) xk − x ) . ( x − x0 ) ( xk − x ) . ( x − x0 ) ( 1 2 + = − .VAR(Yk ) = σ .∑ 2 − 2. 2 S xx n S S . k =1 n xx xx n
n 2 ( xk − x ) ∑ 1 2. ( x − x ) n ( x − x ) 1 ( x − x0 )2 2 k =1 k 0 2 2 = σ . + = σ . − + ( x − x0 ) . .∑ . n S xx n S xx 2 S n k =1 xx Z těchto rovností tedy vyplývá , že
1 ( x0 − x )2 2 a + b.x0 ∼ N α + β .x0 , σ . + n S xx
(12.28)
Pro intervalový odhad skutečných hodnot využijeme vztahů (12.24) a (12.28) . Z těchto vztahů vyplývá a + b.x0 − α − β .x0 1 ( x0 − x ) + . n S xx 2
Sr ( n − 2)
∼ tn − 2
(12.29)
Pomocí tohoto vztahu je možno odhadovat chování a + b.x v daném bodě x0 ( odhadujeme chování parametrů !). Na druhou stranu je někdy důležité vyšetřit chování hodnoty Y v konkrétním bodě x0 . V tomto případě stanovujeme budoucí hodnotu náhodné veličiny Y. Opět platí nám známé vztahy v bodě x0 Y ∼ N (α + β .x0 , σ 2 )
(12.30)
a vztah (12.29). Z nich můžeme stanovit rozdělení následujícího rozdílu 1 ( x0 − x )2 2 Y − a − b.x0 ∼ N 0, σ . 1 + + nebo ekvivalentně n S xx
Y − a − b.x0 n + 1 ( x0 − x ) σ . + n S xx
2
∼ N ( 0,1)
(12.31)
Využijeme – li dále vztahu (12.24) získáváme
Y − a − b.x0 n + 1 ( x0 − x ) + . n S xx 2
Sr ( n − 2)
∼ tn − 2
(12.32)
Tento vztah je rozhodující pro stanovení odhadu budoucí hodnoty Y v bodě x0 .
