0

我有 96 个氨基酸序列，用 MAFFT 比对并手动修剪（FASTA 格式），用 ProtTest 选择氨基酸替换模型（LG+I+G 模型），用 MEGAX 进行系统发育重建（ML 方法，bootstrap 测试 1000复制，Newick 格式的树）和用 PAML 进行的祖先重建，总共有 664 个最终氨基酸位置。但是，我的对齐有插入缺失。我用字母（A 到 T）命名每个插入缺失，并且各自的酰胺酸位置范围：A:89-92、B:66-67、C:181-186、D:208-208、E:214-219 , F:244-250, G:237-296, H:278-280, I:295-295, J:329-334, K:345-349, L:371-375, M:390-425, N ：432-433，O：440-443，P：480-480，Q：500-500，R：541-544，S：600-600。序列的初始部分和结尾部分都是非常可变的，因此从位置 0 到 34（初始）和 600 到 664（最终），

我想知道，在每个祖先节点上，每个 indel 出现在祖先序列中的概率是多少。有人告诉我，包“ape - 系统发育和进化分析”上的 R-studio“ace”功能可以执行此任务。我已经安装了“ape”和“ggtree”。我检查了这个网页https://www.rdocumentation.org/packages/ape/versions/3.0-1/topics/ace，但是，我不知道如何构建脚本。我是生物学家，也是 R 的新手。

有人可以帮忙吗？将不胜感激，谢谢。

4

1 回答 1

0

很难从示例中准确确定您需要什么，但以下内容可能符合总体思路：

1 - 加载你的树R

对于这一步，您可以使用函数read.tree或read.nexus取决于您的树格式：即您的系统发育软件是输出 NEXUS 文件（通常这些文件中的第一行是#NEXUS，最后一行是end;or END;）还是 newick 输出（通常是第一行直接以系统发育开始，((my_species...以 ) 结束;。您可以找到此文件，然后在 R 中使用以下命令读取它：

## Loading the package
library(ape)
## Reading the tree
my_tree <- read.tree("<the_path_to_your_file>")

2 - 加载你的特征数据R

然后，您需要将特征数据（例如您上面列出的 indels 位置）加载为 amatrix或 a data.frame。最简单的方法是将它们设置为一种格式（“逗号分隔值”），然后您可以使用以下函数.csv读取：Rread.csv

## Reading the variables as a matrix
my_variables <- read.csv("<the_path_to_your_file>")

3 - 运行祖先字符估计

最后，您可以使用ace包中的函数为每个变量运行祖先字符估计ape：

## Selecting the variable of interest (e.g. the first column of the dataset)
one_variable <- my_variables[, 1]
## Running the ancestral character estimation for this variable
my_ace <- ace(x = one_variable, phy = my_tree, type = "discrete")
## Looking at the results
my_ace

当然还有更多，但希望这可以让你开始。

于 2021-07-09T11:07:15.227 回答