same old question - lme4 and p-values
Hi Doug:
Perhaps I should clarify. The summary of a fitted lmer model does not provide p-values because I don't know how to calculate them in an acceptable way, not because I am philosophically opposed to them. The estimates and the approximate standard errors can be readily calculated as can their ratio. The problem is determining the appropriate reference distribution for that ratio from which to calculate a p-value. In fixed-effects models (under the "usual" assumptions) that ratio is distributed as a T with a certain number of degrees of freedom. For mixed models it is not clear exactly what distribution it has - except in certain cases of completely balanced data sets (i.e. the sort of data sets that occur in text books). At one time I used a T distribution and an upper bound on the degrees of freedom but I was persuaded that providing p-values that could be strongly "anti-conservative" is worse than not providing any.
Now I understand the situation better and am in agreement that this is clearly the right solution at this point.
The approach that I feel is most likely to be successful in summarizing these models is first to obtain the REML or ML estimates of the parameters then to run a Markov chain Monte Carlo sampler to assess the variability in the parameters (or, if you prefer, the variability in the parameter estimators). (Note: I am not advocating using MCMC to obtain the estimates, I suggest MCMC for assessing the variability.)
I'm a little confused as to what is the Monte Carlo part of this scenario? If you perform REML or ML, theoretically it should always converge to the REML/ML estimates (unless you have a flat or multimodal likelihood which each produce other problems). I understand you are fixing the parameter estimates of something at the REML/ML estimates, but what is the random component? Of course, I could always stop being lazy and just look at the source... ;^)
The current version of the mcmcsamp function suffers from the practical problem that it gets stuck at near-zero values of variance components. There are some approaches to dealing with that. Over the weekend I thought that I had a devastatingly simple way of dealing with such cases until I reflected on it a bit more and realized that it would require a division by zero. Other than that, it was a good idea.
At least the variance estimates were not negative... ;^) Thanks!! Dave H -- David Henderson, Ph.D. Director of Community REvolution Computing 1100 Dexter Avenue North, Suite 250 206-577-4778 x3203 DNADave at Revolution-Computing.Com http://www.revolution-computing.com