So in order to better acquaint myself with Gibbs sampling, I've been working on a fairly simple linear model, written in Python/R. Basically, I have 2-dimensional input data (the xi) and a scalar output vector (the yi). I am looking to fit a beta vector, i.e. βT * xi = yi + εi (εi is the noise).
So I decided I'd use a Gaussian prior for the betas, plus an inverse gamma prior for their variance, giving me a posterior distribution that is Gaussian (thanks conjugacy!). And so I used the formula from the Wikipedia page on conjugate distributions to get the likelihood function I need to generate random samples of my beta coefficients and their variances, given the data (by their notation, xi's):
So, it would seem in order to do Gibbs sampling I simply iterate through β1, γ1, β2, and γ2, generating samples from the distributions listed on that wikipedia page. My confusion is, what exactly are my data points (the "samples" I'm using to calculate the distribution parameters, which I need to sample my own betas, i.e. the xi's on the Wikipedia page)??? The way I see it, it could be one of two things. For discussion purposes, let's discuss β1 for now.
First thought is, I'll just iterate over each i, generating my data points using (yi - β2*xi,2) / xi,1 — basically, I subtract out the influence of the second factor and divide the "leftovers" by the first factor, to obtain the impact of the first factor on the response variable. Then I can just use that population to get my distribution parameters, and finally can sample my β1 and γ1.
Second thought is, I'll just go through each i, and generate my data points as (yi / xi), meaning I do NOT subtract out the influence of the other factors on the response variable; it treats everything independently.
So far I've been using #1, but even with artificial data sets, I'm finding that it is implying enormous beta variances that don't seem to make any sense — BUT, they are darn close to the sample variances in the populations I'm generating (my "data points" I use to calc those distribution parameters). What am I doing wrong? Is #1 or #2 the right way to do a linear model? Or am I missing something entirely?
Thanks in advance for your help! Please let me know if any info would help.