multcomp::glht 中具有交互作用的混合效應模型 (lme4) 的事後測試

R我正在對（lme4包）中的線性混合效應模型執行事後測試。我正在使用multcomp包（glht()函數）來執行事後測試。

我的實驗設計是重複測量，具有隨機塊效應。模型指定為：

mymod <- lmer(variable ~ treatment * time + (1|block), data = mydata, REML = TRUE)

我沒有在此處附加我的數據，而是使用包中調用的warpbreaks數據multcomp。

data <- warpbreaks
warpbreaks$rand <- NA

我添加了一個額外的隨機變量來模仿我的“塊”效果：

warpbreaks$rand <- rep(c("foo", "bar", "bee"), nrow(warpbreaks)/3)

這模仿了我的模型：

mod <- lmer(breaks ~ tension * wool + (1|rand), data = warpbreaks) 

我知道“ Additional Multcomp Examples - 2 Way Anova”中的例子。這個例子引導你比較wool.

如果我想反其道而行之——比較 的 水平之wool內的水平tension怎麼辦？（在我的情況下，這將是在時間水平（三 - 六月、七月、八月）內比較治療水平（二 - 0、1）。

我已經提出了以下代碼來執行此操作，但它似乎不起作用（請參閱下面的錯誤消息）。

首先，從示例中（有wool和tension交換的地方）：

tmp <- expand.grid(wool = unique(warpbreaks$wool), tension = unique(warpbreaks$tension))
X <- model.matrix(~ tension * wool, data = tmp)
glht(mod, linfct = X)

Tukey <- contrMat(table(warpbreaks$wool), "Tukey")

K1 <- cbind(Tukey, matrix(0, nrow = nrow(Tukey), ncol = ncol(Tukey)))
rownames(K1) <- paste(levels(warpbreaks$tension)[1], rownames(K1), sep = ":")

K2 <- cbind(matrix(0, nrow = nrow(Tukey), ncol = ncol(Tukey)), Tukey)
rownames(K2) <- paste(levels(warpbreaks$tension)[2], rownames(K2), sep = ":")

從這裡到底部，我自己的代碼：

K3 <- cbind(matrix(0, nrow = nrow(Tukey), ncol = ncol(Tukey)), Tukey)
rownames(K2) <- paste(levels(warpbreaks$tension)[3], rownames(K3), sep = ":")

K <- rbind(K1, K2, K3)
colnames(K) <- c(colnames(Tukey), colnames(Tukey))

> summary(glht(mod, linfct = K %*% X))
Error in summary(glht(mod, linfct = K %*% X)) : 
 error in evaluating the argument 'object' in selecting a method for function 'summary': Error in K %*% X : non-conformable arguments

使用 lsmeans 包要容易得多

library(lsmeans)
lsmeans(mod, pairwise ~ tension | wool)
lsmeans(mod, pairwise ~ wool | tension)

引用自：https://stats.stackexchange.com/questions/145765