(Re: [R] lme in R-2.0.0: Problem with lmeControl) and parameter specification
Pavel Khomski
pkhomski at wiwi.uni-bielefeld.de
Wed Dec 1 17:49:40 CET 2004
Hello!
Thanks a lot to Douglas Bates for your advice.
concerning the lme(...) function i wanted to put four other questions.
1. in the specification of initial values in the pdMat-constructor i
probably define a standard deviation (sigma_b) and not a variance
(sigma_b^2). For instance
in the Rail example in
Pinheiro/Bates on p.81 if i specify a random parameter as
random=pdIdent(value=lambda<-diag(1000,1),form=~1), (in S-plus),
then the call to the lme(...)
with just 0 iterations produces:
>
fm1Rail.lme<-lme(travel~1,data=Rail,random=pdIdent(value=lambda<-diag(1000,1),form=~1),control=list(msMaxIter=0,msVerbose=TRUE,niterEM=0))
Iteration: 0 , 1 function calls, F= 66.37359
Parameters:
[1] -3.453878
Warning messages:
ITERATION LIMIT REACHED WITHOUT OTHER CONVERGENCE in: ms( ~ -
logLik(lmeSt,
lmePars), start = list(lmePars = c(coef(lmeSt))), ....
if i now print out an estimated std.dev. for sigma_b i get:
> (fm1Rail.lme$sigma)*exp(unlist(fm1Rail.lme$modelStruct))
reStruct.Rail
107.6767
> lambda
[,1]
[1,] 1000
so that the estimated variance would be 107.6767 ^2 = 11594.27 what is
much grater then 1000. But hier we know that with iterations the
value of variance will reduce ( and at the convergence the StdDev is
24.80547 )
so i think that lambda=1000 is specified equal to sigma_b as initial value.
2. What is the meaning for 0-Iteration?
3. are the parameters fixed=beta and random=sigma being
calculated (just on time) only after all iterations have run, or they
also be updated
at every iteration with new value of teta ? if the latter, how
can i get them for each run ?
4. Is it in principle possible to hold a variance components
parameter, say sigma_b, as in Rail-example, fixed (on specified
value) through all the
iteration steps (without changing it) and only optimize for
teta=log(sigma_b/sigma_epsilon) with fixed known value of sigma_b ?
how can it be done ?
Thank you for replay
More information about the R-help
mailing list