estimability (version 1.5.1)

epredict: Estimability Enhancements for lm and Relatives

Description

These functions call the corresponding S3 predict methods in the stats package, but with a check for estimability of new predictions, and with appropriate actions for non-estimable cases.

Usage

# S3 method for lm
epredict(object, newdata, ..., 
    type = c("response", "terms", "matrix", "estimability"), 
    nonest.tol = 1e-8, nbasis = object$nonest)
# S3 method for glm
epredict(object, newdata, ..., 
    type = c("link", "response", "terms", "matrix", "estimability"), 
    nonest.tol = 1e-8, nbasis = object$nonest)
# S3 method for mlm
epredict(object, newdata, ..., 
    type = c("response", "matrix", "estimability"), 
    nonest.tol = 1e-8, nbasis = object$nonest)
    
eupdate(object, ...)

Value

The same as the result of a call to the predict method in the stats package, except rows or elements corresponding to non-estimable predictor combinations are set to NA. The value for type is "matrix" or "estimability" is explained under details.

Arguments

object

An object inheriting from lm

newdata

A data.frame containing predictor combinations for new predictions

...

Arguments passed to predict or update

nonest.tol

Tolerance used by is.estble to check estimability of new predictions

type

Character string specifying the desired result. See Details.

nbasis

Basis for the null space, e.g., a result of a call to nonest.basis. If nbasis is NULL, a basis is constructed from object.

Author

Russell V. Lenth <[e-mail]>

Details

If newdata is missing or object is not rank-deficient, this method passes its arguments directly to the same method in the stats library. In rank-deficient cases with newdata provided, each row of newdata is tested for estimability against the null basis provided in nbasis. Any non-estimable cases found are replaced with NAs.

The type argument is passed to predict when it is one of "response", "link", or "terms". With newdata present and type = "matrix", the model matrix for newdata is returned, with an attribute "estble" that is a logical vector of length nrow(newdata) indicating whether each row is estimable. With type = "estimability", just the logical vector is returned.

If you anticipate making several epredict calls with new data, it improves efficiency to either obtain the null basis and provide it in the call, or add it to object with the name "nonest" (perhaps via a call to eupdate).

eupdate is an S3 generic function with a method provided for "lm" objects. It updates the object according to any arguments in ..., then obtains the updated object's nonestimable basis and returns it in object$nonest.

See Also

predict.lm in the stats package; nonest.basis.

Examples

require("estimability")

# Fake data where x3 and x4 depend on x1, x2, and intercept
x1 <- -4:4
x2 <- c(-2,1,-1,2,0,2,-1,1,-2)
x3 <- 3*x1 - 2*x2
x4 <- x2 - x1 + 4
y <- 1 + x1 + x2 + x3 + x4 + c(-.5,.5,.5,-.5,0,.5,-.5,-.5,.5)

# Different orderings of predictors produce different solutions
mod1234 <- lm(y ~ x1 + x2 + x3 + x4)
mod4321 <- eupdate(lm(y ~ x4 + x3 + x2 + x1))
# (Estimability checking with mod4321 will be more efficient because
#  it will not need to recreate the basis)
mod4321$nonest

 
# test data:
testset <- data.frame( 
              x1 = c(3,  6,  6,  0,  0,  1), 
              x2 = c(1,  2,  2,  0,  0,  2), 
              x3 = c(7, 14, 14,  0,  0,  3), 
              x4 = c(2,  4,  0,  4,  0,  4))

# Look at predictions when we don't check estimability
suppressWarnings( # Disable the warning from stats::predict.lm
    rbind(p1234 = predict(mod1234, newdata = testset),
          p4321 = predict(mod4321, newdata = testset)))

# Compare with results when we do check:
rbind(p1234 = epredict(mod1234, newdata = testset),
      p4321 = epredict(mod4321, newdata = testset))

# stats::predict has same capability for lm objects starting in version 4.3.0:
if((R.Version()$major >= 4) && (R.Version()$minor >= 3))
  stats::predict(mod1234, newdata = testset, rankdeficient = "NA")

# Note that estimable cases have the same predictions

# change mod1234 and include nonest basis 
mod134 <- eupdate(mod1234, . ~ . - x2, subset = -c(3, 7))
mod134$nonest

# When row spaces are the same, bases are interchangeable
# so long as you account for the ordering of parameters:
epredict(mod4321, newdata = testset, type = "estimability",
    nbasis = nonest.basis(mod1234)[c(1,5:2), ])
    
# Comparison with predict.lm in R >= 4.3.0
if((R.Version()$major >= 4) && (R.Version()$minor >= 3))
  stats::predict(mod4321, newdata = testset, rankdeficient = "NA")

if (FALSE) {
### Additional illustration
example(nonest.basis)  ## creates model objects warp.lm1 and warp.lm2

# The two models have different contrast specs. But the empty cell
# is correctly identified in both:
fac.cmb <- expand.grid(wool = c("A", "B"), tension = c("L", "M", "H"))
cbind(fac.cmb, 
      pred1 = epredict(warp.lm1, newdata = fac.cmb), 
      pred2 = epredict(warp.lm2, newdata = fac.cmb))
}