0

我正在尝试将包segmented与混合模型一起使用。我遇到了各种问题，具体取决于我是否尝试使用lmer、lme或glmmTMB. 最后的方法需要是glmmTMB，因为在我的完整模型中，我使用的是负二项分布加上自相关结构。

segmented.default包中功能的帮助segmented说

要最小化的目标函数是 logLik 函数提取的（减）值，或者它可以通过 seg.control 中的 fn.obj 参数传递

所以我做了一个小辅助函数，然后传递给segmented.default. 但是，我在使用任何混合模型时都没有运气。欢迎提出任何建议，并为该glmmTMB块的解决方案提供额外的奖励积分。

可重现的例子：

library(ggplot2)
library(segmented)

# make up data
set.seed(12)
x <- 1:100
z <- runif(100)
y <- 2+1.5*pmax(x + 50-35,0)+15*pmax(z-.5,0)+rnorm(100,0,2)
y[y < 0] <- 0
y1 <- rpois(100, y)

dati <- data.frame(x = x, y = y, y1 = y1, z = z)
# add a subject column
dati$Subject <- sample(c("A", "B", "C"), 100, replace = TRUE) 
# look at the data
ggplot(dati) + geom_point(aes(x = x, y = log(y1), colour = Subject))

# helper function
f = function(x){
    as.numeric(-logLik(x))
                }

# try a glm - this works fine
out.glm.poisson <- glm(y1 ~ x, data=dati, family = "poisson")
o1.glm.poisson <- segmented.default(out.glm.poisson,seg.Z=~x,psi=list(x=c(45)))
confint.segmented(o1.glm.poisson)

# try the same glm, but with the helper function - the results are identical to before
out.glm.poisson <- glm(y1 ~ x, data=dati, family = "poisson")
o1.glm.poisson <- segmented.default(out.glm.poisson,seg.Z=~x,psi=list(x=c(45)), control=seg.control(fn.obj="f(x)"))
confint.segmented(o1.glm.poisson)

现在到混合部分——lmer、lme 和 glmmTMB 都以不同的方式失败。

out.lmer <- lme4:::lmer(y1~x + (1|Subject), data=dati)
o1.lmer <- segmented.default(out.lmer,seg.Z=~x,psi=list(x=c(35)),
           control=seg.control(fn.obj="f(x)")) # Error: $ operator not defined for this S4 class

out.lme <- nlme:::lme(y1~x, random = ~1|Subject, data=dati) # Error: The first fit with U variables does not work..

out.glmmTMB.poisson <- glmmTMB:::glmmTMB(y1~x + (1|Subject), data=dati, family = "poisson")
o1.glmmTMB.poisson <- segmented.default(out.glmmTMB.poisson,seg.Z=~x,psi=list(x=c(30)),
                       control=seg.control(fn.obj="f(x)")) # Error in Z - PSI : non-conformable arrays
4

0 回答 0