BGweights.RdCompute empirical weights based on out of sample forecast variances, following Bates and Granger (1969).
BGWeights(object, ..., data, force.update = FALSE)two or more fitted glm objects, or a
list of such, or an "averaging"=model.avg object.
a data frame containing the variables in the model.
if TRUE, the much less efficient method of
updating glm function will be used rather than directly via
glm.fit. This only applies to glms, in
case of other model types update is always used.
A numeric vector of model weights.
Bates-Granger model weights are calculated using prediction covariance. To
get the estimate of prediction covariance, the models are fitted to
randomly selected half of data and prediction is done on the
remaining half.
These predictions are then used to compute the variance-covariance between
models, \(\Sigma\). Model weights are then calculated as
w_BG = (1' ^-1 1)^-1 1 ^-1
w_BG = (1' ^-1 1)^-1 1 ^-1
w_BG = (1' Sigma^-1 1)^-1 1 \ Sigma^-1,
where \(1\) a vector of 1-s.
Bates-Granger model weights may be outside of the \([0,1]\) range, which may cause the averaged variances to be negative. Apparently this method works best when data is large.
For matrix inversion, MASS::ginv()MASS:ginv is more stable near singularities
than solve. It will be used as a fallback if solve fails and
MASS is available.
Bates, J. M. and Granger, C. W. J. 1969 The combination of forecasts. Journal of the Operational Research Society 20, 451-468.
Dormann, C. et al. (2018) Model averaging in ecology: a review of Bayesian, information-theoretic, and tactical approaches for predictive inference. Ecological Monographs 88, 485–504.
Weights, model.avg
Other model weights: bootWeights, cos2Weights, jackknifeWeights, stackingWeights
fm <- glm(Prop ~ mortality + dose, family = binomial, Beetle, na.action = na.fail)
models <- lapply(dredge(fm, evaluate = FALSE), eval)
#> Fixed term is "(Intercept)"
ma <- model.avg(models)
# this produces warnings because of negative variances:
set.seed(78)
Weights(ma) <- BGWeights(ma, data = Beetle)
coefTable(ma, full = TRUE)
#> Estimate Std. Error
#> (Intercept) -6.956333 6.1472
#> dose 0.095898 0.1183
#> mortality 2.813398 2.8082
# SE for prediction is not reliable if some or none of coefficient's SE
# are available
predict(ma, data = test.data, se.fit = TRUE)
#> $fit
#> [1] -1.9139710 -1.2557521 -0.6828942 0.2847985 1.5891400 2.1629570 2.7639489
#> [8] 3.1970959
#>
#> $se.fit
#> [1] 0.9254979 0.6544394 0.4801025 0.1670266 0.5464026 0.6472613 0.8764116
#> [8] 1.1546316
#>
coefTable(ma, full = TRUE)
#> Estimate Std. Error
#> (Intercept) -6.956333 6.1472
#> dose 0.095898 0.1183
#> mortality 2.813398 2.8082