[R] effects package --- add abline to plot

Tue Apr 28 17:33:36 CEST 2009

Dear David,

> -----Original Message-----
> From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org] On
> Behalf Of David Winsemius
> Sent: April-28-09 10:12 AM
> To: Prew, Paul
> Cc: r-help at r-project.org
> Subject: Re: [R] effects package --- add abline to plot
> 
> 
> On Apr 28, 2009, at 12:00 AM, Prew, Paul wrote:
> 
> > Hello,  I am not having success in a simple task.  Using the effects
> > package, I would like to add reference lines at probability values
> > of 0.1 – 0.6 on a plot of the effects.
> 
> I have concerns that you are considering these probabilities. They are
> not going to be probabilities. They are effects.
> 
> >  The plot command works, but following up with an abline command
> > produces the message “plot .new has not been called yet”, and of
> > course the reference lines were not added.
> >
> > Looking through past R help lists, there was a similar request for
> > help --- trying to add an abline but “got the error "plot.new  has
> > not been called yet".
> >
> > The help list reply was
> >
> > “ ?abline: "This function adds one or more straight lines through
> > the
> > current plot.", i.e. the already existing *current plot*.
> >
> > So plot your data (e.g. with plot(x, y)) before adding a regression
> > line.”
> >
> > I interpreted the above to suggest the following ---
> >
> > plot(allEffects(Clean.label),ask=FALSE, alternating = TRUE,
> > 	ylab="Probability of Rating", xlab="City",main="Cleanliness Ratings
> > by City",
> > 	factor.names=FALSE, ticks=c(0.1,0.2,0.3,0.4,0.5,0.6))
> > abline(h=c(0.1,0.2,0.3,0.4,0.5,0.6))
> 
> I do not know why that is happening and you have not provided a
> minimal executable example. The vectorized use of abline does succeed
> in a simper example:
> 
>  > plot(.5,.5)
>  > abline(h=c(0.1,0.2,0.3,0.4,0.5,0.6))
> 
>   .... so the problem may lie in how the effects package completes its
> plot function for this particular object. You ought to provide at a
> minimum the results of str on that object. Perhaps it executes a
> device call and then turns off the device? However I loaded the
> effects package and ran that abline call after the example:
> 
>  > mod.cowles <- glm(volunteer ~ sex + neuroticism*extraversion,
> +     data=Cowles, family=binomial)
>  > eff.cowles <- allEffects(mod.cowles, xlevels=list(neuroticism=0:24,
> +     extraversion=seq(0, 24, 6)))
>  > eff.cowles
> 
> 
> I did not get what I expected, which would have been a single
> horizontal line at 0.4 but rather got four lines roughly at 0.351,
> 0.378, 0.408, 0.439. Even then, I would have expected one more line
> before the upper limits of that plot, which makes me think these four
> lines were the results of arguments 0.3 ,0.4, 0.5, 0.6.   Most R
> plotting is done in the coordinate system rather than with absolute
> coordinates, but perhaps the mixture of base graphics with lattice
> graphis is ht eproblem

The plot() methods in the effects package make lattice graphs which in most instances will have more than one panel. For a binomial GLM, the default is to plot on the scale of the linear predictor (e.g., the logit scale) but to label the response axis on the scale of the response (i.e., the probability scale). To draw a line on the graph, even if you could do it, would require that you translate to the scale of the linear predictor [e.g., for a logit model, log(p/(1 - p))].

> >
> >
> > Error in int_abline(a = a, b = b, h = h, v = v, untf = untf, ...) :
> >  plot.new has not been called yet
> >
> > Less bothersome is the fact that the tick marks weren’t modified to
> > 0.1, 0.2, 0.3, etc.

>From ?plot.eff

"ticks: a two-item list controlling the placement of tick marks on the vertical axis, with elements at and n. If at=NULL (the default), the program attempts to find `nice' locations for the ticks, and the value of n (default, 5) gives the approximate number of tick marks desired; if at is non-NULL, then the value of n is ignored."

> >
> > Further searching brought the panel.abline command to light, but
> > that didn’t produce any results, not even an error message.
> >
> >> plot(allEffects(Clean.label),ask=FALSE, alternating = TRUE,
> > + ylab="Probability of Rating", xlab="City",main="Cleanliness
> > Ratings by City",
> > +  factor.names=FALSE, ticks=c(0.1,0.2,0.3,0.4,0.5,0.6))
> >
> >> panel.abline(h=c(0.1,0.2,0.3,0.4,0.5,0.6))

You can't just modify a lattice graph on the screen like that. plot() invisibly returns the lattice object; I suppose that you could try to modify that, but I think that my original suggestion -- to make a custom plot from the object returned by effect() -- is likely simpler.

John

> 
> >>
> >
> >
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> ;;;;;;;;;;;;;;
> >> sessionInfo()
> > R version 2.9.0 RC (2009-04-10 r48318)
> > i386-pc-mingw32
> >
> > locale:
> > LC_COLLATE=English_United States.1252;LC_CTYPE=English_United States.
> > 1252;LC_MONETARY=English_United States.
> > 1252;LC_NUMERIC=C;LC_TIME=English_United States.1252
> >
> > attached base packages:
> > [1] tcltk     grid      stats     graphics  grDevices utils
> > datasets  methods
> > [9] base
> >
> > other attached packages:
> > [1] relimp_1.0-1     Rcmdr_1.4-9      car_1.2-13       effects_2.0-4
> > [5] colorspace_1.0-0 nnet_7.2-46      MASS_7.2-46      lattice_0.17-22
> >
> > loaded via a namespace (and not attached):
> > [1] tools_2.9.0
> > Thank you for any advice.
> > Paul
> >
> > Paul Prew   ▪  Statistician
> > 651-795-5942   ▪   fax 651-204-7504
> > Ecolab Research Center   ▪  Mail Stop ESC-F4412-A
> > 655 Lone Oak Drive   ▪   Eagan, MN 55121-1560
> >
> >
> 
> David Winsemius, MD
> Heritage Laboratories
> West Hartford, CT
> 
> ______________________________________________
> 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.