Pbinary: Bayesian analysis for a parametric Bernoulli regression model
In DPpackage: Bayesian Nonparametric Modeling in R

Description Usage Arguments Details Value Author(s) References Examples

This function generates a posterior density sample for a parametric binary regression model.

1
2
3

Pbinary(formula,link="logit",prior,mcmc,state,status,misc=NULL,
         data=sys.frame(sys.parent()),na.action=na.fail)

`formula`	a two-sided linear formula object describing the model fit, with the response on the left of a `~` operator and the terms, separated by `+` operators, on the right.
`link`	a description of the link function to be used in the model. The links considered by `Pbinary` so far are logit (default), probit, cloglog, and cauchy.
`prior`	a list giving the prior information. The list includes the following parameters: beta0 and Sbeta0 giving the hyperparameters of the normal prior distribution for the regression coefficients.
`mcmc`	a list giving the MCMC parameters. The list must include the following integers: nburn giving the number of burn-in scans, nskip giving the thinning interval, nsave giving the total number of scans to be saved, ndisplay giving the number of saved scans to be displayed on the screen (the function reports on the screen when every ndisplay iterations have been carried out), and tune giving the Metropolis tuning parameter.
`state`	a list giving the current value of the parameters. This list is used if the current analysis is the continuation of a previous analysis.
`status`	a logical variable indicating whether this run is new (`TRUE`) or the continuation of a previous analysis (`FALSE`). In the latter case the current value of the parameters must be specified in the object state.
`misc`	misclassification information. When used, this list must include two objects, sens and spec, giving the sensitivity and specificity, respectively. Both can be a vector or a scalar. This information is used to correct for misclassification in the conditional bernoulli model.
`data`	data frame.
`na.action`	a function that indicates what should happen when the data contain `NA`s. The default action (`na.fail`) causes `Pbinary` to print an error message and terminate if there are any incomplete observations.

Pbinary simulates from the posterior density of a parametric Bernoulli regression model,

yi ~ Bernoulli(pii)

where pii = F(Xi beta) and F is a distribution function on the real line known as the inverse of the link function in the context of generalized linear models. The links considered by Pbinary so far are logit (default), probit, cloglog, and cauchy.

To complete the model specification, the following prior distribution is assumed,

beta | beta0, Sbeta0 ~ N(beta0,Sbeta0)

A Metropolis-Hastings step is used to sample the posterior distribution of the regression coefficients. The Metropolis proposal distribution is centered at its current value and the variance-covariance matrix correspond to the variance-covariance matrix of the MLEs times the tunning parameter, tune, specified in the mcmc list.

When the model considers correction for misclassification, a modified link function is used. The modified link is a function of the sensitivity and specificity of the classification (see, e.g., Jara, Garcia-Zattera and Lesaffre, 2006).

An object of class Pbinary representing the parametric regression model fit. Generic functions such as print, plot, summary, predict, and anova have methods to show the results of the fit. The results include only the regression coefficients, beta.

The MCMC samples of the parameters are stored in the object thetasave. This object is included in the list save.state and is a matrix which can be analyzed directly by functions provided by the coda package.

The list state in the output object contains the current value of the parameters necessary to restart the analysis. If you want to specify different starting values to run multiple chains set status=TRUE and create the list state based on this starting values. In this case the list state must include the following objects:

beta

giving the value of the regression coefficients.

Alejandro Jara <atjara@uc.cl>

Jara, A., Garcia-Zattera, M.J., Lesaffre, E. (2006) Semiparametric Bayesian Analysis of Misclassified Binary Data. XXIII International Biometric Conference, July 16-21, Montreal, Canada.

## Not run: 
    # Bioassay Data Example
    # Cox, D.R. and Snell, E.J. (1989). Analysis of Binary Data. 2nd ed. 
    # Chapman and Hall. p. 7  
    # In this example there are 150 subjects at 5 different stimulus levels, 
    # 30 at each level.

      y<-c(rep(0,30-2),rep(1,2),
           rep(0,30-8),rep(1,8),
           rep(0,30-15),rep(1,15),
           rep(0,30-23),rep(1,23),
           rep(0,30-27),rep(1,27))

      x<-c(rep(0,30),
           rep(1,30),
           rep(2,30),
           rep(3,30),
           rep(4,30))

    # Initial state
      state <- NULL

    # MCMC parameters
      nburn<-5000
      nsave<-5000
      nskip<-10
      ndisplay<-1000
      mcmc <- list(nburn=nburn,nsave=nsave,nskip=nskip,ndisplay=ndisplay,
                   tune=1.1)


    # Prior distribution
      prior <- list(beta0=rep(0,2), Sbeta0=diag(10000,2))

    # Fit a logistic regression model
      fit1 <- Pbinary(y~x,link="logit",prior=prior,
                      mcmc=mcmc,state=state,status=TRUE) 
      fit1

    # Fit a probit regression model
      fit2 <- Pbinary(y~x,link="probit",prior=prior, 
                      mcmc=mcmc,state=state,status=TRUE) 
      fit2

    # Fit a cloglog regression model
      fit3 <- Pbinary(y~x,link="cloglog",prior=prior,
                      mcmc=mcmc,state=state,status=TRUE) 
      fit3

    # Fit a cauchy regression model
      fit4 <- Pbinary(y~x,link="cauchy",prior=prior,
                      mcmc=mcmc,state=state,status=TRUE) 
      fit4

## End(Not run)