Skip to contents

Confidence intervals around predicted values

Usage

get_predicted_ci(x, predictions = NULL, ...)

# S3 method for default
get_predicted_ci(
  x,
  predictions = NULL,
  data = NULL,
  se = NULL,
  ci = 0.95,
  ci_type = "confidence",
  ci_method = NULL,
  dispersion_method = "sd",
  vcov = NULL,
  vcov_args = NULL,
  ...
)

Arguments

x

A statistical model (can also be a data.frame, in which case the second argument has to be a model).

predictions

A vector of predicted values (as obtained by stats::fitted(), stats::predict() or get_predicted()).

...

Other argument to be passed, for instance to get_predicted_ci().

data

An optional data frame in which to look for variables with which to predict. If omitted, the data used to fit the model is used. Visualization matrices can be generated using get_datagrid().

se

Numeric vector of standard error of predicted values. If NULL, standard errors are calculated based on the variance-covariance matrix.

ci

The interval level. Default is NULL, to be fast even for larger models. Set the interval level to an explicit value, e.g. 0.95, for 95% CI).

ci_type

Can be "prediction" or "confidence". Prediction intervals show the range that likely contains the value of a new observation (in what range it would fall), whereas confidence intervals reflect the uncertainty around the estimated parameters (and gives the range of the link; for instance of the regression line in a linear regressions). Prediction intervals account for both the uncertainty in the model's parameters, plus the random variation of the individual values. Thus, prediction intervals are always wider than confidence intervals. Moreover, prediction intervals will not necessarily become narrower as the sample size increases (as they do not reflect only the quality of the fit). This applies mostly for "simple" linear models (like lm), as for other models (e.g., glm), prediction intervals are somewhat useless (for instance, for a binomial model for which the dependent variable is a vector of 1s and 0s, the prediction interval is... [0, 1]).

ci_method

The method for computing p values and confidence intervals. Possible values depend on model type.

  • NULL uses the default method, which varies based on the model type.

  • Most frequentist models: "normal" (default).

  • Bayesian models: "quantile" (default), "hdi", "eti", and "spi".

  • Mixed effects lme4 models: "normal" (default), "satterthwaite", and "kenward".

dispersion_method

Bootstrap dispersion and Bayesian posterior summary: "sd" or "mad".

vcov

Variance-covariance matrix used to compute uncertainty estimates (e.g., for robust standard errors). This argument accepts a covariance matrix, a function which returns a covariance matrix, or a string which identifies the function to be used to compute the covariance matrix.

  • A covariance matrix

  • A function which returns a covariance matrix (e.g., stats::vcov())

  • A string which indicates the kind of uncertainty estimates to return.

    • Heteroskedasticity-consistent: "vcovHC", "HC", "HC0", "HC1", "HC2", "HC3", "HC4", "HC4m", "HC5". See ?sandwich::vcovHC

    • Cluster-robust: "vcovCR", "CR0", "CR1", "CR1p", "CR1S", "CR2", "CR3". See ?clubSandwich::vcovCR()

    • Bootstrap: "vcovBS", "xy", "residual", "wild", "mammen", "webb". See ?sandwich::vcovBS

    • Other sandwich package functions: "vcovHAC", "vcovPC", "vcovCL", "vcovPL".

vcov_args

List of arguments to be passed to the function identified by the vcov argument. This function is typically supplied by the sandwich or clubSandwich packages. Please refer to their documentation (e.g., ?sandwich::vcovHAC) to see the list of available arguments.

Examples

# Confidence Intervals for Model Predictions
# ------------------------------------------

data(mtcars)

# Linear model
# ------------
x <- lm(mpg ~ cyl + hp, data = mtcars)
predictions <- predict(x)
ci_vals <- get_predicted_ci(x, predictions, ci_type = "prediction")
head(ci_vals)
#>         SE    CI_low  CI_high
#> 1 3.255505 14.558527 27.87504
#> 2 3.255505 14.558527 27.87504
#> 3 3.305931 19.309850 32.83263
#> 4 3.255505 14.558527 27.87504
#> 5 3.303717  8.687625 22.20134
#> 6 3.266957 14.630713 27.99407
ci_vals <- get_predicted_ci(x, predictions, ci_type = "confidence")
head(ci_vals)
#>          SE   CI_low  CI_high
#> 1 0.7281647 19.72752 22.70605
#> 2 0.7281647 19.72752 22.70605
#> 3 0.9279509 24.17337 27.96911
#> 4 0.7281647 19.72752 22.70605
#> 5 0.9200310 13.56281 17.32616
#> 6 0.7777664 19.72168 22.90310
ci_vals <- get_predicted_ci(x, predictions, ci = c(0.8, 0.9, 0.95))
head(ci_vals)
#>                          SE CI_low_0.8 CI_high_0.8 CI_low_0.9 CI_high_0.9
#> Mazda RX4         0.7281647   20.26184    22.17172   19.97954    22.45403
#> Mazda RX4 Wag     0.7281647   20.26184    22.17172   19.97954    22.45403
#> Datsun 710        0.9279509   24.85429    27.28818   24.49453    27.64794
#> Hornet 4 Drive    0.7281647   20.26184    22.17172   19.97954    22.45403
#> Hornet Sportabout 0.9200310   14.23793    16.65104   13.88124    17.00773
#> Valiant           0.7777664   20.29240    22.33238   19.99087    22.63391
#>                   CI_low_0.95 CI_high_0.95
#> Mazda RX4            19.72752     22.70605
#> Mazda RX4 Wag        19.72752     22.70605
#> Datsun 710           24.17337     27.96911
#> Hornet 4 Drive       19.72752     22.70605
#> Hornet Sportabout    13.56281     17.32616
#> Valiant              19.72168     22.90310

# Bootstrapped
# ------------
if (require("boot")) {
  predictions <- get_predicted(x, iterations = 500)
  get_predicted_ci(x, predictions)
}
#>           SE   CI_low  CI_high
#> 1  0.6778386 20.06367 22.73426
#> 2  0.6778386 20.06367 22.73426
#> 3  1.1213280 24.07299 28.19569
#> 4  0.6778386 20.06367 22.73426
#> 5  0.7137549 13.96361 16.74841
#> 6  0.7196399 20.11708 22.92205
#> 7  0.9616939 11.54478 15.32004
#> 8  1.1099003 24.70988 29.07236
#> 9  1.1281547 24.00743 28.15526
#> 10 0.5968210 19.90476 22.18073
#> 11 0.5968210 19.90476 22.18073
#> 12 0.6867313 13.94831 16.57341
#> 13 0.6867313 13.94831 16.57341
#> 14 0.6867313 13.94831 16.57341
#> 15 0.6619519 13.42214 15.93349
#> 16 0.7053541 13.08691 15.77129
#> 17 0.8151358 12.40637 15.55107
#> 18 1.1012410 24.63705 28.94998
#> 19 1.1440413 24.87906 29.34752
#> 20 1.1031298 24.65526 28.99487
#> 21 1.1356722 23.94017 28.10764
#> 22 0.9243855 14.17452 17.75612
#> 23 0.9243855 14.17452 17.75612
#> 24 0.9616939 11.54478 15.32004
#> 25 0.7137549 13.96361 16.74841
#> 26 1.1012410 24.63705 28.94998
#> 27 1.1152049 24.14757 28.24896
#> 28 1.2189054 23.41192 27.86701
#> 29 1.1785211 10.35989 15.04183
#> 30 0.8033095 18.06171 21.24781
#> 31 2.1096285  6.00274 14.50050
#> 32 1.1944721 23.58991 27.93128

if (require("datawizard") && require("bayestestR")) {
  ci_vals <- get_predicted_ci(x, predictions, ci = c(0.80, 0.95))
  head(ci_vals)
  datawizard::reshape_ci(ci_vals)

  ci_vals <- get_predicted_ci(x,
    predictions,
    dispersion_method = "MAD",
    ci_method = "HDI"
  )
  head(ci_vals)
}
#>          SE   CI_low  CI_high
#> 1 0.6511853 20.08775 22.75861
#> 2 0.6511853 20.08775 22.75861
#> 3 1.1886593 24.12011 28.22016
#> 4 0.6511853 20.08775 22.75861
#> 5 0.6922916 13.92620 16.69426
#> 6 0.7042060 20.11663 22.92984


# Logistic model
# --------------
x <- glm(vs ~ wt, data = mtcars, family = "binomial")
predictions <- predict(x, type = "link")
ci_vals <- get_predicted_ci(x, predictions, ci_type = "prediction")
head(ci_vals)
#>                   CI_low CI_high
#> Mazda RX4           -Inf     Inf
#> Mazda RX4 Wag       -Inf     Inf
#> Datsun 710          -Inf     Inf
#> Hornet 4 Drive      -Inf     Inf
#> Hornet Sportabout   -Inf     Inf
#> Valiant             -Inf     Inf
ci_vals <- get_predicted_ci(x, predictions, ci_type = "confidence")
head(ci_vals)
#>          SE     CI_low   CI_high
#> 1 0.5623444 -0.3931282 1.8112213
#> 2 0.4690190 -0.6974034 1.1411172
#> 3 0.7195076 -0.1279982 2.6924199
#> 4 0.4459072 -1.3016913 0.4462326
#> 5 0.5021936 -1.8418839 0.1266787
#> 6 0.5094490 -1.8943152 0.1026881