Skip to contents

Predictions from gMAP analyses

Source: R/predict_gMAP.R
predict.gMAP.Rd

Produces a sample of the predictive distribution.

Usage

# S3 method for gMAP
predict(
  object,
  newdata,
  type = c("response", "link"),
  probs = c(0.025, 0.5, 0.975),
  na.action = na.pass,
  thin,
  ...
)

# S3 method for gMAPpred
print(x, digits = 3, ...)

# S3 method for gMAPpred
summary(object, ...)

# S3 method for gMAPpred
as.matrix(x, ...)

Arguments

newdata

data.frame which must contain the same columns as input into the gMAP analysis. If left out, then a posterior prediction for the fitted data entries from the gMAP object is performed (shrinkage estimates).

type

sets reported scale (response (default) or link).

probs

defines quantiles to be reported.

na.action

how to handle missings.

thin

thinning applied is derived from the gMAP object.

...

ignored.

x, object

gMAP analysis object for which predictions are performed

digits

number of displayed significant digits.

Details

Predictions are made using the \(\tau\) prediction stratum of the gMAP object. For details on the syntax, please refer to predict.glm and the example below.

See also

gMAP, predict.glm

Examples

## Setting up dummy sampling for fast execution of example
## Please use 4 chains and 20x more warmup & iter in practice
.user_mc_options <- options(RBesT.MC.warmup=50, RBesT.MC.iter=100,
                            RBesT.MC.chains=2, RBesT.MC.thin=1)

# create a fake data set with a covariate
trans_cov <- transform(transplant, country=cut(1:11, c(0,5,8,Inf), c("CH", "US", "DE")))
set.seed(34246)
map <- gMAP(cbind(r, n-r) ~ 1 + country | study,
            data=trans_cov,
            tau.dist="HalfNormal",
            tau.prior=1,
            # Note on priors: we make the overall intercept weakly-informative
            # and the regression coefficients must have tighter sd as these are
            # deviations in the default contrast parametrization
            beta.prior=rbind(c(0,2), c(0,1), c(0,1)),
            family=binomial,
            ## ensure fast example runtime
            thin=1, chains=1)
#> Warning: The largest R-hat is 1.1, indicating chains have not mixed.
#> Running the chains for more iterations may help. See
#> https://mc-stan.org/misc/warnings.html#r-hat
#> Warning: Bulk Effective Samples Size (ESS) is too low, indicating posterior means and medians may be unreliable.
#> Running the chains for more iterations may help. See
#> https://mc-stan.org/misc/warnings.html#bulk-ess
#> Warning: Tail Effective Samples Size (ESS) is too low, indicating posterior variances and tail quantiles may be unreliable.
#> Running the chains for more iterations may help. See
#> https://mc-stan.org/misc/warnings.html#tail-ess
#> Final MCMC sample equivalent to less than 1000 independent draws.
#> Please consider increasing the MCMC simulation size.

# posterior predictive distribution for each input data item (shrinkage estimates)
pred_cov <- predict(map)
pred_cov
#> Meta-Analytic-Predictive Prior Predictions
#> Scale: response 
#> 
#> Summary:
#>     mean     sd  2.5%   50% 97.5%
#> 1  0.205 0.0510 0.123 0.199 0.291
#> 2  0.195 0.0406 0.120 0.195 0.263
#> 3  0.216 0.0412 0.154 0.209 0.302
#> 4  0.253 0.0427 0.193 0.250 0.333
#> 5  0.197 0.0276 0.154 0.198 0.252
#> 6  0.184 0.0458 0.103 0.181 0.273
#> 7  0.238 0.0516 0.149 0.233 0.344
#> 8  0.168 0.0317 0.120 0.166 0.239
#> 9  0.241 0.0521 0.154 0.239 0.345
#> 10 0.174 0.0313 0.121 0.171 0.246
#> 11 0.243 0.0231 0.206 0.246 0.279

# extract sample as matrix
samp <- as.matrix(pred_cov)

# predictive distribution for each input data item (if the input studies were new ones)
pred_cov_pred <- predict(map, trans_cov)
pred_cov_pred
#> Meta-Analytic-Predictive Prior Predictions
#> Scale: response 
#> 
#> Summary:
#>     mean     sd   2.5%   50% 97.5%
#> 1  0.210 0.0650 0.0979 0.198 0.323
#> 2  0.221 0.0749 0.1210 0.220 0.379
#> 3  0.220 0.0748 0.1130 0.212 0.418
#> 4  0.214 0.0720 0.1110 0.206 0.377
#> 5  0.203 0.0697 0.0791 0.196 0.348
#> 6  0.198 0.0747 0.1030 0.180 0.360
#> 7  0.203 0.0638 0.1070 0.200 0.333
#> 8  0.191 0.0617 0.0827 0.188 0.299
#> 9  0.237 0.0985 0.0998 0.217 0.446
#> 10 0.238 0.0885 0.1200 0.234 0.406
#> 11 0.234 0.0792 0.0838 0.232 0.418


# a summary function returns the results as matrix
summary(pred_cov)
#>         mean         sd      2.5%       50%     97.5%
#> 1  0.2046888 0.05103156 0.1233163 0.1989446 0.2908793
#> 2  0.1950595 0.04060381 0.1204606 0.1949910 0.2630753
#> 3  0.2160028 0.04120147 0.1538030 0.2089321 0.3022593
#> 4  0.2529856 0.04270467 0.1929258 0.2503303 0.3333643
#> 5  0.1969653 0.02764795 0.1544012 0.1980138 0.2518814
#> 6  0.1843000 0.04581607 0.1028620 0.1814753 0.2733072
#> 7  0.2383059 0.05162711 0.1487544 0.2326939 0.3438535
#> 8  0.1680299 0.03170941 0.1198606 0.1659361 0.2389757
#> 9  0.2410675 0.05211712 0.1538566 0.2392947 0.3450478
#> 10 0.1742228 0.03131481 0.1213062 0.1707167 0.2455262
#> 11 0.2426294 0.02309088 0.2061749 0.2463644 0.2792656

# obtain a prediction for new data with specific covariates
pred_new <- predict(map, data.frame(country="CH", study=12))
pred_new
#> Meta-Analytic-Predictive Prior Predictions
#> Scale: response 
#> 
#> Summary:
#>    mean     sd   2.5%   50% 97.5%
#> 1 0.217 0.0793 0.0926 0.209 0.404

## Recover user set sampling defaults
options(.user_mc_options)