首页 文章

在R中运行glmer时出现错误消息

提问于
浏览
0

我试图在R中运行两个类似的广义线性混合模型 . 两个模型对于预测变量,协变量和随机因子具有相同的输入变量,但是,响应变量不同 . 型号需要 lme4 包 . Ben Bolker解决了第二个模型的问题 .

在第一个模型中,响应变量是生物量重量和 family = gaussian .

global.model <- lmer(ex.drywght ~ forestloss562*forestloss17*roaddenssec*nearestroadprim + 
                            elevation + soilPC1 + soilPC2 +
                            (1|block/fragment), 
                          data = RespPredComb, 
                          family = "gaussian") 

Predictors have the following units: 
    forestloss562 = %,
    forestloss17 = %, 
    roaddenssec = (km/km2) and
    nearestroadprim = (m).

执行此模型会显示以下警告消息:

警告信息:1:在glmer中(ex.drywght~forestloss562 * forestloss17 * roaddenssec *:使用family = gaussian(身份链接)调用glmer()作为lmer()的快捷方式已弃用;请直接调用lmer()2:一些预测变量的尺度非常不同:考虑重新缩放

然后我执行这些后续步骤(遵循Grueber等人(2011)中描述的步骤顺序:

我标准化预测因子,

stdz.model <- standardize(global.model, standardize.y = FALSE)

(需要包 arm

使用自动模型选择和所提供的“全局”模型的子集

model.set <- dredge(stdz.model)

需要包裹( MuMIn

在这里,我收到以下警告消息:

Warning message:
In dredge(stdz.model2) : comparing models fitted by REML

找到前2个AIC型号和

top.models <- get.models(model.set, subset = delta < 2)

进行模型平均

model.avg(model.set, subset = delta < 2)

在这里,我收到此错误消息:

Error in apply(apply(z, 2L, is.na), 2, all) : 
  dim(X) must have a positive length

任何关于如何可能修复此错误的建议都将非常感激 .

在第二个模型中,响应变量是丰富度,家庭是泊松 .

global.model <- glmer(ex.richness ~ forestloss562*forestloss17*roaddenssec*nearestroadprim + 
                        elevation + soilPC1 + soilPC2 +
                       (1|block/fragment), 
                     data = mydata, 
                     family = "poisson")

当我执行上面的命令时,我收到以下错误和警告消息:

错误:(maxstephalfit)PIRLS步骤减半未能减少pwrssUpdate中的偏差另外:警告消息:1:一些预测变量的尺度非常不同:考虑重新调整2:在pwrssUpdate中(pp,resp,tolPwrss,GQmat,compDev, fac,verbose):文件中的Cholmod警告'不肯定':../ Cholesky / t_cholmod_rowfac.c,第431行3:在pwrssUpdate(pp,resp,tolPwrss,GQmat,compDev,fac,verbose):Cholmod警告'不肯定的'在档案:../ Cholesky / t_cholmod_rowfac.c,第431行

请在下面找到我的数据的可重现子集:

structure(list(plot.code = structure(c(1L, 3L, 2L, 4L, 5L, 6L, 
7L), .Label = c("a100m56r", "b1m177r", "c100m56r", "d1f1r", "e1m177r", 
"f1m17r", "lf10m56r"), class = "factor"), site.code = structure(c(1L, 
3L, 2L, 4L, 5L, 6L, 7L), .Label = c("a100m56", "b1m177", "c100m56", 
"d1f1", "e1m177", "f1m17", "lf10m56"), class = "factor"), block = structure(c(1L, 
3L, 2L, 4L, 5L, 6L, 7L), .Label = c("a", "b", "c", "d", "e", 
"f", "lf"), class = "factor"), fragment = structure(c(1L, 3L, 
2L, 4L, 5L, 6L, 7L), .Label = c("a100", "b1", "c100", "d1", "e1", 
"f1", "lf10"), class = "factor"), elevation = c(309L, 342L, 435L, 
495L, 443L, 465L, 421L), forestloss562 = c(25.9, 56.77, 5.32, 
27.4, 24.25, 3.09, 8.06), forestloss17 = c(7.47, 51.93, 79.76, 
70.41, 80.55, 0, 0), roaddenssec = c(2.99, 3.92, 2.61, 1.58, 
1.49, 1.12, 1.16), nearestroadprim = c(438L, 237L, 2637L, 327L, 
655L, 528L, 2473L), soilPC1 = c(0.31, -0.08, 1.67, 2.39, -1.33, 
-1.84, -0.25), soilPC2 = c(0.4, 0.41, -0.16, 0.15, 0.03, -0.73, 
0.51), ex.richness = c(0L, 0L, 1L, 7L, 0L, 0L, 1L), ex.drywght = c(0, 
0, 1.255, 200.2825, 0, 0, 0.04)), .Names = c("plot.code", "site.code", 
"block", "fragment", "elevation", "forestloss562", "forestloss17", 
"roaddenssec", "nearestroadprim", "soilPC1", "soilPC2", "ex.richness", 
"ex.drywght"), class = "data.frame", row.names = c(NA, -7L))

1 回答

  • 8

    tl;dr 您需要在适合模型之前标准化变量,以获得更高的数值稳定性 . 我也有一些关于你最终保存它们的可行性的评论......

    source("SO_glmer_26904580_data.R")
    library("arm")
    library("lme4")
    library("MuMIn")
    

    尝试第一次适合:

    pmod <- glmer(ex.richness ~
                  forestloss562*forestloss17*roaddenssec*nearestroadprim + 
                      elevation + soilPC1 + soilPC2 +
                          (1|block/fragment), 
                  data = dat, 
                  family = "poisson")
    

    如上所述,这失败了 . 但是,我收到上面没有报告的警告:

    ## 1: Some predictor variables are on very different scales: consider rescaling
    

    这提供了一个线索 .

    缩放数字参数:

    pvars <- c("forestloss562","forestloss17",
               "roaddenssec","nearestroadprim",
               "elevation","soilPC1","soilPC2")
    datsc <- dat
    datsc[pvars] <- lapply(datsc[pvars],scale)
    

    再试一次:

    pmod <- glmer(ex.richness ~
                  forestloss562*forestloss17*roaddenssec*nearestroadprim + 
                      elevation + soilPC1 + soilPC2 +
                          (1|block/fragment), 
                  data = datsc, 
                  family = "poisson",
                  na.action="na.fail")
    

    这是有效的,虽然我们得到一个关于太大梯度的警告信息 - 我认为这实际上是可以忽略的(我们仍在努力使这些误差敏感度阈值正确) .

    据我所知,以下几行似乎正在起作用:

    stdz.model <- standardize(pmod, standardize.y = FALSE)
    ## increases max gradient -- larger warning
    model.set <- dredge(stdz.model)  ## slow, but running ...
    

    以下是我对可行性的评论:

    • 即使对随机效应参数进行计数,您的观测值只有预测变量的8倍 . 这是推动它(根据经验,每个参数应该有10-20个观测值) .
    nrow(datsc) ## 159
     ncol(getME(pmod,"X"))  ## 19
    
    • 疏浚/多模型 - 对有和没有相互作用的模型进行平均可能是危险的 - 至少,为了使其可解释,必须使连续变量居中 . (我不知道 dredge 在这种情况下是否做了任何可以理解的事情 . )

    我也尝试了 glmmLasso 这个问题 - 它最终缩小了所有固定效果条款......

    library("glmmLasso")
    datsc$bf <- interaction(datsc$block,datsc$fragment)
    glmmLasso(ex.richness ~
              forestloss562+forestloss17+roaddenssec+nearestroadprim + 
                  elevation + soilPC1 + soilPC2,
              rnd=list(block=~1,bf=~1),
                  data = datsc, 
                  family = poisson(),
              lambda=500)
    

相关问题