1

I would like to plot the emmeans of a glmmTMB model using plot(). When my glmmTMB model takes in log transformed data, plot(emmeans(glmmTMB_model)) works just fine. However, when I attempt to plot the emmeans of the glmmTMB model of non-transformed data, I get the following error: Error in linkinv(summ$the.emmean) : could not find function "linkinv".

See my code below:

Site <- c(4,4,5,5)
Treatment <- c("Burnt_Intact-Vegetation","Burnt_Intact-Vegetation", "Burnt_Cleared", "Burnt_Cleared")
pH <- c(5.94, 5.91, 5.44, 5.49)

pH_EC_data <- data_frame(Site, Treatment, pH)

pH_model <- glmmTMB(pH~Treatment+(1|Site), data = pH_EC_data)
log10pH_model<- glmmTMB(log10(`pH`)~Treatment+(1|Site), data = pH_EC_data)


log10pH_analysis <- emmeans(log10pH_model, pairwise~Treatment, type = "response")
plot(log10pH_analysis)
##This plot works just fine.

pH_analysis <- emmeans(pH_model, pairwise~Treatment, type = "response")
plot(pH_analysis)
##This code results in the following error: Error in linkinv(summ$the.emmean) : could not find function "linkinv"

Note, log10pH_analysis and pH_analysis differ by one column. Emmeans of glmmTMB of logged data creates a "response" column whereas the same manipulation of non-tranformed data resulted in an "emmeans" column. See below:

log10pH_analysis

$emmeans
 Treatment               response     SE df lower.CL upper.CL
 Burnt_Cleared               5.47 0.0167  2     5.40     5.54
 Burnt_Intact-Vegetation     5.90 0.0180  2     5.82     5.98

Confidence level used: 0.95 
Intervals are back-transformed from the log10 scale 

$contrasts
 contrast                                  estimate      SE df t.ratio p.value
 Burnt_Cleared - (Burnt_Intact-Vegetation)  -0.0329 0.00188  2 -17.520 0.0032 

Note: contrasts are still on the log10 scale

pH_analysis

$emmeans
 Treatment               emmean     SE df lower.CL upper.CL
 Burnt_Cleared             5.47 0.0176  2     5.39     5.55
 Burnt_Intact-Vegetation   5.90 0.0176  2     5.82     5.98

Confidence level used: 0.95 

$contrasts
 contrast                                  estimate     SE df t.ratio p.value
 Burnt_Cleared - (Burnt_Intact-Vegetation)    -0.43 0.0249  2 -17.238 0.0033

Thank you.

4

2 回答 2

0

似乎 type = "response" 导致了未转换数据的问题。

pH_analysis <- emmeans(pH_model, pairwise~Treatment)
plot(pH_analysis)
##Plots beautifully; problem resolved
于 2020-11-23T03:28:23.433 回答
0

我可以用一个类似的例子来复制它。由于编码错误，我将在包的下一次更新中更正。

问题与您的pH_analysis对象实际上是两个emmGrid对象的列表有关 - 一个用于估计的边际均值，另一个用于它们的成对比较。如果你这样做

plot(pH_analysis, which = 1)    # or plot(pH_analysis[[1]])

它会工作得很好。在尝试绘制第二个时也会出现错误，并且编码错误与成对比较不是在变换尺度上但它认为它们是这样的事实有关。

我建议您不要plot()在调用结果中使用emmeans()双面公式或规格列表；而是选择你真正想要绘制的部分。另外，我认为您必须拥有相当旧的emmeans版本，因为plot()最近版本的默认值是which = 1.

感谢您报告此错误。

于 2020-11-24T22:38:49.610 回答