Phylogenetic model with repeated measurements without using species mean

Interfaces brms
1

I am working on a phylogenetic regression model in brms with multiple measurements per species. I have read over the excellent brms_phylogenetics vignette several times, but still have a question about repeated measurements. I will use the sample data variables and scripts from the vignette to ask my question.

Here is the head of the dataset:

       phen  cofactor species phylo spec_mean_cf
1 107.41919 11.223724    sp_1  sp_1    10.309588
2 109.16403  9.805934    sp_1  sp_1    10.309588
3  91.88672 10.308423    sp_1  sp_1    10.309588
4 121.54341  8.355349    sp_1  sp_1    10.309588
5 105.31638 11.854510    sp_1  sp_1    10.309588
6  64.99859  4.314015    sp_2  sp_2     3.673914

spec_mean_cf is the mean value of the cofactor for the species.
There are several examples in the vignette of modeling phen ~ spec_mean_cf, including a way to incorporate the variability of the cofactor for each species by using the distance from the mean of each individual (a separate variable called within_spec_cf in the vignette).

phen ~ spec_mean_cf + (1 | phylo) + (1 | species) + within_spec_cf

However, the mean of the species is always used as a predictor variable.

In brms, is there a way to use the individual values, and not the mean, as the predictor while still accounting for species measurement variability? In MCMCglmm it would be written as:

MCMCglmm( fixed = phen ~ cofactor,  random = ~ species + phylo)

I have tried a few variations with my own data, such as:

phen ~ cofactor + (1 | phylo) +  within_spec_cf

and

phen ~ cofactor + (1 | phylo) +  (1 | within_spec_cf )

but I don’t think this is correct, as I get an estimate of zero for cofactor for each model (not possible in this case)

I apologize if this is very naive. I am trying to learn-by-doing.

Thanks!

Please also provide the following information in addition to your question:

  • Operating System: OSX 10.11.6
  • brms Version: 2.3.0
2

Go for phen ~ cofactor + (1 | phylo) + (1 | species)

In the vignette I just split up cofactor into spec_mean_cf and within_spec_cf following the book from which I have the data.

3

Thanks Paul!

I am still getting strange results, but that is an issue for another post.

4

@paul.buerkner @jonnations did you ever post about your strange results? I am getting strange results with this same example… if I run the model on the species means:

# model_repeat1 <- brm(
#   phen ~ spec_mean_cf + (1|gr(phylo, cov = A)) + (1|species), 
#   data = data_repeat, 
#   family = gaussian(), 
#   data2 = list(A = A),
#   prior = c(
#     prior(normal(0,10), "b"),
#     prior(normal(0,200), "Intercept"),
#     prior(student_t(3,0,20), "sd"),
#     prior(student_t(3,0,20), "sigma")
#   ),
#   sample_prior = TRUE, chains = 2, cores = 2, 
#   iter = 3000, warmup = 1000 #set this way higher for final runs (see website)
# )

The output matches the raw data:

WithSpecMeanCF
WithSpecMeanCF834×464 11.4 KB

whereas if I run the model on cofactor:

# model_repeat1 <- brm(
#   phen ~ cofactor + (1|gr(phylo, cov = A)) + (1|species), 
#   data = data_repeat, 
#   family = gaussian(), 
#   data2 = list(A = A),
#   prior = c(
#     prior(normal(0,10), "b"),
#     prior(normal(0,200), "Intercept"),
#     prior(student_t(3,0,20), "sd"),
#     prior(student_t(3,0,20), "sigma")
#   ),
#   sample_prior = TRUE, chains = 2, cores = 2, 
#   iter = 3000, warmup = 1000 #set this way higher for final runs (see website)
# )

It does not.

WithCofactor
WithCofactor834×464 11.6 KB

What am I missing here?

Thanks!

5

Hi @roxannebeltran, I’m no expert, but I think that because you use the species’ values and not the means, the measurement error is included in the model which causes the effect size to go down. It is less confident in the effect of cofactor on phen because there is a lot more (real world) error in the model. IMHO, this is why we should work to incorporate measurement error (or repeated measurements) in models rather than rely on species’ means.

Also, if phylogenetic signal is high, then we shouldn’t always expect the regression line to match the data. Check out the Figure 7 from this classic : Felsenstein, J. 1985. Phylogenies and the Comparative Method. The American Naturalist 125:1-15.

Maybe that helps ? :)