Успех вложения (биномиальный glmm) в r

Я запускаю GLMM, используя glmer() в R:

glmer(survive ~ fyear + site + fyear * site.x + (1|fyear),
family = binomial(link = logexp(shaffer.sub$exposure)),
data = shaffer.sub)

где выживание равно 0 или 1 в зависимости от того, было ли гнездо успешным или нет. Здесь вы можете увидеть, как выглядят данные:

structure(list(id = structure(1:7, .Label = c("1", "2", "3", 
"4", "5", "6", "7"), class = "factor"), year.x = structure(c(1L, 
1L, 2L, 3L, 3L, 3L, 3L), .Label = c("1994", "1995", "1999"), class = "factor"), 
    survive = structure(c(1L, 2L, 2L, 2L, 2L, 2L, 1L), .Label = c("0", 
    "1"), class = "factor"), fyear = structure(c(1L, 1L, 2L, 
    3L, 3L, 3L, 3L), .Label = c("1994", "1995", "1999"), class = "factor"), 
    site.x = structure(c(1L, 2L, 1L, 1L, 1L, 2L, 1L), .Label = c("N", 
    "S"), class = "factor")), .Names = c("id", "year.x", "survive", 
"fyear", "site.x"), row.names = c(NA, -7L), class = "data.frame")

но я получаю это предупреждающее сообщение:

*Warning messages:
1: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv,  :
  Model failed to converge with max|grad| = 0.0299425 (tol = 0.001, component 12)
2: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv,  :
  Model is nearly unidentifiable: large eigenvalue ratio
 - Rescale variables?*

Мне сказали, что я не должен использовать тот же случайный фактор в качестве фиксированного эффекта для той же модели.

В конце я хотел бы получить результат, в котором я могу видеть год, сайт и год взаимодействия: эффекты сайта. Как в таблице ANOVA (возможно ли это? Я пытался использовать summary(aov(model)), но это не работает; anova(model) тоже не работает.

Я получаю эту ошибку для команды aov():

*Error in summary`(aov(syearXsite))` : 
  error in evaluating the argument 'object' in selecting a method for function 'summary': Error in if (fixed.only) { : argument is not interpretable as logical*

Как я могу увидеть влияние этих переменных на выживаемость?


person ri_qc    schedule 09.04.2015    source источник
comment
спасибо Alexforrence!   -  person ri_qc    schedule 10.04.2015


Ответы (1)


Тот, кто сказал вам не использовать категориальную входную переменную (fyear) как фиксированный и случайный эффект, был правильным. Трудно понять, что именно порекомендовать, это зависит от количества лет и сайтов в вашем наборе данных (это данные, которые вы связали со всеми своими данными (надеюсь, что нет), или просто первые несколько строк? Сколько лет и сколько участков и сколько всего наблюдений у вас есть?)

Если вы хотите рассматривать год как случайный, а сайт как фиксированный (что было бы разумно, если у вас есть только два сайта (N против S, как видно из ваших данных) и довольно много лет, например, более 5), тогда вы можете подойти:

g1 <- glmer(survive~site.x+(site.x|fyear),
      family=binomial(link=logexp(shaffer.sub$exposure)),
      data=shaffer.sub)

Я не знаю, что такое site против site.x: я вижу только site.x в вашем фрагменте данных.

Чтобы получить информацию, попробуйте summary(g1). (Это даст вам дисперсии только для случайных эффектов, а не для фиксированных эффектов; GLM-модели не работают в том же режиме «объяснения дисперсии», как ANOVA, в частности потому, что дисперсии, объясняемые разными терминами, обычно не суммируем с общей дисперсией.)

person Ben Bolker    schedule 09.04.2015
comment
Привет, Бен, у меня есть 11 лет и сотни наблюдений за несколько лет. Всего мой набор данных составляет 4742 наблюдения и 2 участка (север и юг). Меня также интересует год как фиксированный эффект. Правда, это site = site.x, он исходит из вывода, когда я связал выживание вывода. Есть ли способ получить эффект года, сайта и года взаимодействия * сайта из этой модели? Большое спасибо! - person ri_qc; 10.04.2015
comment
Модель работает отлично! но, думаю, тогда невозможно получить взаимодействие yearXsite. Насколько я понимаю, точка пересечения - это год (p-значение = ‹2e-16 ***), поэтому, если он значительный, успех вложения меняется с годами (т.е. различается), и для участка (p-значение = 0,505) гнездование успех одинаков для обоих сайтов? а что я могу сказать об оценках? Оценка Std. Ошибка z значение Pr (›| z |) (Intercept) 4.6105 0.5134 8.981‹ 2e-16 *** site.xS -0.2855 0.4281 -0.667 0.505 Спасибо! - person ri_qc; 10.04.2015