effects package --- add abline to plot

Dear Paul,
-----Original Message-----
From: Prew, Paul [mailto:Paul.Prew at ecolab.com]
Sent: April-28-09 2:19 PM
To: John Fox; David Winsemius
Cc: r-help at r-project.org
Subject: RE: [R] effects package --- add abline to plot

Dear John and David,  thank you for your help.  I apologize for not defining
the analysis as an ordinal regression, or including a structure --- could
have taken some of the guesswork out for you.

John --- for the ticks, I would still like to make this work for future
analyses, but still not sure what specifically needs changing.  Before
initially posting, I did read ?effect, and several other searches around
"tick" and "at", but couldn't find a workable description or example for how
to use "at".  The one example I did find I thought I had copied pretty
closely with my command below.

plot(..., ticks=c(0.1,0.2,0.3,0.4,0.5,0.6))

I tried other combinations that probably look pretty silly...
plot(..., ticks(at=c(0.1,0.2,0.3,0.4,0.5,0.6)))
plot(..., ticks=c(0.1:0.6/0.1))

Just don't know how to properly populate this list ---
"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."
You need to specify ticks as a *list*; in your case, you just need the first element, so ticks=list(at=c(0.1,0.2,0.3,0.4,0.5,0.6)) or more compactly ticks=list(at=0.1*1:6) should do the trick. 

Can you tell me a bit more about why you want to draw horizontal lines on the plot? If this seems like a generally useful idea, I can put it my to-do list.

I hope this helps,
 John
Thanks again.  Effects is a terrific package.
Paul

**************  Model Specification ****************
Clean.label <- polr(Clean.lbl ~ City.Abbr, method="logistic", data=Safeway,
  Hess=TRUE)

str(Clean.label)
List of 18
 $ coefficients : Named num [1:8] -1.0887 -1.1449 0.6923 0.0894 -0.8229 ...
  ..- attr(*, "names")= chr [1:8] "City.Abbr[T.Dublin]" "City.Abbr[T.Englwd]"
"City.Abbr[T.Fairfax]" "City.Abbr[T.Falls Ch]" ...
 $ zeta         : Named num [1:4] -2.529 0.894 3.447 6.406
  ..- attr(*, "names")= chr [1:4] "Excellent|V.Good" "V.Good|Good"
"Good|Fair" "Fair|Poor"
 $ deviance     : num 2327
 $ fitted.values: num [1:1248, 1:5] 0.0568 0.0568 0.0568 0.0568 0.0568 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:1248] "1" "2" "3" "4" ...
  .. ..$ : chr [1:5] "Excellent" "V.Good" "Good" "Fair" ...
 $ lev          : chr [1:5] "Excellent" "V.Good" "Good" "Fair" ...
 $ terms        :Classes 'terms', 'formula' length 3 Clean.lbl ~ City.Abbr
  .. ..- attr(*, "variables")= language list(Clean.lbl, City.Abbr)
  .. ..- attr(*, "factors")= int [1:2, 1] 0 1
  .. .. ..- attr(*, "dimnames")=List of 2
  .. .. .. ..$ : chr [1:2] "Clean.lbl" "City.Abbr"
  .. .. .. ..$ : chr "City.Abbr"
  .. ..- attr(*, "term.labels")= chr "City.Abbr"
  .. ..- attr(*, "order")= int 1
  .. ..- attr(*, "intercept")= int 1
  .. ..- attr(*, "response")= int 1
  .. ..- attr(*, ".Environment")=<environment: R_GlobalEnv>
  .. ..- attr(*, "predvars")= language list(Clean.lbl, City.Abbr)
  .. ..- attr(*, "dataClasses")= Named chr [1:2] "ordered" "factor"
  .. .. ..- attr(*, "names")= chr [1:2] "Clean.lbl" "City.Abbr"
 $ df.residual  : num 1236
 $ edf          : num 12
 $ n            : num 1248
 $ nobs         : num 1248
 $ call         : language polr(formula = Clean.lbl ~ City.Abbr, data =
Safeway, Hess = TRUE,      method = "logistic")
 $ method       : chr "logistic"
 $ convergence  : int 0
 $ niter        : Named int [1:2] 62 22
  ..- attr(*, "names")= chr [1:2] "f.evals.function" "g.evals.gradient"
 $ Hessian      : num [1:12, 1:12] 31.9 0 0 0 0 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:12] "City.Abbr[T.Dublin]" "City.Abbr[T.Englwd]"
"City.Abbr[T.Fairfax]" "City.Abbr[T.Falls Ch]" ...
  .. ..$ : chr [1:12] "City.Abbr[T.Dublin]" "City.Abbr[T.Englwd]"
"City.Abbr[T.Fairfax]" "City.Abbr[T.Falls Ch]" ...
 $ model        :'data.frame':	1248 obs. of  2 variables:
  ..$ Clean.lbl: Ord.factor w/ 5 levels "Excellent"<"V.Good"<..: 1 2 3 3 3 3
3 3 2 2 ...
  ..$ City.Abbr: Factor w/ 9 levels "DC","Dublin",..: 9 9 9 9 9 9 9 9 9 9 ...
  ..- attr(*, "terms")=Classes 'terms', 'formula' length 3 Clean.lbl ~
City.Abbr
  .. .. ..- attr(*, "variables")= language list(Clean.lbl, City.Abbr)
  .. .. ..- attr(*, "factors")= int [1:2, 1] 0 1
  .. .. .. ..- attr(*, "dimnames")=List of 2
  .. .. .. .. ..$ : chr [1:2] "Clean.lbl" "City.Abbr"
  .. .. .. .. ..$ : chr "City.Abbr"
  .. .. ..- attr(*, "term.labels")= chr "City.Abbr"
  .. .. ..- attr(*, "order")= int 1
  .. .. ..- attr(*, "intercept")= int 1
  .. .. ..- attr(*, "response")= int 1
  .. .. ..- attr(*, ".Environment")=<environment: R_GlobalEnv>
  .. .. ..- attr(*, "predvars")= language list(Clean.lbl, City.Abbr)
  .. .. ..- attr(*, "dataClasses")= Named chr [1:2] "ordered" "factor"
  .. .. .. ..- attr(*, "names")= chr [1:2] "Clean.lbl" "City.Abbr"
 $ contrasts    :List of 1
  ..$ City.Abbr: chr "contr.Treatment"
 $ xlevels      :List of 1
  ..$ City.Abbr: chr [1:9] "DC" "Dublin" "Englwd" "Fairfax" ...
 - attr(*, "class")= chr "polr"

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

-----Original Message-----
From: John Fox [mailto:jfox at mcmaster.ca]
Sent: Tuesday, April 28, 2009 10:34 AM
To: 'David Winsemius'
Cc: r-help at r-project.org; Prew, Paul
Subject: RE: [R] effects package --- add abline to plot

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.

CONFIDENTIALITY NOTICE:
This e-mail communication and any attachments may contain proprietary and
privileged information for the use of the designated recipients named above.
Any unauthorized review, use, disclosure or distribution is prohibited.
If you are not the intended recipient, please contact the sender by reply e-
mail and destroy all copies of the original message.

effects package --- add abline to plot

Thread (7 messages)