post.pred.check {CoMiRe}

Posterior predictive check plot

Description

A plot for an object of classCoMiRe class. The plot is a goodness-of-fit assessment of CoMiRe model. Since Version 0.8 if z is provided into the fit object, an error message is returned. If family = 'continuous', a smoothed empirical estimate of F(a|x) = pr(y < a | x) is computed from the observed data (black line) and from some of the data sets simulated from the posterior predictive distribution in the fit object (grey lines). If family = 'binary', a smoothed empirical estimate of the proportion of events (black line) and of the smoothed empirical proportion of data simulated from the posterior predictive distribution in the fit object (grey lines). In the x axis are reported the observed exposures.

Usage

post.pred.check(y, x, fit, mcmc, J=10, H=10, a, max.x=max(x), 
xlim=c(0, max(x)), bandwidth = 20, oneevery = 20)

Arguments

y

numeric vector for the response used in comire.gibbs

x

numeric vector for the covariate relative to the dose of exposure used in comire.gibbs

fit

the output of comire.gibbs opportunely trasformed in classCoMiRe class

mcmc

a list giving the MCMC parameters

J

parameter controlling the number of elements of the I-spline basis

H

total number of components in the mixture at x_0

a

threshold of clinical interest to compute the F(a|x,z)

max.x

maximum value allowed for x

xlim

numeric vectors of length 2, giving the x coordinates ranges for the plot

bandwidth

the kernel bandwidth smoothing parameter

oneevery

integer number representing how many MCMC draws to plot in the posterior predictive check. It draws one sample every oneevery.

Author(s)

Antonio Canale, Arianna Falcioni

Examples

{
data(CPP)
attach(CPP)

n <- NROW(CPP)
J <- H <- 10

premature <- as.numeric(gestage<=37)

mcmc <- list(nrep=5000, nb=2000, thin=5, ndisplay=4)

## too few iterations to be meaningful. see below for safer and more comprehensive results

mcmc <- list(nrep=10, nb=2, thin=1, ndisplay=4) 

prior <- list(mu.theta=mean(gestage), k.theta=10, eta=rep(1, J)/J, 
              alpha=rep(1,H)/H, a=2, b=2, J=J, H=H)
              
fit.dummy <- comire.gibbs(gestage, dde, family="continuous", 
                     mcmc=mcmc, prior=prior, seed=1, max.x=180)
                     
post.pred.check(y = gestage, x = dde, fit = fit.dummy, mcmc = mcmc, J = 10, H = 10, a = 37, 
                max.x = max(dde), xlim = c(0,150), oneevery = 4)
                

## safer procedure with more iterations (it may take some time)

mcmc <- list(nrep=5000, nb=2000, thin=5, ndisplay=4)

## Fit the model for continuous y 

prior <- list(mu.theta=mean(gestage), k.theta=10, eta=rep(1, J)/J, 
              alpha=rep(1,H)/H, a=2, b=2, J=J, H=H)
              
fit1 <- comire.gibbs(gestage, dde, family="continuous", 
                     mcmc=mcmc, prior=prior, seed=5, max.x=180)

post.pred.check(y = gestage, x = dde, fit = fit1, mcmc = mcmc, J = 10, H = 10, a = 37, 
                max.x = max(dde), xlim = c(0,150))


}