[R] lme questions

[Christian Hennig]

Dear list,

the following should fit the model
log(PEPC.Wert)=fKTemp+fHerk+interaction(fKTemp,fHerk)+fMub+error,
where fKTemp, fHerk are fixed effect factors and fMub is a random effects
factor nested in fHerk (values are different fopr different values of
fHerk).

logpepcr1 <- lme(log(PEPC.Wert) ~ fKTemp*fHerk, random= ~ 1 | fMub, na.action=na.omit)
logpepcr0 <- lme(log(PEPC.Wert) ~ fKTemp*fHerk,na.action=na.omit)

Error in getGroups.data.frame(dataMix, groups) :
Invalid formula for groups

logpepcr0 <- lm(log(PEPC.Wert) ~ fKTemp*fHerk,na.action=na.omit)

...works.

Is the following correct for testing the significance of the
random factor?

anova.lme(logpepcr0,logpepcr1)

Model df AIC BIC logLik Test L.Ratio p-value
logpepcr0 1 10 3977.354 4030.318 -1978.677
logpepcr1 2 11 3885.392 3943.653 -1931.696 1 vs 2 93.96144 <.0001

Where can I get an estimator for the variance of the random factor fMub?

Is there a way to specify an interaction between fMub and fKTemp?

And why does this happen?

plot(logpepcr1)

Error in do.call(plotFun, as.list(args)) :
couldn't find function "xyplot"

(plot(logpepcr0) works; I use R 1.3.0 on Solaris.)

Thank you for help.
Christian

***********************************************************************
Christian Hennig
University of Hamburg, Faculty of Mathematics - SPST/ZMS
(Schwerpunkt Mathematische Statistik und Stochastische Prozesse,
Zentrum fuer Modellierung und Simulation)
Bundesstrasse 55, D-20146 Hamburg, Germany
Tel: x40/42838 4907, privat x40/631 62 79
hennig at math.uni-hamburg.de, http://www.math.uni-hamburg.de/home/hennig/
#######################################################################
ich empfehle www.boag.de


-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._

[Prof Brian Ripley]

Post by Christian Hennig
Dear list,
the following should fit the model
log(PEPC.Wert)=fKTemp+fHerk+interaction(fKTemp,fHerk)+fMub+error,
where fKTemp, fHerk are fixed effect factors and fMub is a random effects
factor nested in fHerk (values are different fopr different values of
fHerk).

logpepcr1 <- lme(log(PEPC.Wert) ~ fKTemp*fHerk, random= ~ 1 | fMub, na.action=na.omit)
logpepcr0 <- lme(log(PEPC.Wert) ~ fKTemp*fHerk,na.action=na.omit)

Invalid formula for groups

logpepcr0 <- lm(log(PEPC.Wert) ~ fKTemp*fHerk,na.action=na.omit)

...works.
Is the following correct for testing the significance of the
random factor?

anova.lme(logpepcr0,logpepcr1)

Model df AIC BIC logLik Test L.Ratio p-value
logpepcr0 1 10 3977.354 4030.318 -1978.677
logpepcr1 2 11 3885.392 3943.653 -1931.696 1 vs 2 93.96144 <.0001

No. You need to fit the lme model by maximum likelihood, not the default
REML, to be comparable with the lm fit.

Post by Christian Hennig
Where can I get an estimator for the variance of the random factor fMub?

The standard error is in the summary() output.

Post by Christian Hennig
Is there a way to specify an interaction between fMub and fKTemp?
And why does this happen?

plot(logpepcr1)

couldn't find function "xyplot"

That's a legacy of the S origins: xyplot is in package lattice in the
Devel section on CRAN (needs grid too) and I don't know if it
yet suffices for plot.lme.

Post by Christian Hennig
(plot(logpepcr0) works; I use R 1.3.0 on Solaris.)

That's using plot.lm.
--
Brian D. Ripley, ripley at stats.ox.ac.uk
Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel: [phone] (self)
1 South Parks Road, [phone] (secr)
Oxford OX1 3TG, UK Fax: [phone]

-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._

[Saikat DebRoy]

Post by Christian Hennig
And why does this happen?

plot(logpepcr1)

couldn't find function "xyplot"

Brian> That's a legacy of the S origins: xyplot is in package lattice in the
Brian> Devel section on CRAN (needs grid too) and I don't know if it
Brian> yet suffices for plot.lme.

It should work with nlme 3.1-17 (I think it is already on
CRAN). Deepayan Sarkar made some changes in the nlme plot methods to
make it work with lattice.

Saikat
--
Department of Statistics Email: saikat at stat.wisc.edu
University of Wisconsin - Madison Phone: [phone]
1210 West Dayton Street Fax: [phone]
Madison, WI 53706-1685
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._