[R-sig-ME] prior specification for multinomial model MCMCglmm

Vincenzo Ellis vincenzoaellis at gmail.com
Thu Apr 6 21:44:06 CEST 2017


Dear list members,

I'm wondering if someone could give me some guidance on how to set up the
priors, specifically regarding random effects, for a multinomial model in
MCMCglmm. The response variable has three levels. The model has an
intercept, but no explanatory variables and it has three nested random
effects. The point of the model is to partition variance in the response
variable among the three random effect terms.

Following the course notes for the package, I think I have the R structure
right for the prior, but I'm finding very little information online to help
with specifying the G structure--both what kind of prior is appropriate in
this case and how to code it. There is at least one worked example (
https://hlplab.wordpress.com/2009/05/07/multinomial-random-effects-models-in-r/)
but it cites a tutorial for most of the explanations that I cannot find.

Any advice would be much appreciated. Data and code (with the prior
specification incomplete) follow.

Thanks!

Vincenzo

## load data
library(gsheet)
dat <- gsheet2tbl('
https://docs.google.com/spreadsheets/d/1qVABAjNSBIlsoG3FYjREejBKyGkJvjaPiIJ_T_bzO_Y/edit?usp=sharing
')

## group taxonomic categories into unique categories for MCMCglmm for use
as nested random effects
dat$Family.Genus <- paste0(dat$Family, dat$Genus)
dat$Family.Genus.Species <- paste0(dat$Family, dat$Genus, dat$Species)

## make variables factors...they come out as characters in the download
dat <- as.data.frame(unclass(dat))

## set up model
library(MCMCglmm)

## k = number of categories in response variable
k <- length(levels(dat$Nest.Type))

## I and J are matrices that will set up constraints on the residuals of
the model
I <- diag(k-1)
J <- matrix(1, k-1, k-1)

## set up prior; how to set up the G structure?
prior1 <- list(R = list(V = (1/k)*(I + J), fix = 1))

## model
mod <- MCMCglmm(Nest.Type ~ trait - 1, rcov = ~us(trait):units,
                random = ~us(trait):Family + us(trait):Family.Genus +
us(trait):Family.Genus.Species,
                prior = prior1, data = dat, family = "categorical")

	[[alternative HTML version deleted]]



More information about the R-sig-mixed-models mailing list