Convex Mixture Regression

Description

Posterior inference under the convex mixture regression (CoMiRe) models introduced by Canale, Durante, and Dunson (2018) <doi:10.1111/biom.12917>.

Details

The CoMiRe package implements the convex mixture regresion approach of Canale, Durante, and Dunson (2018) and some extensions to deal with binary response variables or to account for the presence of continuous and categorical confunders. Estimation is conducted via Gibbs sampler. Posterior plots for inference and goodness-of-fit tests are also available.

Author(s)

Antonio Canale [aut, cre], Daniele Durante [ctb], Arianna Falcioni [aut], Luisa Galtarossa [aut], Tommaso Rigon [ctb] Maintainer: Antonio Canale <antonio.canale@unipd.it>

References

Canale, A., Durante, D., and Dunson, D. (2018), Convex Mixture Regression for Quantitative Risk Assessment, Biometrics, 74, 1331-1340

Benchmark dose

Description

Benchmark dose associated to a particular risk

Usage

BMD(level, risk, x, alpha=0.05)

Arguments

Value

A dataframe containing as variables:

• q the dose level of interest.

• BMD the benchmark dose.

• low lower credible limit.

• upp upper credible limit.

• BMDL a more conservative benchmark dose.

Author(s)

Antonio Canale

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))
bmd.data <- BMD(seq(0,.20, length=50), risk.data$mcmc.risk, 
x=seq(0,max(dde), length=100), alpha=0.05)
bmd.plot(bmd.data)       
                     

## 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)
              
fit <- comire.gibbs(gestage, dde, family="continuous", 
                     mcmc=mcmc, prior=prior, seed=5, max.x=180)


risk.data <- add.risk(y = gestage, x = dde, fit = fit, mcmc = mcmc,
a = 37, x.grid = seq(0, max(dde), length = 100))
bmd.data <- BMD(seq(0,.20, length=50), risk.data$mcmc.risk, 
x=seq(0,max(dde), length=100), alpha=0.05)
bmd.plot(bmd.data)       


}

Collaborative Perinatal Project data

Description

A subset of the Collaborative Perinatal Project data set (Klebanoff, 2009) focusing on studying the effect of DDE exposure on pregnancies (Longnecker et al., 2001). The dataset contains the following variables for each pregnant women enrolled in the study:

• hosp, factor denoting the hospital where the woman was hospitalized;

• dde, Dichlorodiphenyldichloroethylene (DDE) concentration in maternal serum;

• gestage, gestational age (in weeks);

• mage, age of the moter (in years);

• mweight, pre pregnancy weigth of the mother (in lbs);

• mbmi, pre pregnancy BMI of the mother;

• sei, socio economic index of the mother;

• smoke, factor. It takes value 2 if the woman is a smoker, 1 otherwise;

Usage

data(CPP)

Format

A data.frame

References

Klebanoff M. A. (2009) The collaborative perinatal project: a 50-year retrospective. Paediatric and perinatal epidemiology, 23, 2.

Longnecker, M. P., Klebanof, M. A., Zhou, H., Brock, J. (2001) Association between maternal serum concentration of the DDT metabolite DDE and preterm and small-for-gestational-age babies at birth. The Lancet, 358, 110-114.

Examples

data(CPP)
str(CPP)

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

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))


}

classCoMiRe class constructor

Description

A constructor for the classCoMiRe class. The class classCoMiRe is a named list containing the output of the posterior estimation of CoMiRe model implemented in comire.gibbs

Usage

as.classCoMiRe(call = NULL, out = NULL, z = NULL, z.val = NULL, f0 = NULL, f1 = NULL, 
nrep, nb, bin = FALSE, univariate = TRUE)

Arguments

Author(s)

Antonio Canale, Arianna Falcioni

\beta(x) plot

Description

Posterior mean (continuous lines) and pointwise credible bands (shaded areas) for \beta(x).

Usage

betaplot(x, fit, x.grid = NULL, xlim = c(0, max(x)), xlab = "x")

Arguments

Author(s)

Antonio Canale

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)
                     
betaplot(x=dde, fit=fit.dummy, x.grid=seq(0,180, length=100), 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)

                         
betaplot(x=dde, fit=fit1, x.grid=seq(0,180, length=100), xlim=c(0,150))


}

Benchmark dose plot

Description

Posterior mean (continuous lines) and pointwise credible bands (shaded areas) for the benchmark dose in function of the increase in risk.

Usage

bmd.plot(bmd.data)

Arguments

Author(s)

Antonio Canale

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))
bmd.data <- BMD(seq(0,.20, length=50), risk.data$mcmc.risk, 
x=seq(0,max(dde), length=100), alpha=0.05)
bmd.plot(bmd.data)       


## 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)
              
fit <- comire.gibbs(gestage, dde, family="continuous", 
                     mcmc=mcmc, prior=prior, seed=5, max.x=180)
                     
risk.data <- add.risk(y = gestage, x = dde, fit = fit, mcmc = mcmc, 
a = 37, x.grid = seq(0, max(dde), length = 100))
bmd.data <- BMD(seq(0,.20, length=50), risk.data$mcmc.risk, 
x=seq(0,max(dde), length=100), alpha=0.05)
bmd.plot(bmd.data)       


}

Gibbs sampler for CoMiRe model

Description

Posterior inference via Gibbs sampler for CoMiRe model

Usage

comire.gibbs(y, x, z = NULL, family = 'continuous', 
       grid = NULL, mcmc, prior, 
       state = NULL, seed, max.x = max(x), z.val = NULL, verbose = TRUE)

Arguments

Details

The function fit a convex mixture regression (CoMiRe) model (Canale, Durante, Dunson, 2018) via Gibbs sampler. For continuous outcome y \in \mathcal{Y}, adverse esposure level x \in \mathcal{X} and no confunding variables, one can set family = 'continuous' and z = NULL and fit model
f_x(y) = \{1-\beta(x)\} \sum_{h=1}^{H}\nu_{0h} \phi(y; \theta_{0h}, \tau_{0h}^{-1}) + \beta(x) \phi(y; \theta_{\infty}, \tau_{\infty}^{-1}) ;

where \beta(x) = \sum_{j=1}^{J} \omega_j \psi_j(x), x\ge0, is a a monotone nondecreasing interpolation function, constrained between 0 and 1 and \psi_1,...,\psi_J are monotone nondecreasing I-splines basis.
If p \ge 1 confounding covariates z \in \mathcal{Z} are available, passing the argument z the function fits model

f(y; x,z) = \{1-\beta(x)\} f_0(y;z) + \beta(x) f_\infty(y;z) ;

where:
f_0(y;z)= \sum_{h=1}^{H} \nu_{0h} \phi(y;\theta_{0h}+z^\mathsf{T}\gamma,\tau_{0h}^{-1}), and f_\infty(y;z)= \phi(y;\theta_\infty+ z^\mathsf{T}\gamma,\tau_{\infty}^{-1}).

Finally, if y is a binary response, one can set family = 'binary' and fit model

p_x(y) = (\pi_x)^y (1 - \pi_x)^{1-y} ;

where \pi_x = P(Y=1 | x) is \pi_x = \{1-\beta(x)\} \pi_0 + \beta(x) \pi_\infty.

Value

An object of the class classCoMiRe, i.e. a list of arguments for generating posterior output. It contains:

• callthe model formula

• post.means a list containing the posterior mean density beta over the grid, of all the mixture parameters and, if family = "continuous" and z = NULL, of f_0 and f_{inf} over the y.grid.

• ci a list containing the 95% credible intervals for all the quantities stored in post.means.

• mcmc a list containing all the MCMC chains.

• z the same of the input

• z.val the same of the input

• f0,f1 MCMC replicates of the density in the two extremes (only if family = 'continuous')

• nrep,nb the same values of the list mcmc in the input arguments

• bin logical, equal to TRUE if family = 'binary'

• univariate logical, equal to TRUE if z is null or a vector

Author(s)

Antonio Canale [aut, cre], Daniele Durante [ctb], Arianna Falcioni [aut], Luisa Galtarossa [aut], Tommaso Rigon [ctb]

References

Canale, A., Durante, D., and Dunson, D. (2018), Convex Mixture Regression for Quantitative Risk Assessment, Biometrics, 74, 1331-1340

Galtarossa, L., Canale, A., (2019), A Convex Mixture Model for Binomial Regression, Book of Short Papers SIS 2019

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)
                     
summary(fit.dummy)
 


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

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

## 1. binary case ##

prior <- list(pi0=mean(gestage), eta=rep(1, J)/J, 
             a.pi0=27, b.pi0=360, J=J)
             
fit_binary<- comire.gibbs(premature, dde, family="binary", 
                          mcmc=mcmc, prior=prior, seed=5, max.x=180)
                          
                          
## 2. continuous case ##

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)


## 2.2 One confunder ##

mage_std <- scale(mage, center = TRUE, scale = TRUE) 

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


## 2.3 More confunders ##

Z <- cbind(mage, mbmi, sei)
Z <- scale(Z, center = TRUE, scale = TRUE)
Z <- as.matrix(cbind(Z, CPP$smoke))
colnames(Z) <- c("age", "BMI", "sei", "smoke")

mod <- lm(gestage ~ dde + Z)
prior <- list(mu.theta = mod$coefficients[1], k.theta = 10,
              mu.gamma = mod$coefficients[-c(1, 2)], sigma.gamma = diag(rep(10, 4)),
              eta = rep(1, J)/J, alpha = 1/H, a = 2, b = 2, H = H, J = J)
              
fit3 <- comire.gibbs(y = gestage, x = dde, z = Z, family = "continuous", mcmc = mcmc, 
                     prior = prior, seed = 5)

 
}

Internal Functions of CoMiRe package

Description

Internal Functions of CoMiRe package

Usage

.labelling_b_C(w, phi, f0i, f1i)

.labelling_bb_C(w, phi, P0, P1)

.labelling_c_C(y, logpi, mu, tau)

.mixdensity_C(y, pi, mu, tau)

.pssq_gaussian(index, data, cluster, locations)

.labelling_b_uni(i, w, phi, f0i, f1i)

.labelling_c_uni(i, y, z, nu, theta, tau, ga)

.mixdensity_uni(i, y, z, nu, theta, tau, ga)

.pssq_uni(index, y, z, cluster, theta, gamma)

.psdp_uni(index, y, z, cluster)

.labelling_b_multi(i, w, phi, f0i, f1i)

.labelling_c_multi(i, y, z, nu, theta, tau, ga)

.mixdensity_multi(i, y, z, nu, theta, tau, ga)

.pssq_multi(y, times, z, gamma, theta)

.psdp_multi(index, y, z, cluster)

Posterior mean density plot for dose intervals

Description

Pointwise posterior mean (continuous blue lines), and credible bands (shaded blue areas) for f (y | x, z) calculated in x.val under the the model fitted in fit.

Usage

fit.pdf.mcmc(y, x, fit, mcmc, J=10, H = 10, alpha = 0.05, 
max.x = max(x), x.val, y.grid = NULL, xlim = c(0, max(x)), 
ylim = c(0, 1), xlab = NULL)

Arguments

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)
                     
fit.pdf.mcmc(y = gestage, x = dde, fit = fit.dummy, mcmc = mcmc, J = 10, H = 10, 
                         alpha = 0.05, max.x = max(dde), x.val = 125, 
                         xlim = c(25,48), ylim = c(0,0.25),
                         xlab = "Gest. age. for DDE = 125")
                         

## 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)
          
fit.pdf.mcmc(y = gestage, x = dde, fit = fit1, mcmc = mcmc, J = 10, H = 10, 
                         alpha = 0.05, max.x = max(dde), x.val = 125, 
                         xlim = c(25,48), ylim = c(0,0.25),
                         xlab = "Gest. age. for DDE = 125")
                         

}

CoMiRe plot

Description

An S3 plot method for an object of classCoMiRe class.

Usage

## S3 method for class 'classCoMiRe'
plot(
  x,
  y,
  xobs,
  mcmc,
  J = 10,
  H = 10,
  a = NULL,
  max.x = max(xobs),
  bandwidth = 20,
  x.grid = NULL,
  xlim = c(0, max(xobs)),
  ylim = c(0, 1),
  xlab = "x",
  alpha = 0.05,
  risk = TRUE,
  bmd = TRUE,
  level,
  oneevery = 20,
  ...
)

Arguments

Details

The output is a list of ggplot2 plots containing the result of the betaplot function and, if the threshold a is provided, of post.pred.check, riskplot, bmd.plot.

Value

If a=NULL returns only betaplot otherwise, if risk=FALSE and bmd=FALSE returns a list containing betaplot (which is automatically plotted) and post.pred.check plot. Finally, if a is provided, risk=TRUE and bmd=TRUE returns a list with betaplot, post.pred.check, riskplot and bmd.plot.

Author(s)

Antonio Canale, Arianna Falcioni

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

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))


}

comire.gibbs for different fixed values of z

Description

This function computes the predicted values of the density al low dose f_0 and of the density at high dose f_{inf}, for fixed values of the confounders z.

Usage

predict_new_z(fit, y, z.val)

Arguments

Value

An object of class classCoMiRe.

Author(s)

Antonio Canale, Arianna Falcioni

CoMiRe print

Description

The print.classCoMiRe method prints the type of a classCoMiRe object.

Usage

## S3 method for class 'classCoMiRe'
print(x, ...)

Arguments

Author(s)

Antonio Canale, Arianna Falcioni

Additional risk function plot

Description

Posterior mean (continuous lines) and pointwise credible bands (shaded areas) for Ra(x, a).

Usage

riskplot(risk.data, xlab = NULL, x = NULL, ylim=c(0,1), xlim=c(0, max(x)))

Arguments

Author(s)

Antonio Canale

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)
              
fit <- comire.gibbs(gestage, dde, family="continuous", 
                     mcmc=mcmc, prior=prior, seed=5, max.x=180)
 
risk.data <- add.risk(y = gestage, x = dde, fit = fit, 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))


}

CoMiRe summary

Description

The summary.classCoMiRe method provides summary information on classCoMiRe objects.

Usage

## S3 method for class 'classCoMiRe'
summary(object, ...)

Arguments

Author(s)

Antonio Canale Arianna Falcioni