每个物种使用多个条目的系统发育模型

Mar*_*mba 4 mixed-models phylogeny

我对系统发育回归模型比较陌生。过去,当我的树中每个物种只有 1 个条目时,我使用了 PGLS。现在我有一个包含总共 9 个物种的数千条记录的数据集,我想运行一个系统发育模型。我阅读了最常见软件包(例如 caper)的教程,但我不确定如何构建模型。

当我尝试为雀跃创建对象时,即使用:

obj <- comparative.data(phy = Tree, data = Data, names.col = species, vcv = TRUE, na.omit = FALSE, warn.dropped = TRUE)
Run Code Online (Sandbox Code Playgroud)

我收到消息:

row.names<-.data.frame( *tmp*, value = value) 中的错误:不允许重复的“row.names”此外:警告消息:设置“row.names”时的非唯一值:“Species1”、“Species2”、“Species3”、“Species4” '、'物种 5'、'物种 6'、'物种 7'、'物种 8'、'物种 9'

我知道我可以通过应用 MCMCglmm 模型来解决这个问题,但我不熟悉贝叶斯模型。

在此先感谢您的帮助。

Tho*_*rme 5

这确实不适用于简单的 PGLS from,caper因为它无法将个体作为随机效应处理。我建议您使用MCMCglmm它来理解并不会复杂得多,并且可以让您将个人作为随机效果。你可以找到包的作者出色的文档在这里在这里或更处理的包(即树的不确定性)的某些具体方面的替代文件在这里

真的很简单,让你开始:

## Your comparative data
comp_data <- comparative.data(phy = my_tree, data =my_data,
      names.col = species, vcv = TRUE)
Run Code Online (Sandbox Code Playgroud)

请注意,您可以拥有一个如下所示的样本列:

   taxa        var1 var2 specimen
1     A  0.08730689    a    spec1
2     B  0.47092692    a    spec1
3     C -0.26302706    b    spec1
4     D  0.95807782    b    spec1
5     E  2.71590217    b    spec1
6     A -0.40752058    a    spec2
7     B -1.37192856    a    spec2
8     C  0.30634567    b    spec2
9     D -0.49828379    b    spec2
10    E  1.42722363    b    spec2
Run Code Online (Sandbox Code Playgroud)

然后你可以设置你的公式(类似于一个简单的lm公式):

## Your formula
my_formula <- variable1 ~ variable2
Run Code Online (Sandbox Code Playgroud)

以及您的 MCMC 设置:

## Setting the prior list (see the MCMCglmm course notes for details)
prior <- list(R = list(V=1, nu=0.002),
              G = list(G1 = list(V=1, nu=0.002)))

## Setting the MCMC parameters
## Number of interactions
nitt <- 12000

## Length of burnin
burnin <- 2000

## Amount of thinning
thin <- 5
Run Code Online (Sandbox Code Playgroud)

然后您应该能够运行默认值MCMCglmm

## Extracting the comparative data
mcmc_data <- comp_data$data

## As MCMCglmm requires a column named animal for it to identify it as a phylo
## model we include an extra column with the species names in it.
mcmc_data <- cbind(animal = rownames(mcmc_data), mcmc_data)
mcmc_tree <- comp_data$phy

## The MCMCglmmm
mod_mcmc <- MCMCglmm(fixed = my_formula, 
                     random = ~ animal + specimen, 
                     family = "gaussian",
                     pedigree = mcmc_tree, 
                     data = mcmc_data,
                     nitt = nitt,
                     burnin = burnin,
                     thin = thin,
                     prior = prior)
Run Code Online (Sandbox Code Playgroud)