[R] converting a BMDP 8V mixed model to R / nlme

[Peter Dalgaard BSA]

"Wilhelm B. Kloke" <wb at vestein.arb-phys.uni-dortmund.de> writes:

> latm <- read.table("latm.dat" )
> latm.aov <- aov( llat ~ person + retest * direct * flicker * ground
> 	+ Error ( person : ( retest * direct * flicker * ground )))
> 
> summary(latm.aov) at least reports the same F-values as 8V does. Though
> the
> R computations thrash my system (FreeBSD, K6-233, 64MB). I have to
> triple the swapspace before starting, and I wait about an hour for
> either rebooting the system or having a lucky outcome. There seems to be
> something fundamentally wrong.
> 
> It seems I have to move to nmle to get variance components. I did not
> manage to translate my problem to mle language yet.

(That's nlme/lme - (non-)linear mixed effect)

aov() does take rather a long time on problems of this size (2268
observations, 16 error strata). I tried it with simulated data and saw
the the memory footprint go to about 300MB for a while... 

The algorithm in R is not as good at taking advantage of balancedness
between error terms as programs like Genstat and does the orthogonal
decomposition by brute force, leading to linear problems with a huge
number of predictors. Of course, given sufficient knowledge of the
particular model, you can optimize very heavily and everything
reduces to means and sums of squares thereof. That is probably what
BMDP 8V does.

Variance components can be obtained from the estimated stratum error
terms by inverting the linear relation between them, which is a bit
tedious to work out, but I suspect you'll find it even more tedious to
set up the model in lme()... The latter involves explicitly setting up
the block diagonal covariance matrix of random effects, with each
block proportional to the identity matrix. Something like

random=list(person=pdBlocked(
   list(pdIdent(~retest-1),
        ...., 
        pdIdent(~retest:direct:flicker:ground-1))))


-- 
   O__  ---- Peter Dalgaard             Blegdamsvej 3  
  c/ /'_ --- Dept. of Biostatistics     2200 Cph. N   
 (*) \(*) -- University of Copenhagen   Denmark      Ph: (+45) 35327918
~~~~~~~~~~ - (p.dalgaard at biostat.ku.dk)             FAX: (+45) 35327907
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
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
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._