这两天看用vcf文件做单倍型网络的内容,找到了一篇plos one上的论文
论文题目是
A workflow with R: Phylogenetic analyses and visualizations using mitochondrial cytochrome b gene sequences
image.png
论文提供了完整的R语言代码和示例数据
里面一小部分内容是关于进化树的可视化展示并且关联多序列比对的结果的。记录下这个代码
我自己的数据是vcf文件,论文中提供的fasta格式的文件
读取vcf文件
代码语言:javascript复制library(vcfR)
vcf.example<-read.vcfR("popgenome/KiwifruitPathogenFiltered.recode.vcf")
df<-vcfR2DNAbin(vcf.example)
做进化树并使用ggtree可视化展示
代码语言:javascript复制library(ggtree)
library(ape)
dnbin<-dist.dna(df,model = "K80")
dnbin
tree<-nj(dnbin)
tree
ggtree(tree,branch.length = "none")
geom_tiplab()
#theme_tree2()
xlim(0,10)
image.png
关联fasta序列内容
这里使用到的是msa这个R包
首先是安装
代码语言:javascript复制BiocManager::install("msa")
library(msa)
help(package="msa")
可视化展示
代码语言:javascript复制ggtree(tree)
xlim(0,0.1)
geom_tiplab(align = T) -> p1
p1
njmsaplot<-msaplot(p1, df,
offset = 0.02,
width=1,
height = 0.5,
color = c(rep("rosybrown", 1),
rep("sienna1", 1),
rep("lightgoldenrod1", 1),
rep("lightskyblue1", 1),
"red"))
image.png
msa这个包是第一次接触,还没有学会其中函数的用法,先知道有这个功能,等到用到的时候再来学习吧