I'm confused on how to evaluate the posterior predictive distribution for Bayesian linear regression, past the basic case described [here](http://www.biostat.umn.edu/~sudiptob/ph5470/Lecture6.pdf) on page 3, and copied below.

$$ p(\tilde y \mid y) = \int p(\tilde y \mid \beta, \sigma^2) p(\beta, \sigma^2 \mid y) $$

The basic case is this linear regression model:

$$ y = X \beta + \epsilon, \hspace{10mm} y \sim N(X \beta, \sigma^2)$$

If we use either a uniform prior on $\beta$, with a scale-Inv$\chi^2$ prior on $\sigma^2$, OR the normal-inverse-gamma prior (see [here](http://www.biostat.umn.edu/~brad/7440/BayesianLinearModelGoryDetails.pdf)) the posterior predictive distribution is analytic and is student t. 

What about for this model?
$$ y = X \beta + \epsilon, \hspace{10mm} y \sim N(X \beta, \Sigma)$$

When $y \sim N(X\beta, \Sigma)$, but $\Sigma$ is known, the posterior predictive distribution is multivariate Gaussian. Usually, you don't know $\Sigma$, but have to estimate it. Maybe you say its diagonal and make the diagonal a function of the covariates in some way. This is discussed in the linear regression chapter of Gelman's *Bayesian Data Analysis*.

**Is there an analytic form for the posterior predictive distribution in this case? Can I just plug my estimate of it into a multivariate student t?** If you estimate more than one variance, is the distribution still multivariate student t?

I am asking because say I have some $\tilde y$ already on hand. I want to know whether it is more likely to have been predicted by e.g. linear regression A, linear regression B