[R] error applying user-defined link function to lmer
Kathleen Vancleef
Kathleen.Vancleef at ppw.kuleuven.be
Fri Aug 10 15:01:09 CEST 2012
Dear R users,
I'm struggling with applying a user-defined link function in lmer. For analyzing data of a 2AFC psychophysical experiment, I would like to model my binary data with a logistic function with a lower limit at 0.5 instead of 0. In a previous question this has been described as a halflogit function. To do so I wrote my own link function and would like to submit it to lmer, however this results in error messages. Below I've listed some of my attempts to solve this problem based on previous postings on this list. I've also added the corresponding problems I encountered with each of these attempts.
1) I've used the link function as specified at https://stat.ethz.ch/pipermail/r-help/2007-August/138436.html
# define link function
halflogit=function(){
half.logit=function(mu) qlogis(2*mu-1)
half.logit.inv=function(eta) .5*plogis(eta)+.5
half.logit.deriv=function(eta) .5*(exp(eta/2)+exp(-eta/2))^-2
half.logit.inv.indicator=function(eta) TRUE
half.logit.indicator=function(mu) mu>.5 & mu<1
link <- "half.logit"
structure(list(linkfun = half.logit, linkinv = half.logit.inv,
mu.eta = half.logit.deriv, validmu =
half.logit.indicator ,valideta = half.logit.inv.indicator, name = link),
class = "link-glm")
}
# submit halflogit function to lmer
rcNrPSrcNxPSxP<-lmer(Correct~-1+cNoise*f_Pos+f_Shape_con+f_Shape_con:f_Pos+(-1+cNoise*f_Pos|f_Subject), data=dens_data, family=binomial(halflogit()))
#resulting error message
Error in if (!(validmu(mu) && valideta(eta))) stop("cannot find valid starting values: please specify some", :
missing value where TRUE/FALSE needed
In addition: Warning message:
In qlogis(p, location, scale, lower.tail, log.p) : NaNs produced
2) Also making the variable global seems not to help (see also http://permalink.gmane.org/gmane.comp.lang.r.lme4.devel/4501)
mafc.logit <- function (.m = 2)
{
.m <<- as.integer(.m) # .m goes global
if (.m < 2)
stop(".m must be an integer > 1")
linkfun <- function(mu) {
mu <- pmax(mu, 1/.m + .Machine$double.eps)
qlogis((.m * mu - 1)/(.m - 1))
}
linkinv <- function(eta) {
1/.m + (.m - 1)/.m * .Call("logit_linkinv", eta, PACKAGE = "stats")
}
mu.eta <- function(eta) ((.m - 1)/.m) * .Call("logit_mu_eta",
eta, PACKAGE = "stats")
valideta <- function(eta) TRUE
link <- paste("mafc.logit(", .m, ")", sep = "")
structure(list(linkfun = linkfun, linkinv = linkinv, mu.eta = mu.eta,
valideta = valideta, name = link), class = "link-glm")
}
lmer(Correct~-1+cNoise*f_Pos+f_Shape_con+f_Shape_con:f_Pos+(-1+cNoise*f_Pos|f_Subject), data=dens_data, family=binomial(mafc.logit(2)))
Error in famType(glmFit$family) : unknown link: 'mafc.logit(2)'
3) The link function does work with glm so the problem is in lmer, I guess, not in the newly defined link function. Unfortunately I would like to model random effect, so I prefer to use the lmer function.
4) I've also tried with other versions of the package: lme4a and lme4.0, but lme4a seems not to work under the current R version, and lme4.0 gives the same errors as lme4.
5) I've tried to transform the data myself with my own function and fit a linear model on the transformed data. However transforming results in Inf, -Inf and NaN for proportion of 1, 0.5 and <0.5 respectively and a linear model cannot handle Inf, -Inf and NaN. In addition, I'm not sure if this is an appropriate way to do my analyses because the variance on the data would be estimated as an additional parameter in a linear model while in logistic regression the variance is prop(1-prop) and is not estimated as an additional parameter.
6) Last I tried using glmmPQL but this does not reach convergence and I'm not sure if it's suitable for my analyses.
I would appreciate any help on applying a halflogit function to lmer.
Kind regards,
Kathleen
---------------------------------------------------------------------------
Kathleen Vancleef
PhD student
Laboratory of Experimental Psychology
University of Leuven
Tiensestraat 102 box 3711
B-3000 Leuven
Belgium
Tel: +32 (0)16/32.62.83
Email: Kathleen.Vancleef at ppw.kuleuven.be
http://www.gestaltrevision.be
More information about the R-help
mailing list