Skip to content

using predict.lm() within a function

2 messages · Michael Friendly, Dennis Murphy

Michael Friendly
# 
I've written a simple function to draw a regression line in a plot and 
annotate the line showing the slope
with a label.  It works, as I'm using it, when the horizontal variable 
is 'x', but gives incorrect results otherwise.
What's wrong?

# simple function to show the slope of a line
show.beta <- function(model, x="x", x1, x2, label, col="black", ...) {
     abline(model, col=col, lwd=2)
#    x <- deparse(substitute(x))      # doesn't help
     xs <- data.frame(x=c(x1, x2, x2))
     ys <- predict(model, xs)
     lines(cbind(xs,ys[c(1,1,2)]), col=col)
     text(x2, mean(ys[1:2]), label, col=col, ...)
}

x=rnorm(10)
DF <- data.frame(x, y=.25*x+rnorm(10))

# OK
with(DF, {
     plot(y ~ x)
     mod <- lm(y ~ x)
     show.beta(mod, "x", -0.5, 0, "b", pos=4)
     })

# not OK
xx=rnorm(10)
DF2 <- data.frame(xx, y=.25*x+rnorm(10))

with(DF2, {
     plot(y ~ xx)
     mod <- lm(y ~ xx)
     show.beta(mod, "xx", -0.5, 0, "b", pos=4)
     })

 From the latter, I get:

Warning message:
'newdata' had 3 rows but variable(s) found have 10 rows
 >
Dennis Murphy
# 
Hi Michael:

Try this:

show.beta <- function(model, x = 'x',
                      x1, x2, label, col="black", ...) {
   abline(model, col=col, lwd=2)
   xs <- data.frame(c(x1, x2, x2))
   names(xs) <- attr(model$coefficients, 'names')[2]
   ys <- predict(model, xs)
   lines(cbind(xs,ys[c(1,1,2)]), col=col)
   text(x2, mean(ys[1:2]), label, col=col, ...)
}

Both of your examples worked for me. The trick is to match the
variable name in the model with the variable name in the prediction
data frame. There are several ways to access that, but I chose the
first opportunity I saw:
List of 12
 $ coefficients : Named num [1:2] 0.501 0.351
  ..- attr(*, "names")= chr [1:2] "(Intercept)" "xx"

The line where I assign the names to xs is where the attr() function
comes into play. The first argument is the choice of list component
from the model object, while the second is the name of the attribute
to extract. This one returns a vector whose second component is what
we want.

HTH,
Dennis
On Mon, Oct 24, 2011 at 2:30 PM, Michael Friendly <friendly at yorku.ca> wrote: