Marginal Student posterior
Corollary taught in Bayesian linear regression, part of Bayesian statistics
Marginal Student posterior
| Type | Corollary |
| Taught in | Bayesian linear regression |
| Part | Bayesian statistics |
| Worked examples | 2 |
\boxed{ {\boldsymbol{\beta}}\mid{\boldsymbol{y}},{\color{blue}\boldsymbol{X}} \sim t_{n-p}\!\left( \widehat{\boldsymbol{\beta}}, s^2({\color{blue}\boldsymbol{X}}^\top{\color{blue}\boldsymbol{X}})^{-1}\right),} where the displayed matrix is the multivariate t scale matrix. In particular, \beta_j\mid{\boldsymbol{y}},{\color{blue}\boldsymbol{X}} \sim t_{n-p}\!\left(\widehat\beta_j,s\sqrt{c_{jj}}\right), so an equal-tailed 100(1-\alpha)\% credible interval is \widehat\beta_j\ \pm\ t_{1-\alpha/2,n-p}s\sqrt{c_{jj}}. Its endpoints coincide with the usual classical interval under this reference prior. The interpretations do not coincide: the Bayesian statement concerns posterior uncertainty about \beta_j for the observed data, whereas the frequentist statement concerns repeated-sample coverage by a random interval.
Examples
Q24. Bayesian regression and prediction
Question page: Q24. Bayesian regression and prediction
[12 marks]
Consider the regression model Y_i=\beta_0+\beta_1x_i+\varepsilon_i, \qquad \varepsilon_i\overset{\mathrm{iid}}{\sim}N(0,\sigma^2), with reference prior p(\beta_0,\beta_1,\sigma^2)\propto1/\sigma^2 for (\beta_0,\beta_1)\in\mathbb{R}^2 and \sigma^2>0. For n=10 observations, the design matrix has rows (1,x_i), and \widehat{\boldsymbol\beta}=\begin{pmatrix}1\\2\end{pmatrix}, \qquad (X^\top X)^{-1}=\begin{pmatrix}1/10&0\\0&1/40\end{pmatrix}, \qquad \operatorname{SSE}=16. The observed predictor values range from -3 to 3.
You may use the standard regression posterior formulas under this prior. Use the inverse-Gamma convention f(v)\propto v^{-a-1}e^{-b/v}. Write t_\nu(m,s) for a Student-t distribution with location m and scale s, and use t_{0.975,8}=2.306.
(a) [2 marks] State the marginal posterior distribution of \sigma^2 and calculate \mathbb{E}(\sigma^2\mid\mathbf y,X).
(b) [3 marks] State the marginal posterior distribution of the slope \beta_1, and calculate its equal-tailed 95% credible interval.
(c) [4 marks] At x_*=0, derive the posterior distributions and equal-tailed 95% intervals for:
(i) the unknown mean response \mu_*=\beta_0+\beta_1x_*;
(ii) a new observation \widetilde Y.
Explain precisely why the second interval is wider.
(d) [3 marks] Calculate h(x_*)=(1,x_*)(X^\top X)^{-1}(1,x_*)^\top at x_*=0 and x_*=4. Compare the posterior predictive scales at these points and explain how the calculation reflects uncertainty from extrapolation.
Definitions. Write D=(\mathbf y,X), v=\sigma^2, and \operatorname{SSE}=\|\mathbf y-X\widehat{\boldsymbol\beta}\|^2. The inverse-Gamma kernel is v^{-a-1}e^{-b/v} on v>0. t_\nu(m,s) denotes location m, scale s, and \nu degrees of freedom: (T-m)/s\sim t_\nu. A scale is not a variance; t_{q,\nu} denotes the q quantile of standard t_\nu.
Assumptions. Use independent Normal errors, fixed covariates, full column rank X, n>p, positive SSE, and the reference prior p(\boldsymbol\beta,v)\propto1/v with respect to d\boldsymbol\beta\,dv. These rank and residual conditions ensure a proper posterior despite the improper prior.
There are p=2 regression coefficients, including the intercept. Define \nu=n-p=8,\qquad s^2=\frac{\operatorname{SSE}}{n-p}=\frac{16}{8}=2, \qquad C=(X^\top X)^{-1}. Here s^2 is the residual variance estimate; it is not the posterior mean of \sigma^2.
Why the regression posterior has this form. The likelihood times the reference prior gives p(\boldsymbol\beta,v\mid D) \propto v^{-n/2-1} \exp\{-\|\mathbf y-X\boldsymbol\beta\|^2/(2v)\}, \qquad v=\sigma^2>0. Expand about the least-squares estimate. The cross term vanishes because X^\top(\mathbf y-X\widehat{\boldsymbol\beta})=\mathbf0, giving \|\mathbf y-X\boldsymbol\beta\|^2 =\operatorname{SSE} +(\boldsymbol\beta-\widehat{\boldsymbol\beta})^\top X^\top X(\boldsymbol\beta-\widehat{\boldsymbol\beta}). The second term is a Normal quadratic. At fixed v it identifies \boldsymbol\beta\mid v,D\sim N_p(\widehat{\boldsymbol\beta},vC). Its unnormalised exponential integrates to (2\pi v)^{p/2}|X^\top X|^{-1/2}, which explains the v^{p/2} factor in the variance marginal below.
(a) Marginal variance posterior. Under the stated reference prior, \sigma^2\mid D\sim\operatorname{InvGamma}\left(\frac{n-p}{2}, \frac{\operatorname{SSE}}2\right). Thus \boxed{\sigma^2\mid D\sim\operatorname{InvGamma}(4,8),\qquad\sigma^2>0.} An \operatorname{InvGamma}(a,b) variable has mean b/(a-1) when a>1, so \boxed{\mathbb{E}(\sigma^2\mid D)=\frac8{4-1}=\frac83\approx2.6667.} For completeness, integrating the coefficient kernel out of the joint posterior contributes a factor (\sigma^2)^{p/2}. Therefore p(\sigma^2\mid D) \propto (\sigma^2)^{-(n/2+1)}(\sigma^2)^{p/2} \exp\!\left(-\frac{\operatorname{SSE}}{2\sigma^2}\right) =(\sigma^2)^{-5}e^{-8/\sigma^2}, confirming the shape 4, not n/2=5.
(b) Slope posterior and credible interval. The conditional coefficient posterior is \boldsymbol\beta\mid\sigma^2,D \sim N_2(\widehat{\boldsymbol\beta},\sigma^2 C). The slope entry of C is 1/40, so \beta_1\mid\sigma^2,D\sim N\!\left(2,\frac{\sigma^2}{40}\right). Integrating over \sigma^2\mid D gives \boxed{\beta_1\mid D\sim t_8\!\left(2,\sqrt{\frac2{40}}\right) =t_8(2,\sqrt{0.05}).} Equivalently, (\beta_1-2)/\sqrt{0.05}\mid D\sim t_8. The equal-tailed 95% credible interval is therefore 2\pm2.306\sqrt{0.05}=2\pm0.5156, \qquad\boxed{(1.484,\,2.516).} The second argument of t_8 is a scale, not a variance. There is no additional division by \sqrt n: sample-size information is already in C.
(c) Mean response versus a new observation at x_*=0. For a general predictor value, set \mathbf z_*=(1,x_*)^\top and write m_*=\mathbf z_*^\top\widehat{\boldsymbol\beta}, \qquad h(x_*)=\mathbf z_*^\top C\mathbf z_*. The linear-transformation rule for a normal vector gives \mu_*\mid\sigma^2,D \sim N\bigl(m_*,\sigma^2 h(x_*)\bigr). For a new response \widetilde Y=\mu_*+\widetilde\varepsilon, the new error is independent of the coefficient uncertainty conditional on \sigma^2,D, and has variance \sigma^2. Thus \widetilde Y\mid\sigma^2,D,x_* \sim N\bigl(m_*,\sigma^2\{1+h(x_*)\}\bigr). Integrating over the inverse-Gamma variance posterior gives \mu_*\mid D\sim t_8\bigl(m_*,s\sqrt{h(x_*)}\bigr), \widetilde Y\mid D,x_*\sim t_8\bigl(m_*,s\sqrt{1+h(x_*)}\bigr). At x_*=0, m_*=(1,0)\begin{pmatrix}1\\2\end{pmatrix}=1, \qquad h(0)=\frac1{10}. Therefore the unknown mean response has posterior \boxed{\mu_*\mid D\sim t_8(1,\sqrt{0.2}),} with equal-tailed 95% credible interval 1\pm2.306\sqrt{0.2}=1\pm1.0313, \qquad\boxed{(-0.031,\,2.031).} The new observation has posterior predictive distribution \boxed{\widetilde Y\mid D,x_*=0\sim t_8(1,\sqrt{2.2}),} with equal-tailed 95% predictive interval 1\pm2.306\sqrt{2.2}=1\pm3.4204, \qquad\boxed{(-2.420,\,4.420).} Both intervals have the same centre and use the same t_8 quantile. The predictive interval is wider because a new observation includes both uncertainty about the regression mean and fresh observational noise. That additional noise produces the extra 1 in 1+h(x_*).
(d) Extrapolation and predictive scale. Using the supplied diagonal matrix, \begin{aligned} h(x_*)&=(1,x_*) \begin{pmatrix}1/10&0\\0&1/40\end{pmatrix} \begin{pmatrix}1\\x_*\end{pmatrix}\\ &=\frac1{10}+\frac{x_*^2}{40}. \end{aligned} Hence \boxed{h(0)=0.1,\qquad h(4)=0.1+\frac{16}{40}=0.5.} The corresponding predictive scales are \boxed{s_{\mathrm{pred}}(0)=\sqrt{2(1+0.1)}=\sqrt{2.2}\approx1.4832,} \boxed{s_{\mathrm{pred}}(4)=\sqrt{2(1+0.5)}=\sqrt3\approx1.7321.} At x_*=4, the predicted centre is 1+2(4)=9, so the predictive distribution is t_8(9,\sqrt3). Its interval is wider than at x_*=0 because it has the same degrees of freedom but a larger scale.
The value 4 lies outside the observed range [-3,3]. Uncertainty about the slope has a greater effect on the fitted mean farther from the centred design, increasing h(x_*). The future-error contribution is unchanged; the increased width comes from greater parameter uncertainty in the prediction. This calculation assumes the linear model remains appropriate beyond the observed range; it does not account for possible failure of that model there.
Key result. The mean-response scale is s\sqrt h; the new-observation scale is s\sqrt{1+h}. The extra 1 represents new observational noise, while h measures uncertainty about the fitted mean.
Q43. Unknown Normal mean and variance: conditional versus marginal
[14 marks]
Suppose X_1,\ldots,X_n\mid\mu,v \overset{\mathrm{iid}}{\sim}N(\mu,v), \qquad p(\mu,v)\propto\frac1v, where \mu\in\mathbb R, v>0, n\ge2, and S=\sum_{i=1}^n(x_i-\bar x)^2>0. Use the inverse-Gamma convention f(v)\propto v^{-a-1}e^{-b/v}.
(a) Establish the identity \sum_{i=1}^n(x_i-\mu)^2=S+n(\mu-\bar x)^2. [2 marks]
(b) Derive the full conditional distributions \mu\mid v,{\boldsymbol{x}} and v\mid\mu,{\boldsymbol{x}}. [4 marks]
(c) Integrate out \mu to derive p(v\mid{\boldsymbol{x}}). Explain why it differs from p(v\mid\mu,{\boldsymbol{x}}). [4 marks]
(d) Describe a Gibbs sampler using part (b). A student instead samples v from its marginal posterior, then \mu conditional on that draw, repeating with fresh random numbers. Is this also a valid way to sample the joint posterior? Are its resulting pairs independent? [4 marks]
Definitions. Write v=\sigma^2, \bar x=n^{-1}\sum_i x_i, and S=\sum_i(x_i-\bar x)^2, the unnormalised sum of squares. We use N(m,v) with variance v and \operatorname{InvGamma}(a,b) with kernel v^{-a-1}e^{-b/v} on v>0. A full conditional fixes the other parameter; a marginal posterior integrates it out.
Assumptions. The Normal observations are conditionally iid, with prior p(\mu,v)\propto1/v relative to d\mu\,dv. Because this prior is improper, require n\ge2 and S>0 for the posterior below to be proper.
(a) Sum-of-squares identity. Write x_i-\mu=(x_i-\bar x)+(\bar x-\mu) and expand: \begin{aligned} \sum_{i=1}^n(x_i-\mu)^2 &=\sum_{i=1}^n(x_i-\bar x)^2 +2(\bar x-\mu)\sum_{i=1}^n(x_i-\bar x) +n(\bar x-\mu)^2\\ &=\boxed{S+n(\mu-\bar x)^2}, \end{aligned} because \sum_i(x_i-\bar x)=0.
(b) Full conditional distributions. The likelihood and joint posterior kernels are L(\mu,v;{\boldsymbol{x}})\propto v^{-n/2}\exp\left\{-\frac{S+n(\mu-\bar x)^2}{2v}\right\}, p(\mu,v\mid{\boldsymbol{x}})\propto v^{-n/2-1}\exp\left\{-\frac{S+n(\mu-\bar x)^2}{2v}\right\}, \qquad \mu\in\mathbb R,\ v>0. For fixed v, retain only factors involving \mu: p(\mu\mid v,{\boldsymbol{x}})\propto \exp\left\{-\frac{n(\mu-\bar x)^2}{2v}\right\}. Recognising a Normal kernel, \boxed{\mu\mid v,{\boldsymbol{x}}\sim N\!\left(\bar x,\frac vn\right).} For fixed \mu, retain all the v-dependent factors: p(v\mid\mu,{\boldsymbol{x}})\propto v^{-n/2-1}\exp\left\{-\frac{S+n(\mu-\bar x)^2}{2v}\right\}. Under the stated inverse-Gamma convention, \boxed{ v\mid\mu,{\boldsymbol{x}}\sim \operatorname{InvGamma}\!\left(\frac n2,\frac{S+n(\mu-\bar x)^2}{2}\right). } Its support is v>0. In this conditional distribution, \mu is held fixed; it has not been replaced by \bar x or integrated out.
(c) Marginal posterior for v. Integrate the joint posterior kernel over all \mu\in\mathbb R: \begin{aligned} p(v\mid{\boldsymbol{x}}) &\propto v^{-n/2-1}e^{-S/(2v)} \int_{-\infty}^{\infty} \exp\left\{-\frac{n(\mu-\bar x)^2}{2v}\right\}\,d\mu\\ &=v^{-n/2-1}e^{-S/(2v)}\sqrt{\frac{2\pi v}{n}}\\ &\propto v^{-(n-1)/2-1}e^{-S/(2v)}. \end{aligned} Therefore \boxed{v\mid{\boldsymbol{x}}\sim\operatorname{InvGamma}\!\left(\frac{n-1}{2},\frac S2\right), \qquad v>0.} The assumptions n\ge2 and S>0 make both parameters positive, so the posterior is proper despite the improper prior.
The marginal distribution accounts for uncertainty in \mu by integrating over it. The integral contributes a factor v^{1/2}, changing the shape from n/2 to (n-1)/2. Even setting \mu=\bar x in the full conditional does not produce the marginal posterior, because it misses this factor.
(d) Gibbs sampling and direct joint sampling. Start from any \mu^{(0)}\in\mathbb R and v^{(0)}>0. At iteration t=1,\ldots,T, draw \mu^{(t)}\sim N\!\left(\bar x,\frac{v^{(t-1)}}n\right), v^{(t)}\sim\operatorname{InvGamma}\!\left( \frac n2,\frac{S+n(\mu^{(t)}-\bar x)^2}{2}\right). Each update uses the latest available value of the other parameter. The resulting Markov chain has the joint posterior as its equilibrium distribution; successive pairs are generally dependent, and convergence from the chosen starting value must be assessed.
The student’s alternative is also valid. Independently for each repetition, draw v^{(t)}\sim\operatorname{InvGamma}\!\left(\frac{n-1}{2},\frac S2\right), \qquad \mu^{(t)}\mid v^{(t)}\sim N\!\left(\bar x,\frac{v^{(t)}}n\right). This samples the factorisation p(\mu,v\mid{\boldsymbol{x}})=p(v\mid{\boldsymbol{x}})p(\mu\mid v,{\boldsymbol{x}}) directly. With fresh independent random numbers for each repetition, the pairs are iid joint-posterior draws, with no burn-in required. The two components within a pair are not independent: the conditional variance of \mu depends on v.
Used in
- Taught in Bayesian linear regression, where the derivation and examples are.
- All 8 questions on Bayesian linear regression
See also
Bayesian credible interval · Gaussian linear model · Error and residual · Ordinary least squares · Sampling distribution of the OLS estimator · Studentised coefficient pivot · Least-squares decomposition · Conditional coefficient posterior