# 1
d.spe<-read.csv2("doubs_spe.csv")
d.env<-read.csv2("doubs_env.csv")

library(vegan)
spRich<-specnumber(d.spe)

lm.full<-lm(spRich~., data=d.env[,-1])
summary(lm.full)
lm.s<-step(lm.full)
summary(lm.s)
lm.s<-update(lm.s, .~.-bdo-har-dfs)

summary(lm.s)
# Call:
# lm(formula = spRich ~ slo + pH + nit + amm, data = d.env[, -1])
# 
# Residuals:
#     Min      1Q  Median      3Q     Max 
# -7.0782 -2.8613  0.7445  2.3079  6.9526 
# 
# Coefficients:
#              Estimate Std. Error t value Pr(>|t|)    
# (Intercept) 106.66097   40.07022   2.662  0.01339 *  
# slo          -3.40756    1.00436  -3.393  0.00231 ** 
# pH           -1.10215    0.47952  -2.298  0.03017 *  
# nit           0.05893    0.01125   5.237 2.02e-05 ***
# amm          -0.22872    0.03511  -6.514 8.00e-07 ***
# ---
# Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# 
# Residual standard error: 4.141 on 25 degrees of freedom
# Multiple R-squared:  0.7888,	Adjusted R-squared:  0.755 
# F-statistic: 23.34 on 4 and 25 DF,  p-value: 3.928e-08

#RDA

rda.full<-rda(spRich~., data=d.env[,-1])
rda.s<-ordistep(rda.full, permutations = how(nperm=999))
anova(rda.s, by="margin")
# Permutation test for rda under reduced model
# Marginal effects of terms
# Permutation: free
# Number of permutations: 999
# 
# Model: rda(formula = spRich ~ slo + pH + nit + amm, data = d.env[, -1])
#          Df Variance       F Pr(>F)    
# slo       1   6.8049 11.5109  0.004 ** 
# pH        1   3.1231  5.2829  0.037 *  
# nit       1  16.2109 27.4216  0.001 ***
# amm       1  25.0881 42.4378  0.001 ***
# Residual 25  14.7793                   


RsquareAdj(rda.s)
# $r.squared
# [1] 0.7888148
# 
# $adj.r.squared
# [1] 0.7550252

###2
spe<-read.delim2("clipboard", header = T)
env<-read.delim2("clipboard", header=T)
summary(env)

sep<-spe[,colSums(spe>0)>1]# Odstraneni druhu, ktere se vyskytly jen jednou

env$Year<-env$Year-2012
env$Block<-as.character(env$Block)
env$Plot<-as.character(env$Plot)


names(env)

rda.1<-rda(sqrt(spe)~(Mown.twice+Rhinanthus.sown)*Year
           +Condition(Block), data=env)
anova(rda.1, by="terms", permutations=how(
  blocks=env$Block, 
  plots=Plots(strata=env$Plot)))

rda.1<-update(rda.1, .~.-Mown.twice:Year)
anova(rda.1, by="terms", permutations=how(
  blocks=env$Block, 
  plots=Plots(strata=env$Plot)))
#                      Df Variance       F Pr(>F)   
# Mown.twice            1   0.5677  2.8964  0.005 **
# Rhinanthus.sown       1   1.7326  8.8396  0.005 **
# Year                  1   3.5007 17.8606  0.005 **
# Rhinanthus.sown:Year  1   0.7316  3.7324  0.005 **
# Residual             86  16.8564 

### 3
bru<-read.delim2("clipboard")
bru.6<-bru[bru$instar==6,]
names(bru.6)

library(vegan)
pca.6<-rda(bru.6[,4:11], scale=T)

rda.6<-rda(bru.6[,4:11]~species, data=bru.6, scale=T)
anova(rda.6)
# Permutation test for rda under reduced model
# Permutation: free
# Number of permutations: 999
# 
# Model: rda(formula = bru.6[, 4:11] ~ species, data = bru.6, scale = T)
#           Df Variance      F Pr(>F)    
# Model      5    6.352 87.881  0.001 ***
# Residual 114    1.648                  
# ---
# Signif. codes:  
# 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

library(fastDummies)
sp.d<-dummy_cols(bru.6$species)
cca.full<-cca(sp.d[,-1]~., data=bru.6[,4:11] )
cca.s<-ordistep(cca.full, permutations = how(nperm=999))

library(MASS)

lda.1<-lda(bru.6[,4:11], bru.6$species)

aa<-ordiplot(lda.1, display=c("si", "sp"), type="n")
ordispider(aa, groups=bru.6$species)
text(aa, what="sp", arrows = T, length=0.03, col=3)
text(aa, what="sp")