Bayesian Robust Regression with the Horseshoe+ Estimator

Abstract

The horseshoe$+$ estimator for Gaussian linear regression models is a novel extension of the horseshoe estimator that enjoys many favourable theoretical properties. We develop the first efficient Gibbs sampling algorithm for the horseshoe$+$ estimator for linear and logistic regression models. Importantly, our sampling algorithm incorporates robust data models that naturally handle non-Gaussian data and are less sensitive to outliers. The resulting software implementation provides a powerful, flexible and robust tool for building prediction and classification models from potentially high-dimensional data and represents the state-of-the-art in Bayesian machine learning techniques.

A Appendix

1.1 A.1 Inverse Gamma Distribution

The inverse gamma probability density function is given by

$$\begin{aligned} p(x | \alpha , \beta ) = \frac{\beta ^\alpha }{\varGamma (\alpha )} x^{-\alpha - 1} \exp \left( - \frac{\beta }{x} \right) , \quad (x > 0), \end{aligned}$$

with shape parameter $(\alpha >0)$ and scale parameter $(\beta > 0)$. The first two moments are

$$\begin{aligned} \mathrm{E} \left( x \right) = \frac{\beta }{\alpha - 1}, \quad \mathrm{Var} \left( x \right) = \frac{\beta ^2}{(\alpha -1)^2 (\alpha - 2)}, \end{aligned}$$

where the mean and variance only exist for $(\alpha >1)$ and $(\alpha >2)$ respectively.

1.2 A.2 Inverse Gaussian Distribution

The inverse Gaussian probability density function is given by

$$\begin{aligned} p(x | \mu , \lambda ) = \left( \frac{\lambda }{2 \pi x^3} \right) ^{\frac{1}{2}} \exp \left( -\frac{\lambda (x - \mu )^2}{2 \mu ^2 x}\right) \!\!, \end{aligned}$$

for $(x > 0)$, where $(\mu > 0)$ is the mean and $(\lambda > 0)$ is the shape parameter. The first two moments are

$$\begin{aligned} \mathrm{E} \left( x \right) = \mu , \quad \mathrm{Var} \left( x \right) = \frac{\mu ^3}{\lambda }. \end{aligned}$$

1.3 A.3 Student-t Distribution

The Student-t distribution probability density function is given by

$$\begin{aligned} p(x | \mu , \sigma ^2, \nu ) = \frac{\varGamma \left( \frac{\nu + 1}{2} \right) }{\varGamma \left( \frac{\nu }{2} \right) \sqrt{\pi \nu \sigma ^2}} \left( 1 + \frac{1}{\nu }\frac{(x - \mu )^2}{\sigma ^2} \right) ^{-\frac{\nu +1}{2}}, \end{aligned}$$

where $(x \in \mathbb {R})$, $(\mu \in \mathbb {R})$, $(\sigma ^2 > 0)$ and the degrees of freedom $(\nu > 0)$. The first two moments are

$$\begin{aligned} \mathrm{E} \left( x \right) = \mu , \quad (\nu > 1), \quad \mathrm{Var} \left( x \right) = \sigma ^2 \left( \frac{\nu }{\nu - 2}\right) \!\!, \end{aligned}$$

(13)

where the mean and variance only exist for $(\nu > 1)$ and $(\nu > 2)$ respectively.

1.4 A.4 Laplace Distribution

The probability density function of the Laplace distribution is

$$\begin{aligned} p(x | \mu , b) = \frac{1}{2b} \exp \left( - \frac{|x - \mu |}{b} \right) , \end{aligned}$$

where $(x \in \mathbb {R})$, $(\mu \in \mathbb {R})$ is the location parameter and $(b > 0)$ is the scale parameter. The first two moments are

$$\begin{aligned} \mathrm{E} \left( x \right) = \mu , \quad \mathrm{Var} \left( x \right) = 2 b^2. \end{aligned}$$

1.5 A.5 Pólya-Gamma Distribution

A random variable x follows a Pólya-gamma distribution [12], $x \sim \mathrm{PG}(b, c)$, if

$$\begin{aligned} x \mathop {=}\limits ^{D} \frac{1}{2\pi ^2} \sum _{k=1}^\infty \frac{g_k}{(k-1/2)^2 + c^2/(4\pi ^2)}, \end{aligned}$$

where $g_k \sim \mathrm{Ga}(b,1)$ are independent gamma random variables, $(b > 0)$ and $(c \in \mathbb {R})$ are the parameters and $\mathop {=}\limits ^{D}$ denotes equality in distribution. The first two moments of x are

$$\begin{aligned} \mathrm{E} \left( x \right) = \frac{b}{2c} \mathrm{tanh} \left( \frac{c}{2}\right) \!\!, \quad \mathrm{Var} \left( x \right) = \frac{b}{4c^3} \left( \mathrm{sinh}(c) - c\right) \mathrm{sech}^2 \left( \frac{c}{2} \right) \!\!. \end{aligned}$$