Past exam of the mathematics course of the University of Cambridge / 2009 / ii / Paper 4 / 13I / Solution

Put $M=X^{T}X$, which is positive definite because $X$ has full column rank. The log likelihood differs from $-\Vert Y-X\beta {\Vert }^2/(2{\sigma }^2)$ by a term independent of $\beta$. Differentiation gives the normal equations $M\hat{\beta }=X^{T}Y$, henceThe distribution follows by writing $\hat{\beta }=\beta +M^{-1}{X}^T\epsilon$ and applying the linear transformation rule for a multivariate normal distribution.

Let $s={\sigma }^2$ and $R=\Vert Y-X\hat{\beta }{\Vert }^2$. Orthogonality of the residual to the columns of $X$ gives the exact decompositionThe prior is density proportional to $s^{-1}$ with respect to $d\beta ds$. Thus the joint posterior distribution isCompleting the quadratic identifiesThe sampling distribution treats the true $\beta$ as fixed and varies the data, whereas the conditional posterior fixes the observed data and distributes $\beta$. Their covariance formulas agree, but their centres and their interpretations differ.

For a new row vector $x^{\ast }$, let $h^{\ast }=x^{\ast }{M}^{-1}(x^{\ast }{)}^T$ and suppose its new error is independent. With known ${\sigma }^2$, the maximum-likelihood plug-in prediction uses mean $x^{\ast }\hat{\beta }$ and noise variance ${\sigma }^2$, while a calibrated frequentist prediction must account for the estimator: $y^{\ast }-x^{\ast }\hat{\beta }\sim N(0,{\sigma }^2(1+h^{\ast }))$. If variance is estimated, the maximum-likelihood estimator is $R/n$, whereas the exact frequentist prediction pivot uses $\hat{s}=R/(n-p)$ and has Student t-distribution with $n-p$ degrees of freedom.

In the Bayesian calculation, integrating $\beta$ first gives $y^{\ast }\mid s,X,Y\sim N(x^{\ast }\hat{\beta },s(1+h^{\ast }))$. Integrating $\beta$ out of the joint posterior gives $s\mid X,Y\sim \mathrm{InvGamma}((n-p)/2,R/2)$. For $\nu =n-p$, integration of the predictive normal against this inverse-gamma distribution has density proportional toThe substitution $t=[R+(y^{\ast }-x^{\ast }\hat{\beta }{)}^2/(1+h^{\ast })]/(2s)$ reduces the integral to a Gamma function. Normalizing identifies the Student t-distribution, givingwhere the second parameter is a scale, not its variance. These Bayesian intervals coincide with the usual exact Gaussian-regression prediction intervals for this prior, although their probability interpretations differ. Posterior propriety requires $R>0$ and $n>p$; the former holds almost surely under the nondegenerate Gaussian model. Prediction of the noiseless conditional mean omits the added one in $1+h^{\ast }$.