1

如果为n 个核酸序列构建了一个词频表（序列 ATG 对应于长度为 2 的两个词，AT 和 TG），则可以使用该表（直接或通过 PCA 降维后）来计算距离这些序列的矩阵，然后可以聚集成系统发育树（doi：10.1007/s00285-002-0185-3）：

library(sequinr)
Bat1 <- read.fasta(file="bat1.FASTA")
Bat1.seq <- Bat1[[1]]
Bat1.count <- as.vector(count(Bat1.seq, 2)) # count word frequencies, k < log4(Sequence length)
...
Counts <- rbind(Bat1.count, ...)
rownames(Counts) <- c("Bat1", ...)
colnames(Counts) <- c(rownames(count(Bat1.seq, 2)))
RowCounts <- rowSums(Counts)     
Counts.norm <- Counts/RowCounts  # normalise word counts for different sequence length
distance <- dist(Counts.norm, method = "euclidian")
hc <- hclust(distance, method = "average")
plot(hc)

几种病毒序列的系统发育树 几种病毒序列的系统发育树

这工作得非常好，输出看起来类似于用 ClustalX 进行多序列比对获得的树，但计算时间是几秒而不是几小时。

问题：如何衡量这些树的质量，选择最佳字长k或（如果使用 PCA）最佳分量q数量，以及距离和聚类方法？最好不要使用随机序列进行冗长的引导；-)。

4

1 回答 1

0

这棵树最重要的特点是这不是系统发育！

在系统发育中，边缘反映了进化过程，我们询问两个分类群是否有共同的祖先以及这种可能性有多大。相反，OP 图像中的树状图代表了分类群之间的 DNA 序列组成相似性，因此是一棵表型树。了解系统发育树和表型树之间的区别对于决定是否使用建议的方法至关重要。如果测试的目的是推断病毒之间的进化关系，那么这种方法是不合适的。

由于这棵树不是系统发育树，因此无需在进化历史意义上测试这些关系。

于 2021-01-21T10:57:14.473 回答