Доверительный интервал случайных эффектов с lmer

Я использую lmer от lme4 пакет для расчета доверительного интервала для дисперсионного компонента.

Когда я подгоняю модель есть предупреждающие сообщения:

fit <- lmer(Y~X+Z+X:Z+(X|group),data=sim_data)
Warning messages:
1: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv,  :
  unable to evaluate scaled gradient
2: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv,  :
  Model failed to converge: degenerate  Hessian with 1 negative eigenvalues

Я много искал, чтобы понять, почему возникает ошибка, и, наконец, пришел к решению, что there is difference between error and warning in the R world,

Я хочу вычислить доверительный интервал для параметров модели и запустить код, который показывает ошибку:

 confint.merMod(fit,oldNames=FALSE)
 Computing profile confidence intervals ...
 Error in if (all(x[iu <- upper.tri(x)] == 0)) t(x[!iu]) else t(x)[!iu] : 
 missing value where TRUE/FALSE needed

Есть ли другой способ получения КИ случайных эффектов с помощью lmer?

РЕДАКТИРОВАТЬ:

simfun <- function(J,n_j,g00,g10,g01,g11,sig2_0,sig01,sig2_1){
    N <- sum(rep(n_j,J))  
    x <- rnorm(N)         
    z <- rnorm(J)        

    mu <- c(0,0)
    sig <- matrix(c(sig2_0,sig01,sig01,sig2_1),ncol=2)
    u   <- rmvnorm(J,mean=mu,sigma=sig)

    b_0j <- g00 + g01*z + u[,1]
    b_1j <- g10 + g11*z + u[,2]

    y <- rep(b_0j,each=n_j)+rep(b_1j,each=n_j)*x + rnorm(N,0,0.5)
    data <- data.frame(Y=y,X=x,Z=rep(z,each=n_j),group=rep(1:J,each=n_j))
  } 

noncoverage <- function(J,n_j,g00,g10,g01,g11,sig2_0,sig01,sig2_1){
    dat <- simfun(J,n_j,g00,g10,g01,g11,sig2_0,sig01,sig2_1)
    fit <- lmer(Y~X+Z+X:Z+(X|group),data=dat)
 }

comb1 = replicate(1000,noncoverage(10,5,1,.3,.3,.3,(1/18),0,(1/18)))
comb26 = replicate(1000,noncoverage(100,50,1,.3,.3,.3,(1/8),0,(1/8)))
Источник

2 ответа

Это зависит от того, что вы ищете именно от доверительных интервалов, но функция sim в arm пакет обеспечивает отличный способ получить повторные образцы из задней части lmer или же glmer Объект, чтобы получить представление об изменчивости коэффициентов как фиксированных, так и случайных членов.

в merTools В пакете мы написали оболочку, которая упрощает процесс извлечения этих значений и взаимодействия с ними:

library(merTools)
randomSims <- REsim(fit, n.sims = 500)
# and to plot it
plotREsim(REsim(fit, n.sims = 500))

Существует ряд других инструментов для их изучения в merTools, Если вы хотите получить фактические результаты моделирования, вам лучше использовать arm::sim,

Источник

Модель Лмера показывает результаты в матрице. Вы можете получить доступ к оценкам и стандартным ошибкам модели для расчета доверительных интервалов.

Поскольку первая строка - это оценка пересечения модели, вторая строка - это оценка вашей управляемой переменной. Первый столбец - это оценка эффекта, а второй столбец - стандартная ошибка.

Так меняется MyModel для названия вашей модели и округления числа до двух, round(coef(summary(MyModel))[2,1],2)-round(coef(summary(MyModel))[2,2],2)*2 даст вам нижний предел доверительного интервала, в то время как простое изменение вычитания для сложения в предыдущей формуле даст вам верхний предел оценки: round(coef(summary(MyModel))[2,1],2)+round(coef(summary(MyModel))[2,2],2)*2

Источник

Вы можете использовать пакет параметров, который предлагает вам методы для вычисления CI для фиксированных эффектов, включая приближение Кенварда-Роджера для небольших размеров выборки, или стандартные ошибки для фиксированных и случайных эффектов, или полные сводки моделей. См. Примеры ниже.

library(parameters)
library(lme4)
#> Loading required package: Matrix
model <- lmer(mpg ~ wt + (1 | gear), data = mtcars)

ci(model)
#>     Parameter CI   CI_low   CI_high
#> 1 (Intercept) 95 31.89749 40.482616
#> 2          wt 95 -6.29877 -3.791242

ci_kenward(model)
#>     Parameter CI    CI_low   CI_high
#> 1 (Intercept) 95 31.208589 41.171513
#> 2          wt 95 -6.573142 -3.516869

standard_error(model)
#>     Parameter        SE
#> 1 (Intercept) 2.1901246
#> 2          wt 0.6396873

standard_error(model, effects = "random")
#> $gear
#>   (Intercept)
#> 3   0.6457169
#> 4   0.6994964
#> 5   0.9067223

model_parameters(model, df_method = "satterthwaite")
#> Parameter   | Coefficient |   SE |         95% CI |     t |    df |      p
#> --------------------------------------------------------------------------
#> (Intercept) |       36.19 | 2.19 | [31.90, 40.48] | 16.52 | 13.85 | < .001
#> wt          |       -5.05 | 0.64 | [-6.30, -3.79] | -7.89 | 21.92 | < .001

Создано 07.02.2020 с помощью пакета REPEX (v0.3.0)

Источник

Проблема здесь в том, что ваш случайный эффект (X | группа), и я предполагаю, что это должна быть (1 | группа), модель случайного перехвата или включать обе как (1+X| группа). Наличие только (X|group) вызовет проблемы, поскольку вы допускаете разные уклоны, но не позволяете изменяться перехвату.

Источник
Другие вопросы по тегам