R-beta: more model.matrix
Douglas Bates
bates at stat.wisc.edu
Fri Oct 17 19:25:48 CEST 1997
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.
R> fm1 <- lm(y ~ x.1 + x.2 + x.3 + x.4 + x.5, data = e3.7, model = T)
R> summary(fm1)
Call:
lm(formula = y ~ x.1 + x.2 + x.3 + x.4 + x.5, data = e3.7, model = T)
Residuals:
Min 1Q Median 3Q Max
-0.3944667 -0.1184741 0.0005345 0.0831318 0.5623167
Coefficients:
Estimate Std.Error t Value Pr(>|t|)
(Intercept) -2.1561 0.9135 -2.3603 0.0333
x.1 0.0000 0.0005 -0.0174 0.9864
x.2 0.0013 0.0013 1.0415 0.3153
x.3 0.0001 0.0001 1.6618 0.1188
x.4 0.0079 0.0140 0.5642 0.5815
x.5 0.0001 0.0001 1.9208 0.0754
Residual standard error: 0.2618 on 14 degrees of freedom
Multiple R-Squared: 0.8107, Adjusted R-squared: 0.743
F-statistic: 11.99 on 5 and 14 degrees of freedom, p-value: 0.0001184
R> fm1Sampl <- bsCoefSample(fm1, 500)
Error: Object "y" not found
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.
--
Douglas Bates bates at stat.wisc.edu
Statistics Department 608/262-2598
University of Wisconsin - Madison http://www.stat.wisc.edu/~bates/
=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=
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