### Name: blasso
### Title: Bayesian Lasso/NG and Ridge Regression
### Aliases: blasso bridge
### Keywords: regression

### ** Examples

## following the lars diabetes example
data(diabetes)
attach(diabetes)

## Ordinary Least Squares regression
reg.ols <- regress(x, y)

## Lasso regression
reg.las <- regress(x, y, method="lasso")

## Bayesian Lasso regression
reg.blas <- blasso(x, y)

## summarize the beta (regression coefficients) estimates
plot(reg.blas, burnin=200)
points(drop(reg.las$b), col=2, pch=20)
points(drop(reg.ols$b), col=3, pch=18)
legend("topleft", c("blasso-map", "lasso", "lsr"),
       col=c(2,2,3), pch=c(21,20,18))

## plot the size of different models visited
plot(reg.blas, burnin=200, which="m")

## get the summary
s <- summary(reg.blas, burnin=200)

## calculate the probability that each beta coef != zero
s$bn0

## summarize s2
plot(reg.blas, burnin=200, which="s2")
s$s2

## summarize lambda2
plot(reg.blas, burnin=200, which="lambda2")
s$lambda2

## fit with Student-t errors
## (~400-times slower due to automatic thinning level)
regt.blas <- blasso(x, y, theta=0.1)

## plotting some information about nu, and quantiles
plot(regt.blas, "nu", burnin=200)
quantile(regt.blas$nu[-(1:200)], c(0.05, 0.95))

## Bayes Factor shows strong evidence for Student-t model
mean(exp(regt.blas$llik[-(1:200)] - regt.blas$llik.norm[-(1:200)]))

## clean up
detach(diabetes)

##
## a big-p small-n example
##

n <- 25; m <- 51
xmuS <- randmvn(n, m)
X <- xmuS$x[,1:(m-1)]
Y <- drop(xmuS$x[,m])
obl <- blasso(X, Y, verb=0)

## plot summary of the model order
plot(obl, burnin=10, which="m")

## fit a standard lasso model
oml <- regress(X, Y, method="lasso")

## compare via RMSE, most often blasso will win
beta <- xmuS$S[m,-m] %*% solve(xmuS$S[-m,-m])
sqrt(mean((apply(obl$beta, 2, mean) - beta)^2))
sqrt(mean((oml$b[-1] - beta)^2))

## now try both Bayesian & ML ridge regression
obr <- bridge(X, Y, verb=0)
omr <- regress(X, Y, method="ridge")
sqrt(mean((apply(obr$beta[-c(1:200),], 2, mean) - beta)^2))
sqrt(mean((omr$b[-1] - beta)^2))



