glmmTMB中变量间的对比

Question

作为一个可重复的例子，让我们使用下一个无意义的例子：

> library(glmmTMB)
> summary(glmmTMB(am ~ disp + hp + (1|carb), data = mtcars))
 Family: gaussian  ( identity )
Formula:          am ~ disp + hp + (1 | carb)
Data: mtcars

     AIC      BIC   logLik deviance df.resid 
    34.1     41.5    -12.1     24.1       27 

Random effects:

Conditional model:
 Groups   Name        Variance  Std.Dev. 
 carb     (Intercept) 2.011e-11 4.485e-06
 Residual             1.244e-01 3.528e-01
Number of obs: 32, groups:  carb, 6

Dispersion estimate for gaussian family (sigma^2): 0.124 

Conditional model:
              Estimate Std. Error z value Pr(>|z|)    
(Intercept)  0.7559286  0.1502385   5.032 4.87e-07 ***
disp        -0.0042892  0.0008355  -5.134 2.84e-07 ***
hp           0.0043626  0.0015103   2.889  0.00387 ** 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

事实上，我真正的模型家庭是nbinom2。我想在disp和hp之间做一个对比测试。所以，我试着：

> glht(glmmTMB(am ~ disp + hp + (1|carb), data = mtcars), linfct = matrix(c(0,1,-1)))
Error in glht.matrix(glmmTMB(am ~ disp + hp + (1 | carb), data = mtcars),  : 
  ‘ncol(linfct)’ is not equal to ‘length(coef(model))’

如何避免此错误？

谢谢!

Answer 1

问题实际上相当简单：linfct需要是一个矩阵，其列数等于参数数。您指定了matrix(c(0,1,-1))，而没有指定行数或列数，因此R默认生成列矩阵。添加nrow=1似乎有效。

library(glmmTMB)
library(multcomp)
m1<- glmmTMB(am ~ disp + hp + (1|carb), data = mtcars)
modelparm.glmmTMB <- function (model, coef. = function(x) fixef(x)[[component]],
                               vcov. = function(x) vcov(x)[[component]],
                               df = NULL, component="cond", ...) {
    multcomp:::modelparm.default(model, coef. = coef., vcov. = vcov.,
                        df = df, ...)
}        
glht(m1, linfct = matrix(c(0,1,-1),nrow=1))