[R-sig-ME] Likelihood Ratio
Iasonas Lamprianou
lamprianou at yahoo.com
Sun Apr 4 17:49:06 CEST 2010
Dear all,
the good old advice of using
sm1 <- mcmcsamp(m3_06, 10)
HPDinterval(sm1)
to evaluate the confidence intervals of the random variance do not seem to work any more. What is the new acceptable advice, other than running simulations with the likelihood ratio test? Or maybe there is no additional/new advice on this?
jason
Dr. Iasonas Lamprianou
Assistant Professor (Educational Research and Evaluation)
Department of Education Sciences
European University-Cyprus
P.O. Box 22006
1516 Nicosia
Cyprus
Tel.: +357-22-713178
Fax: +357-22-590539
Honorary Research Fellow
Department of Education
The University of Manchester
Oxford Road, Manchester M13 9PL, UK
Tel. 0044 161 275 3485
iasonas.lamprianou at manchester.ac.uk
--- On Sun, 4/4/10, r-sig-mixed-models-request at r-project.org <r-sig-mixed-models-request at r-project.org> wrote:
> From: r-sig-mixed-models-request at r-project.org <r-sig-mixed-models-request at r-project.org>
> Subject: R-sig-mixed-models Digest, Vol 40, Issue 6
> To: r-sig-mixed-models at r-project.org
> Date: Sunday, 4 April, 2010, 11:00
> Send R-sig-mixed-models mailing list
> submissions to
> r-sig-mixed-models at r-project.org
>
> To subscribe or unsubscribe via the World Wide Web, visit
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
> or, via email, send a message with subject or body 'help'
> to
> r-sig-mixed-models-request at r-project.org
>
> You can reach the person managing the list at
> r-sig-mixed-models-owner at r-project.org
>
> When replying, please edit your Subject line so it is more
> specific
> than "Re: Contents of R-sig-mixed-models digest..."
>
>
> Today's Topics:
>
> 1. Re: lme and prediction intervals (D
> Chaws)
>
>
> ----------------------------------------------------------------------
>
> Message: 1
> Date: Sun, 4 Apr 2010 01:12:52 -0400
> From: D Chaws <cat.dev.urandom at gmail.com>
> To: r-sig-mixed-models at r-project.org
> Subject: Re: [R-sig-ME] lme and prediction intervals
> Message-ID:
> <t2ydcf23fb81004032212p656253e4v42234d8a5b90a669 at mail.gmail.com>
> Content-Type: text/plain; charset=ISO-8859-1
>
> Ok, issue solved for the most straightforward random
> effects cases.
> Not sure about nested random effects or more complex
> cases.
>
> vcov(fm1) will give you the appropriate covariance matrix
> that can
> then be applied to the design matrix to get SEs.
>
> Brian S. Yandell has a nice version of lsmean for his pda
> package that
> handles prediction SEs
> well, and which appear match very closely to those computed
> from SAS.
>
> The link is:
>
> http://r-forge.r-project.org/plugins/scmsvn/viewcvs.php/*checkout*/pkg/R/lsmean.R?rev=2&root=pda
>
> -- DC
>
> On Fri, Feb 26, 2010 at 12:32 AM, D Chaws <cat.dev.urandom at gmail.com>
> wrote:
> > And the saga continues...
> >
> > I checked on the SAS website to see how they compute
> standard errors
> > for predictions in LSMEANS, just for giggles.
> According to
> > http://support.sas.com/documentation/cdl/en/statug/63033/HTML/default/statug_mixed_sect014.htm
> >
> > SE = L(X'V-1X)-L'
> >
> > or in R terms
> >
> > SE = Designmat %*%
> solve(formXtViX(extract.lme.cov2(fm1, Orthodont),
> > extract.lmeDesign(fm1)$X), t(Designmat))
> >
> > which is, yes...you guessed it, the exact same thing
> as Designmat %*%
> > fm1$varFix %*% t(Designmat).
> >
> > And of course, the SE's don't match those in
> LSMEANS. Any thoughts?
> > I can't believe a solution to this has never come
> up. I don't care
> > about replicating LSMEANS, I just want to be confident
> that the SE's I
> > present are accurate and meaningful. Since they
> currently leave out
> > random-effects variance, I suspect something important
> is missing.
> >
> > -- DC
> >
> > On Thu, Feb 18, 2010 at 12:25 PM, D Chaws <cat.dev.urandom at gmail.com>
> wrote:
> >> Dear lme(r) users,
> >>
> >> I have seen this issue come up many times over the
> years, but haven't
> >> come across an answer as of yet.
> >> Take for example the growth model:
> >>
> >> fm1 <- lme(distance ~ age*Sex, random = ~ 1 +
> age | Subject, data = Orthodont)
> >>
> >> Distance increases with age, but more so in
> Males.
> >>
> >> I want to obtain the model predicted values of
> distance at each age
> >> (c(8,10,12,14)) for males and female separately to
> explore this
> >> interaction.
> >>
> >> So,
> >>
> >> newdat <- expand.grid(age=c(8,10,12,14),
> Sex=c("Male","Female"))
> >> newdat$pred <- predict(fm1, newdat, level = 0)
> >>
> >> R# newdat
> >> age Sex pred
> >> 1 8 Male 22.62
> >> 2 10 Male 24.18
> >> 3 12 Male 25.75
> >> 4 14 Male 27.32
> >> 5 8 Female 21.21
> >> 6 10 Female 22.17
> >> 7 12 Female 23.13
> >> 8 14 Female 24.09
> >>
> >> Yup, males have a steeper increase with age than
> females.
> >>
> >> The question is, how to go about getting
> prediction intervals around
> >> these predictions. It seems reasonable to
> need to know
> >> the precision of these predictions, and of course
> most journals
> >> require the reporting of error bars etc...
> However, predict.lme
> >> doesn't
> >> have a se.fit or intervals argument.
> >>
> >> The only answer I have found at the moment is to
> use the design matrix
> >> and $varFix from the model.
> >>
> >> So,
> >>
> >> Designmat <-
> model.matrix(eval(eval(fm1$call$fixed)[-2]), newdat[-3])
> >> R# Designmat
> >> (Intercept) age SexFemale age:SexFemale
> >> 1
> 1 8
> 0
> 0
> >> 2
> 1 10
> 0
> 0
> >> 3
> 1 12
> 0
> 0
> >> 4
> 1 14
> 0
> 0
> >> 5
> 1 8
> 1
> 8
> >> 6
> 1 10
> 1
> 10
> >> 7
> 1 12
> 1
> 12
> >> 8
> 1 14
> 1
> 14
> >>
> >> newdat$SE <- sqrt(diag(Designmat %*% fm1$varFix
> %*% t(Designmat)))
> >> R# newdat
> >> age Sex pred
> SE
> >> 1 8 Male 22.62
> 0.5265
> >> 2 10 Male 24.18 0.4848
> >> 3 12 Male 25.75 0.5021
> >> 4 14 Male 27.32 0.5730
> >> 5 8 Female 21.21 0.6350
> >> 6 10 Female 22.17 0.5847
> >> 7 12 Female 23.13 0.6056
> >> 8 14 Female 24.09 0.6910
> >>
> >> Are these true pointwise prediction intervals?
> >>
> >> Any help would be greatly appreciated. I
> refuse to use SAS for this!
> >>
> >> -- D. Chaws
> >>
> >
>
>
>
> ------------------------------
>
> _______________________________________________
> R-sig-mixed-models mailing list
> R-sig-mixed-models at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
>
> End of R-sig-mixed-models Digest, Vol 40, Issue 6
> *************************************************
>
More information about the R-sig-mixed-models
mailing list