### Cvika 2
doubs.spe<-read.csv2("doubs_spe.csv", row.names = 1)
doubs.env<-read.csv2("doubs_env.csv",row.names = 1)

summary(doubs.spe)


library(vegan)


sh<-diversity(doubs.spe)# Shannon index
sp.rich<-specnumber(doubs.spe)# Pocet druhu
eve<-sh/log(sp.rich) # Vyrovnanost (Pielou eveness)

df.div<-data.frame(sh, sp.rich, eve) #Vytvoreni data frame div, ktery obsahuje druhovou bohatost, Shannon index a Pielou eveness
summary(df.div)

install.packages("psych") #Instalace
library(psych) # a aktivace balicku psych

df.div<-df.div[-c(1,8),]# Odstraneni 1. a 8. radku - vzorky s jednim a zadnym druhem
pairs.panels(df.div) # Parovy graf ukazujici vztahy mezi promennymi a histogramy jednotlivych promennych

### Regrese

lm.0<-lm(df.div$sp.rich~+1, data=doubs.env[-c(1,8),])# Nulovy (intercept-only) model 

names(doubs.env)# Vygenerovani jmen promennych z df - pro manualni edit a nasledne pouziti v add1

pdf("pairs.pdf", 10, 10)
pairs.panels(doubs.env[-c(1,8),])
dev.off() # parovy graf pro prediktory - rovnou ulozeny do pdf


corrs<-add1(lm.0, .~.+dfs+slo+flo+pH+har+phos+nit+amm+oxy+bdo, 
     test="F") # Test jednoduchych korelaci (simple effects) mezi jednotlivymi prediktory a poctem druhu - ulozeni do objekdu corrs
corrs
#        Df Sum of Sq     RSS    AIC F value    Pr(>F)    
# <none>              1720.43 117.31                      
# dfs     1    773.36  947.06 102.59 21.2314 9.467e-05 ***
# slo     1    609.64 1110.79 107.06 14.2698 0.0008329 ***
# flo     1    702.05 1018.38 104.62 17.9238 0.0002536 ***
# pH      1      6.98 1713.45 119.19  0.1060 0.7474029    
# har     1    246.68 1473.75 114.97  4.3520 0.0469192 *  
# phos    1     11.48 1708.94 119.12  0.1747 0.6793751    
# nit     1    263.85 1456.57 114.65  4.7098 0.0393114 *  
# amm     1     20.63 1699.80 118.97  0.3155 0.5791072    
# oxy     1     23.72 1696.71 118.92  0.3635 0.5517733    
# bdo     1     38.00 1682.43 118.68  0.5872 0.4504065  

p.adjust(corrs$`Pr(>F)`, method = "holm")# Uprava p-hodnot Holmovou korekci

lm.sat<-update(lm.0, .~.+alt+slo+flo+pH+har+phos+nit+amm+oxy+bdo)# saturovany model vytvoreny pomoci update nuloveho modelu

summary(lm.sat)# Spodni cast vystupu - ukazuje test signifikance pro saturovany model

# Residual standard error: 4.387 on 17 degrees of freedom
# Multiple R-squared:  0.8098,	Adjusted R-squared:  0.698 
# F-statistic:  7.24 on 10 and 17 DF,  p-value: 0.0002118

### Vystavba modelu pomoci forward selection
lm.1<-update(lm.0, .~.+dfs)# pridani prefiktoru dfs do modelu
add1(lm.1, .~.+slo+flo+pH+har+phos+nit+amm+oxy+bdo, 
     test="F")# test dalsich prediktoru
lm.2<-update(lm.1, .~.+bdo) #Pridani prediktoru bdo do modelu
add1(lm.2, .~.+slo+flo+pH+har+phos+nit+amm+oxy+bdo, 
     test="F")# test dalsich prediktoru - uz nebylo nic prukazne

plot(lm.2)# Regresni diagnostika modelu lm.2


summary(lm.2)# Summary lm.2
# Residuals:
#     Min      1Q  Median      3Q     Max 
# -7.3250 -3.5428 -0.6032  3.5232  8.1681 
# 
# Coefficients:
#               Estimate Std. Error t value Pr(>|t|)    
# (Intercept)  8.4186686  1.7325421   4.859 5.37e-05 ***
# dfs          0.0052006  0.0007305   7.119 1.84e-07 ***
# bdo         -0.1067380  0.0254420  -4.195    3e-04 ***
# ---
# Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# 
# Residual standard error: 4.715 on 25 degrees of freedom
# Multiple R-squared:  0.677,	Adjusted R-squared:  0.6511 
# F-statistic: 26.19 on 2 and 25 DF,  p-value: 7.342e-07