[Rd] model.response() as a generic function?

[Olivia Lau]

Hi,

We were thinking of the following:  R has a great user interface -- formula 
+ data, followed by terms, model.frame(), model.matrix(), 
model.response() -- for single equation models (e.g., univariate response 
models). We're trying to extend that user interface to multiple equation 
models so that multivariate response models, hierarchical models, 
multi-level models, etc. will have a similar, easy to use UI and programmer 
interface.

An example:  In a case like SURM, having Y as a matrix makes sense and we 
can implement that using the current model.response(model.frame(cbind(y1, 
y2) ~ X, data)).  In the case of multi-level models, however, you may want 
to manipulate the responses in each level separately and need a way to 
identify which response is associated with which formula.  For example, 
let's say that you are trying to write a user interface for a multi-level 
model that takes counties (or some other sub-unit) at the first level, and 
US states at the second level -- at the first level, there will be n_k 
units; at the second level, there will be k units -- clearly, cbind(y1, y2) 
will not work in this case.

One solution is to simply write a function that takes myfn(formula1, data1, 
formula2, data2, ...) for however many levels that are required for your 
model.  Using that approach, however, means that you need to either cobble 
together the RHS and LHS of the various formulas into one long formula (so 
that you can use model.matrix.default and model.frame.default) *or* write 
your own version of model.frame and model.matrix.  These are non-trivial 
tasks if your model accepts constraints across equations (e.g., users can 
specify that parameter effects be constrained to be equal for x1 in formula1 
and x1 in formula2).  You can ask users to input a constraint matrix, but 
this is (in practice) difficult for the average R user.  What we've been 
working on is a way to use formula plus some new specials that we've written 
to construct constraint matrices on the fly (e.g., in the formula).

To do this, we wrote a version of model.matrix.ourclass() that takes all the 
arguments of model.matrix.default() and an additional eqn argument, which 
specifies the equations (say mu, which corresponds to the sub-unit level, 
and gamma, which corresponds to the unit level) for which you would like to 
construct a model.matrix().  We would like to do the same to 
model.response(), which currently takes a model frame.  We are writing a 
version of model.response.ourclass() that takes a model frame and an 
additional argument that specifies the equation.  Since we are just adding 
an argument to the current model.response() function, the dispatch method 
would require a model frame as the first input, with dots being passed to 
other methods -- all code that currently works with model.response() will 
still work since the current model.response() would become 
model.response.default().

Let "MM" be a multi-level model with levels "mu" corresponding to the unit 
and "gamma" corresponding to the sub-unit, and let ourclass = "multiple". 
Then our proposed procedure is as follows:

fml <- parse.formula(formula, model = "MM")
            # checks formula input to make sure it corresponds to
            # what your model needs, where where fml is of class "multiple"
D <- model.frame(fml, data = mydata)
           # since fml is of class "multiple", this actually goes to 
model.frame.multiple()
           # D is a data.frame with classes c("multiple", "data.frame")
Xgamma <- model.matrix(fml, data = D, eqn = "gamma")
Ygamma <- model.response(D, eqn = "gamma")
            # Goes to model.response.multiple() since the first class of D 
is  "multiple"

The versions of model.frame.multiple() and model.matrix.multiple() that 
we've written also take terms as the first input rather than fml, and if 
data is null, looks in the parent environment for the appropriate data frame 
given in terms.  (e.g., works just as the .default() functions work, with 
the added eqn argument)

I've omitted a few details (esp re the parse.formula() function), but this 
is the general gist.  Please let me know if you have further questions!

Thanks,

Olivia Lau

----- Original Message ----- 
From: "Prof Brian Ripley" <ripley at stats.ox.ac.uk>
To: "Kosuke Imai" <kimai at princeton.edu>
Cc: <r-devel at r-project.org>; "Olivia Lau" <olau at fas.harvard.edu>; "Gary 
King" <king at harvard.edu>
Sent: Saturday, October 29, 2005 2:17 AM
Subject: Re: [Rd] model.response() as a generic function?


> There is a considerable difference: model.response is documented to work 
> on a model frame.  Why does model.response(model.frame(object))) not work 
> in the generality you need?
>
> Please give us some examples of why you are trying to do, and how you 
> would envisage a generic model.response being documented.  (It is 
> currently documented as equivalent to model.extract(,"response").)
>
> On Fri, 28 Oct 2005, Kosuke Imai wrote:
>
>> Dear R contributors,
>>  Gary King, Olivia Lau (both at Harvard) and I are working on an R
>> package where we are trying to write functions equivalent to
>> model.frame(), model.matrix(), and model.response() for multiple 
>> equations
>> models (i.e., models that require the specification of multiple 
>> formulae).
>> However, we noticed that while model.frame() and model.matrix() are
>> generic functions, model.response() is not. It would be nice if
>> model.response() is also a generic function so that we can make use of 
>> the
>> same function name for a large class of models we are working on. Is 
>> there
>> any possibility that this change can be made in the future version of R?
>> Thank you very much for your consideration in advance,
>> Kosuke
>
> -- 
> Brian D. Ripley,                  ripley at stats.ox.ac.uk
> Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
> University of Oxford,             Tel:  +44 1865 272861 (self)
> 1 South Parks Road,                     +44 1865 272866 (PA)
> Oxford OX1 3TG, UK                Fax:  +44 1865 272595
>