Package 'car'

Description

Calculates type-II or type-III analysis-of-variance tables for model objects produced by lm, glm, multinom (in the nnet package), polr (in the MASS package), coxph (in the survival package), coxme (in the coxme pckage), svyglm and svycoxph (in the survey package), rlm (in the MASS package), lmer (in the lme4 package), lme (in the nlme package), clm and clmm (in the ordinal package), and (by the default method) for most models with a linear predictor and asymptotically normal coefficients (see details below). For linear models, F-tests are calculated; for generalized linear models, likelihood-ratio chisquare, Wald chisquare, or F-tests are calculated; for multinomial logit and proportional-odds logit models, likelihood-ratio tests are calculated. Various test statistics are provided for multivariate linear models produced by lm or manova. Partial-likelihood-ratio tests or Wald tests are provided for Cox models. Wald chi-square tests are provided for fixed effects in linear and generalized linear mixed-effects models. Wald chi-square or F tests are provided in the default case.

Usage

Anova(mod, ...)

Manova(mod, ...)

## S3 method for class 'lm'
Anova(mod, error, type=c("II","III", 2, 3),
	white.adjust=c(FALSE, TRUE, "hc3", "hc0", "hc1", "hc2", "hc4"),
	vcov.=NULL, singular.ok, ...)

## S3 method for class 'aov'
Anova(mod, ...)

## S3 method for class 'glm'
Anova(mod, type=c("II","III", 2, 3),
    test.statistic=c("LR", "Wald", "F"),
    error, error.estimate=c("pearson", "dispersion", "deviance"),
   vcov.=vcov(mod, complete=TRUE),  singular.ok, ...)

## S3 method for class 'multinom'
Anova(mod, type = c("II","III", 2, 3), ...)

## S3 method for class 'polr'
Anova(mod, type = c("II","III", 2, 3), ...)

## S3 method for class 'mlm'
Anova(mod, type=c("II","III", 2, 3), SSPE, error.df,
    idata, idesign, icontrasts=c("contr.sum", "contr.poly"), imatrix,
    test.statistic=c("Pillai", "Wilks", "Hotelling-Lawley", "Roy"),...)

## S3 method for class 'manova'
Anova(mod, ...)

## S3 method for class 'mlm'
Manova(mod, ...)

## S3 method for class 'Anova.mlm'
print(x, ...)

## S3 method for class 'Anova.mlm'
summary(object, test.statistic, univariate=object$repeated,
    multivariate=TRUE, p.adjust.method, ...)

## S3 method for class 'summary.Anova.mlm'
print(x, digits = getOption("digits"),
    SSP=TRUE, SSPE=SSP, ... )

## S3 method for class 'univaov'
print(x, digits = max(getOption("digits") - 2L, 3L),
                          style=c("wide", "long"),
                          by=c("response", "term"),
                          ...)

## S3 method for class 'univaov'
as.data.frame(x, row.names, optional, by=c("response", "term"), ...)

## S3 method for class 'coxph'
Anova(mod, type=c("II", "III", 2, 3),
	test.statistic=c("LR", "Wald"), ...)

## S3 method for class 'coxme'
Anova(mod, type=c("II", "III", 2, 3),
    test.statistic=c("Wald", "LR"), ...)

## S3 method for class 'lme'
Anova(mod, type=c("II","III", 2, 3),
		vcov.=vcov(mod, complete=FALSE), singular.ok, ...)

## S3 method for class 'mer'
Anova(mod, type=c("II", "III", 2, 3),
	test.statistic=c("Chisq", "F"), vcov.=vcov(mod, complete=FALSE), singular.ok, ...)

## S3 method for class 'merMod'
Anova(mod, type=c("II", "III", 2, 3),
    test.statistic=c("Chisq", "F"), vcov.=vcov(mod, complete=FALSE), singular.ok, ...)

## S3 method for class 'svyglm'
Anova(mod, ...)

## S3 method for class 'svycoxph'
Anova(mod, type=c("II", "III", 2, 3),
  test.statistic="Wald", ...)

## S3 method for class 'rlm'
Anova(mod, ...)

## S3 method for class 'clm'
Anova(mod, ...)

## S3 method for class 'clmm'
Anova(mod, ...)

## Default S3 method:
Anova(mod, type=c("II", "III", 2, 3),
	test.statistic=c("Chisq", "F"), vcov.=vcov(mod, complete=FALSE),
	singular.ok, error.df, ...)

Arguments

Details

The designations "type-II" and "type-III" are borrowed from SAS, but the definitions used here do not correspond precisely to those employed by SAS. Type-II tests are calculated according to the principle of marginality, testing each term after all others, except ignoring the term's higher-order relatives; so-called type-III tests violate marginality, testing each term in the model after all of the others. This definition of Type-II tests corresponds to the tests produced by SAS for analysis-of-variance models, where all of the predictors are factors, but not more generally (i.e., when there are quantitative predictors). Be very careful in formulating the model for type-III tests, or the hypotheses tested will not make sense.

As implemented here, type-II Wald tests are a generalization of the linear hypotheses used to generate these tests in linear models.

For tests for linear models, multivariate linear models, and Wald tests for generalized linear models, Cox models, mixed-effects models, generalized linear models fit to survey data, and in the default case, Anova finds the test statistics without refitting the model. The svyglm method simply calls the default method and therefore can take the same arguments.

The standard R anova function calculates sequential ("type-I") tests. These rarely test interesting hypotheses in unbalanced designs.

A MANOVA for a multivariate linear model (i.e., an object of class "mlm" or "manova") can optionally include an intra-subject repeated-measures design. If the intra-subject design is absent (the default), the multivariate tests concern all of the response variables. To specify a repeated-measures design, a data frame is provided defining the repeated-measures factor or factors via idata, with default contrasts given by the icontrasts argument. An intra-subject model-matrix is generated from the formula specified by the idesign argument; columns of the model matrix corresponding to different terms in the intra-subject model must be orthogonal (as is insured by the default contrasts). Note that the contrasts given in icontrasts can be overridden by assigning specific contrasts to the factors in idata. As an alternative, the within-subjects model matrix can be specified directly via the imatrix argument. Manova is essentially a synonym for Anova for multivariate linear models.

If univariate tests are requested for the summary of a multivariate linear model, the object returned contains a univaov component of "univaov"; print and as.data.frame methods are provided for the "univaov" class.

For the default method to work, the model object must contain a standard terms element, and must respond to the vcov, coef, and model.matrix functions. If any of these requirements is missing, then it may be possible to supply it reasonably simply (e.g., by writing a missing vcov method for the class of the model object).

Value

An object of class "anova", or "Anova.mlm", which usually is printed. For objects of class "Anova.mlm", there is also a summary method, which provides much more detail than the print method about the MANOVA, including traditional mixed-model univariate F-tests with Greenhouse-Geisser and Huynh-Feldt corrections.

Warning

Be careful of type-III tests: For a traditional multifactor ANOVA model with interactions, for example, these tests will normally only be sensible when using contrasts that, for different terms, are orthogonal in the row-basis of the model, such as those produced by contr.sum, contr.poly, or contr.helmert, but not by the default contr.treatment. In a model that contains factors, numeric covariates, and interactions, main-effect tests for factors will be for differences over the origin. In contrast (pun intended), type-II tests are invariant with respect to (full-rank) contrast coding. If you don't understand this issue, then you probably shouldn't use Anova for type-III tests.

Author(s)

John Fox [email protected]; the code for the Mauchly test and Greenhouse-Geisser and Huynh-Feldt corrections for non-spericity in repeated-measures ANOVA are adapted from the functions stats:::stats:::mauchly.test.SSD and stats:::sphericity by R Core; summary.Anova.mlm and print.summary.Anova.mlm incorporates code contributed by Gabriel Baud-Bovy.

References

Fox, J. (2016) Applied Regression Analysis and Generalized Linear Models, Third Edition. Sage.

Fox, J. and Weisberg, S. (2019) An R Companion to Applied Regression, Third Edition, Sage.

Hand, D. J., and Taylor, C. C. (1987) Multivariate Analysis of Variance and Repeated Measures: A Practical Approach for Behavioural Scientists. Chapman and Hall.

O'Brien, R. G., and Kaiser, M. K. (1985) MANOVA method for analyzing repeated measures designs: An extensive primer. Psychological Bulletin 97, 316–333.

See Also

linearHypothesis, anova anova.lm, anova.glm, anova.mlm, anova.coxph, svyglm.

Examples

## Two-Way Anova

mod <- lm(conformity ~ fcategory*partner.status, data=Moore,
  contrasts=list(fcategory=contr.sum, partner.status=contr.sum))
Anova(mod)
Anova(mod, type=3)  # note use of contr.sum in call to lm()

## use of vcov.; the following are equivalent:

Anova(mod, white.adjust = TRUE)
Anova(mod, vcov. = hccm) # vcov. is a function, type = "hc3" is the default
Anova(mod, vcov. = hccm(mod, type = "hc3")) # vcov. is a matrix
Anova(mod, vcov. = function(m) hccm(m, type = "hc3")) # passing type as an argument

## One-Way MANOVA
## See ?Pottery for a description of the data set used in this example.

summary(Anova(lm(cbind(Al, Fe, Mg, Ca, Na) ~ Site, data=Pottery)))

## MANOVA for a randomized block design (example courtesy of Michael Friendly:
##  See ?Soils for description of the data set)

soils.mod <- lm(cbind(pH,N,Dens,P,Ca,Mg,K,Na,Conduc) ~ Block + Contour*Depth,
    data=Soils)
Manova(soils.mod)
summary(Anova(soils.mod), univariate=TRUE, multivariate=FALSE,
    p.adjust.method=TRUE)

## a multivariate linear model for repeated-measures data
## See ?OBrienKaiser for a description of the data set used in this example.

phase <- factor(rep(c("pretest", "posttest", "followup"), c(5, 5, 5)),
    levels=c("pretest", "posttest", "followup"))
hour <- ordered(rep(1:5, 3))
idata <- data.frame(phase, hour)
idata

mod.ok <- lm(cbind(pre.1, pre.2, pre.3, pre.4, pre.5,
                     post.1, post.2, post.3, post.4, post.5,
                     fup.1, fup.2, fup.3, fup.4, fup.5) ~  treatment*gender,
                data=OBrienKaiser)
(av.ok <- Anova(mod.ok, idata=idata, idesign=~phase*hour))

summary(av.ok, multivariate=FALSE)

## A "doubly multivariate" design with two  distinct repeated-measures variables
## (example courtesy of Michael Friendly)
## See ?WeightLoss for a description of the dataset.

imatrix <- matrix(c(
	1,0,-1, 1, 0, 0,
	1,0, 0,-2, 0, 0,
	1,0, 1, 1, 0, 0,
	0,1, 0, 0,-1, 1,
	0,1, 0, 0, 0,-2,
	0,1, 0, 0, 1, 1), 6, 6, byrow=TRUE)
colnames(imatrix) <- c("WL", "SE", "WL.L", "WL.Q", "SE.L", "SE.Q")
rownames(imatrix) <- colnames(WeightLoss)[-1]
(imatrix <- list(measure=imatrix[,1:2], month=imatrix[,3:6]))
contrasts(WeightLoss$group) <- matrix(c(-2,1,1, 0,-1,1), ncol=2)
(wl.mod<-lm(cbind(wl1, wl2, wl3, se1, se2, se3)~group, data=WeightLoss))
Anova(wl.mod, imatrix=imatrix, test="Roy")

## mixed-effects models examples:

## Not run:  # loads nlme package
	library(nlme)
	example(lme)
	Anova(fm2)

## End(Not run)

## Not run:  # loads lme4 package
	library(lme4)
	example(glmer)
	Anova(gm1)

## End(Not run)

Description

These functions construct added-variable, also called partial-regression, plots for linear and generalized linear models.

Usage

avPlots(model, ...)

## Default S3 method:
avPlots(model, terms=~., intercept=FALSE, 
  layout=NULL, ask, main, ...)

avp(...)

avPlot(model, ...)

## S3 method for class 'lm'
avPlot(model, variable,
	id=TRUE, col = carPalette()[1], col.lines = carPalette()[2],
	xlab, ylab, pch = 1, lwd = 2, cex = par("cex"), pt.wts = FALSE,
	main=paste("Added-Variable Plot:", variable),
	grid=TRUE,
	ellipse=FALSE,
  marginal.scale=FALSE, ...)

## S3 method for class 'glm'
avPlot(model, variable, 
	id=TRUE,
	col = carPalette()[1], col.lines = carPalette()[2],
	xlab, ylab, pch = 1, lwd = 2, cex = par("cex"), pt.wts = FALSE,
	type=c("Wang", "Weisberg"), 
	main=paste("Added-Variable Plot:", variable), grid=TRUE,
	ellipse=FALSE, ...)
	
avPlot3d(model, coef1, coef2, id=TRUE, ...)

## S3 method for class 'lm'
avPlot3d(model, coef1, coef2, id=TRUE, fit="linear", ...)

## S3 method for class 'glm'
avPlot3d(model, coef1, coef2, id=TRUE, type=c("Wang", "Weisberg"), 
  fit="linear", ...)

Arguments

Details

The functions intended for direct use are avPlots (for which avp is an abbreviation) and avPlot3d.

Value

These functions are used for their side effect of producing plots, but also invisibly return the coordinates of the plotted points.

Author(s)

John Fox [email protected], Sanford Weisberg [email protected]

References

Cook, R. D. and Weisberg, S. (1999) Applied Regression, Including Computing and Graphics. Wiley.

Fox, J. (2016) Applied Regression Analysis and Generalized Linear Models, Third Edition. Sage.

Fox, J. and Weisberg, S. (2019) An R Companion to Applied Regression, Third Edition, Sage.

Wang, P C. (1985) Adding a variable in generalized linear models. Technometrics 27, 273–276.

Weisberg, S. (2014) Applied Linear Regression, Fourth Edition, Wiley.

See Also

residualPlots, crPlots, ceresPlots, link{dataEllipse}, showLabels, dataEllipse.

Examples

avPlots(lm(prestige ~ income + education + type, data=Duncan))

avPlots(glm(partic != "not.work" ~ hincome + children, 
  data=Womenlf, family=binomial), id=FALSE, pt.wts=TRUE)
  
m1 <- lm(partic ~ tfr + menwage + womwage + debt + parttime, Bfox)
par(mfrow=c(1,3))
# marginal plot, ignoring other predictors:
with(Bfox, dataEllipse(womwage, partic, levels=0.5)) 
abline(lm(partic ~ womwage, Bfox), col="red", lwd=2)
# AV plot, adjusting for others:
avPlots(m1, ~ womwage, ellipse=list(levels=0.5)) 
# AV plot, adjusting and scaling as in marginal plot
avPlots(m1, ~ womwage, marginal.scale=TRUE, ellipse=list(levels=0.5)) 

# 3D AV plot, requires the rgl package
if (interactive() && require("rgl")){
avPlot3d(lm(prestige ~ income + education + type, data=Duncan), 
  "income", "education")
}

Description

Transform the elements of a vector or columns of a matrix using, the Box-Cox, Box-Cox with negatives allowed, Yeo-Johnson, or simple power transformations.

Usage

bcPower(U, lambda, jacobian.adjusted=FALSE, gamma=NULL)

bcnPower(U, lambda, jacobian.adjusted = FALSE, gamma)

bcnPowerInverse(z, lambda, gamma)

yjPower(U, lambda, jacobian.adjusted = FALSE)

basicPower(U,lambda, gamma=NULL)

Arguments

.

Details

The Box-Cox family of scaled power transformations equals $(x^{\lambda}-1)/\lambda$ for $\lambda \neq 0$, and $\log(x)$ if $\lambda =0$. The bcPower function computes the scaled power transformation of $x = U + \gamma$, where $\gamma$ is set by the user so $U+\gamma$ is strictly positive for these transformations to make sense.

The Box-Cox family with negatives allowed was proposed by Hawkins and Weisberg (2017). It is the Box-Cox power transformation of

$z = .5 (U + \sqrt{U^2 + \gamma^2)})$

where for this family $\gamma$ is either user selected or is estimated. gamma must be positive if $U$ includes negative values and non-negative otherwise, ensuring that $z$ is always positive. The bcnPower transformations behave similarly to the bcPower transformations, and introduce less bias than is introduced by setting the parameter $\gamma$ to be non-zero in the Box-Cox family.

The function bcnPowerInverse computes the inverse of the bcnPower function, so U = bcnPowerInverse(bcnPower(U, lambda=lam, jacobian.adjusted=FALSE, gamma=gam), lambda=lam, gamma=gam) is true for any permitted value of gam and lam.

If family="yeo.johnson" then the Yeo-Johnson transformations are used. This is the Box-Cox transformation of $U+1$ for nonnegative values, and of $|U|+1$ with parameter $2-\lambda$ for $U$ negative.

The basic power transformation returns $U^{\lambda}$ if $\lambda$ is not 0, and $\log(\lambda)$ otherwise for $U$ strictly positive.

If jacobian.adjusted is TRUE, then the scaled transformations are divided by the Jacobian, which is a function of the geometric mean of $U$ for skewPower and yjPower and of $U + gamma$ for bcPower. With this adjustment, the Jacobian of the transformation is always equal to 1. Jacobian adjustment facilitates computing the Box-Cox estimates of the transformation parameters.

Missing values are permitted, and return NA where ever U is equal to NA.

Value

Returns a vector or matrix of transformed values.

Author(s)

Sanford Weisberg, <[email protected]>

References

Fox, J. and Weisberg, S. (2019) An R Companion to Applied Regression, Third Edition, Sage.

Hawkins, D. and Weisberg, S. (2017) Combining the Box-Cox Power and Generalized Log Transformations to Accomodate Nonpositive Responses In Linear and Mixed-Effects Linear Models South African Statistics Journal, 51, 317-328.

Weisberg, S. (2014) Applied Linear Regression, Fourth Edition, Wiley Wiley, Chapter 7.

Yeo, In-Kwon and Johnson, Richard (2000) A new family of power transformations to improve normality or symmetry. Biometrika, 87, 954-959.

See Also

powerTransform, testTransform

Examples

U <- c(NA, (-3:3))
## Not run: bcPower(U, 0)  # produces an error as U has negative values
bcPower(U, 0, gamma=4)
bcPower(U, .5, jacobian.adjusted=TRUE, gamma=4)
bcnPower(U, 0, gamma=2)
basicPower(U, lambda = 0, gamma=4)
yjPower(U, 0)
V <- matrix(1:10, ncol=2)
bcPower(V, c(0, 2))
basicPower(V, c(0,1))

Description

This function provides a simple front-end to the boot function in the boot package that is tailored to bootstrapping based on regression models. Whereas boot is very general and therefore has many arguments, the Boot function has very few arguments.

Usage

Boot(object, f=coef, labels=names(f(object)), R=999,
  method=c("case", "residual"), ncores=1, ...)

## Default S3 method:
Boot(object, f=coef, labels=names(f(object)),
  R=999, method=c("case", "residual"), ncores=1,
  start = FALSE, ...)

## S3 method for class 'lm'
Boot(object, f=coef, labels=names(f(object)),
  R=999, method=c("case", "residual"), ncores=1, ...)

## S3 method for class 'glm'
Boot(object, f=coef, labels=names(f(object)),
  R=999, method=c("case", "residual"), ncores=1, ...)

## S3 method for class 'nls'
Boot(object, f=coef, labels=names(f(object)),
  R=999, method=c("case", "residual"), ncores=1, ...)

Arguments

Details

Boot uses a regression object and the choice of method, and creates a function that is passed as the statistic argument to the boot function in the boot package. The argument R is also passed to boot. If ncores is greater than 1, then the parallel and ncpus arguments to boot are set appropriately to use multiple codes, if available, on your computer. All other arguments to boot are kept at their default values unless you pass values for them.

The methods available for lm and nls objects are “case” and “residual”. The case bootstrap resamples from the joint distribution of the terms in the model and the response. The residual bootstrap fixes the fitted values from the original data, and creates bootstraps by adding a bootstrap sample of the residuals to the fitted values to get a bootstrap response. It is an implementation of Algorithm 6.3, page 271, of Davison and Hinkley (1997). For nls objects ordinary residuals are used in the resampling rather than the standardized residuals used in the lm method. The residual bootstrap for generalized linear models has several competing approaches, but none are without problems. If you want to do a residual bootstrap for a glm, you will need to write your own call to boot.

For the default object to work with other types of regression models, the model must have methods for the the following generic functions: residuals(object, type="pearson") must return Pearson residuals; fitted(object) must return fitted values; hatvalues(object) should return the leverages, or perhaps the value 1 which will effectively ignore setting the hatvalues. In addition, the data argument should contain no missing values among the columns actually used in fitting the model, as the resampling may incorrectly attempt to include cases with missing values. For lm, glm and nls, missing values cause the return of an error message.

An attempt to fit using a bootstrap sample may fail. In a lm or glm fit, the bootstrap sample could have a different rank from the original fit. In an nls fit, convergence may not be obtained for some bootstraps. In either case, NA are returned for the value of the function f. The summary methods handle the NAs appropriately.

Fox and Weisberg (2017) cited below discusses this function and provides more examples.

Value

See boot for the returned value of the structure returned by this function.

Warning

C=A call like car::Boot(object, method="residual") will fail for the residual method if not preceded by library(car). If method="case" the library(car) command is not required.

Author(s)

Sanford Weisberg, [email protected]. Achim Zeileis added multicore support, and also fixed the default method to work for many more regression models.

References

Davison, A, and Hinkley, D. (1997) Bootstrap Methods and their Applications. Oxford: Oxford University Press.

Fox, J. and Weisberg, S. (2019) Companion to Applied Regression, Third Edition. Thousand Oaks: Sage.

Fox, J. and Weisberg, S. (2019) Bootstrapping Regression Models in R, https://www.john-fox.ca/Companion/appendices/Appendix-Bootstrapping.pdf.

Weisberg, S. (2014) Applied Linear Regression, Fourth Edition, Wiley Wiley, Chapters 4 and 11.

See Also

Functions that work with boot objects from the boot package are boot.array, boot.ci, plot.boot and empinf. Additional functions in the car package are summary.boot, confint.boot, and hist.boot.

Examples

m1 <- lm(Fertility ~ ., swiss)
betahat.boot <- Boot(m1, R=199) # 199 bootstrap samples--too small to be useful
summary(betahat.boot)  # default summary
confint(betahat.boot)
hist(betahat.boot)
# Bootstrap for the estimated residual standard deviation:
sigmahat.boot <- Boot(m1, R=199, f=sigmaHat, labels="sigmaHat")
summary(sigmahat.boot)
confint(sigmahat.boot)

Description

Computes and optionally plots profile log-likelihoods for the parameter of the Box-Cox power family, the Yeo-Johnson power family, or for either of the parameters in a bcnPower family. This is a slight generalization of the boxcox function in the MASS package that allows for families of transformations other than the Box-Cox power family. the boxCox2d function produces a contour plot of the two-dimensional likelihood profile for the bcnPower family.

Usage

boxCox(object, ...)

## Default S3 method:
boxCox(object,
        lambda = seq(-2, 2, 1/10), plotit = TRUE,
        interp = plotit, eps = 1/50,
        xlab=NULL, ylab=NULL, main= "Profile Log-likelihood",
        family="bcPower",
        param=c("lambda", "gamma"), gamma=NULL,
        grid=TRUE, ...)

## S3 method for class 'formula'
boxCox(object, lambda = seq(-2, 2, 1/10), plotit = TRUE, family = "bcPower",
    param = c("lambda", "gamma"), gamma = NULL, grid = TRUE,
    ...)

## S3 method for class 'lm'
boxCox(object, lambda = seq(-2, 2, 1/10), plotit = TRUE, ...)

boxCox2d(x, ksds = 4, levels = c(0.5, 0.95, 0.99, 0.999),
                 main = "bcnPower Log-likelihood", grid=TRUE, ...)

Arguments

Details

The boxCox function is an elaboration of the boxcox function in the MASS package. The first 7 arguments are the same as in boxcox, and if the argument family="bcPower" is used, the result is essentially identical to the function in MASS. Two additional families are the yjPower and bcnPower families that allow a few values of the response to be non-positive. The bcnPower family has two parameters: a power $\lambda$ and a start or location parameter $\gamma$, and the boxCox function can be used to obtain a profile log-likelihood for either parameter with $\lambda$ as the default. Alternatively, the boxCox2d function can be used to get a contour plot of the profile log-likelihood.

Value

Both functions ae designed for their side effects of drawing a graph. The boxCox function returns a list of the lambda (or possibly, gamma) vector and the computed profile log-likelihood vector, invisibly if the result is plotted. If plotit=TRUE plots log-likelihood vs lambda and indicates a 95% confidence interval about the maximum observed value of lambda. If interp=TRUE, spline interpolation is used to give a smoother plot.

Author(s)

Sanford Weisberg, <[email protected]>

References

Box, G. E. P. and Cox, D. R. (1964) An analysis of transformations. Journal of the Royal Statisistical Society, Series B. 26 211-46.

Cook, R. D. and Weisberg, S. (1999) Applied Regression Including Computing and Graphics. Wiley.

Fox, J. (2016) Applied Regression Analysis and Generalized Linear Models, Third Edition. Sage.

Fox, J. and Weisberg, S. (2019) An R Companion to Applied Regression, Third Edition, Sage.

Hawkins, D. and Weisberg, S. (2017) Combining the Box-Cox Power and Generalized Log Transformations to Accomodate Nonpositive Responses In Linear and Mixed-Effects Linear Models South African Statistics Journal, 51, 317-328.

Weisberg, S. (2014) Applied Linear Regression, Fourth Edition, Wiley.

Yeo, I. and Johnson, R. (2000) A new family of power transformations to improve normality or symmetry. Biometrika, 87, 954-959.

See Also

boxcox, yjPower, bcPower, bcnPower, powerTransform, contour

Examples

with(trees, boxCox(Volume ~ log(Height) + log(Girth), data = trees,
         lambda = seq(-0.25, 0.25, length = 10)))

  data("quine", package = "MASS")
  with(quine, boxCox(Days ~ Eth*Sex*Age*Lrn, 
         lambda = seq(-0.05, 0.45, len = 20), family="yjPower"))

Description

Computes a constructed variable for the Box-Cox transformation of the response variable in a linear model.

Usage

boxCoxVariable(y)

Arguments

Details

The constructed variable is defined as $y[\log(y/\widetilde{y}) - 1]$, where $\widetilde{y}$ is the geometric mean of y.

The constructed variable is meant to be added to the right-hand-side of the linear model. The t-test for the coefficient of the constructed variable is an approximate score test for whether a transformation is required.

If $b$ is the coefficient of the constructed variable, then an estimate of the normalizing power transformation based on the score statistic is $1 - b$. An added-variable plot for the constructed variable shows leverage and influence on the decision to transform y.

Value

a numeric vector of the same length as y.

Author(s)

John Fox [email protected]

References

Atkinson, A. C. (1985) Plots, Transformations, and Regression. Oxford.

Box, G. E. P. and Cox, D. R. (1964) An analysis of transformations. JRSS B 26 211–246.

Fox, J. (2016) Applied Regression Analysis and Generalized Linear Models, Third Edition. Sage.

Fox, J. and Weisberg, S. (2019) An R Companion to Applied Regression, Third Edition, Sage.

See Also

boxcox, powerTransform, bcPower

Examples

mod <- lm(interlocks + 1 ~ assets, data=Ornstein)
mod.aux <- update(mod, . ~ . + boxCoxVariable(interlocks + 1))
summary(mod.aux)
# avPlots(mod.aux, "boxCoxVariable(interlocks + 1)")

Description

Boxplot is a wrapper for the standard R boxplot function, providing point identification, axis labels, and a formula interface for boxplots without a grouping variable.

Usage

Boxplot(y, ...)

## Default S3 method:
Boxplot(y, g, id=TRUE, xlab, ylab, ...)

## S3 method for class 'formula'
Boxplot(formula, data=NULL, subset, na.action=NULL, 
    id=TRUE, xlab, ylab, ...)

## S3 method for class 'list'
Boxplot(y, xlab="", ylab="", ...)

## S3 method for class 'data.frame'
Boxplot(y, id=TRUE, ...)

## S3 method for class 'matrix'
Boxplot(y, ...)

Arguments

Author(s)

John Fox [email protected], with a contribution from Steve Ellison to handle at argument (see boxplot).

References

Fox, J. and Weisberg, S. (2019) An R Companion to Applied Regression, Third Edition, Sage.

See Also

boxplot

Examples

Boxplot(~income, data=Prestige, id=list(n=Inf)) # identify all outliers
Boxplot(income ~ type, data=Prestige)
Boxplot(income ~ type, data=Prestige, at=c(1, 3, 2))
Boxplot(k5 + k618 ~ lfp*wc, data=Mroz)
with(Prestige, Boxplot(income, id=list(labels=rownames(Prestige))))
with(Prestige, Boxplot(income, type, id=list(labels=rownames(Prestige))))
Boxplot(scale(Prestige[, 1:4]))

Description

Computes the Box-Tidwell power transformations of the predictors in a linear model.

Usage

boxTidwell(y, ...)

## S3 method for class 'formula'
boxTidwell(formula, other.x=NULL, data=NULL, subset, 
  na.action=getOption("na.action"), verbose=FALSE, tol=0.001, 
  max.iter=25, ...)

## Default S3 method:
boxTidwell(y, x1, x2=NULL, max.iter=25, tol=0.001, 
  verbose=FALSE, ...)
  
## S3 method for class 'boxTidwell'
print(x, digits=getOption("digits") - 2, ...)

Arguments

Details

The maximum-likelihood estimates of the transformation parameters are computed by Box and Tidwell's (1962) method, which is usually more efficient than using a general nonlinear least-squares routine for this problem. Score tests for the transformations are also reported.

Value

an object of class boxTidwell, which is normally just printed.

Author(s)

John Fox [email protected]

References

Box, G. E. P. and Tidwell, P. W. (1962) Transformation of the independent variables. Technometrics 4, 531-550.

Fox, J. (2016) Applied Regression Analysis and Generalized Linear Models, Third Edition. Sage.

Fox, J. and Weisberg, S. (2019) An R Companion to Applied Regression, Third Edition, Sage.

Examples

boxTidwell(prestige ~ income + education, ~ type + poly(women, 2), data=Prestige)

Description

Print data objects and statistical model summaries in abbreviated form.

Usage

brief(object, ...)

## S3 method for class 'data.frame'
brief(object, rows = if (nr <= 10) c(nr, 0) else c(3, 2),
    cols, head=FALSE, tail=FALSE, elided = TRUE,
    classes = inherits(object, "data.frame"), ...)
## S3 method for class 'tbl'
brief(object, ...)
## S3 method for class 'matrix'
brief(object, rows = if (nr <= 10) c(nr, 0) else c(3, 2), ...)

## S3 method for class 'numeric'
brief(object, rows = c(2, 1), elided = TRUE, ...)
## S3 method for class 'integer'
brief(object, rows = c(2, 1), elided = TRUE, ...)
## S3 method for class 'character'
brief(object, rows = c(2, 1), elided = TRUE, ...)
## S3 method for class 'factor'
brief(object, rows=c(2, 1), elided=TRUE, ...)

## S3 method for class 'list'
brief(object, rows = c(2, 1), elided = TRUE, ...)

## S3 method for class 'function'
brief(object, rows = c(5, 3), elided = TRUE, ...)

## S3 method for class 'lm'
brief(object,  terms = ~ .,
    intercept=missing(terms), pvalues=FALSE,
    digits=3, horizontal=TRUE, vcov., ...)

## S3 method for class 'glm'
brief(object, terms = ~ .,
    intercept=missing(terms), pvalues=FALSE,
    digits=3, horizontal=TRUE, vcov., dispersion, exponentiate, ...)

## S3 method for class 'multinom'
brief(object, terms = ~ .,
    intercept=missing(terms), pvalues=FALSE,
    digits=3, horizontal=TRUE, exponentiate=TRUE, ...)

## S3 method for class 'polr'
brief(object, terms = ~ .,
    intercept, pvalues=FALSE,
    digits=3, horizontal=TRUE, exponentiate=TRUE, ...)

## Default S3 method:
brief(object, terms = ~ .,
    intercept=missing(terms), pvalues=FALSE,
    digits=3, horizontal=TRUE, ...)

Arguments

Value

Invisibly returns object for a data object, or summary for a model object.

Note

The method brief.matrix calls brief.data.frame, and brief.tbl (for tibbles) calls print.

Author(s)

John Fox [email protected]

References

Fox, J. and Weisberg, S. (2019) An R Companion to Applied Regression, Third Edition, Sage.

See Also

S

Examples

brief(rnorm(100))
brief(Duncan)
brief(OBrienKaiser, elided=TRUE)
brief(matrix(1:500, 10, 50))
brief(lm)

mod.prestige <- lm(prestige ~ education + income + type, Prestige)
brief(mod.prestige, pvalues=TRUE)
brief(mod.prestige, ~ type)
mod.mroz <- glm(lfp ~ ., data=Mroz, family=binomial)
brief(mod.mroz)

Description

These functions are were deprecated in 2009 and are now defunct.

Usage

av.plot(...)
av.plots(...)
box.cox(...)
bc(...)
box.cox.powers(...)
box.cox.var(...)
box.tidwell(...)
cookd(...)
confidence.ellipse(...)
ceres.plot(...)
ceres.plots(...)
cr.plot(...)
cr.plots(...)
data.ellipse(...)
durbin.watson(...)
levene.test(...)
leverage.plot(...)
leverage.plots(...)
linear.hypothesis(...)
ncv.test(...)
outlier.test(...)
qq.plot(...)
skewPower(...)
spread.level.plot(...)

Arguments

Details

av.plot and av.plots are replaced by avPlot and avPlots functions.

box.cox and bc are now replaced by bcPower.

box.cox.powers is replaced by powerTransform.

box.cox.var is replaced by boxCoxVariable.

box.tidwell is replaced by boxTidwell.

cookd is replaced by cooks.distance in the stats package.

confidence.ellipse is replaced by confidenceEllipse.

ceres.plot and ceres.plots are now replaced by the ceresPlot and ceresPlots functions.

cr.plot and cr.plots are now replaced by the crPlot and crPlots functions.

data.ellipse is replaced by dataEllipse.

durbin.watson is replaced by durbinWatsonTest.

levene.test is replaced by leveneTest function.

leverage.plot and leverage.plots are now replaced by the leveragePlot and leveragePlots functions.

linear.hypothesis is replaced by the linearHypothesis function.

ncv.test is replaced by ncvTest.

outlier.test is replaced by outlierTest.

qq.plot is replaced by qqPlot.

skewPower is replaced by bcnPower.

spread.level.plot is replaced by spreadLevelPlot.

Description

These functions are provided for compatibility with older versions of the car package only, and may be removed eventually. Commands that worked in versions of the car package prior to version 3.0-0 will not necessarily work in version 3.0-0 and beyond, or may not work in the same manner.

Usage

bootCase(...)

nextBoot(...)

Arguments

Details

These functions are replaced by Boot.

See Also

See Also Boot

Description

These objects (currently only the .carEnv environment) are exported for technical reasons and are not for direct use.

Author(s)

John Fox [email protected]

Description

Open the official hex sticker for the car package in your browser

Usage

carHexsticker()

Value

Used for its side effect of openning the hex sticker for the car package in your browser.

Author(s)

John Fox [email protected]

Examples

## Not run:  # intended for interactive use
carHexsticker()

## End(Not run)

Description

This function is used to set or retrieve colors to be used in car package graphics functions.

Usage

carPalette(palette)

Arguments

Details

This function sets or returns the value of options(carPalette=pallete) that will be use in car graphics functions to determine colors. The default is c("black", "blue", "magenta", "cyan", "orange", "gray", "green3", "red")), which is nearly a permutation of the colors returned by the standard palette function that minimizes the use of red and green in the same graph, and that substitutes orange for the often hard to see yellow.

Value

Invisibly returns the previous value of the car palette.

Author(s)

Sanford Weisberg and John Fox

References

Fox, J. and Weisberg, S. (2019) An R Companion to Applied Regression, Third Edition, Sage.

See Also

palette, colors

Examples

# Standard color palette
palette()
# car standard color palette
carPalette()
# set colors to all black
carPalette(rep("black", 8))
# Use a custom color palette with 12 distinct colors
carPalette(sample(colors(distinct=TRUE), 12))
# restore default
carPalette("default")

Description

This function will access the website for An R Companion to Applied Regression.

Usage

carWeb(...)

Arguments

Value

Displays the web page for An R Companion to Applied Regression, Third Edition, Sage

Description

These functions draw Ceres plots for linear and generalized linear models.

Usage

ceresPlots(model, ...)

## Default S3 method:
ceresPlots(model, terms = ~., layout = NULL, ask, main,
    ...)

ceresPlot(model, ...)

## S3 method for class 'lm'
ceresPlot(model, variable, id=FALSE,
  line=TRUE,  smooth=TRUE, col=carPalette()[1], col.lines=carPalette()[-1],
  xlab, ylab, pch=1, lwd=2,  grid=TRUE, ...)

## S3 method for class 'glm'
ceresPlot(model, ...)

Arguments

Details

Ceres plots are a generalization of component+residual (partial residual) plots that are less prone to leakage of nonlinearity among the predictors.

The function intended for direct use is ceresPlots.

The model cannot contain interactions, but can contain factors. Factors may be present in the model, but Ceres plots cannot be drawn for them.

Value

NULL. These functions are used for their side effect: producing plots.

Author(s)

John Fox [email protected]

References

Cook, R. D. and Weisberg, S. (1999) Applied Regression, Including Computing and Graphics. Wiley.

Fox, J. (2016) Applied Regression Analysis and Generalized Linear Models, Third Edition. Sage.

Fox, J. and Weisberg, S. (2019) An R Companion to Applied Regression, Third Edition, Sage.

See Also

crPlots, avPlots, showLabels

Examples

ceresPlots(lm(prestige~income+education+type, data=Prestige), terms= ~ . - type)

Description

This function extracts estimates of regression parameters and their standard errors from one or more models and prints them in a table.

Usage

compareCoefs(..., se = TRUE, zvals = FALSE, pvals = FALSE, vcov.,
    print = TRUE, digits = 3)

Arguments

Value

This function is mainly used for its side-effect of printing the result. It also invisibly returns a matrix of estimates, standard errors, Wald statistics, and p-values.

Author(s)

John Fox [email protected]

References

Fox, J. and Weisberg, S. (2019) An R Companion to Applied Regression, Third Edition, Sage.

Examples

mod1 <- lm(prestige ~ income + education, data=Duncan)
mod2 <- update(mod1, subset=-c(6,16))
mod3 <- update(mod1, . ~ . + type)
mod4 <- update(mod1, . ~ . + I(income + education)) # aliased coef.
compareCoefs(mod1)
compareCoefs(mod1, mod2, mod4)
compareCoefs(mod1, mod2, mod3, zvals=TRUE, pvals=TRUE)
compareCoefs(mod1, mod2, se=FALSE)
compareCoefs(mod1, mod1, vcov.=list(vcov, hccm))

Description

These are substitutes for similarly named functions in the stats package (note the uppercase letter starting the second word in each function name). The only difference is that the contrast functions from the car package produce easier-to-read names for the contrasts when they are used in statistical models.

The functions and this documentation are adapted from the stats package.

Usage

contr.Treatment(n, base = 1, contrasts = TRUE)

contr.Sum(n, contrasts = TRUE)

contr.Helmert(n, contrasts = TRUE)

Arguments

Details

These functions are used for creating contrast matrices for use in fitting analysis of variance and regression models. The columns of the resulting matrices contain contrasts which can be used for coding a factor with n levels. The returned value contains the computed contrasts. If the argument contrasts is FALSE then a square matrix is returned.

Several aspects of these contrast functions are controlled by options set via the options command:

decorate.contrasts

This option should be set to a 2-element character vector containing the prefix and suffix characters to surround contrast names. If the option is not set, then c("[", "]") is used. For example, setting options(decorate.contrasts=c(".", "")) produces contrast names that are separated from factor names by a period. Setting options( decorate.contrasts=c("", "")) reproduces the behaviour of the R base contrast functions.

decorate.contr.Treatment

A character string to be appended to contrast names to signify treatment contrasts; if the option is unset, then "T." is used.

decorate.contr.Sum

Similar to the above, with default "S.".

decorate.contr.Helmert

Similar to the above, with default "H.".

contr.Sum.show.levels

Logical value: if TRUE (the default if unset), then level names are used for contrasts; if FALSE, then numbers are used, as in contr.sum in the base package.

Note that there is no replacement for contr.poly in the base package (which produces orthogonal-polynomial contrasts) since this function already constructs easy-to-read contrast names.

Value

A matrix with n rows and k columns, with k = n - 1 if contrasts is TRUE and k = n if contrasts is FALSE.

Author(s)

John Fox [email protected]

References

Fox, J. and Weisberg, S. (2019) An R Companion to Applied Regression, Third Edition, Sage.

See Also

contr.treatment, contr.sum, contr.helmert, contr.poly

Examples

# contr.Treatment vs. contr.treatment in the base package:

lm(prestige ~ (income + education)*type, data=Prestige, 
    contrasts=list(type="contr.Treatment"))

##  Call:
##  lm(formula = prestige ~ (income + education) * type, data = Prestige,
##      contrasts = list(type = "contr.Treatment"))
##  
##  Coefficients:
##          (Intercept)                  income               education  
##              2.275753                0.003522                1.713275  
##          type[T.prof]              type[T.wc]     income:type[T.prof]  
##              15.351896              -33.536652               -0.002903  
##      income:type[T.wc]  education:type[T.prof]    education:type[T.wc]  
##              -0.002072                1.387809                4.290875  

lm(prestige ~ (income + education)*type, data=Prestige, 
    contrasts=list(type="contr.treatment"))    

##  Call:
##  lm(formula = prestige ~ (income + education) * type, data = Prestige,
##      contrasts = list(type = "contr.treatment"))
##  
##  Coefficients:
##      (Intercept)              income           education  
##          2.275753            0.003522            1.713275  
##          typeprof              typewc     income:typeprof  
##          15.351896          -33.536652           -0.002903  
##      income:typewc  education:typeprof    education:typewc  
##          -0.002072            1.387809            4.290875

Description

These functions construct component+residual plots, also called partial-residual plots, for linear and generalized linear models.

Usage

crPlots(model, ...)

## Default S3 method:
crPlots(model, terms = ~., layout = NULL, ask, main, 
    ...)

crp(...)

crPlot(model, ...)

## S3 method for class 'lm'
crPlot(model, variable, id=FALSE,
    order=1, line=TRUE, smooth=TRUE, 
    col=carPalette()[1], col.lines=carPalette()[-1],
    xlab, ylab, pch=1, lwd=2, grid=TRUE, ...)
    
crPlot3d(model, var1, var2, ...)

## S3 method for class 'lm'
crPlot3d(model, var1, var2,
    xlab = var1,
    ylab = paste0("C+R(", eff$response, ")"), zlab = var2,
    axis.scales = TRUE, axis.ticks = FALSE, revolutions = 0,
    bg.col = c("white", "black"),
    axis.col = if (bg.col == "white") c("darkmagenta", "black", "darkcyan") 
        else c("darkmagenta", "white", "darkcyan"),
    surface.col = carPalette()[2:3], surface.alpha = 0.5,
    point.col = "yellow", text.col = axis.col,
    grid.col = if (bg.col ==  "white") "black" else "gray",
    fogtype = c("exp2", "linear", "exp", "none"),
    fill = TRUE, grid = TRUE, grid.lines = 26,
    smoother = c("loess", "mgcv", "none"), df.mgcv = NULL, loess.args = NULL,
    sphere.size = 1, radius = 1, threshold = 0.01, speed = 1, fov = 60,
    ellipsoid = FALSE, level = 0.5, ellipsoid.alpha = 0.1,
    id = FALSE, 
    mouseMode=c(none="none", left="polar", right="zoom", middle="fov", 
	               wheel="pull"),
    ...)

Arguments

Details

The functions intended for direct use are crPlots, for which crp is an abbreviation, and, for 3D C+R plots, crPlot3d.

For 2D plots, the model cannot contain interactions, but can contain factors. Parallel boxplots of the partial residuals are drawn for the levels of a factor. crPlot3d can handle models with two-way interactions.

For 2D C+R plots, the fit is represented by a broken blue line and a smooth of the partial residuals by a solid magenta line. For 3D C+R plots, the fit is represented by a blue surface and a smooth of the partial residuals by a magenta surface.

Value

These functions are used for their side effect of producing plots, but also invisibly return the coordinates of the plotted points.

Author(s)

John Fox [email protected]

References

Cook, R. D. and Weisberg, S. (1999) Applied Regression, Including Computing and Graphics. Wiley.

Fox, J. (2016) Applied Regression Analysis and Generalized Linear Models, Third Edition. Sage.

Fox, J. and Weisberg, S. (2019) An R Companion to Applied Regression, Third Edition, Sage.

See Also

ceresPlots, avPlots

Examples

crPlots(m<-lm(prestige ~ income + education, data=Prestige)) 

crPlots(m, terms=~ . - education) # get only one plot

crPlots(lm(prestige ~ log2(income) + education + poly(women,2), data=Prestige))

crPlots(glm(partic != "not.work" ~ hincome + children, 
  data=Womenlf, family=binomial), smooth=list(span=0.75))

# 3D C+R plot, requires the rgl, effects, and mgcv packages
if (interactive() && require(rgl) && require(effects) && require(mgcv)){
crPlot3d(lm(prestige ~ income*education + women, data=Prestige), 
    "income", "education") 
}

Description

deltaMethod is a generic function that uses the delta method to get a first-order approximate standard error for a nonlinear function of a vector of random variables with known or estimated covariance matrix.

Usage

deltaMethod(object, ...)

## Default S3 method:
deltaMethod(object, g., vcov., func=g., constants, level=0.95, 
  rhs, singular.ok=FALSE, ..., envir=parent.frame())
## S3 method for class 'lm'
 deltaMethod(object, g., vcov.=vcov(object, complete=FALSE), 
           parameterNames=names(coef(object)), ..., envir=parent.frame())
## S3 method for class 'nls'
deltaMethod(object, g., vcov.=vcov(object, complete=FALSE), ..., envir=parent.frame())
## S3 method for class 'multinom'
 deltaMethod(object, g., vcov. = vcov(object, complete=FALSE), 
           parameterNames = if (is.matrix(coef(object))) 
           colnames(coef(object)) else names(coef(object)), ..., envir=parent.frame()) 
## S3 method for class 'polr'
 deltaMethod(object, g., vcov.=vcov(object, complete=FALSE), ..., envir=parent.frame())
## S3 method for class 'survreg'
 deltaMethod(object, g., vcov. = vcov(object, complete=FALSE), 
           parameterNames = names(coef(object)), ..., envir=parent.frame())
## S3 method for class 'coxph'
 deltaMethod(object, g., vcov. = vcov(object, complete=FALSE), 
           parameterNames = names(coef(object)), ..., envir=parent.frame())
## S3 method for class 'mer'
 deltaMethod(object, g., vcov. = vcov(object, complete=FALSE),
           parameterNames = names(fixef(object)), ..., envir=parent.frame())
## S3 method for class 'merMod'
 deltaMethod(object, g., vcov. = vcov(object, complete=FALSE),
           parameterNames = names(fixef(object)), ..., envir=parent.frame())
## S3 method for class 'lme'
 deltaMethod(object, g., vcov. = vcov(object, complete=FALSE),
           parameterNames = names(fixef(object)), ..., envir=parent.frame())
## S3 method for class 'lmList'
 deltaMethod(object, g.,  ..., envir=parent.frame())

Arguments

Details

Suppose $x$ is a random vector of length $p$ that is at least approximately normally distributed with mean $\beta$ and estimated covariance matrix $C$. Then any function $g(\beta)$ of $\beta$, is estimated by $g(x)$, which is in large samples normally distributed with mean $g(\beta)$ and estimated variance $h'Ch$, where $h$ is the first derivative of $g(\beta)$ with respect to $\beta$ evaluated at $x$. This function returns both $g(x)$ and its standard error, the square root of the estimated variance.

The default method requires that you provide $x$ in the argument object, $C$ in the argument vcov., and a text expression in argument g. that when evaluated gives the function $g$. The call names(object) must return the names of the elements of x that are used in the expression g..

Since the delta method is often applied to functions of regression parameter estimates, the argument object may be the name of a regression object from which the estimates and their estimated variance matrix can be extracted. In most regression models, estimates are returned by the coef(object) and the variance matrix from vcov(object). You can provide an alternative function for computing the sample variance matrix, for example to use a sandwich estimator.

For mixed models using lme4 or nlme, the coefficient estimates are returned by the fixef function, while for multinom, lmList and nlsList coefficient estimates are returned by coef as a matrix. Methods for these models are provided to get the correct estimates and variance matrix.

The argument g. must be a quoted character string that gives the function of interest. For example, if you set m2 <- lm(Y ~ X1 + X2 + X1:X2), then deltaMethod(m2,"X1/X2") applies the delta method to the ratio of the coefficient estimates for X1 and X2. The argument g. can consist of constants and names associated with the elements of the vector of coefficient estimates.

In some cases the names may include characters such as the colon : used in interactions, or mathematical symbols like + or - signs that would confuse the function that computes numerical derivatives, and for this case you can replace the names of the estimates with the parameterNames argument. For example, the ratio of the X2 main effect to the interaction term could be computed using deltaMethod(m2, "b1/b3", parameterNames=c("b0", "b1", "b2", "b3")). The name “(Intercept)” used for the intercept in linear and generalized linear models is an exception, and it will be correctly interpreted by deltaMethod. Another option is to use back-ticks to quote nonstandard R names, as in deltaMethod(m2,"X1/`X1:X2`").

For multinom objects, the coef function returns a matrix of coefficients, with each row giving the estimates for comparisons of one category to the baseline. The deltaMethod function applies the delta method to each row of this matrix. Similarly, for lmList and nlsList objects, the delta method is computed for each element of the list of models fit.

For nonlinear regression objects produced by the nls function, the call coef(object) returns the estimated coefficient vectors with names corresponding to parameter names. For example, m2 <- nls(y ~ theta/(1 + gamma * x), start = list(theta=2, gamma=3)) will have parameters named c("theta", "gamma"). In many other familiar regression models, such as those produced by lm and glm, the names of the coefficient estimates are the corresponding regressor names, not parameter names.

For mixed-effects models fit with lmer and glmer from the lme4 package or lme and nlme from the nlme package, only fixed-effect coefficients are considered.

For regression models for which methods are not provided, you can extract the named vector of coefficient estimates and and estimate of its covariance matrix and then apply the default deltaMethod function.

Note: Earlier versions of deltaMethod included an argument parameterPrefix that implemented the same functionality as the parameterNames argument, but which caused several problems that were not easily fixed without the change in syntax.

Value

An object of class "deltaMethod", inheriting from "data.frame", for which a print method is provided. The object contains columns named Estimate for the estimate, SE for its standard error, and columns for confidence limits and possibly a hypothesis test. The value of g. is given as a row label.

Author(s)

Sanford Weisberg, [email protected], John Fox [email protected], and Pavel Krivitsky.

References

Fox, J. (2016) Applied Regression Analysis and Generalized Linear Models, Third Edition. Sage.

Fox, J. and Weisberg, S. (2019) An R Companion to Applied Regression, Third Edition, Sage.

Weisberg, S. (2014) Applied Linear Regression, Fourth Edition, Wiley, Section 6.1.2.

See Also

First derivatives of g. are computed using symbolic differentiation by the function D.

Examples

m1 <- lm(time ~ t1 + t2, data = Transact) 
deltaMethod(m1, "b1/b2", parameterNames= paste("b", 0:2, sep="")) 
deltaMethod(m1, "t1/t2", rhs=1) # use names of preds. rather than coefs.
deltaMethod(m1, "t1/t2", vcov=hccm) # use hccm function to est. vars.
deltaMethod(m1, "1/(Intercept)")
# The next example invokes the default method by extracting the
# vector of estimates and covariance matrix explicitly
deltaMethod(coef(m1), "t1/t2", vcov.=vcov(m1))

Description

densityPlot contructs and graphs nonparametric density estimates, possibly conditioned on a factor, using the standard R density function or by default adaptiveKernel, which computes an adaptive kernel density estimate. depan provides the Epanechnikov kernel and dbiwt provides the biweight kernel.

Usage

densityPlot(x, ...)

## Default S3 method:
densityPlot(x, g, method=c("adaptive", "kernel"),
    bw=if (method == "adaptive") bw.nrd0 else "SJ", adjust=1,
    kernel, xlim, ylim,
    normalize=FALSE, xlab=deparse(substitute(x)), ylab="Density", main="",
    col=carPalette(), lty=seq_along(col), lwd=2, grid=TRUE,
    legend=TRUE, show.bw=FALSE, rug=TRUE, ...)

## S3 method for class 'formula'
densityPlot(formula, data=NULL, subset,
    na.action=NULL, xlab, ylab, main="", legend=TRUE, ...)

adaptiveKernel(x, kernel=dnorm, bw=bw.nrd0, adjust=1.0, n=500,
    from, to, cut=3, na.rm=TRUE)
    
depan(x)
dbiwt(x)

Arguments

Details

If you use a different kernel function than the default dnorm that has a standard deviation different from 1 along with an automatic rule like the default function bw.nrd0, you can attach an attribute to the kernel function named "scale" that gives its standard deviation. This is true for the two supplied kernels, depan and dbiwt

Value

densityPlot invisibly returns the "density" object computed (or list of "density" objects) and draws a graph. adaptiveKernel returns an object of class "density" (see density).

Author(s)

John Fox [email protected]

References

Fox, J. and Weisberg, S. (2019) An R Companion to Applied Regression, Third Edition, Sage.

W. N. Venables and B. D. Ripley (2002) Modern Applied Statistics with S. New York: Springer.

B.W. Silverman (1986) Density Estimation for Statistics and Data Analysis. London: Chapman and Hall.

See Also

density, bw.SJ, plot.density

Examples

densityPlot(~ income, show.bw=TRUE, method="kernel", data=Prestige)
densityPlot(~ income, show.bw=TRUE, data=Prestige)
densityPlot(~ income, from=0, normalize=TRUE, show.bw=TRUE, data=Prestige)

densityPlot(income ~ type, data=Prestige)
densityPlot(~ income, show.bw=TRUE, method="kernel", data=Prestige)
densityPlot(~ income, show.bw=TRUE, data=Prestige)
densityPlot(~ income, from=0, normalize=TRUE, show.bw=TRUE, data=Prestige)

densityPlot(income ~ type, kernel=depan, data=Prestige)
densityPlot(income ~ type, kernel=depan, legend=list(location="top"), data=Prestige)

plot(adaptiveKernel(UN$infantMortality, from=0, adjust=0.75), col="magenta")
lines(density(na.omit(UN$infantMortality), from=0, adjust=0.75), col="blue")
rug(UN$infantMortality, col="cyan")
legend("topright", col=c("magenta", "blue"), lty=1,
       legend=c("adaptive kernel", "kernel"), inset=0.02)


plot(adaptiveKernel(UN$infantMortality, from=0, adjust=0.75), col="magenta")
lines(density(na.omit(UN$infantMortality), from=0, adjust=0.75), col="blue")
rug(UN$infantMortality, col="cyan")
legend("topright", col=c("magenta", "blue"), lty=1,
       legend=c("adaptive kernel", "kernel"), inset=0.02)

Description

These functions display index plots of dfbeta (effect on coefficients of deleting each observation in turn) and dfbetas (effect on coefficients of deleting each observation in turn, standardized by a deleted estimate of the coefficient standard error). In the plot of dfbeta, horizontal lines are drawn at 0 and +/- one standard error; in the plot of dfbetas, horizontal lines are drawn and 0 and +/- 1.

Usage

dfbetaPlots(model, ...)

dfbetasPlots(model, ...)

## S3 method for class 'lm'
dfbetaPlots(model, terms= ~ ., intercept=FALSE, layout=NULL, ask, 
	    main, xlab, ylab, labels=rownames(dfbeta), 
        id.method="y",  
        id.n=if(id.method[1]=="identify") Inf else 0, id.cex=1,
        id.col=carPalette()[1], id.location="lr", col=carPalette()[1], grid=TRUE, ...)

## S3 method for class 'lm'
dfbetasPlots(model, terms=~., intercept=FALSE, layout=NULL, ask, 
	    main, xlab, ylab,
        labels=rownames(dfbetas), id.method="y",
        id.n=if(id.method[1]=="identify") Inf else 0, id.cex=1, 
        id.col=carPalette()[1], id.location="lr", col=carPalette()[1], grid=TRUE, ...)

Arguments

.

Value

NULL. These functions are used for their side effect: producing plots.

Author(s)

John Fox [email protected]

References

Fox, J. (2016) Applied Regression Analysis and Generalized Linear Models, Third Edition. Sage.

Fox, J. and Weisberg, S. (2019) An R Companion to Applied Regression, Third Edition, Sage.

See Also

dfbeta ,dfbetas

Examples

dfbetaPlots(lm(prestige ~ income + education + type, data=Duncan))

dfbetasPlots(glm(partic != "not.work" ~ hincome + children, 
  data=Womenlf, family=binomial))

Description

Computes residual autocorrelations and generalized Durbin-Watson statistics and their bootstrapped p-values. dwt is an abbreviation for durbinWatsonTest.

Usage

durbinWatsonTest(model, ...)

dwt(...)

## S3 method for class 'lm'
durbinWatsonTest(model, max.lag=1, simulate=TRUE, reps=1000,
    method=c("resample","normal"),
    alternative=c("two.sided", "positive", "negative"), ...)

## Default S3 method:
durbinWatsonTest(model, max.lag=1, ...)

## S3 method for class 'durbinWatsonTest'
print(x, ...)

Arguments

Value

Returns an object of type "durbinWatsonTest".

Note

p-values are available only from the lm method.

Author(s)

John Fox [email protected]

References

Fox, J. (2016) Applied Regression Analysis and Generalized Linear Models, Third Edition. Sage.

Examples

durbinWatsonTest(lm(fconvict ~ tfr + partic + degrees + mconvict, data=Hartnagel))

Description

These functions draw ellipses, including data ellipses, and confidence ellipses for linear, generalized linear, and possibly other models.

Usage

ellipse(center, shape, radius, log="", center.pch=19, center.cex=1.5, 
  segments=51, draw=TRUE, add=draw, xlab="", ylab="", 
  col=carPalette()[2], lwd=2, fill=FALSE, fill.alpha=0.3, grid=TRUE, ...)

dataEllipse(x, ...)

## S3 method for class 'formula'
dataEllipse(formula, data, subset, weights, xlab, ylab, id=FALSE, ...)

## Default S3 method:
dataEllipse(x, y, groups, 
  group.labels=group.levels, ellipse.label, 
  weights, log="", levels=c(0.5, 0.95), center.pch=19, 
  center.cex=1.5, draw=TRUE, plot.points=draw, add=!plot.points, segments=51, 
  robust=FALSE, 
  xlab=deparse(substitute(x)), 
  ylab=deparse(substitute(y)), 
  col=if (missing(groups)) carPalette()[1:2] else carPalette()[1:length(group.levels)],
  pch=if (missing(groups)) 1 else seq(group.levels),
  lwd=2, fill=FALSE, fill.alpha=0.3, grid=TRUE, id=FALSE, 
  label.pos=NULL, label.cex=1.25, label.xpd=FALSE, ...)

confidenceEllipse(model, ...)

## S3 method for class 'lm'
confidenceEllipse(model, which.coef, vcov.=vcov, 
  L, levels=0.95, Scheffe=FALSE, dfn,
  center.pch=19, center.cex=1.5, segments=51, xlab, ylab, 
  col=carPalette()[2], lwd=2, fill=FALSE, fill.alpha=0.3, draw=TRUE, add=!draw,
  grid=TRUE, ...)
  
## S3 method for class 'glm'
confidenceEllipse(model, chisq, ...)

## S3 method for class 'mlm'
confidenceEllipse(model, xlab, ylab, which.coef=1:2, ...)

## Default S3 method:
confidenceEllipse(model, which.coef, vcov.=vcov,
  L, levels=0.95, Scheffe=FALSE,  dfn,
  center.pch=19, center.cex=1.5, segments=51, xlab, ylab, 
  col=carPalette()[2], lwd=2, fill=FALSE, fill.alpha=0.3, draw=TRUE, add=!draw, 
  grid=TRUE, ...)
  
confidenceEllipses(model, ...)

## Default S3 method:
confidenceEllipses(model, coefnames,  main, grid=TRUE, ...)

## S3 method for class 'mlm'
confidenceEllipses(model, coefnames, main, ...)

Arguments

.

Details

The ellipse is computed by suitably transforming a unit circle.

dataEllipse superimposes the normal-probability contours over a scatterplot of the data.

confidenceEllipses plots a matrix of all pairwise confidence ellipses; each panel of the matrix is created by confidenceEllipse.

Value

These functions are mainly used for their side effect of producing plots. For greater flexibility (e.g., adding plot annotations), however, ellipse returns invisibly the (x, y) coordinates of the calculated ellipse. dataEllipse and confidenceEllipse return invisibly the coordinates of one or more ellipses, in the latter instance a list named by levels; confidenceEllipses invisibly returns NULL.

Author(s)

Georges Monette, John Fox [email protected], and Michael Friendly.

References

Fox, J. (2016) Applied Regression Analysis and Generalized Linear Models, Third Edition. Sage.

Fox, J. and Weisberg, S. (2019) An R Companion to Applied Regression, Third Edition, Sage.

Monette, G. (1990) Geometry of multiple regression and 3D graphics. In Fox, J. and Long, J. S. (Eds.) Modern Methods of Data Analysis. Sage.

See Also

cov.trob, cov.wt, linearHypothesis.

Examples

dataEllipse(Duncan$income, Duncan$education, levels=0.1*1:9, 
    ellipse.label=0.1*1:9, lty=2, fill=TRUE, fill.alpha=0.1)

# formula method
dataEllipse(education ~ income | type, data=Duncan,
            levels = 0.68)

dataEllipse(education ~ income | type, data=Duncan, 
            levels=0.68, label.xpd=TRUE)

# vary label positions
dataEllipse(education ~ income | type, data=Duncan,
            levels = 0.68,
            fill = TRUE, fill.alpha = 0.1,
            label.pos = c(3, 2, 4),
            xlim = c(0, 100),
            ylim = c(2, 120))

    
confidenceEllipse(lm(prestige ~ income + education, data=Duncan), Scheffe=TRUE)

confidenceEllipse(lm(prestige ~ income + education, data=Duncan), vcov.=hccm)

confidenceEllipse(lm(prestige ~ income + education, data=Duncan), 
	L=c("income + education", "income - education"))
	
confidenceEllipses(lm(prestige ~ income + education + type, data=Duncan),
  fill=TRUE)
cov2cor(vcov(lm(prestige ~ income + education + type, 
  data=Duncan))) # correlations among coefficients

wts <- rep(1, nrow(Duncan))
wts[c(6, 16)] <- 0 # delete Minister, Conductor
with(Duncan, {
	dataEllipse(income, prestige, levels=0.68)
	dataEllipse(income, prestige, levels=0.68, robust=TRUE, 
	    plot.points=FALSE, col="green3")
	dataEllipse(income, prestige, weights=wts, levels=0.68, 
	    plot.points=FALSE, col="brown")
	dataEllipse(income, prestige, weights=wts, robust=TRUE, levels=0.68, 
		plot.points=FALSE, col="blue")
	})
    
with(Prestige, dataEllipse(income, education, type, 
    id=list(n=2, labels=rownames(Prestige)), pch=15:17,
    xlim=c(0, 25000), center.pch="+",
    group.labels=c("Blue Collar", "Professional", "White Collar"),
    ylim=c(5, 20), level=.95, fill=TRUE, fill.alpha=0.1))