Fits a regression model using a polynomial basis expansion of the input variables, with penalization via the adaptive LASSO or SCAD to provide oracle variable selection.
1 2 3 4 5 6 7 8
polywog(formula, data, subset, weights, na.action, degree = 3, family = c("gaussian", "binomial"), method = c("alasso", "scad"), penwt.method = c("lm", "glm"), unpenalized = character(0), .parallel = FALSE, boot = 0, control.boot = control.bp(.parallel = .parallel), lambda = NULL, nlambda = 100, lambda.min.ratio = 1e-04, nfolds = 10, foldid = NULL, thresh = ifelse(method == "alasso", 1e-07, 0.001), maxit = ifelse(method == "alasso", 1e+05, 5000), model = TRUE, X = FALSE, y = FALSE)
model formula specifying the response and input variables. See "Details" for more information.
a data frame, list or environment containing the variables specified in the model formula.
an optional vector specifying a subset of observations to be used in fitting.
an optional vector specifying weights for each observation to be used in fitting.
a function specifying what to do with observations
integer specifying the degree of the polynomial expansion of the input variables.
variable selection method:
estimator for obtaining first-stage estimates in
logistic models when
names of model terms to be exempt from the adaptive
penalty (only available when
logical: whether to perform k-fold cross-validation in
parallel (only available when
number of bootstrap iterations (0 for no bootstrapping).
list of arguments to be passed to
a vector of values from which the penalty factor is to be
selected via k-fold cross-validation.
number of values of the penalty factor to examine via
ratio of the lowest value to the highest in the
generated sequence of values of the penalty factor if
number of folds to use in cross-validation to select the penalization factor.
optional vector manually assigning fold numbers to each
observation used for fitting (only available when
convergence threshold, passed as the
maximum number of iterations to allow in adaptive LASSO or SCAD fitting.
logical: whether to include the model frame in the returned object.
logical: whether to include the raw design matrix (i.e., the matrix of input variables prior to taking their polynomial expansion) in the returned object.
logical: whether to include the response variable in the returned object.
The design matrix for the regression is a polynomial basis expansion of the
matrix of raw input variables. This includes all powers and interactions of
the input variables up to the specified
degree. For example, the
following terms will be included in
polywog(y ~ x1 + x2, degree = 3,
terms of degree 0: intercept
terms of degree 1:
terms of degree 2:
terms of degree 3:
To exclude certain terms from the basis expansion, use a model formula like
y ~ x1 + x2 | z1 + z2. Only the degree 1 terms of
z2 will be included.
It is possible that the "raw" basis expansion will be rank-deficient, such
as if there are binary input variables (in which case x_i = x_i^n for
all n > 0). The procedure detects collinearity via
removes extraneous columns before fitting.
For both the adaptive LASSO and SCAD, the penalization factor λ
is chosen by k-fold cross-validation. The selected value minimizes the
average mean squared error of out-of-sample fits. (To select both
λ and the polynomial degree simultaneously via cross-validation,
The cross-validation process may be run in parallel via
foreach by registering an appropriate backend and specifying
.parallel = TRUE. The appropriate backend is system-specific; see
foreach for information on selecting and registering a
backend. The bootstrap iterations may also be run in parallel by
control.boot = control.bp(.parallel = TRUE).
An object of class
"polywog", a list containing:
the estimated coefficients.
value of the penalty factor λ used to fit the final model.
a list containing the results of the cross-validation procedure used to select the penalty factor:
values of the penalty factor tested in cross-validation.
out-of-fold prediction error corresponding to
each value of
lambda with the minimal
minimized value of the cross-validation error.
the fitted mean values for each observation used in fitting.
coefficients from an unpenalized least-squares regression of the response variable on the polynomial expansion of the input variables.
adaptive weight given to each term in the LASSO
NULL for models fit via SCAD).
model formula, as a
degree of the polynomial basis expansion.
observation weights if specified.
the specified regularization method.
the specified method for calculating
the adaptive LASSO weights (
NULL for models fit via SCAD).
logical vector indicating which terms were not included in the LASSO penalty.
convergence threshold used in fitting.
iteration limit used in fitting.
terms object used to construct the
a matrix indicating the power of each raw input term (columns) in each term of the polynomial expansion used in fitting (rows).
the number of observations used to fit the model.
information on how
NA values in the input
data were handled.
levels of factor variables used in fitting.
names of the raw input variables included in the model formula.
the original function call.
model = TRUE, the model frame used in
X = TRUE, the raw model matrix (i.e., prior to
taking the polynomial expansion); otherwise
NULL. For calculating
the expanded model matrix, see
y = TRUE, the response variable used in
boot > 0, a sparse matrix of class
"dgCMatrix" where each column is the estimate from a
bootstrap replicate. See
bootPolywog for more information
Brenton Kenkel and Curtis S. Signorino
Brenton Kenkel and Curtis S. Signorino. 2012. "A Method for Flexible Functional Form Estimation: Bootstrapped Basis Regression with Variable Selection." Typescript, University of Rochester.
To estimate variation via the bootstrap, see
bootPolywog. To generate fitted values, see
predVals (and the underlying method
predict.polywog). For plots, see
The polynomial degree may be selected via cross-validation using
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27
## Using occupational prestige data data(Prestige, package = "carData") Prestige <- transform(Prestige, income = income / 1000) ## Fit a polywog model with bootstrap iterations ## (note: using low convergence threshold to shorten computation time of the ## example, *not* recommended in practice!) set.seed(22) fit1 <- polywog(prestige ~ education + income + type, data = Prestige, degree = 2, boot = 5, thresh = 1e-4) ## Basic information print(fit1) summary(fit1) ## See how fitted values change with education holding all else fixed predVals(fit1, "education", n = 10) ## Plot univariate relationships plot(fit1) ## Use SCAD instead of adaptive LASSO fit2 <- update(fit1, method = "scad", thresh = 1e-3) cbind(coef(fit1), coef(fit2))
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.