Bayesian Analysis of Complete Contingency Tables
These functions are the workhorses behind bcct, bcctu, bcctsubset, and bcctsubsetu.
bcct.fit(priornum, maximal.mod, IP, eta.hat, ini.index, ini.beta, ini.sig, iters, save, name, null.move.prob, a, b, progress) bcctsubset.fit(priornum, subset.index, maximal.mod, IP, eta.hat, ini.index, ini.beta, ini.sig, iters, save, name, null.move.prob, a, b, progress)
priornum |
A numeric scalar indicating which prior is to be used: 1 = |
subset.index |
A matrix where each row gives the index for each model in the subset of models under consideration. |
maximal.mod |
An object of class |
IP |
A p by p matrix giving the inverse of the prior scale matrix for the maximal model. |
eta.hat |
A vector of length n (number of cells) giving the posterior mode of the linear predictor under the maximal model. |
ini.index |
A binary vector, of the same length as the number of log-linear parameters in the maximal model, indicating which parameters are present in the initial model. |
ini.beta |
A numeric vector giving the starting values of the log-linear parameters for the MCMC algorithm. |
ini.sig |
A numeric scalar giving the starting value of sigma^2 for the MCMC algorithm. |
iters |
The number of iterations of the MCMC algorithm to peform. |
save |
If positive, the function will save the MCMC output to external text
files every |
name |
A prefix to the external files saved if the argument |
null.move.prob |
A scalar argument giving the probability of performing a null move, i.e. proposing a move to the current model. |
a |
The shape hyperparameter of the Sabanes-Bove & Held prior, see Overstall & King (2014). |
b |
The scale hyperparameter of the Sabanes-Bove & Held prior, see Overstall & King (2014). |
progress |
Logical argument. If |
The function will return a list with the following components:
BETA |
An |
MODEL |
A vector of length |
SIG |
A vector of length |
rj_acc |
A binary vector of the same length as the number of reversible jump moves attempted. A 0 indicates that the proposal was rejected, and a 1 that the proposal was accepted. |
mh_acc |
A binary vector of the same length as the number of Metropolis-Hastings moves attempted. A 0 indicates that the proposal was rejected, and a 1 that the proposal was accepted. |
This function will not typically be called by the user.
Antony M. Overstall A.M.Overstall@soton.ac.uk.
Overstall, A.M. & King, R. (2014) conting: An R package for Bayesian analysis of complete and incomplete contingency tables. Journal of Statistical Software, 58 (7), 1–27. http://www.jstatsoft.org/v58/i07/
data(AOH) ## Load the AOH data. maximal.mod<-glm(formula=y~(alc+hyp+obe)^3,data=AOH,x=TRUE,y=TRUE, contrasts=list(alc="contr.sum",hyp="contr.sum",obe="contr.sum")) ## Set up the maximal model which in this case is the saturated ## model. curr.index<-formula2index(big.X=maximal.mod$x,formula=y~alc+hyp+obe+hyp:obe,data=AOH) ## Set up the binary vector for the model containing all main effects and the ## hyp:obe interaction. IP<-t(maximal.mod$x)%*%maximal.mod$x/length(maximal.mod$y) IP[,1]<-0 IP[1,]<-0 ## Set up the inverse scale matrix for the prior distribution under ## the maximal model. bmod<-beta_mode(X=maximal.mod$x,prior="UIP",y=maximal.mod$y,IP=IP) ## Find the posterior mode under the maximal model eta.hat<-as.vector(maximal.mod$x%*%bmod) ## Find the posterior mode of the linear predictor ## under the maximal model. set.seed(1) ## Set seed for reproducibility test1<-bcct.fit(priornum=1, maximal.mod=maximal.mod, IP=IP, eta.hat=eta.hat, ini.index=curr.index, ini.beta=bmod[curr.index==1], ini.sig=1, iters=5, save=0, name=NULL,null.move.prob=0.5, a=0.001, b=0.001, progress=TRUE) ## Run for 5 iterations starting at model defined by curr.index. test1$MODEL ## Look at sampled model indicators. Should be: ## [1] "fe00c0" "fe0000" "fe0000" "fe0000" "fe0000" model2index(test1$MODEL,dig=24) ## Convert these to binary indicators of the log-linear parameters. ## Will get: # [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12] [,13] #fe00c0 1 1 1 1 1 1 1 0 0 0 0 0 0 #fe0000 1 1 1 1 1 1 1 0 0 0 0 0 0 #fe0000 1 1 1 1 1 1 1 0 0 0 0 0 0 #fe0000 1 1 1 1 1 1 1 0 0 0 0 0 0 #fe0000 1 1 1 1 1 1 1 0 0 0 0 0 0 # [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21] [,22] [,23] [,24] #fe00c0 0 0 0 1 1 0 0 0 0 0 0 #fe0000 0 0 0 0 0 0 0 0 0 0 0 #fe0000 0 0 0 0 0 0 0 0 0 0 0 #fe0000 0 0 0 0 0 0 0 0 0 0 0 #fe0000 0 0 0 0 0 0 0 0 0 0 0 ## See how the hyp:obe interactions in columns 17 and 18 gets dropped after ## the 1st iteration.
Please choose more modern alternatives, such as Google Chrome or Mozilla Firefox.