[R-sig-ME] 3 bugs in lmer
Douglas Bates
bates at stat.wisc.edu
Fri Sep 26 23:56:42 CEST 2008
On Thu, Sep 25, 2008 at 12:51 PM, Terry Therneau <therneau at mayo.edu> wrote:
> These are of the better-user-interface variety; features that are common to all
> the other model functions that I know of, but seem to be missing in lmer. One
> could argue that "shortfalls" would be a better label than "bugs".
>> version
> _
> platform sparc-sun-solaris2.10
> arch sparc
> os solaris2.10
> system sparc, solaris2.10
> status
> major 2
> minor 7.1
> year 2008
> month 06
> day 23
> svn rev 45970
> language R
> version.string R version 2.7.1 (2008-06-23)
> #1 - subsets
> fit1c <- lmer(il2 ~ (1|plate/id), data=allsets, subset=(il2>0) )
> fit1d <- lmer(il2 ~ (1|plate/id), data=allsets, subset= il2>0)
> I always put subset arguments in parenthesis for readability. (In the actual
> code that sparked this message, the subset expression was longer, and it helps
> more). Both models fit ok, but the first one doesn't print.
Turns out that was a problem with the utility asOneSidedFormula in the
"stats" package. I'm not sure why I used in the print method but the
problem is fixed in lme4_0.999375-27 which I will release over the
weekend if the tests run on R-forge don't show further problems. At
least the problem is fixed in the examples that I concocted. You did
not provide a reproducible example.
>> all.equal(coef(fit1c), coef(fit1d)
> [1] TRUE
>> print(fit1c)
> Linear mixed model fit by REML
> Formula: il2 ~ (1 | plate/id)
> Data: allsets
> Error in sprintf(gettext(fmt, domain = domain), ...) :
> object "x" not found
>
> #2 coefficients
> Like everyone else, I almost never type the full "coefficients" method and use
> the "coef" abbreviation. But I happened to do so today:
I am somewhat reluctant to add a coefficients method as a synonym for
coef because the coef method does not return what is typically meant
by "the coefficients" in a model based on a linear predictor. The
fixef method returns the fixed effects estimates and the ranef method
returns the conditional modes of the random effects. The current coef
method is a sort of a hodge podge that returns "subject-specific"
coefficients. In a way I wish that I hadn't created that method
because, although those values may be meaningful for specific kinds of
mixed models, it is not clear what the result should be in general
models. I would be more amenable to the suggestion that I remove the
coef method altogether than to the suggestion that I add
"coefficients" as a synonym. I don't feel that the coefficients <->
coef equivalence should hold in this case because the coef method
isn't what it means in other model forms.
>> coefficients(fit1d)
> Error in object$coefficients : $ operator not defined for this S4 class
> #3 character variables
> I almost never, ever, turn a character variable into a factor. We can argue
> about the wisdom of my choice another day, but you'll likely not sway me.
> Character variables work fine in model statements, being treated as factors; for
> models other than lmer.
> Here is a fit on my original data set, before I modified it for the benefit of
> lmer:
The description of the formula argument in the man page for lmer states:
\item{formula}{a two-sided linear formula object describing the
fixed-effects part of the model, with the response on the left of a
\code{~} operator and the terms, separated by \code{+} operators, on
the right. The vertical bar character \code{"|"} separates an
expression for a model matrix and a grouping factor.}
Notice how the description ends - the expression on the right hand
side of the vertical bar should be a factor. I am having a little
difficulty responding in a collegial manner to a bug report of the
form "I don't like factors so you have to rewrite your code to handle
grouping factors that aren't factors."
It appears that the problem here is in the expansion of the / but
again I can't really tell because you didn't include a reproducible
example.
Do you encounter problems if you write the model as
lmer(il2 ~ (1|plate) + (1|plate:id), data=data2)
which is an equivalent specification? If so then the underlying
problem is with the interpretation of the ':' operator. Within a
model formula it should be interpreted as an interaction - i.e. the
definition that applies when the arguments are factors. The error
messages you quote indicate that it is being used as the sequence
operator and that would come from numeric values, not character
values.
I suppose that I could, grudgingly, add code that would coerce the
expression on the right hand side of the | to a factor but I must
admit that I am very tempted instead to add code that would throw an
error with the error message
"The grouping factor must be a factor - and this means you too, Terry!"
Now that I have gotten that off my chest, :-) could you check if the
underlying problem is the ":" operator?
>> fit1b <- lmer(il2 ~ (1|plate/id), data=data2)
> Error in id:plate : NA/NaN argument
> In addition: Warning messages:
> 1: In id:plate :
> numerical expression has 3184 elements: only the first used
> 2: In id:plate :
> numerical expression has 3184 elements: only the first used
> 3: In inherits(x, "factor") : NAs introduced by coercion
>
> By the way, there are no odd or missing values in plate or id.
>
> Oops, I just realized the above is wrong: plate is an integer. So try again:
>
>> table(data2$plate)
>
> 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26
> 80 80 70 80 80 80 46 78 80 80 80 80 80 64 80 80 80 80 80 80 76 80 80 80 80 80
> 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41
> 80 70 80 80 80 80 80 80 80 80 64 80 80 80 80
The output from str(data2) would have been more informative.
>> fit1b <- lmer(il2 ~ (1|factor(plate)/id), data=data2)
> Error in inherits(x, "factor") : object "plate" not found
That one is a problem and it is a nontrivial problem related to what
Thomas Lumley has called "the standard, non-standard evaluation
model". The call to model.frame in lmer is somehow not getting the
right environment in the formula under certain circumstances. It is
very difficult to debug that one.
>> data2$plate2 <- factor(data2$plate)
>> fit1b <- lmer(il2 ~ (1|plate2/id), data=data2)
>
> Error in id:plate2 : NA/NaN argument
> In addition: Warning messages:
> 1: In id:plate2 :
> numerical expression has 3184 elements: only the first used
> 2: In id:plate2 :
> numerical expression has 3184 elements: only the first used
> 3: In inherits(x, "factor") : NAs introduced by coercion
>
>> data2$id2 <- factor(data2$id)
>> fit1b <- lmer(il2 ~ (1|plate2/id2), data=data2)
>> # Finally works
>
>
> Terry Therneau
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
More information about the R-sig-mixed-models
mailing list