#################################################################################
# Bayesian model selection and statistical modeling 
# Chapman & Hall/CRC Taylor and Francis Group 
#
# Chapter6: Section6.7.1
#
# Bayesian analysis of the probit model
#
# Author : Tomohiro Ando
#
################################################################################
library(MCMCpack)

# Dutch voting behavior data

D <- matrix(scan("Defaultdata.txt"),ncol=4,byrow=T)
D <- as.data.frame(D)
names(D) <- c("ROA","logSales","CF-Sales ratio","Ratings")
n <- nrow(D)


posterior <- MCMCprobit(Default~., data = D, burnin = 1000, mcmc = 4000, b0=0, B0 =1)

Beta <- as.vector(rep(1,len=nrow(posterior))%*%posterior)/nrow(posterior)
X <- as.matrix(cbind(1,D[,1:3]))
O <- X%*%Beta

p1 <- pnorm(O)
p2 <- 1-pnorm(O)

Loglike <- sum( log(p1)*D[,4]+log(p2)*(1-D[,4]) )

Prior <- 4/2*log(2*pi)-sum( (Beta)^2 )/2

V <- as.matrix(var(posterior))

Post <- 0

for(j in 1:nrow(posterior)){

B <- as.vector(posterior[j,])

A <- 4/2*log(2*pi)-log(det(V))/2-sum( t(B-Beta)%*%solve(V)%*%(B-Beta) )/2

Post <- Post + A

}

Post <- Post/nrow(posterior)

LogMarginallikelihood <- Loglike+Prior-Post


print(LogMarginallikelihood)