[Rd] model.matrix evaluation challenges

Ben Bolker bolker at ufl.edu
Mon Aug 10 02:30:13 CEST 2009


  I am having difficulty with evaluation/environment construction
for a formula to be evaluated by model.matrix().  Basically, I
want to construct a model matrix that first looks in "newdata"
for the values of the model parameters, then in "object at data".
Here's what I've tried:

  1. model.matrix(~f,data=c(newdata,object at data)) -- fails because
something (terms()?) within model.matrix() tries to coerce "data"
to a data frame -- but newdata and object at data have different
numbers of rows.

  2. with(c(newdata,object at data),model.matrix(~f))   works --
BUT not when ~f is assigned to another object (as it will be
when it is being dug out of a previously fit model).

  3. If "f" exists as a column within newdata and/or object at data,
but the formula is assigned to another object (see below), then
with(c(newdata,object at data),model.matrix(z))  fails -- because
model.matrix() is explicitly evaluating the variables in the formula
in the environment of z (i.e., ignoring the first argument of "with" ...)


  Any advice on how to solve this without making a bigger mess?

  sincerely
    Ben Bolker


## set up a data frame for prediction

set.seed(1001)
f = factor(rep(letters[1:4],each=20))
x = runif(80)
u = rnorm(4)
y = rnorm(80,mean=2+x*(3+u[f]),sd=0.1)
dat = data.frame(f,x,y)

## fit a model ... could easily do by lm() but want to
##   demonstrate the problem

library(bbmle)
m1 = mle2(y~dnorm(a+b*x,sd=exp(logs)),parameters=list(b~f),data=dat,
  start=list(a=0,b=2,logs=-3))

## data frame for prediction
pp0 = expand.grid(x=seq(0,1,length=11),
  f=levels(dat$f))

## combine frame and model data: have to keep the model data
##  around, because it contain other information needed for
##  prediction.

cc = c(pp0,m1 at data)
try(mm <- model.matrix(~f,cc))  ## fails -- different number of rows
  ## (at some point R tries to coerce cc to a data frame)
nrow(with(cc,model.matrix(~f))) ## works ... but ...

## here's what happens within predict.mle2()
## extract information about parameter models
params <- eval(m1 at call$parameters)
params <- params[[1]]
params[[2]] <- NULL

nrow(with(cc,model.matrix(params)))
rm(f)
## object "f" not found -- because model.matrix
##  evaluates in the parent environment of params ...
try(nrow(with(cc,model.matrix(params))))



-- 
Ben Bolker
Associate professor, Biology Dep't, Univ. of Florida
bolker at ufl.edu / www.zoology.ufl.edu/bolker
GPG key: www.zoology.ufl.edu/bolker/benbolker-publickey.asc

-------------- next part --------------
A non-text attachment was scrubbed...
Name: signature.asc
Type: application/pgp-signature
Size: 261 bytes
Desc: OpenPGP digital signature
URL: <https://stat.ethz.ch/pipermail/r-devel/attachments/20090809/5c009c7b/attachment.bin>


More information about the R-devel mailing list