A Gibbs Sampler for Multivariate Linear Regression

TL;DR

This paper extends Kelly's Gibbs sampling algorithm for linear regression to handle multiple responses and models the covariate prior with a Dirichlet process, enabling flexible, data-driven modeling of complex astrophysical relations.

Contribution

The paper introduces a generalized Gibbs sampler for multivariate regression with Dirichlet process priors, enhancing modeling flexibility for complex data.

Findings

Successfully modeled galaxy cluster properties with the extended algorithm

Demonstrated the ability to learn the number of mixture components from data

Provided an R implementation of the Gibbs sampler

Abstract

Kelly (2007, hereafter K07) described an efficient algorithm, using Gibbs sampling, for performing linear regression in the fairly general case where non-zero measurement errors exist for both the covariates and response variables, where these measurements may be correlated (for the same data point), where the response variable is affected by intrinsic scatter in addition to measurement error, and where the prior distribution of covariates is modeled by a flexible mixture of Gaussians rather than assumed to be uniform. Here I extend the K07 algorithm in two ways. First, the procedure is generalized to the case of multiple response variables. Second, I describe how to model the prior distribution of covariates using a Dirichlet process, which can be thought of as a Gaussian mixture where the number of mixture components is learned from the data. I present an example of multivariate regression using the extended algorithm, namely fitting scaling relations of the gas mass, temperature, and luminosity of dynamically relaxed galaxy clusters as a function of their mass and redshift. An implementation of the Gibbs sampler in the R language, called LRGS, is provided.