Lineární modely (stručné shrnutí)
(Obecné) lineární modely
• Jedna odpověď, jeden a více prediktorů
- Prediktory kvantitativní i kategoriální
- ANOVAs prediktorem o n úrovních je analogická mnohonásobné regresi o n-1 prediktorech
• Rovnice: y = a + bx1 + cx2 + ...+£
- y: odpověď
- a: intercept
- b,c: regresní koeficienty prediktorů
- z: reziduály (residua): předpokládá se, že všechny pocházejí ze stejného normálního rozložení N(0, a)
Modely s více prediktory
• Dvoucestná (Vícecestná) ANOVA
- odpověď ~ faktor. 1 + faktor.2 + ...
• Mnohonásobná regrese
- odpověď - lin.prediktor.1 + lin.prediktor.2 + ...
• Lineárni model (Analýza kovariance): kategoriální i lineárni prediktory
• Aditivní efekty prediktorů vs. interakce - odpověď na faktor. 1 (prediktor.1) závisí na hodnotě faktoru.2 (prediktoru.2)
- aditivitu lze statisticky testovat a případně zamítnout ve prospěch interakce
Interakce mezi prediktory
•Průkazná interakce znamená vzájemné ovlivňování vlivu prediktorů na odpověď
-y = a + bx1 + cx2 + dx1x2 + £ •Test interakce HO: d = 0 -d > 0: pozitivní int., vyšší hodnoty než aditivita -d < 0 : negativní int., nižší hodnoty než aditivita
•dfint = dfx1*dfx2 •Jak napsat
-x (Alt + 0215) formálne v odborných textech -: v R pouze interakční člen -*vR znamená additivitu dohromady s interakcí -rovnici nahoře v R napíšeme takto
y ~ x1*x2
•Interaction plot Znázornění interakce
-Pokud int. není přítomna, budou čáry paralelní Interakce neznamená korelaci prediktorů!
—i
LU LU
m
E o
to Cl
Interaction plot
o o
cn
O
O O
O Li")
O -
Control
x.factor
data
biomass
Reakce růstu rostlin na zalévání a
hnojení
Experimentální design: 24 květináčů s rostlinami rozdělenými náhodně do 4 skupin s faktoriálním uspořádáním hnojení a velikosti zálivky
Otázka: Jaký je vliv hnojení a zalévání na produkci nadzemní biomasy rostlin?
aov.l<-aov(biomass-watering*fertil,  data=plants) summary(aov.1)
watering	fertil
65 Control	Control
58 Control	Control
74 Control	Control
65 Control	Control
81 Control	Control
78 Control	Control
92 Watered	Control
86 Watered	Control
94 Watered	Control
100 Watered	Control
89 Watered	Control
90 Watered	Control
110 Control	Fertilized
118 Control	Fertilized
128 Control	Fertilized
121 Control	Fertilized
119 Control	Fertilized
116 Control	Fertilized
185 Watered	Fertilized
196 Watered	Fertilized
201 Watered	Fertilized
195 Watered	Fertilized
193 Watered	Fertilized
191 Watered	Fertilized
	Df	Sum Sq	Mean Sq	F value	Pr(>F)	
watering	1	13968	13968	336.38	5.604e-	]_ 4 * * *
fertil	1	33825	33825	814.57	< 2.2e-	]_ g -k-k-k
watering:fertil	1	4240	4240	102.11	2.654e-	Q g    -k -k -k
Residuals	20	831	42			
Interaction plot
cl
Control
Fertilized
x.f a cto r
data
biomass
Reakce růstu rostlin na zalévání a
hnojení
Experimentální design: 24 květináčů s rostlinami rozdělenými náhodně do 4 skupin s faktoriálním uspořídáním hnojení a velikosti zálivky
Otázka: Jaký je vliv hnojení a zalévání na produkci nadzemní biomasy rostlin?
aov.l<-aov(biomass-watering*fertil,  data=plants) summary(aov.1)
Df Sum Sq Mean Sq F value Pr(>F)
336.38 5.604e-14 *** 814.57 < 2.2e-16 *** 102.11 2.654e-09 ***
Interaction plot
Závěr: Hnojení i zalévání mají průkazný pozitivní vliv růst rostlin. Mezi prediktory je navíc pozitivní interakce.
watering	fertil
65 Control	Control
58 Control	Control
74 Control	Control
65 Control	Control
81 Control	Control
78 Control	Control
92 Watered	Control
86 Watered	Control
94 Watered	Control
100 Watered	Control
89 Watered	Control
90 Watered	Control
110 Control	Fertilized
118 Control	Fertilized
128 Control	Fertilized
121 Control	Fertilized
119 Control	Fertilized
116 Control	Fertilized
185 Watered	Fertilized
196 Watered	Fertilized
201 Watered	Fertilized
195 Watered	Fertilized
193 Watered	Fertilized
191 Watered	Fertilized
watering	1	13968	13968
fertil	1	33825	33825
watering:fertil	1	4240	4240
Residuals	20	831	42
E ^
Q_
Control
Fertilized
x.f actor
Testování vlivu jednotlivých prediktorů / výběr prediktorů do
modelu
• Cíle
- ukázat, které prediktory mají průkazný vliv
• lze testovat jednoduché (marginální) efekty, tj. párové korelace mezi jednotlivými prediktory a odpovědí
- sestavit minimální adekvátní model, který bude zahrnovat pouze průkazné prediktory
• kondicionální (parciální) efekty prediktorů
• Statistické testy tohle moc neumí
- umí otestovat model, porovnat kvalitu více modelů
• test signifikance
• AIC
Akaike information Criterion
• Kvantifikuje množství informaci v odpovědi, kterou vysvětluje model
- umožňuje porovnání dvou modelů, které se liší počtem stupňu volnosti (modelu)
- nižší AIC značí lepší fit než vyšší AIC (absolutní hodnoty nejsou důležité)
• AIC = 2k - 2log(L), (log = přirozený log), k je počet parametrů (tedy df model v lin. modelu)
- v lin. mod. AIC = 2k - 2 log (n/RSS) + C, kde RSS je reziduálni suma čtverců, C je konstanta (lze ignorovat)
• Různé názory na možnost kombinovat s F-testem průkaznosti
Výběr prediktorů do modelu
• Postupný
- Forward selection: k nulovému modelu přidávám postupně prediktory
• vhodnější pro observační data
- Backward selection: ze saturovaného modelu odebírám nevýznamné prediktory
• vhodnější pro experimentální data
- Oba směry: zvažují v každém kroku přidání I odebrání prediktorů
Problém s mnohonásobným paralelním testováním
• Provedu-li více testů s pravděpodobností chyby I. druhu 0.05, pravděpodobnost, že udělám chybu aspoň jednou velmi vzrůstá.
- pro dva testy p = 0.05 + 0.05 - 0.05 x 0.05 = 0.0975
• Řešení: různé korekce (Holm, Bonferroni, false detection rate -FDR) upravují p-hodnoty nahom, čímž kontrolují/redukují riziko chyby
• Alternativa: "protected multiple testing"
- spočtu test saturovaného modelu se všemi prediktory. Je-li průkazný, další korekci kvůli multiple testing už neřeším
- Je-li ovšem neprůkazný, tak s testováním končím se závěrem, že regrese není průkazná.
• Skutečný problém hlavně v analýzách "velkých" dat z databází.
Lineární model pro experimentální
data
Sample
seedlings
1
2
3
4
5
6
7
8
9
10
11
treatment
7 control 5 control
12 control
8 control
13 control 12 control
7 control
5 control 7 control
6 control 19 grazing
productivity temperature
714	7.2
518	4.5
379	7.4
686	5.4
703	5
775	6
651	6.2
630	7.6
470	8.4
557	4.7
394	6.8
to be continued
Q: Jaký je vliv obhospodařování louky (kosení, pastva) na počet semenáčků, které se objeví na jaře? Ovlivňují semenáčky ještě další faktory prostředí?
Design: 30 experimentálních ploch (10 na každý typ obhospodařování) náhodně rozmístěné v krajině (různé lokality). Control = opuštěná nekosená louka. Zaznamenána byla i produktivita a průměrná teplota lokality.
lm.full<-lm(seedlings-treatment*productivity*temperature, data=seedl] anova(lm.full)
# Analysis of Variance Table
#
Response: seedlings
#
#
#
#
#
#
#
#
#
#
#
#
Df
treatment 2
productivity 1
temperature 1
treatment:productivity 2
treatment:temperature 2
productivity:temperature 1
treatment:productivity:temperature 2
Residuals 18
Signif. codes
0
0 . 001
Sum Sq Mean Sq F value
618.07 309.033 31.5710
290.61 290.605 29.6883 1.88 1.882 0.1922 5.02 2.509 0.2563 6.87 3.434 0.3508 0.25 0.250 0.0256 6.48 3.242 0.3312 176.19 9.789
0.01 0.05 0.1 *
Pr (>F) , 301e-06 , 553e-05
~k ~k ~k ~k ~k ~k
0 0 0 0 0
6663 7767 7088 8747 7223
Tuto tabulku lze použít jako výstup testu jednotlivých prediktorů a jejich interakcí
Lineární model pro experimentální
data
Další krok: odstraňte neprůkazné prediktory z modelu (pomocí backward selekce; ponechejte I neprůkazné main efekty pokud je průkazná interakce)
lm.final<-lm(seedlings-treatment+productivity,  data=seedl) summary(lm.final) Call:
lm(formula = seedlings ~ treatment + productivity,  data = seedl)
Residuals
Min -4.698 -1
1Q Median !40 -0.315
3Q 1. 975
Max 4 .741
Coefficients:
(Intercept) treatmentgrazing treatmentmowing productivity
Signif.  codes: 0
Estimate Std.  Error t value Pr(>|t|)
22.315470 2.437909 9.154 1.29e-09 ***
9.581024 1.251180 7.658 3.98e-08 ***
2.173362 1.268726 1.713 0.0986 .
-0.024764 0.003996 -6.198 1.48e-06 ***
y -k-k-k r
0.001
y -k-k r
0.01
y -k r
0 . 05
0.1
Residual standard error:  2.75 on 26 degrees of freedom Multiple R-squared:  0.8221,     Adjusted R-squared: 0.8015 F-statistic:  40.04 on 3 and 26 DF,    p-value: 6.857e-10
Koeficienty prediktorů (nebo kontrasty pro faktory)
Celkový test modelu
Lineární model pro experimentální
data
Další krok: odstraňte neprůkazné prediktory z modelu (pomocí backward selekce; ponechejte I neprůkazné main efekty pokud je průkazná interakce)
lm.final<-lm(seedlings-treatment+productivity,  data=seedl) summary(lm.final) Call:
lm(formula
seedlings ~ treatment + productivity,  data = seedl
Residuals
Min -4.698 -1
1Q Median :40 -0.315
3Q 1. 975
Max 4 .741
Coefficients:
(Intercept) treatmentgrazing treatmentmowing productivity
Estimate Std.  Error t value Pr(>|t|)
22.315470 2.437909 9.154 1.29e-09 ***
9.581024 1.251180 7.658 3.98e-08 ***
2.173362 1.268726 1.713 0.0986 .
-0.024764 0.003996 -6.198 1.48e-06 ***
Koeficienty prediktorů (nebo kontrasty pro faktory)
Signif.  codes:    0 0.001 0.01 0.05
0.1
Residual standard error:  2.75 on 26 degrees of freedom Multiple R-sguared:  0.8221,     Adjusted R-sguared: 0.8015 F-statistic:  40.04 on 3 and 26 DF,    p-value: 6.857e-10
Celkový test modelu
Závěr: Průkazný vliv na počet semenáčků má způsob obhospodařování a produktivita mají. Jejich efekty jsou aditivní. Produktivita snižuje počet semenáčků. Obhospodařivání jejich počet zvyšuje. Průkazné zvýšení je ale způsobené pouze pastvou.
Lineární model pro observační data
• Typicky mnoho potenciálních prediktorů
- mnohdy nemožnost fitovat saturovaný model (potenciálních prediktorů více než pozorování)
• Forward selection s korekcí - akceptovatelná varianta
• Korelované potenciální prediktory: těžko řešitelný problém
- např. v ČR: nadmořská výška, průměrná teplota a úhrn srážek
- Možnost prezentovat jednoduché korelace spolu s finální lineárním modelem
Výběr prediktorů do lineárního
modelu
• Spousta možností jak to udělat v zásadě dobře
- postupný výběr s korekcí (ať už jakoukoliv)
- využití "test protection"
• Jednoznačně špatně - Statistical fishing
- Forward selekce bez korekce pro observační data
• speciálně při použití mnoha potenciálních prediktorů
- Úpravy testů "aby to vyšlo průkazně"