add.risk {CoMiRe}

Additional risk function

Description

Additional risk function estimated from the object fit

Usage

add.risk(y, x, fit, mcmc, a, alpha=0.05, 
x.grid=NULL, y.grid=NULL)

Arguments

y

optional numeric vector for the response used in comire.gibbs. If y is missing, y.grid must be provided.

x

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

fit

the output of comire.gibbs. an object of the class classCoMiRe.

mcmc

a list giving the MCMC parameters.

a

threshold of clinical interest for the response variable

alpha

level of the credible bands.

x.grid

optional numerical vector giving the actual values of the grid for x for plotting the additional risk function. If x.gird is not provided, standard grids are automatically used.

y.grid

optional numerical vector giving the actual values of the grid for y for plotting the additional risk function. If y.gird is not provided, standard grids are automatically used.

Value

A list of arguments for generating posterior output. It contains:

• mcmc.risk a matrix containing in the lines the MCMC chains, after thinning, of the additional risk function over x.grid, in the columns.

• summary.risk a data frame with four variables: the posterior means of the additional risk function over x.grid, the respective \alpha/2 and 1-\alpha/2 quantiles, and x.grid.

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)
                     
risk.data <- add.risk(y = gestage, x = dde, fit = fit.dummy, mcmc = mcmc, 
    a = 37, x.grid = seq(0, max(dde), length = 100))
riskplot(risk.data$summary.risk, xlab="DDE", x = dde, xlim = c(0,150))
                    

## 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)
 
risk.data <- add.risk(y = gestage, x = dde, fit = fit1, mcmc = mcmc, 
a = 37, x.grid = seq(0, max(dde), length = 100))
riskplot(risk.data$summary.risk, xlab="DDE", x = dde, xlim = c(0,150))


}