An embedded and charset-unspecified text was scrubbed... Name: not available URL: <https://stat.ethz.ch/pipermail/r-help/attachments/20090223/eb26d61f/attachment-0001.pl>
lmer, estimation of p-values and mcmcsamp
3 messages · Maureen Ryan, Mark Difford, levyofi
Hi Maureen, There is a specialist mailing list for mixed-models, where you are likely to get a well-informed answer to this, quite possibly from Prof. Bates himself (the author of lmer). https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models Regards, Mark. PS: R Hompage -> Mailing Lists for other specialist list.
Maureen Ryan wrote:
(To the list moderator: I just subscribed to the list. Apologies for not
having done so longer before trying to post.)
Hi all,
I am currently using lmer to analyze data from an experiment with a
single fixed factor (treatment, 6 levels) and a single random factor
(block). I've been trying to follow the online guidance for estimating
p-values for parameter estimates on these and other help threads:
https://stat.ethz.ch/pipermail/r-help/2006-May/094765.html
http://pidgin.ucsd.edu/pipermail/r-lang/2007-August/000057.html
but have not yet been successful. This is my code:
lnmass <- lmer(log.mass ~ treatment + (1|block), data=exp1)
summary(lnmass)
samp <- rnorm(n=10000)
mcmcpvalue <- function(samp)
{std <- backsolve(chol(var(samp)),
cbind(0,t(samp)) - colMeans(samp),
transpose = TRUE)
sqdist <- colSums(std*std)
sum(sqdist[-1] > sqdist[1]/nrow(samp) }
markov1 <- mcmcsamp(lnmass, 10000)
HPDinterval(markov1)
mcmcpvalue(as.matrix(markov1[,1]))
mcmcpvalue(as.matrix(markov1[,2]))
mcmcpvalue(as.matrix(markov1[,3]))
mcmcpvalue(as.matrix(markov1[,4]))
mcmcpvalue(as.matrix(markov1[,5]))
mcmcpvalue(as.matrix(markov1[,6]))
The first time I tried it HPDinterval generated CIs for each treatment
level, etc. but after each of the last 6 lines of code I received the
following error messages:
Error in markov[,1] : object is not subsettable
Error in as.matrix(markov1[,1]) : error in evaluating the argument 'x' in
selecting a method for function 'as.matrix'
Error in chol(val(samp)) : error in evaluating the argument 'x' in
selecting
a method for function 'chol'
When I came back to it today HPDinterval(markov1) also no longer worked
and
gave me the following error message:
Error in UseMethod("HPDinterval"): no applicable method for "HPDinterval"
I imagine at least some of this has to do with the fact that I don't fully
understand what exactly the mcmcpvalue function is doing. Any help on
either problem would be very much appreciated.
Thanks, Maureen
*****
Maureen Ryan
Department of Evolution and Ecology
University of California Davis
One Shields Avenue
Davis, CA 95616
(530)304-2266
[[alternative HTML version deleted]]
______________________________________________ R-help at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
View this message in context: http://www.nabble.com/lmer%2C-estimation-of-p-values-and-mcmcsamp-tp22177307p22178938.html Sent from the R help mailing list archive at Nabble.com.
For getting the p-values from a lmer object you need to use the pvals.fnc() function. do the following: library(lme4) library(languageR) lnmass <- lmer(log.mass ~ treatment + (1|block), data=exp1) pvals.fnc(lnmass, nsim=10000) See the documentation for futher details. Best regards, Ofir.
Maureen Ryan wrote:
(To the list moderator: I just subscribed to the list. Apologies for not
having done so longer before trying to post.)
Hi all,
I am currently using lmer to analyze data from an experiment with a
single fixed factor (treatment, 6 levels) and a single random factor
(block). I've been trying to follow the online guidance for estimating
p-values for parameter estimates on these and other help threads:
https://stat.ethz.ch/pipermail/r-help/2006-May/094765.html
http://pidgin.ucsd.edu/pipermail/r-lang/2007-August/000057.html
but have not yet been successful. This is my code:
lnmass <- lmer(log.mass ~ treatment + (1|block), data=exp1)
summary(lnmass)
samp <- rnorm(n=10000)
mcmcpvalue <- function(samp)
{std <- backsolve(chol(var(samp)),
cbind(0,t(samp)) - colMeans(samp),
transpose = TRUE)
sqdist <- colSums(std*std)
sum(sqdist[-1] > sqdist[1]/nrow(samp) }
markov1 <- mcmcsamp(lnmass, 10000)
HPDinterval(markov1)
mcmcpvalue(as.matrix(markov1[,1]))
mcmcpvalue(as.matrix(markov1[,2]))
mcmcpvalue(as.matrix(markov1[,3]))
mcmcpvalue(as.matrix(markov1[,4]))
mcmcpvalue(as.matrix(markov1[,5]))
mcmcpvalue(as.matrix(markov1[,6]))
The first time I tried it HPDinterval generated CIs for each treatment
level, etc. but after each of the last 6 lines of code I received the
following error messages:
Error in markov[,1] : object is not subsettable
Error in as.matrix(markov1[,1]) : error in evaluating the argument 'x' in
selecting a method for function 'as.matrix'
Error in chol(val(samp)) : error in evaluating the argument 'x' in
selecting
a method for function 'chol'
When I came back to it today HPDinterval(markov1) also no longer worked
and
gave me the following error message:
Error in UseMethod("HPDinterval"): no applicable method for "HPDinterval"
I imagine at least some of this has to do with the fact that I don't fully
understand what exactly the mcmcpvalue function is doing. Any help on
either problem would be very much appreciated.
Thanks, Maureen
*****
Maureen Ryan
Department of Evolution and Ecology
University of California Davis
One Shields Avenue
Davis, CA 95616
(530)304-2266
[[alternative HTML version deleted]]
______________________________________________ R-help at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
View this message in context: http://www.nabble.com/lmer%2C-estimation-of-p-values-and-mcmcsamp-tp22177307p22180956.html Sent from the R help mailing list archive at Nabble.com.