Prediction intervals for future observations based on linear random effects models
Prediction intervals for future observations based on linear random effects models
lmer_pi_futvec() calculates a bootstrap calibrated prediction interval for one or more future observation(s) based on linear random effects models. With this approach, the experimental design of the future data is taken into account (see below).
futvec: an integer vector that defines the structure of the future data based on the row numbers of the historical data. If length(futvec) is one, a PI for one future observation is computed
newdat: a data.frame with the same column names as the historical data on which model depends
alternative: either "both", "upper" or "lower". alternative specifies if a prediction interval or an upper or a lower prediction limit should be computed
alpha: defines the level of confidence (1-alpha)
nboot: number of bootstraps
delta_min: lower start value for bisection
delta_max: upper start value for bisection
tolerance: tolerance for the coverage probability in the bisection
traceplot: if TRUE: Plot for visualization of the bisection process
n_bisec: maximal number of bisection steps
algorithm: either "MS22" or "MS22mod" (see details)
Returns
lmer_pi_futvec() returns an object of class c("predint", "normalPI")
with prediction intervals or limits in the first entry ($prediction).
Details
This function returns bootstrap-calibrated prediction intervals as well as lower or upper prediction limits.
If algorithm is set to "MS22", both limits of the prediction interval are calibrated simultaneously using the algorithm described in Menssen and Schaarschmidt (2022), section 3.2.4. The calibrated prediction interval is given as
[l,u]=μ^±qcalibvar(μ^)+c=1∑C+1σ^c2
with μ^ as the expected future observation (historical mean) and σ^c2 as the c=1,2,...,C variance components and σ^C+12
as the residual variance obtained from the random effects model fitted with lme4::lmer() and qcalib as the as the bootstrap-calibrated coefficient used for interval calculation.
If algorithm is set to "MS22mod", both limits of the prediction interval are calibrated independently from each other. The resulting prediction interval is given by
Please note, that this modification does not affect the calibration procedure, if only prediction limits are of interest.
Be aware that the sampling structure of the historical data must contain the structure of the future data. This means that the observations per random factor must be less or equal in the future data compared to the historical data.
This function is an implementation of the PI given in Menssen and Schaarschmidt 2022 section 3.2.4 except that the bootstrap calibration values are drawn from bootstrap samples that mimic the future data.
Examples
# loading lme4library(lme4)# Fitting a random effects model based on c2_dat1fit <- lmer(y_ijk~(1|a)+(1|b)+(1|a:b), c2_dat1)summary(fit)#----------------------------------------------------------------------------### Prediction interval using c2_dat3 as future data# without printing c2_dat3 in the output# Row numbers of the historical data c2_dat1 that define the structure of# the future data c2_dat3futvec <- c(1,2,4,5,10,11,13,14)# Calculating the PIpred_int <- lmer_pi_futvec(model=fit, futvec=futvec, nboot=100)summary(pred_int)#----------------------------------------------------------------------------### Calculating the PI with c2_dat3 printed in the outputpred_int_new <- lmer_pi_futvec(model=fit, futvec=futvec, newdat=c2_dat3, nboot=100)summary(pred_int_new)#----------------------------------------------------------------------------### Upper prediction limit for m=1 future observationpred_u <- lmer_pi_futvec(model=fit, futvec=1, alternative="upper", nboot=100)summary(pred_u)#----------------------------------------------------------------------------# Please note that nboot was set to 100 in order to decrease computing time# of the example. For a valid analysis set nboot=10000.
References
Menssen and Schaarschmidt (2022): Prediction intervals for all of M future observations based on linear random effects models. Statistica Neerlandica, tools:::Rd_expr_doi("10.1111/stan.12260")