Nejdříve načteme data:
ohrazeni.spe <- read.csv2 ('http://www.davidzeleny.net/anadat-r/data-download/ohrazeni-spe.csv', row.names = 1, dec = '.') ohrazeni.env <- read.csv2 ('http://www.davidzeleny.net/anadat-r/data-download/ohrazeni-env.csv', dec = '.')
Z původních dat vyštípneme jenom řádky (vzorky) z třetího roku pozorování, a z druhových dat odstraníme první sloupeček (molicaer):
spe <- sqrt (ohrazeni.spe[ohrazeni.env$YEAR == 3,-1]) env <- ohrazeni.env[ohrazeni.env$YEAR == 3,] library (vegan)
Který z faktorů má největší vliv na vegetaci po odčerpání vlivu ostatních?
rda.mow.cond <- rda (spe ~ MOWING + Condition (FERTIL + REMOV), data = env) rda.fer.cond <- rda (spe ~ FERTIL + Condition (MOWING + REMOV), data = env) rda.rem.cond <- rda (spe ~ REMOV + Condition (MOWING + FERTIL), data = env) RsquareAdj (rda.mow.cond) #$r.squared #[1] 0.09939733 #$adj.r.squared #[1] 0.07165358 RsquareAdj (rda.fer.cond) #$r.squared #[1] 0.1878872 #$adj.r.squared #[1] 0.168571 RsquareAdj (rda.rem.cond) #$r.squared #[1] 0.03322556 #$adj.r.squared #[1] -0.0008202652
A test jednotlivých vlivů:
anova (rda.mow.cond) #Permutation test for rda under reduced model #Model: rda(formula = spe ~ MOWING + Condition(FERTIL + REMOV), data = env) # Df Var F N.Perm Pr(>F) #Model 1 5.599 2.9256 199 0.005 ** #Residual 20 38.273 #--- #Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 anova (rda.fer.cond) #Permutation test for rda under reduced model #Model: rda(formula = spe ~ FERTIL + Condition(MOWING + REMOV), data = env) # Df Var F N.Perm Pr(>F) #Model 1 10.583 5.5302 199 0.005 ** #Residual 20 38.273 #--- #Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 anova (rda.rem.cond) #Permutation test for rda under reduced model #Model: rda(formula = spe ~ FERTIL + Condition(MOWING + REMOV), data = env) # Df Var F N.Perm Pr(>F) #Model 1 10.583 5.5302 199 0.005 ** #Residual 20 38.273 #--- #Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Největší vliv na vegetaci má tedy hnojení (R2adj = 16.9% vysvětlené variability, p < 0.005), druhé je kosení (R2adj = 7.2%, p < 0.005). Naopak odstranění dominanty nemělo výrazný vliv - jeho adjustovaný R2 je dokonce záporný, což značí že proměnná vysvětlí míň variability než by vysvětlila proměnná náhodná, a test je neprůkazný (p = 0.51).
Použití funkce varpart
:
varpart.ohr <- varpart (spe, ~MOWING, ~ FERTIL, ~ REMOV, data = env) plot (varpart.ohr)