Skip to content

problem running mixed model including a phylogeny suggested in vignette "Covariance structure"

4 messages · Ben Bolker, Pablo Inchausti

#
Hi,
I have been trying to run the mixed model incorporating a phylogeny that is
suggested in
https://cran.r-project.org/web/packages/glmmTMB/vignettes/covstruct.html
(section Proptional). I have encountered errors that I cannot resolve and
hence this post.

This is what I did:
require(ade4); require(ape); require(glmmTMB)

data(carni70) # loads the data

carnidat <- data.frame(species = rownames(carni70$tab), carni70$tab)

phylo <- read.tree(text=carni70$tre) # reads the phylogeny

A <- vcv(phylo) # uses package ape to calculate the var-covar matrix from
the phylogeny

# verifies whether the species names of  rows and columns of var-covar
phylogenetic matrix are the same as those in the dataframe

setdiff(carnidat$species, rownames(A)) # No!. I then correct them.

carnidat$species=gsub("_", ".", carnidat$species) # corrects all species
names

setdiff(colnames(A), rownames(A)) # final check. All is fine.

setdiff(carnidat$species, rownames(A)) # another final check. All is fine.

carnidat$dummy <- factor(1) # a dummy grouping variable must be added to
the dataset

Then I ran the model:

fit_phylo=glmmTMB(log(range)~log(size)+propto(0+species|dummy, A),
data= carnidat)

And I obtain the following error message:

Error: column or row names of the matrix do not match the terms.
Expecting names:
?speciesAtelocynus.microtis??speciesBassariscus.alleni??
speciesBassariscus.astutus??
speciesBassariscus.beddardi??speciesBassariscus.gabbii??speciesBassariscus.lasius??speciesBassariscus.pauli??
speciesBassariscus.sumichrasti??speciesCanis.latrans??speciesCanis.lupus??speciesCerdocyon.

(This is much longer, but one can get the idea)

Just in case, the same model does run using the package brms.

I would appreciate very much any suggestion that you might have.
Thanks in advance.
Cheers
Pablo
Platform: x86_64-pc-linux-gnu
Running under: Ubuntu 22.04.5 LTS


Attached base packages:

[1] stats     graphics  grDevices utils     datasets  methods   base

other attached packages:
[1] glmmTMB_1.1.10 ape_5.8

loaded via a namespace (and not attached):
 [1] nlme_3.1-166        TH.data_1.1-2       estimability_1.5.1
reformulas_0.3.0    xtable_1.8-4
 [6] minqa_1.2.8         zoo_1.8-12          TMB_1.9.15
lme4_1.1-35.5       grid_4.4.3
[11] MASS_7.3-61         mvtnorm_1.3-1       numDeriv_2016.8-1.1
compiler_4.4.3      multcomp_1.4-26
[16] codetools_0.2-20    sandwich_3.1-1      emmeans_1.10.5
coda_0.19-4.1       Rcpp_1.0.13
[21] mgcv_1.9-1          rstudioapi_0.16.0   lattice_0.22-6
digest_0.6.37       nloptr_2.1.1
[26] Rdpack_2.6.1        parallel_4.4.3      splines_4.4.3
rbibutils_2.3       Matrix_1.7-0
[31] tools_4.4.3         boot_1.3-31         survival_3.7-0

  
    
3 days later
#
The short answer is that the *values* of the names are OK, but the
*order* is wrong (setdiff() only checks that the values are
equivalent, not the order).  I will try to add some better
tests/diagnostics to the code.

On Thu, Mar 13, 2025 at 12:10?PM Pablo Inchausti
<pablo.inchausti.f at gmail.com> wrote:
2 days later
#
Hi,
First, thanks again for your quick response.
I have sorted the dataframe carnidat to match the order of  rows and
columns of species names of the var-covar matrix obtained from the
phylogeny. Nevertheless, there is still a problem that I cannot solve.
--------------------------------
For completeness, I recopy the entire set of commands that I have used.

data(carni70)  #loads the data
carnidat <- data.frame(species = rownames(carni70$tab), carni70$tab) #
generate the data frame
phylo <- read.tree(text=carni70$tre)  # the phylogeny
A <- vcv(phylo) # var-covar matrix from the phylogeny

# verifies if names of species are the same in rows of var-covar
phylogenetic matrix and in the data frame
setdiff(carnidat$species, rownames(A))
carnidat$species=gsub("_", ".", carnidat$species) # corrects all species
names
setdiff(colnames(A), rownames(A)) # final check
setdiff(carnidat$species, rownames(A))  # another final check
row.names(carnidat)=NULL # deletes the row names from the data frame

# checks that the order of levels in carnidat$Species is the same as in
rows and columns of  var-covar matrix from phylogeny
carnidat=carnidat[match(carnidat$species, rownames(A)),] # sorts the
dataframe in the order given by species from the var-covar matrix
[22] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
[43] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
[64] TRUE TRUE TRUE TRUE TRUE TRUE TRUE

carnidat$dummy <- factor(1) # a dummy grouping variable must be added to
the dataset

I show just a few species names in the data frame and in the var-covar
matrix A:
[4] "Leopardus.pardalis"       "Oreailurus.jacobita"
"Oncifelis.colocolo"
 [7] "Oncifelis.guigna"         "Oncifelis.geoffroyi"
"Leopardus.tigrinus"
[10] "Lynx.rufus"               "Lynx.canadensis"          "Panthera.onca"
[4] "Leopardus.pardalis"       "Oreailurus.jacobita"
"Oncifelis.colocolo"
 [7] "Oncifelis.guigna"         "Oncifelis.geoffroyi"
"Leopardus.tigrinus"
[10] "Lynx.rufus"               "Lynx.canadensis"
"Panthera.onca"


When fitting the model, I obtained the same error message as before:
(the error message is longer but you can get the idea)
------------------

Again, any suggestion would be appreciated.
Cheers

Pablo
On Mon, 17 Mar 2025 at 11:47, Ben Bolker <bbolker at gmail.com> wrote:

            

  
    
#
It is the *levels* of the species factor that need to be ordered properly.

Adding

carnidat$species <- factor(carnidat$species, levels = rownames(A))

before fitting the model should work.

 I've added some more details to the vignette (since this is obviously
not clear enough yet!).
  I think  we should also consider automatically re-ordering in this case ...

On Thu, Mar 13, 2025 at 12:10?PM Pablo Inchausti
<pablo.inchausti.f at gmail.com> wrote: