ridge_regression: Gibbs Sampler for Ridge Regression

In dcbdan/s525: Implementation of common Bayesian models

Description

The ridge regression model coeffecients are given by

\arg\min_{β} \parallel Y - Xβ \parallel^2 + λ \parallel β \parallel^2.

This is the implementation of the bayesian interpretation. Namely,

[y_i|x_i, β, τ] \sim N(x_i^tβ, τ^{-1}), i = 1,...,n ,

[β|η] \sim N(0, η^{-1}I),

τ \sim Gamma(\frac{1}{100}, \frac{1}{100}),

σ_β \sim unif(0, A)

where A = 1000 and τ = \frac{1}{σ_β^2}.

Usage

ridge_regression <- function(
  Y,
  X,
  niter = 10000)

Arguments

Y	n by 1
X	n by p predictor matrix
niter	number of gibbs sampling iterations

Details

This function returns the generated parameters from the gibbs sampling markov chain.

Value

beta	An niter by p matrix
sigma	An niter by 1 matrix
sigma_beta	An niter by 1 matrix

Examples

# Load the data
 prostate.data = "https://web.stanford.edu/~hastie/ElemStatLearn/datasets/prostate.data"
prostate = read.table(file = "prostate.data", sep="", header = TRUE)
# Training data:
prostate_train = prostate[which(prostate$train),-10]
# Testing data:
prostate_test = prostate[which(!prostate$train),-10]
# Response:
y = prostate_train$lpsa
# Center and scale the data:
y = scale(y)
# And the predictors
X = scale(prostate_train[,names(prostate_train) != "lpsa"])

gibbs_ridge <- ridge_regression(y, X, niter=10000)
shrinkage_regression_plot(gibbs_ridge$beta, y, X, main = "Ridge Prior")

dcbdan/s525 documentation built on May 19, 2019, 10:48 p.m.