library(vegan)
#1-2
d.spe<-read.csv2("doubs_spe.csv", row.names = 1)
summary(d.spe)

colSums(d.spe>0)

#3
d.spe.hell<-decostand(d.spe, method="hellinger")

#4
pca.1<-rda(d.spe.hell)
summary(pca.1)#neni moc informativni

pca.1#lepsi vystup ukaze celkovou variabilitu 
# a variability vysvetlene prvnimi osami

# Vypocty vysvetlene variability pro 1.-3. osu
24.9/50.23
6.59/50.23
4.61/50.23

#5
screeplot(pca.1, bstick = T)

abline(h=mean(pca.1$CA$eig))# Kaiser-Guttman kriterium


#7
plot(pca.1)# Rychly ale ne moc pekny graf

#Malovani pomoci ordiplot
#Pro dokumentaci je treba mrknout na 
#?ordiplot i ?plot.cca
#Ordinacni diagram pro 1. a 2. osu
ordiplot.obj<-ordiplot(pca.1, type="n",
                       scaling=2,
xlab="PC1 (49.6%)", ylab="PC2 (13.1%)")
text(ordiplot.obj, what="sp", arrows=T,
     length=0.05)
#Pro 1. a 3. osu
ordiplot.obj<-ordiplot(pca.1, type="n",
                       scaling=2,
xlab="PC1 (49.6%)", ylab="PC3 (9.2%)",
choices = c(1,3))
text(ordiplot.obj, what="sp", arrows=T,
     length=0.05)
# Diagram vzorku pro 1. a 2. osu
ordiplot.obj<-ordiplot(pca.1, type="n",
                       scaling=1,
xlab="PC1 (49.6%)", ylab="PC2 (13.1%)")
points(ordiplot.obj, what="si", arrows=F)
text(ordiplot.obj, what="si", arrows=F)

#8
bk.env<-read.csv2("env.csv", row.names=1)
#9
summary(bk.env)

#10
model.matrix(~as.factor(bk.env$hemiparasite))
### Dummy kodovani je nejlepsi udelat pomoci funkce
#dummy_cols v package fastDummies
install.packages("fastDummies")
library(fastDummies)
hemi_dummy<-dummy_cols(bk.env,"hemiparasite")# Zakoduje dummy promennou
# z faktoru a ulozi ji na konec data frame
bk.env.2<-hemi_dummy[,-2]# vymazani puvodniho faktoru z df
names(bk.env.2)
names(bk.env.2)[14:17]<-c("MNE", "RAL", "RMA", "RMI")# prejmenovani dummy promennych

pca.bk<-rda(bk.env.2, scale=T)
pca.bk
# Call: rda(X = bk.env.2,
# scale = T)
# 
#               Inertia Rank
# Total              17     
# Unconstrained      17   16
# Inertia is correlations 
# 
# Eigenvalues for unconstrained axes:
#   PC1   PC2   PC3   PC4   PC5 
# 4.605 2.519 2.090 1.915 1.320 

#11
screeplot(pca.bk, bstick=T)
abline(h=mean(pca.bk$CA$eig))

# Vyrvoreni ordinacniho diagramu pro
# 1-4 osu a ulozeni do pdf.
pdf("PCA1-4.pdf", 8,5)
par(mfrow=c(1,2))
bk.op<-ordiplot(pca.bk, type="n", 
                scaling=3)
text(bk.op, what="sp", arrows=T, 
     length=0.03, col="darkblue")
points(bk.op, what="si")

bk.op<-ordiplot(pca.bk, type="n",
                choices=c(3,4), 
                scaling=3)
text(bk.op, what="sp", arrows=T, 
     length=0.03, col="darkblue")
points(bk.op, what="si")
dev.off()