重测序群体绘制 NJ Tree 进化树

Guikong
Guikong
发布于 2024-08-08 / 89 阅读
0

重测序群体绘制 NJ Tree 进化树

这里的 NJ 树绘制适用于基于全基因组重测序后,得到多个个体的测序数据的情况。待输入的数据为变异信息文件,这里以 SNP 编译信息为例,因此待输入文件应为 vcf 格式(或者 vcf.gz),假设待输入文件名为 snp.vcf。

需要使用到的工具有 plink 1.9、 R 以及iTol 网站。

1.将vcf文件转换为plink格式

plink --vcf snp.vcf --recode --out snp_plink --allow-extra-chr

这里我们输出文件名为 snp_plink ,如果你的染色体名称为标准格式,参数--allow-extra-chr可以不加。

2.使用plink计算进化树获取dist文件

plink --file snp_plink --distance square --out NJ_distance --allow-extra-chr

这里基于我们刚才生成的 plink 格式的文件,计算 distance 数据,输出文件命名为 NJ_distance。这里我们会得到一个NJ_distance.dist文件以及其他若干文件,注意保存。

3.使用 R 代码获取 iTol 支持的 nwk 格式

首先是第一个版本的代码,这个版本的代码直接将生成的NJ_distance.dist转换为 nwk 文件,转换后导入 iTol 后绘制无个体名称。

# 加载必要的库
if (!requireNamespace("ape", quietly = TRUE)) {
  install.packages("ape")
}
library(ape)
# 读取距离文件
dist.mat <- as.matrix(read.table("NJ_distance.dist"))
# 转换为距离对象
dist.obj <- as.dist(dist.mat)
# 使用Neighbor-Joining算法构建树
tree <- nj(dist.obj)
# 导出Newick格式文件
write.tree(tree, file="NJ_distance_no_individual_name.nwk")

接下来这个版本会将个体名称加上去,导入 iTol 后会有个体名称标注。

# 加载必要的库
if (!requireNamespace("ape", quietly = TRUE)) {
  install.packages("ape")
}
library(ape)

# 读取距离文件
dist.mat <- read.table("NJ_distance.dist")

# 读取个体名称文件
ids <- read.table("NJ_distance.dist.id", header = FALSE)
individual_names <- ids$V1  # 假设个体名称在第一列

# 检查并确保个体名称唯一
unique_individual_names <- make.unique(as.character(individual_names))

# 设置矩阵的行和列名称
rownames(dist.mat) <- unique_individual_names
colnames(dist.mat) <- unique_individual_names
# 转换为距离对象
dist.obj <- as.dist(dist.mat)

tree <- nj(dist.obj)

# 导出Newick格式文件
write.tree(tree, file="NJ_distance_with_individual_name.nwk")

4.导入 iTol 进行绘制

这一步我们就不具体讲了,直接将文件上传或者复制粘贴进去就可以进行美化了,教程非常多。