今天的推文內(nèi)容主要參考
首先是讀入數(shù)據(jù)
今天推文用到的示例數(shù)據(jù)是參考鏈接2中提供的
usflu.fasta盏筐，fasta文件已經(jīng)比對(duì)好顶掉，R語言里讀入fasta格式的數(shù)據(jù)可以使用adegenet包中的fasta2DNAbin函數(shù)
#install.packages("adegenet")
library(adegenet)
dna<-fasta2DNAbin(file = "usflu.fasta")
dna
計(jì)算距離矩陣
library(ape)
dd<-dist.dna(dna)
用到的是
ape包中的dist.dna()函數(shù)
構(gòu)建NJ樹
tree<-nj(dd)
用到的是
ape包中的nj()函數(shù)
ggtree進(jìn)行可視化
library(ggtree)
ggtree(tree)+
geom_tiplab(size=2)
image.png
做bootstrap檢驗(yàn)
bs.tree<-boot.phylo(tree,dna,
function(x)nj(dist.dna(x)),1000,
trees = TRUE)
將得到的bootstrap值合并到tree中
tree$node.label<-bs.tree$BP
這一步不知道對(duì)不對(duì)皆愉，好像是有問題，暫時(shí)還不知道如何驗(yàn)證
結(jié)果里展示bootstrap值
ggtree(tree)+
geom_tiplab(size=2)+
geom_nodelab(hjust=-0.2,size=2)
image.png
歡迎大家關(guān)注我的公眾號(hào)
小明的數(shù)據(jù)分析筆記本
