Филогенетическая модель с использованием нескольких записей для каждого вида

Я относительно новичок в моделях филогенетической регрессии. В прошлом я использовал PGLS, когда у меня было только 1 запись для каждого вида в моем дереве. Теперь у меня есть набор данных с тысячами записей для 9 видов, и я хотел бы запустить филогенетическую модель. Я прочитал руководство по наиболее распространенным пакетам (например, caper), но не знаю, как построить модель.

Когда я пытаюсь создать объект для каперсов, т.е. используя:

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

Я получаю сообщение:

Ошибка в row.names<-.data.frame(*tmp*, значение = значение): дубликаты «row.names» не допускаются. Кроме того: Предупреждающее сообщение: неуникальные значения при настройке «row.names»: «Species1», «Species2», «Species3», «Вид4», «Вид5», «Вид6», «Вид7», «Вид8», «Вид9».

Я понял, что могу решить эту проблему, применив модель MCMCglmm, но я не знаком с байесовскими моделями.

Заранее спасибо за вашу помощь.


mixed-models phylogeny
person Marco Gamba    schedule 02.07.2018    source источник

Ответы (1)


arrow_upward
5
arrow_downward

Это действительно не будет работать с простым PGLS из caper, потому что он не может работать с отдельными людьми как со случайным эффектом. Я предлагаю вам использовать MCMCglmm, который не намного сложнее для понимания и позволит вам иметь людей в качестве случайного эффекта. Вы можете найти отличную документацию от автора пакета здесь или здесь или другую документацию, которая больше касается некоторые специфические аспекты пакета (а именно неопределенность дерева) здесь .

Очень кратко, чтобы вы начали:

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

Обратите внимание, что у вас может быть столбец образцов, который может выглядеть следующим образом:

   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

Затем вы можете настроить свою формулу (похожую на простую формулу lm):

## Your formula
my_formula <- variable1 ~ variable2

И ваши настройки 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

И тогда вы сможете запустить 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)
person Thomas Guillerme    schedule 04.07.2018
comment
Большое спасибо. Попробую и сообщу, справлюсь ли я с этим или нет. - person Marco Gamba; 04.07.2018
comment
В чем выгода от указания индивидуумов в качестве случайного фактора? Насколько я понимаю, все, что он делает, это учитывает чрезмерную дисперсию (как мы поступили бы в обобщенной модели со смешанными эффектами, как это было предложено Зууром и др. в одной из их книг). Я прав? Заранее спасибо! - person Marco Plebani; 14.11.2018
comment
Есть ли функция (возможно, в MCMCglmm), позволяющая собирать сравнительные данные при наличии нескольких наблюдений за видом? В коде, предложенном @thomas-guillerme, comp_data по-прежнему создается с использованием comparative.data из caper, что дает ошибку, описанную @marco-gamba. Ваше здоровье - person Marco Plebani; 14.11.2018
comment
Вы можете использовать функцию as.mulTree из пакета mulTree, которая создает объекты, работающие с несколькими образцами и передаваемые до MCMCglmm через mulTree. См. пример документации mulTree для более подробной информации. - person Thomas Guillerme; 15.11.2018