[R-sig-ME] contrast of predicted slopes from MCMCglmm model

[Lenth, Russell V]

Daniel,

As long as you can get the regression coefficients, covariance of the coefficients, and a suitable grid of linear functions, you can obtain an lsmobj object that can be summarized and further analyzed in the 'lsmeans' package. For your question, I think something like this will work (I'm assuming your fitted mcmcGLMM object is named 'mod'):

# Get coefs and covariance matrix
bhat <- apply(mod$Sol, 2, mean)
V <- cov(mod$Sol)

# Set up desired factor levels (you likely want to change this)
levels <- list(day = 1:4, treatment = c("control", "cold"), species = factor(1:10))

# Obtain a differece quotient of model matrices 
lv1 <- lv2 <- levels
lv1$day = lv1$day - .0001
lv2$day = lv2$day + .0001
grid1 <- do.call(expand.grid, lv1)
grid2 <- do.call(expand.grid, lv2)
lf1 <- model.matrix(fixedForm, data = grid1)
lf2 <- model.matrix(fixedForm, data = grid2)
linfct <- (lf2 - lf1) / .0002

# Create the needed lsmobj
library(lsmeans)
mod.lsm <- lsmobj(bhat = bhat, V = V, levels = levels, 
    linfct = linfct, by = "species")

# Can now use summary, contrast, pairs, etc. on this result


I hope this helps

Russ

Russell V. Lenth  -  Professor Emeritus
Department of Statistics and Actuarial Science   
The University of Iowa  -  Iowa City, IA 52242  USA   
Voice (319)335-0712 (Dept. office)  -  FAX (319)335-3017


-----Original Message-----
Date: Thu, 09 Apr 2015 18:32:11 -0700
From: Daniel Fulop <dfulop.ucd at gmail.com>
To: R_MixedModels <r-sig-mixed-models at r-project.org>
Subject: [R-sig-ME] contrast of predicted slopes from MCMCglmm model
Message-ID: <5527281B.3080607 at gmail.com>
Content-Type: text/plain; charset=ISO-8859-1; format=flowed

Hi mixed modelers,

I would like to contrast average slopes using predictions from an LMM model fit with MCMCglmm.  In essence I want to do something analogous to lsmeans' lstrends() function.  The predicted slopes are growth rates ...more below.

I have 12 days of plant growth data for 10 closely related species at 2 temperatures.  There are several individuals per species at each temperature, and they're fully randomized.  I want to asses whether each species' growth rate differs between temperatures.  What I mean by that is that I would like to contrast the average slopes for each species in control vs. cold temperature; these slopes are the growth rates.

I am modeling the data in MCMCglmm so I can account for phylogeny and also to be able to try multi-response models (to jointly model the stem, plastochrons, root, etc).  I am starting out with the stem length data and after trying different time dependencies settled on this 3rd degree polynomial model:

stemLen ~ poly(day, 3, raw=TRUE) * treatment * species + (1 + day | indiv)

where day is a numeric time variable, treatment a 2 level factor (control and cold), and species a 10 level factor.

I have a slight idea of how I would go about contrasting the average slopes for each species and that I would use the predict.MCMCglmm() function, but I would really appreciate some guidance before I go down the wrong the wrong path.

Thanks for your help!
Dan.

--
Daniel Fulop, Ph.D.
Postdoctoral Scholar
Dept. Plant Biology, UC Davis
Maloof Lab, Rm. 2220
Life Sciences Addition, One Shields Ave.
Davis, CA 95616