[R] plotting multilevel / lme lines

[Deepayan Sarkar]

On 4/7/07, Rense Nieuwenhuis <r.nieuwenhuis at student.ru.nl> wrote:
> Dear Chuck,
>
> thanks for your reply. I think the graphs you showed me looks
> beautifully. The only problem I see with it, is that the regression-
> lines are potentially extended beyond the original data. So, I
> changed some of your suggested code.
>
> The weird thing is, when I try to implement the working code into a
> function, it stops working. I can't figure out why, but it seems like
> the called variables are not transferred properly to the xyplot-
> function within the overall function.
>
> I tried the following:
>
> # using the fm2-estimated model
>
> library(nlme)
> library(lattice)
> fm2 <- lme(distance ~ poly(age, 2) * Sex,
>                        data = Orthodont, random = ~ 1)
>
>
> newdat <- fm2$data
>
> newdat$newpred <- predict(fm2, newdat, level = 0)
>
> xyplot(newpred ~ age, groups=Sex, data = newdat,        panel = function(x,
> y, ...)
>                         {
>                         panel.superpose(x, y, type="l", ...)
>                 }
>             )
>
> #This works properly and is basically a simplified version of what
> you suggested, but it doesn't extent the model beyond the data it was
> based on. The following does not work though:
>
> # now inside function
>
> vis <- function(model, pred, grp, lvl)
>         {
>                 newdat <- model$data
>
>                 newdat$newpred <- predict(model, newdat, level = lvl)
>
>                 xyplot(newpred ~ pred,
>                 groups=grp,
>                 data = newdat,
>                 panel = function(x, y, ...)
>                         {
>                                 panel.superpose(x, y, type="l", ...)
>                         }
>          )
>         }

Formulas don't work that way. You might try this instead:

vis <- function(model, pred, grp, lvl)
{
    newdat <- model$data

    newdat$newpred <- predict(model, newdat, level = lvl)

    xyplot(as.formula(sprintf("newpred ~ %s", pred)),
           groups=newdat[[grp]],
           data = newdat,
           type = "l")
}

-Deepayan


>
>  > vis(fm2,"age","Sex",0)
> Error in !all.na : invalid argument type
> In addition: Warning messages:
> 1: NAs introduced by coercion
> 2: NAs introduced by coercion
>  > vis(fm2,age,Sex,0)
> Error in eval(expr, envir, enclos) : object "Sex" not found
>
>
> I was wondering what I am missing here. Has it something to do with
> the way a lattice function behaves within a function, or am I missing
> something else here?
>
> I hope someone can help me a bit further here.
>
> Greetings,
>
> Rense Nieuwenhuis
>
>
>
> On Apr 6, 2007, at 14:24 , Chuck Cleland wrote:
>
> > Rense Nieuwenhuis wrote:
> >> Dear expeRts,
> >>
> >> I am trying to plot a lme-object {package nlme) in such a way, that
> >> on a selected level the x-axis represents the value on a selected
> >> predictor and the y-axis represents the predicted-outcome variable.
> >> The graphs would than consist of several lines that each represent
> >> one group. I can't find such a plotting function.
> >>
> >> I could write such a function myself, based on ranef() and fixef(),
> >> but it would be a waste of time if such a function would already
> >> exist.
> >>
> >> Does any of you such a function?
> >
> >   I don't know of a single function with an lme object as argument,
> > but
> > for what I think you have in mind, here is how you might go about it:
> >
> > library(nlme)
> >
> > fm2 <- lme(distance ~ poly(age, 2) * Sex,
> >                       data = Orthodont, random = ~ 1)
> >
> > newdat <- expand.grid(age = 8:14, Sex = c("Male","Female"))
> >
> > newdat$PREDDIST <- predict(fm2, newdat, level = 0)
> >
> > library(lattice)
> >
> > xyplot(PREDDIST ~ age, groups=Sex, ylab="Model Predicted Distance",
> >            data = newdat, xlab="Age",
> >            panel = function(x, y, ...){
> >                    panel.grid(h=6,v=6)
> >                    panel.superpose(x, y, type="l", ...)},
> >                    main="Orthodont Growth Model",
> >                    key = simpleKey(levels(newdat$Sex),
> >                                    lines=TRUE, points=FALSE)
> >                    )
> >
> >> Regards,
> >>
> >> Rense Nieuwenhuis
> >>
> >> ______________________________________________
> >> R-help at stat.math.ethz.ch mailing list
> >> https://stat.ethz.ch/mailman/listinfo/r-help
> >> PLEASE do read the posting guide http://www.R-project.org/posting-
> >> guide.html
> >> and provide commented, minimal, self-contained, reproducible code.
> >
> > --
> > Chuck Cleland, Ph.D.
> > NDRI, Inc.
> > 71 West 23rd Street, 8th floor
> > New York, NY 10010
> > tel: (212) 845-4495 (Tu, Th)
> > tel: (732) 512-0171 (M, W, F)
> > fax: (917) 438-0894
> >
>
>
>         [[alternative HTML version deleted]]
>
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>