# hratky se stromy

# workspace phylo.RData
http://is.muni.cz/www/151317/phylo.RData

                    15 13 17


# balicek ape
install.packages("ape")
library(ape)


# Co je phylo?
par(mar=c(1,1,1,1))
plot(strom)

read.tree()
read.nexus()
write.tree()
write.nexus()

#k cemu to? random example:
library(phytools)
library(nlme)

tree<-rtree(20, rooted=T, br=runif)
x<-fastBM(tree, mu=2, sig2=10)
y<-fastBM(tree, mu=1, sig2=0.1)
data<-as.data.frame(x)
data$y<-y
all(row.names(data)==tree$tip.label)

plot.phylo(tree, use.edge.length=F, label.offset=1)
tiplabels(pch=16, cex=data$x/10)
tiplabels(pch=16, cex=data$y, col="blue", adj=1)

plot(data$x ~ data$y)
cor.test(data$x, data$y, method="s") #Spearman
m1<-gls(x~y, data=data, cor=corBrownian(phy=tree)) #Generalized least squares
anova(m1)
abline(m1, col="red")



#zpatky k obrazkum
plot.phylo(strom)
strom
summary(strom)
head(strom)


#priprava
all(row.names(gen)==strom$tip.label)
gen<-gen[strom$tip.label,]
strom<-ladderize(strom, r=F)


#ZOBRAZIT PROMENNE VEDLE JMEN - barevne i textove:

plot.phylo(strom, cex=0.75, label.offset=4, use.edge.length=F)

#standardizace: 1-(hodnota-min/max-min)
gen$stand<-(1-((gen$GS)-min(gen$GS))/(max(gen$GS)-min(gen$GS)))
# nebo
library(vegan)
gen$GS_stand<-1-decostand(gen$GS, method="range")

#exaktni zpusob:
gen$GS_barva<-rgb(red=(gen$stand)*255, 
                  green=(gen$stand)*255, 
                  blue=255, 
                  alpha=255, maxColorValue=255)
tiplabels(pch=22, bg=gen$GS_barva, col="black", adj=1, cex = 1.25)

#veselejsi zpusob:
moje.paleta<-colorRampPalette(c("white","navy"))
barva<-moje.paleta(length(gen$GS))
#gen$barva<-barva[order(order(gen$GS))]
gen$GS_barva2<-barva[rank(gen$GS, ties.method='first')]
tiplabels(pch=22, bg=gen$GS_barva2, col="black", adj=2, cex = 1.25)

#hodnoty pro kontrolu
tiplabels(text=gen$GS, fr="n", cex=0.75, adj=-1, col="red")


# to same pro transposony:
moje.paleta<-colorRampPalette(c("white","green3"))
barva<-moje.paleta(length(gen$transp))
#gen$barva<-barva[order(order(gen$GS))]
gen$t_barva<-barva[rank(gen$transp, ties.method='first')]
tiplabels(pch=22, bg=gen$t_barva, col="black", adj=2, cex = 1.25)




#OBARVIT NODY:
plot.phylo(strom, cex=0.75, use.edge.length=F, edge.color="gray")
nodelabels()
#nodelabels(pch=16, col=as.character(traity$NodeMean_A), cex=1.75) #promenna kodovana red/blue
#nodelabels(pch=16, cex=(traity$CI_A)*5) #promenna primo cislo

nodelabels(pch=16, col=as.character(res$Mean_GS), cex=1.75) #vetsi/mensi nez by mel byt
nodelabels(pch=16, col=as.character(res$SD_GS), cex=1.25) #vic/min diverzifikovany
nodelabels(pch=1, cex=res$CI_GS*5)

#spojenim se ctverecky vznikne:
plot.phylo(strom, cex=0.75, label.offset=1, use.edge.length=F, edge.color="gray")
tiplabels(pch=22, bg=gen$GS_barva2, col="black", adj=1, cex = 1.25)

#mozno nadpis (nutne mit margins v defaultni podobe a ne 1,1,1,1)
title(cex.main=0.75, main="VELIKOST GENOMU z phylocomu:
      - modre snizeni a cervene zvyseni VG v cladu z nodu vychazejicim
      - zelene narust nebo zlute pokles diverzity VG v cladu z nodu vychazejicim")




#OBARVIT VETVE
plot.phylo(strom, edge.color="red", tip.color="blue")

plot.phylo(strom, use.edge.length=F, cex=0.75, 
           edge.col=c("black", "red")[(strom$edge.length>0.01)+1])



#OBARViT ZVOLENY CLADE:
plot(strom, use.edge.length=F, label.offset=1)
strom$edge
nodelabels(fr="n", col="red")
tiplabels(fr="n", col="blue")

which.edge()
which.edge(strom, 1:10) #1:9 jsou cisla tipu!
  #vysledkem je vektor s odpovidajicimi radky strom$edge

#prirazeni barev vybranym cladeum
barvy<-rep("black", Nedge(strom))
barvy[which.edge(strom, 1:10)]<-"red"
barvy[which.edge(strom, 11:16)]<- "green3"
barvy

plot.phylo(strom, edge.col=barvy)


#prirazeni ruznych tyou car vybranum cladeum
carky<-rep(1, Nedge(strom))
carky[which.edge(strom, 1:9)]<-2
carky[which.edge(strom, 11:16)]<-3

plot.phylo(strom, edge.col=barvy, edge.lty=carky)


#spojit barevne 2 taxony
barvy<-rep("black", Nedge(strom))
barvy[which.edge(strom, c("Hieracium_porrifolium", "Hieracium_transsilvanicum"))]<-"blue"
plot.phylo(strom, edge.col=barvy)

#pro zadavani pomoci jmen tipu
library(geiger)
plot.phylo(strom, use.edge.length=F)
node.leaves()
cladeA <- node.leaves(strom, mrca(strom)["Hieracium_alpinum", "Hieracium_laniferum"])
#samozrejme lze zadat rovnou cislo prislusneho mrca nodu naseho cladu:
node.leaves(strom, 24)

which.edge(strom, cladeA) #do argumentu group funkce which.edge nelze zadat funkci, musi byt vektor
barvy<-rep("black", Nedge(strom))
barvy[which.edge(strom, cladeA)]<-"red" 
plot.phylo(strom, edge.col=barvy)