這兩天看用vcf文件做單倍型網(wǎng)絡(luò)的內(nèi)容，找到了一篇plos one上的論文

論文題目是

A workflow with R: Phylogenetic analyses and visualizations using mitochondrial cytochrome b gene sequences

image.png

論文提供了完整的R語言代碼和示例數(shù)據(jù)

里面一小部分內(nèi)容是關(guān)于進(jìn)化樹的可視化展示并且關(guān)聯(lián)多序列比對(duì)的結(jié)果的。記錄下這個(gè)代碼

我自己的數(shù)據(jù)是vcf文件，論文中提供的fasta格式的文件

讀取vcf文件

library(vcfR)

vcf.example<-read.vcfR("popgenome/KiwifruitPathogenFiltered.recode.vcf")
df<-vcfR2DNAbin(vcf.example)

做進(jìn)化樹并使用ggtree可視化展示

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

關(guān)聯(lián)fasta序列內(nèi)容

這里使用到的是msa這個(gè)R包

首先是安裝

BiocManager::install("msa")
library(msa)
help(package="msa")

可視化展示

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這個(gè)包是第一次接觸，還沒有學(xué)會(huì)其中函數(shù)的用法，先知道有這個(gè)功能，等到用到的時(shí)候再來學(xué)習(xí)吧

今天推文的示例數(shù)據(jù)偏大，運(yùn)行時(shí)間可能會(huì)比較長(zhǎng)

歡迎大家關(guān)注我的公眾號(hào)

小明的數(shù)據(jù)分析筆記本

小明的數(shù)據(jù)分析筆記本 公眾號(hào) 主要分享：1、R語言和python做數(shù)據(jù)分析和數(shù)據(jù)可視化的簡(jiǎn)單小例子；2、園藝植物相關(guān)轉(zhuǎn)錄組學(xué)、基因組學(xué)、群體遺傳學(xué)文獻(xiàn)閱讀筆記；3、生物信息學(xué)入門學(xué)習(xí)資料及自己的學(xué)習(xí)筆記！