[R] FW: error predicting values from the LME

Andrej Kveder andrejk at zrc-sazu.si
Wed Oct 1 11:57:36 CEST 2003


Thanks for the hint. And it's a general one I intend to use more often.
I managed to sort out the problem proceeding with the simpler models first
and getting to more komplex ones. Even the one I presented in my question
suddenly worked.
However I think I found the problem, but I'm unable to solve it. With the
problematic prediction I used the prespecified model using the as.formula.
This is my exapmle:

# first way
m.gr<-as.formula("y ~ v11+v21+v22+v23 | inter")
m.fix<-as.formula("y ~ v11+v21+v22+v23+v21:v11+v22:v11+v23:v11")
m.ran<-as.formula("~ v11 | inter")
d.n.gr.1<-grouped.data(m.gr data=d.n)
m.1<-lme(m.fix, random=m.ran, data=d.n.gr.1, na.action=na.omit)
dnNew<-subset(d.n,is.na(d.n$y),select=c(v11,v21,v22,v23,inter))
p.1<-predict(m.1,dnNew)
# second way
d.n.gr.2<-grouped.data(y ~ v11+v21+v22+v23 | inter data=d.n)
m.2<-lme(y ~ v11+v21+v22+v23+v21:v11+v22:v11+v23:v11, random=~ v11 | inter,
data=d.n.gr.2, na.action=na.omit)
p.2<-predict(m.1,dnNew)

Second way works, while the first way doesn't. I get all the model
specifications from the first one including fitted and residuals, everything
that is derived from the lme object itself. However with the predictions I
get the error.

Error in eval(expr, envir, enclos) : 1 argument passed to "$" which requires
2.

I have saved the 3 objects: both models and the test dataset and attached it
to the message in hope to help solve my problem, however my message was
rejected as being to big (the attached file is 80K). If somebody has an idea
how to solve the problem, I can send the data to him/her.

Is there a solution to the problem?

Thanks for the help so far.
I would appreciate any other suggestions.


Andrej



-----Original Message-----
From: Spencer Graves [mailto:spencer.graves at pdf.com]
Sent: Tuesday, September 30, 2003 3:13 PM
To: Andrej Kveder
Cc: R-Help
Subject: Re: [R] FW: error predicting values from the LME


      Could you dumb it down to a toy example with 4 observations for a
model like "y ~ 1 | inter", 2 observations for each of 2 levels of
"inter"?  If that works, then you can play with the example that works
and the example that doesn't;  this is one of the strategies mentioned
in Poly (1971) How to Solve It (Princeton U. Pr.).  With luck, this will
help you figure out what you need to do to get the answers you want.  If
not, it should help you produce a small toy example that doesn't work,
which you can then send us.  Please include a data.frame call, so
someone can copy your example into R and try it in 2 seconds.  That
should increase the chances that you would get a helpful reply.

      Also, have you read Pinhiero and Bates (2000) Mixed-Effect Models
in S and S-Plus (Springer)?  I've found that book to be indispensible
for using "lme".

      hope this helps.  spencer graves

Andrej Kveder wrote:

>HI all,
>
>I might add some more information in order to possibly solve my problem.
I'm
>really stuck and no obvious solutions do the trick.
>I'm using R 1.7.1 on Windows 2000 with the packages regurarly updated.
>I'm using hypothetical data constructed as a pseudo population conforming
to
>a certain Var-Cov structure.
>I might add that just
>
>
>
>>predict(level2)
>>
>>
>
>works. But when I add the new dataset it doesn't. Following a suggestion I
>even tried refactoring of the grouping variable (inter) after I created the
>subset. It didn't work. I have no other factor variables in the model. I
>really have got no clue what could be wrong.
>
>There is a sample from my data:
>
>
>>dnNew
>>
>>
>Grouped Data: y ~ v11 + v21 + v22 + v23 | inter
>         v11             v21          v22         v23    inter
>4 5.55186635 5.6620022 24.18033 5.003409 1
>13 2.03852426 5.6620022 24.18033 5.003409 1
>15 2.19825772 7.5676798 31.03986 4.746891 2
>16 4.51368278 7.5676798 31.03986 4.746891 2
>18 3.35322702 7.5676798 31.03986 4.746891 2
>19 2.46414346 7.5676798 31.03986 4.746891 2
>20 2.66670834 7.5676798 31.03986 4.746891 2
>
>and this is the model:
>
>
>>level2
>>
>>
>Linear mixed-effects model fit by REML
>Data: d.n.gr.2
>Log-restricted-likelihood: -533.0011
>Fixed: model$fixed
>(Intercept) v11 v21 v22 v23 v11:v21
>3.205519074 0.298941539 -0.017743958 0.016007280 -0.410760471 0.002700954
>v11:v22 v11:v23
>-0.003680952 -0.018005717
>Random effects:
>Formula: ~v11 | inter
>Structure: General positive-definite, Log-Cholesky parametrization
>StdDev Corr
>(Intercept) 0.385620605 (Intr)
>v11 0.003147431 -0.048
>Residual 0.450012367
>Number of Observations: 729
>Number of Groups: 50
>If this give you some more insight to my problem.
>
>I would reallly appreciate any suggestion.
>
>Thanks
>
>Andrej
>
>-----Original Message-----
>From: Andrej Kveder [mailto:andrejk at zrc-sazu.si]
>Sent: Monday, September 29, 2003 7:05 PM
>To: R-Help
>Subject: predicting values from the LME
>
>
>Dear listers,
>
>I experinced a problem prdicting the values using the LME with multilevel
>data.
>I have NA's in my dependent variable and the model is fitted only on the
>completed cases.
>I want to estimate the predicted values for the rest of the data (those
>cases with missing dep. variable)
>I extracted a subset from the original file containing the variables used
in
>the model as well as the second level indicator.
>I used the following command
>
>p<-predict(level2,newdata=d.n.new,level=0:1)
>
>where level2 is my LME model.
>But, I get the following error:
>
>Error in eval(expr, envir, enclos) : 1 argument passed to "$" which
requires
>2.
>
>I tried with omitting the level specification (which is 0 by default) and I
>transformed the new data to be groupedData with no luck.
>
>I have tried the example from the Pinheiro,Bates book and it works - mine
>doesn't. Does anybody have an idea what could be wrong?
>
>Thanks for all the suggestions.
>
>Andrej
>
>_________
>Andrej Kveder, M.A.
>researcher
>Institute of Medical Sciences SRS SASA; Novi trg 2, SI-1000 Ljubljana,
>Slovenia
>phone: +386 1 47 06 440   fax: +386 1 42 61 493
>
>	[[alternative HTML version deleted]]
>
>______________________________________________
>R-help at stat.math.ethz.ch mailing list
>https://www.stat.math.ethz.ch/mailman/listinfo/r-help
>
>




More information about the R-help mailing list