#install.packages("ape")
#install.packages("phylotools")
#install.packages("phangorn")
#install.packages("maps")
#install.packages("phytools")
#install.packages("TreeTools")

library(ape)
library(phylotools)
library(phangorn)
library(maps)
library(phytools)
library(TreeTools)

##NGS data
#Spocitani distancni matice dist.p ze souboru *.snps
input <- read.dna("C:/Users/kuba/Plocha/CIR-HYB.snps", format="sequential", as.character = TRUE)
aaPhyDat <- as.phyDat(input)
DMat <- dist.p(aaPhyDat, cost = "polymorphism", ignore.indels = TRUE)
write.nexus.dist(DMat, file="C:/Users/kuba/Plocha/OUTPUT3")

DMat

##Jednoducha prehledova heatmapa
heatmap(x=as.matrix(DMat), Rowv=TRUE, Colv=TRUE, symm=TRUE)

##Pokrocilejsi heatmapa
#vychozi strom
myTreeNGS <- upgma(DMat)
plotTree(myTreeNGS,ftype="i",fsize=0.6,lwd=1)

#preskladani stromu
NewTree3 <- root(myTreeNGS, "DCAN_9.1")
plotTree(NewTree3,ftype="i",fsize=0.6,lwd=1)
NGS_Tree <- SortTree(myTreeNGS, how = "TipLabels", order =c("DCAN_2.1", "DCAN_10.1", "DCAN_6.1", "DCAN_4.1", "DCAN_7.1", "DCAN_9.1", "DCAN_1.1", "DCAN_5.1", "DCAN_3.1", "DCAN_8.1", "CIR-D27-2.1", "CIR-D27-3.1", "CIR-D42-2.1", "CIR-D20-4.1", "CIR-D68-3.1", "CIR-D34-3.1", "CIR-D57-1.1", "CIR-D29-4.1", "CIR-DB9-3.1", "CIR-D26-1.1", "CIR-DB13-3.1", "DPAL_6.1", "CIR-D29-1.1", "DPAL_4.1", "DPAL_11.1", "DPAL_10.1", "DPAL_5.1", "DPAL_9.1", "DPAL_3.1", "DPAL_8.1", "DPAL_2.1", "DPAL_7.1"))
plotTree(NGS_Tree,ftype="i",fsize=0.6,lwd=1)

##heatmapa a strom
#prevod distancni matice na formát matice
DM_NGS <- as.matrix (DMat)

# uprava poradi vzorku dle stromu (vzorky jsou v matici dle abecedniho poradi - je potřeba je preskladat dle stromu)
DM_NGS2 <- DM_NGS[c("DCAN_2.1", "DCAN_10.1", "DCAN_6.1", "DCAN_4.1", "DCAN_7.1", "DCAN_9.1", "DCAN_1.1", "DCAN_5.1", "DCAN_3.1", "DCAN_8.1", "CIR-D27-2.1", "CIR-D27-3.1", "CIR-D42-2.1", "CIR-D20-4.1", "CIR-D68-3.1", "CIR-D34-3.1", "CIR-D57-1.1", "CIR-D29-4.1", "CIR-DB9-3.1", "CIR-D26-1.1", "CIR-DB13-3.1", "DPAL_6.1", "CIR-D29-1.1", "DPAL_4.1", "DPAL_11.1", "DPAL_10.1", "DPAL_5.1", "DPAL_9.1", "DPAL_3.1", "DPAL_8.1", "DPAL_2.1", "DPAL_7.1"), c("DCAN_2.1", "DCAN_10.1", "DCAN_6.1", "DCAN_4.1", "DCAN_7.1", "DCAN_9.1", "DCAN_1.1", "DCAN_5.1", "DCAN_3.1", "DCAN_8.1", "CIR-D27-2.1", "CIR-D27-3.1", "CIR-D42-2.1", "CIR-D20-4.1", "CIR-D68-3.1", "CIR-D34-3.1", "CIR-D57-1.1", "CIR-D29-4.1", "CIR-DB9-3.1", "CIR-D26-1.1", "CIR-DB13-3.1", "DPAL_6.1", "CIR-D29-1.1", "DPAL_4.1", "DPAL_11.1", "DPAL_10.1", "DPAL_5.1", "DPAL_9.1", "DPAL_3.1", "DPAL_8.1", "DPAL_2.1", "DPAL_7.1")]

#heatmapa ML - oprava
phylo.heatmap(NGS_Tree, X=DM_NGS2, split=c(0.7,0.7),fsize=c(0.5,0.5,0.5), colors=NULL, standardize=FALSE, lwd=1, ftype="reg", pts=TRUE, labels=TRUE)
phylo.heatmap(NGS_Tree, X=DM_NGS2, split=c(0.7,0.7),fsize=c(0.5,0.5,0.5), colors=NULL, standardize=FALSE, lwd=1, ftype="reg", pts=TRUE, labels=FALSE)