Ramzi Nahhas
2011-Jul-07 17:09 UTC
[R] How do I overlay two trellis plots of lme fitted lines produced by plot.augPred?
Hello,
I want to use lme to fit two (or more) models, and then compare the fits
on each individual. I know how to write my own code to do this (for each
individual, plot the raw data, followed by lines() to plot each fitted
curve) but I would like to use plot(augPred(... as it produces a nice
trellis plot. I thought I could do this with par(new=T) but it does not
seem to work. trellis.par.set() does not have a component called "new"
so I cannot use that to imitate what I would do with par.
Does anyone know how to overlay two or more plots produced by plot.augPred?
Thank you!
Ramzi
# Here is code that reproduces my problem:
tmpdat structure(list(ptno = c(4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L,
8L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 10L, 10L, 10L, 10L, 10L, 10L,
10L, 10L, 10L, 10L, 20L, 20L, 20L, 20L, 20L, 20L, 20L, 20L, 20L,
20L, 20L, 20L, 20L, 23L, 23L, 23L, 23L, 23L, 23L, 23L, 23L, 23L,
23L, 24L, 24L, 24L, 24L, 24L, 24L, 24L, 24L, 24L, 24L), age =
c(6.9897330595,
9.0075290897, 10.009582478, 12.005475702, 13.013004791, 13.990417522,
16.101300479, 17.062286105, 17.995893224, 20.996577687, 8.0027378508,
8.9965776865, 10.006844627, 15.085557837, 15.994524298, 16.993839836,
17.984941821, 22.431211499, 7.0034223135, 8, 8.9965776865, 11.047227926,
11.997262149, 14.124572211, 15.006160164, 15.986310746, 17.059548255,
18.001368925, 6.9952087611, 8.0027378508, 9.0075290897, 10.001368925,
11.000684463, 11.997262149, 12.988364134, 14.009582478, 14.986995209,
16.002737851, 16.963723477, 17.995893224, 22.702258727, 7.0034223135,
7.9972621492, 9.0047912389, 12.027378508, 12.982888433, 13.995893224,
15.282683094, 16.041067762, 17.029431896, 19.498973306, 7.0061601643,
8.9993155373, 10.015058179, 11.011635866, 12.999315537, 14.006844627,
15.030800821, 16.065708419, 17.015742642, 18.015058179), y =
c(40.632301331,
40.996398926, 44.707000732, 39.288200378, 42.174301147, 50.073101044,
45.052200317, 44.862499237, 53.557498932, 54.667499542, 42.447601318,
43.082500458, 45.340499878, 47.740501404, 50.521499634, 52.097400665,
53.834999084, 53.002498627, 31.815000534, 41.812698364, 45.349998474,
42.350700378, 43.663898468, 43.941699982, 47.007400513, 46.071899414,
48.099998474, 50.319999695, 33.178501129, 38.366100311, 33.740398407,
40.996498108, 45.250499725, 40.964000702, 49.529201508, 46.635799408,
43.848800659, 49.038299561, 59.698799133, 57.349998474, 56.979999542,
37.90530014, 39.726600647, 34.828800201, 37.519298553, 36.681400299,
41.897899628, 48.111301422, 52.746299744, 45.144901276, 53.650001526,
42.632099152, 42.175498962, 46.064498901, 48.142799377, 51.856700897,
55.925800323, 51.002101898, 56.268901825, 61.367401123, 57.844799042
)), .Names = c("ptno", "age", "y"), row.names =
c(1L, 3L, 4L,
5L, 6L, 7L, 8L, 9L, 10L, 11L, 25L, 26L, 27L, 28L, 29L, 30L, 31L,
32L, 44L, 45L, 46L, 47L, 48L, 49L, 50L, 51L, 52L, 53L, 77L, 78L,
79L, 80L, 81L, 82L, 83L, 84L, 85L, 86L, 87L, 88L, 89L, 90L, 91L,
92L, 93L, 94L, 95L, 96L, 97L, 98L, 99L, 100L, 101L, 102L, 103L,
104L, 105L, 106L, 107L, 108L, 109L), class = "data.frame")
library(nlme)
grpdat = groupedData(y ~ age | ptno,
data = tmpdat,
FUN = mean,
labels = list(x="Age", y="Y"))
fit.lme = lme(grpdat)
plot(augPred( fit.lme ), ylim=c(20,65))
plot(augPred(update(fit.lme, y ~ log(age))), ylim=c(20,65))
# I want these two plots to be in the same window so I can visually
compare the two fits
# The following does not work:
plot(augPred( fit.lme ), ylim=c(20,65),
col.line="green")
par(new=T)
plot(augPred(update(fit.lme, y ~ log(age))), ylim=c(20,65),
col.line="red")
par(new=F)
# trellis.par.set() does not have a component called "new" so I cannot
use that to imitate what I would do with par
--
---------------------------------------------------------------------
- Ramzi W. Nahhas, Ph.D.
- Assistant Professor
- Lifespan Health Research Center, Department of Community Health
- Boonshoft School of Medicine, Wright State University
- 3171 Research Boulevard, Dayton OH 45420
- Phone: (937) 775-1421 Cell: (937) 269-6797 Fax: (937) 775-1456
- www.med.wright.edu/lhrc/fac/rn.html
Dennis Murphy
2011-Jul-07 21:10 UTC
[R] How do I overlay two trellis plots of lme fitted lines produced by plot.augPred?
Hi: Here's one way, using the latticeExtra package: library(lattice) library(latticeExtra) p1 <- plot(augPred( fit.lme ), ylim=c(20,65)) p2 <- plot(augPred(update(fit.lme, y ~ log(age))), ylim=c(20,65), col.line = 'red') p1 + p2 HTH, Dennis On Thu, Jul 7, 2011 at 10:09 AM, Ramzi Nahhas <ramzi.nahhas at wright.edu> wrote:> Hello, > > I want to use lme to fit two (or more) models, and then compare the fits on > each individual. I know how to write my own code to do this (for each > individual, plot the raw data, followed by lines() to plot each fitted > curve) but I would like to use plot(augPred(... as it produces a nice > trellis plot. ?I thought I could do this with par(new=T) but it does not > seem to work. trellis.par.set() does not have a component called "new" so I > cannot use that to imitate what I would do with par. > > Does anyone know how to overlay two or more plots produced by plot.augPred? > > Thank you! > > Ramzi > > # Here is code that reproduces my problem: > tmpdat > structure(list(ptno = c(4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, > 8L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 10L, 10L, 10L, 10L, 10L, 10L, > 10L, 10L, 10L, 10L, 20L, 20L, 20L, 20L, 20L, 20L, 20L, 20L, 20L, > 20L, 20L, 20L, 20L, 23L, 23L, 23L, 23L, 23L, 23L, 23L, 23L, 23L, > 23L, 24L, 24L, 24L, 24L, 24L, 24L, 24L, 24L, 24L, 24L), age > c(6.9897330595, > 9.0075290897, 10.009582478, 12.005475702, 13.013004791, 13.990417522, > 16.101300479, 17.062286105, 17.995893224, 20.996577687, 8.0027378508, > 8.9965776865, 10.006844627, 15.085557837, 15.994524298, 16.993839836, > 17.984941821, 22.431211499, 7.0034223135, 8, 8.9965776865, 11.047227926, > 11.997262149, 14.124572211, 15.006160164, 15.986310746, 17.059548255, > 18.001368925, 6.9952087611, 8.0027378508, 9.0075290897, 10.001368925, > 11.000684463, 11.997262149, 12.988364134, 14.009582478, 14.986995209, > 16.002737851, 16.963723477, 17.995893224, 22.702258727, 7.0034223135, > 7.9972621492, 9.0047912389, 12.027378508, 12.982888433, 13.995893224, > 15.282683094, 16.041067762, 17.029431896, 19.498973306, 7.0061601643, > 8.9993155373, 10.015058179, 11.011635866, 12.999315537, 14.006844627, > 15.030800821, 16.065708419, 17.015742642, 18.015058179), y = c(40.632301331, > 40.996398926, 44.707000732, 39.288200378, 42.174301147, 50.073101044, > 45.052200317, 44.862499237, 53.557498932, 54.667499542, 42.447601318, > 43.082500458, 45.340499878, 47.740501404, 50.521499634, 52.097400665, > 53.834999084, 53.002498627, 31.815000534, 41.812698364, 45.349998474, > 42.350700378, 43.663898468, 43.941699982, 47.007400513, 46.071899414, > 48.099998474, 50.319999695, 33.178501129, 38.366100311, 33.740398407, > 40.996498108, 45.250499725, 40.964000702, 49.529201508, 46.635799408, > 43.848800659, 49.038299561, 59.698799133, 57.349998474, 56.979999542, > 37.90530014, 39.726600647, 34.828800201, 37.519298553, 36.681400299, > 41.897899628, 48.111301422, 52.746299744, 45.144901276, 53.650001526, > 42.632099152, 42.175498962, 46.064498901, 48.142799377, 51.856700897, > 55.925800323, 51.002101898, 56.268901825, 61.367401123, 57.844799042 > )), .Names = c("ptno", "age", "y"), row.names = c(1L, 3L, 4L, > 5L, 6L, 7L, 8L, 9L, 10L, 11L, 25L, 26L, 27L, 28L, 29L, 30L, 31L, > 32L, 44L, 45L, 46L, 47L, 48L, 49L, 50L, 51L, 52L, 53L, 77L, 78L, > 79L, 80L, 81L, 82L, 83L, 84L, 85L, 86L, 87L, 88L, 89L, 90L, 91L, > 92L, 93L, 94L, 95L, 96L, 97L, 98L, 99L, 100L, 101L, 102L, 103L, > 104L, 105L, 106L, 107L, 108L, 109L), class = "data.frame") > > library(nlme) > grpdat = groupedData(y ? ? ?~ age | ptno, > ? ? ? ? ? ? ? ? ? ? data ? = tmpdat, > ? ? ? ? ? ? ? ? ? ? FUN ? ?= mean, > ? ? ? ? ? ? ? ? ? ? labels = list(x="Age", y="Y")) > fit.lme = lme(grpdat) > plot(augPred( ? ? ? fit.lme ? ? ? ? ? ? ? ), ylim=c(20,65)) > plot(augPred(update(fit.lme, y ~ log(age))), ylim=c(20,65)) > > # I want these two plots to be in the same window so I can visually compare > the two fits > > # The following does not work: > plot(augPred( ? ? ? fit.lme ? ? ? ? ? ? ? ), ylim=c(20,65), > col.line="green") > par(new=T) > plot(augPred(update(fit.lme, y ~ log(age))), ylim=c(20,65), col.line="red") > par(new=F) > > # trellis.par.set() does not have a component called "new" so I cannot use > that to imitate what I would do with par > > -- > --------------------------------------------------------------------- > - Ramzi W. Nahhas, Ph.D. > - Assistant Professor > - Lifespan Health Research Center, Department of Community Health > - Boonshoft School of Medicine, Wright State University > - 3171 Research Boulevard, Dayton OH 45420 > - Phone: (937) 775-1421 ? Cell: (937) 269-6797 ? Fax: (937) 775-1456 > - www.med.wright.edu/lhrc/fac/rn.html > > ______________________________________________ > R-help at r-project.org 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. >
rwnahhas
2011-Jul-08 12:41 UTC
[R] How do I overlay two trellis plots of lme fitted lines produced by plot.augPred?
Thank you! Exactly what I wanted. Ramzi -- View this message in context: http://r.789695.n4.nabble.com/How-do-I-overlay-two-trellis-plots-of-lme-fitted-lines-produced-by-plot-augPred-tp3652204p3653988.html Sent from the R help mailing list archive at Nabble.com.