[R] User error in calling predict/model.frame

[Joshua Wiley]

On Sat, Jan 29, 2011 at 12:31 PM, David Winsemius
<dwinsemius at comcast.net> wrote:
> Huh?. With the same model and data, they should be the same:

I must be missing something, because that is what I would have
expected too, but it is not what I get (at least when I run the code
as shown below).

>> lm.mod2 <- lm(out ~ scale(xxA)*scale(xxB), data=dat)
>> newdata2 <- expand.grid(xxA=c(-1,0,1),xxB=c(-1,0,1))
>> newdata2$preds <- predict(lm.mod, newdata2)
                                                    ^^^^^ is that lm.mod2?
>> newdata2
>  xxA xxB    preds
> 1  -1  -1 218.6366
> 2   0  -1 232.4330
> 3   1  -1 246.2295
> 4  -1   0 230.9129
> 5   0   0 244.8696
> 6   1   0 258.8263
> 7  -1   1 243.1892
> 8   0   1 257.3062
> 9   1   1 271.4232
>
>> Predict(rms.res,xxA=c(-1,0,1),xxB=c(-1,0,1), conf.int=FALSE)
>  xxA xxB     yhat
> 1  -1  -1 218.6366
> 2   0  -1 232.4330
> 3   1  -1 246.2295
> 4  -1   0 230.9129
> 5   0   0 244.8696
> 6   1   0 258.8263
> 7  -1   1 243.1892
> 8   0   1 257.3062
> 9   1   1 271.4232

This is copied directly from my console in a clean session.  I get
different values for the predicted values from predict(glm or lm) and
"yhat" from Predict(Glm).


> sessionInfo()
R version 2.12.1 (2010-12-16)
Platform: i486-pc-linux-gnu (32-bit)

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C
 [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8
 [5] LC_MONETARY=C              LC_MESSAGES=en_US.UTF-8
 [7] LC_PAPER=en_US.UTF-8       LC_NAME=C
 [9] LC_ADDRESS=C               LC_TELEPHONE=C
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base
> #################################
> require(rms)
Loading required package: rms
Loading required package: Hmisc
Loading required package: survival
Loading required package: splines

Attaching package: 'Hmisc'

The following object(s) are masked from 'package:survival':

    untangle.specials

The following object(s) are masked from 'package:base':

    format.pval, round.POSIXt, trunc.POSIXt, units


Attaching package: 'rms'

The following object(s) are masked from 'package:survival':

    Surv

> set.seed(10)
> dat <- data.frame(xxA = rnorm(20, 10), xxB = rnorm(20, 20))
> dat$out <- with(dat, xxA+xxB+xxA*xxB+rnorm(20,20))
>
> rms.res <- Glm(out ~ scale(xxA) * scale(xxB), data=dat)
> glm.mod <- glm(out ~ scale(xxA) * scale(xxB), data=dat)
> lm.mod2 <- lm(out ~ scale(xxA)*scale(xxB), data=dat)
>
> newdata2 <- expand.grid(xxA = -1:1, xxB = -1:1)
> newdata2$preds1 <- predict(glm.mod, newdata2)
> newdata2$preds2 <- predict(lm.mod2, newdata2)
> newdata2
  xxA xxB   preds1   preds2
1  -1  -1 143.3433 143.3433
2   0  -1 131.6110 131.6110
3   1  -1 119.8787 119.8787
4  -1   0 137.1149 137.1149
5   0   0 126.9708 126.9708
6   1   0 116.8267 116.8267
7  -1   1 130.8865 130.8865
8   0   1 122.3306 122.3306
9   1   1 113.7747 113.7747
>
> Predict(rms.res, xxA= -1:1, xxB= -1:1, conf.int=FALSE)
  xxA xxB     yhat
1  -1  -1 212.2324
2   0  -1 229.6472
3   1  -1 247.0620
4  -1   0 221.7795
5   0   0 240.6413
6   1   0 259.5030
7  -1   1 231.3266
8   0   1 251.6353
9   1   1 271.9441

Response variable (y): X * Beta
> ##################################