R-beta: more model.matrix
Thomas Lumley
thomas at biostat.washington.edu
Fri Oct 17 22:02:21 CEST 1997
On Fri, 17 Oct 1997, Douglas Bates wrote:
> I am trying to show some techniques to my graduate regression class.
> The textbook mentioned using bootstrap samples of regression
> coefficients for assessing variability. I decided to show them
> reasonably effective ways of doing the resampling.
>
> The following is a function I wrote to create bootstrap samples of
> coefficients from a fitted linear regression model.
>
> bsCoefSample <-
> ## Construct a bootstrap sample of coefficients from a
> ## fitted regression model
> function(fittedModel, Nsampl, ...)
> {
> coef(fittedModel) +
> qr.coef(qr(model.matrix(fittedModel)),
> matrix(sample(resid(fittedModel),
> length(resid(fittedModel)) * Nsampl,
> repl = T),
> ncol = Nsampl))
> }
>
> It works ok with S-PLUS but the R version of model.matrix encounters
> difficulties.
> Can anyone suggest a modified function that would work in both S-PLUS
> and R?
>
> It appears that model.matrix is a generic in S-PLUS
> S-PLUS> methods(model.matrix)
> $SHOME/stat/.Functions $SHOME/stat/.Functions
> "model.matrix.default" "model.matrix.lm"
>
> Could this be done in R as well? I am willing to try to do that if
> it seems like a good idea.
We can define
model.matrix.lm<-function(lmobj){
model.matrix(delete.response(terms(lmobj)),model.frame(lmobj))
}
so that an explicit call to model.matrix.lm instead of model.matrix would
work in both R and S. It might well be a good idea to make model.matrix
generic.
Another option is to use update(,qr=T) to refit the model and return the
qr decomposition. This is a bit more tricky, since you need to evaluate
things in the right environment.
Thomas Lumley
------------------------------------------------------+------
Biostatistics : "Never attribute to malice what :
Uni of Washington : can be adequately explained by :
Box 357232 : incompetence" - Hanlon's Razor :
Seattle WA 98195-7232 : :
------------------------------------------------------------
=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch
=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=
More information about the R-help
mailing list