glmulti: назначение функции предсказания для glmer с двумя вложенными случайными величинами - PullRequest
0 голосов
/ 26 октября 2018

Я пытаюсь использовать glmulti с glmer для усреднения модели и для получения усредненных прогнозов по модели.Я следовал за примерами в документации по glmulti («Использование glmulti со статистической моделью любого типа, с примерами», прилагаемой к пакету) и обновлениями, представленными на этом веб-сайте ( glmulti и линейные смешанные модели ) и наБлог сопровождающего пакета (https://vcalcagnoresearch.wordpress.com/package-glmulti/). Мне удалось создать оболочку для функции glmer:

glmer.glmulti <- function (formula, data, family, random = "") {
glmer(paste(deparse(formula), random), data = data, family = binomial)
}

И мне удалось назначить метод getfit для glmer, чтобы я мог получить модельусредненные коэффициенты:

setMethod('getfit', 'merMod', function(object, ...) {
summ=summary(object)$coef
summ1=summ[,1:2]
if (length(dimnames(summ)[[1]])==1) {
summ1=matrix(summ1, nr=1, dimnames=list(c("(Intercept)"),c("Estimate","Std. Error")))
}
cbind(summ1, df=rep(10000,length(fixef(object))))
})

Следующим шагом является назначение функции предсказания для glmer. Это пример, приведенный в документации пакета:

predict.mer=function(objectmer,random=random, newdata, withRandom=F,se.fit=F, ...){
if (missing(newdata) || is.null(newdata)) {
DesignMat <- model.matrix(objectmer) }
else {
DesignMat=model.matrix(delete.response(terms(objectmer)),newdata)
}
output=DesignMat %*% fixef(objectmer)
if(withRandom){
z=unlist(ranef(objectmer))
if (missing(newdata) || is.null(newdata)) {
Zt<- objectmer@Zt
} else {
Zt<-as(as.factor(newdata[,names(ranef(objectmer))]),"sparseMatrix")
}
output = as.matrix(output + t(Zt) %*% z)
}
if(se.fit){
pvar <- diag(DesignMat %*% tcrossprod(vcov(objectmer),DesignMat))
if(withRandom){
pvar <- pvar+ VarCorr(objectmer)[[1]]
}
output=list(fit=output,se.fit=sqrt(pvar))
}
return(output)
}

Затем, чтобы получить усредненную модельпрогнозы (bab - это подобранный объект glmulti в примере):

> predict(bab, se.fit=T, withR=T)
Error in predict.merMod(coffee[[i]], se.fit = se.fit, ...) :
cannot calculate predictions with both standard errors and random effects

Я также пробовал:

> predict(bab, se.fit=T, withR=F)
Error in predict.merMod(coffee[[i]], se.fit = se.fit, ...) :
cannot calculate predictions with both standard errors and random effects

И:

> predict(bab, se.fit=F, withR=T)
Error in waou %*% t(matrix(unlist(preds), nrow = nbpo)) :
non-conformable arguments
In addition: There were 50 or more warnings (use warnings() to see the first 50)

Я не совсем уверен, что не так, хотя это может быть чем-то очевидным. Я понял, что пакет lme4 был обновлен и изменен с момента написания этого примера, так что это может быть связано с этим (?).

Другая возможность заключается в том, что в документации сказано, что эта функция будет обрабатывать только одну случайную переменную.Моя модель имеет две вложенные случайные величины: x ~ y + z + w + (1|u/v).

Мне нужно а) заставить это работать, б) обновить функцию для обработки двух случайных переменных.Любые предложения будут высоко оценены.

...