Normal–inverse-Gamma conjugacy

Theorem taught in Bayesian linear regression, part of Bayesian statistics

Normal–inverse-Gamma conjugacy

Type Theorem
Taught in Bayesian linear regression
Part Bayesian statistics
Worked examples 2
Questions using it 4

{\boldsymbol{\beta}}\mid\sigma^2,{\boldsymbol{y}} \sim\operatorname{N}_p({\boldsymbol{m}}_n,\sigma^2{\color{blue}\boldsymbol{V}}_n), \qquad \sigma^2\mid{\boldsymbol{y}}\sim\operatorname{InvGamma}(a_n,b_n).

Examples

Question 5. Normal–inverse-\chi^2 regression

Question page: Question 5. Normal–inverse-\chi^2 regression

Suppose

y\mid\beta,\sigma^2,X\sim \operatorname{N}_n(X\beta,\sigma^2 I_n),

with prior

\beta\mid\sigma^2\sim \operatorname{N}_p(\beta_0,\sigma^2\Sigma_0), \qquad \sigma^2\sim\operatorname{Inv\text{-}\chi^2}(\mu_0,\sigma_0^2).

Here X is known, \Sigma_0 is positive definite, and \mu_0,\sigma_0^2>0. The inverse-\chi^2 convention is

\frac{\mu_0\sigma_0^2}{\sigma^2}\sim\chi^2_{\mu_0}.

(a) State the inverse-\chi^2 prior kernel and its equivalent inverse-gamma parameters.

(b) State the updated parameters \Sigma_n and \beta_n.

(c) State the updated parameters \mu_n and \sigma_n^2.

(d) Express the joint posterior through a conditional posterior for \beta and a marginal posterior for \sigma^2.

(e) Express the updated scale using a residual sum of squares and a prior penalty. Explain why it is positive.

(f) Is \Sigma_n the conditional covariance matrix of \beta? Is \sigma_n^2 the posterior mean of \sigma^2? Give the correct expressions and the posterior mode of \sigma^2.

(g) State the full conditional posterior of \sigma^2\mid\beta,y,X. How does it differ from its marginal posterior?

(h) Without consulting the updates, complete the square to reconstruct them. Explain why the marginal degrees of freedom are \mu_0+n.

Definitions. There are n responses and p regression coefficients; X is an n\times p known design matrix. Write v=\sigma^2. The matrix \Sigma_0 scales the conditional prior covariance v\Sigma_0. In this question \mu_0 denotes prior degrees of freedom, not a normal mean. Our \operatorname{Inv\text{-}\chi^2}(\nu,s^2) convention means \nu s^2/v\sim\chi^2_\nu, or equivalently v\sim\operatorname{InvGamma}(\nu/2,\nu s^2/2).

(a) The transformation \mu_0\sigma_0^2/\sigma^2\sim\chi^2_{\mu_0} has kernel

\pi(\sigma^2) \propto (\sigma^2)^{-\mu_0/2-1} \exp\Bigl(-\frac{\mu_0\sigma_0^2}{2\sigma^2}\Bigr), \qquad\sigma^2>0.

Matching (\sigma^2)^{-a-1}\exp(-b/\sigma^2) gives a=\mu_0/2 and b=\mu_0\sigma_0^2/2.

For completeness, if U=\mu_0\sigma_0^2/v, then |dU/dv|=\mu_0\sigma_0^2/v^2. Multiplying the \chi^2_{\mu_0} density at U by this Jacobian gives the displayed inverse-gamma kernel. The Jacobian supplies the power v^{-2} and must be retained.

\boxed{ \sigma^2\sim\operatorname{InvGamma}\Bigl(\frac{\mu_0}{2},\frac{\mu_0\sigma_0^2}{2}\Bigr). }

(b) Since \Sigma_0^{-1} is positive definite and X^\top X is positive semidefinite, \Sigma_n^{-1} is positive definite. Precision addition gives

\boxed{\Sigma_n^{-1}=X^\top X+\Sigma_0^{-1},} \qquad \boxed{\beta_n=\Sigma_n(X^\top y+\Sigma_0^{-1}\beta_0).}

The data and prior precision contributions add, and their combined precision determines the posterior centre.

At fixed v, the likelihood contributes X^\top X/v to the precision of \beta, and the prior contributes \Sigma_0^{-1}/v. Their sum is \Sigma_n^{-1}/v. A proper positive-definite prior therefore gives a proper conditional update even if X is not of full column rank; an unregularised least-squares inverse (X^\top X)^{-1} is not needed.

(c)

\boxed{\mu_n=\mu_0+n,} \qquad \boxed{ \sigma_n^2 = \frac{ \mu_0\sigma_0^2+y^\top y+\beta_0^\top\Sigma_0^{-1}\beta_0 -\beta_n^\top\Sigma_n^{-1}\beta_n }{\mu_n}. }

(d)

\boxed{\beta\mid\sigma^2,y,X\sim \operatorname{N}_p(\beta_n,\sigma^2\Sigma_n),} \qquad \boxed{\sigma^2\mid y,X\sim\operatorname{Inv\text{-}\chi^2}(\mu_n,\sigma_n^2).}

Together these specify a normal–inverse-\chi^2 joint posterior.

This is a conditional–marginal factorisation. To draw an independent posterior pair, draw v from the marginal inverse-\chi^2 distribution and then \beta from its normal distribution conditional on that v. The components within a pair are dependent through the conditional covariance, even though different pairs can be generated independently.

(e) Equivalently,

\mu_n\sigma_n^2 = \mu_0\sigma_0^2 +(y-X\beta_n)^\top(y-X\beta_n) +(\beta_n-\beta_0)^\top\Sigma_0^{-1}(\beta_n-\beta_0).

Both quadratic forms are nonnegative and \mu_0\sigma_0^2>0, so \sigma_n^2>0. The algebra is completed in (h).

(f) No. The conditional covariance is \sigma^2\Sigma_n, not \Sigma_n:

\boxed{\operatorname{Var}(\beta\mid\sigma^2,y,X)=\sigma^2\Sigma_n.}

The parameter \sigma_n^2 is a scale. By the inverse-gamma equivalence in (a),

\boxed{ \mathbb{E}(\sigma^2\mid y,X) = \frac{\mu_n\sigma_n^2}{\mu_n-2}, \qquad\mu_n>2, }

and

\boxed{ \operatorname{mode}(\sigma^2\mid y,X) = \frac{\mu_n\sigma_n^2}{\mu_n+2}. }

(g)

\boxed{ \sigma^2\mid\beta,y,X \sim \operatorname{Inv\text{-}\chi^2}\!\left( \mu_n+p,\; \frac{\mu_n\sigma_n^2+(\beta-\beta_n)^\top\Sigma_n^{-1}(\beta-\beta_n)}{\mu_n+p} \right). }

The full conditional fixes \beta. The marginal posterior integrates \beta out, which removes p from the degrees of freedom.

The difference comes from the normal density’s factor v^{-p/2}. It remains when \beta is fixed and therefore adds p/2 to the inverse-gamma shape. Integrating over \beta instead contributes a factor v^{p/2}, cancelling it. Substituting \beta=\beta_n in the full conditional would leave the wrong degrees of freedom.

(h) Expand the likelihood quadratic and the prior quadratic:

\begin{aligned} &(y-X\beta)^\top(y-X\beta) +(\beta-\beta_0)^\top\Sigma_0^{-1}(\beta-\beta_0)\\ &\quad= \beta^\top(X^\top X+\Sigma_0^{-1})\beta -2\beta^\top(X^\top y+\Sigma_0^{-1}\beta_0) +y^\top y+\beta_0^\top\Sigma_0^{-1}\beta_0\\ &\quad= (\beta-\beta_n)^\top\Sigma_n^{-1}(\beta-\beta_n) +y^\top y+\beta_0^\top\Sigma_0^{-1}\beta_0 -\beta_n^\top\Sigma_n^{-1}\beta_n. \end{aligned}

Adding \mu_0\sigma_0^2, the terms free of \beta are the numerator of \mu_n\sigma_n^2 in (c), with \mu_n=\mu_0+n. Thus

\boxed{ Q(\beta) = (\beta-\beta_n)^\top\Sigma_n^{-1}(\beta-\beta_n) +\mu_n\sigma_n^2. }

The joint posterior kernel is

(\sigma^2)^{-(\mu_0+n+p)/2-1} \exp\Bigl(-\frac{Q(\beta)}{2\sigma^2}\Bigr).

The integral over \beta of the normal kernel equals (2\pi\sigma^2)^{p/2}|\Sigma_n|^{1/2}. The remaining power of \sigma^2 is

-\frac{\mu_0+n+p}{2}-1+\frac p2 = -\frac{\mu_0+n}{2}-1.

\therefore \text{the marginal inverse-}\chi^2\text{ degrees of freedom are }\mu_0+n.

Q58. Joint Normal conjugacy and exact sampling

Question page: Q58. Joint Normal conjugacy and exact sampling

[20 marks; approximately 25 minutes]

Let X_1,\ldots,X_n\mid\mu,\sigma^2 \overset{\mathrm{iid}}{\sim}N(\mu,\sigma^2), with both \mu\in\mathbb R and \sigma^2>0 unknown, and n\ge1. Write \theta=(\mu,\sigma^2). Prior densities below are with respect to d\mu\,d\sigma^2. Use the inverse-Gamma convention p(v)\propto v^{-a-1}e^{-b/v} for v>0.

(i) Show the Normal belongs to the exponential family in the form f(x\mid\theta) =\frac{b(x)\exp\{c_1(\theta)T_1(x)+c_2(\theta)T_2(x)\}} {a(\theta)}. Identify b,a,c_1,c_2,T_1,T_2 (signs go in the c_j). Hence write the likelihood of the whole sample in terms of \sum_i x_i and \sum_i x_i^2. [4 marks]

(ii) The conjugate prior for this family is \pi(\theta)\propto \frac{\exp\{c_1(\theta)\alpha_1+c_2(\theta)\alpha_2\}} {a(\theta)^\beta}. Write it out explicitly in \mu and \sigma^2. By completing the square in \mu, show it factorises as \mu\mid\sigma^2\sim N\!\left(m_0,\frac{\sigma^2}{\kappa_0}\right), \qquad\sigma^2\sim\operatorname{InvGamma}(a_0,b_0), and find m_0,\kappa_0,a_0,b_0 in terms of \alpha_1,\alpha_2,\beta. State the conditions on the hyperparameters for a proper distribution. [6 marks]

(iii) Show that the posterior is in the same family, with hyperparameters (\alpha_1+\sum_i x_i,\,\alpha_2+\sum_i x_i^2,\,\beta+n). Then give the posterior m_n,\kappa_n,a_n,b_n in terms of \bar x, s^2=\sum_i(x_i-\bar x)^2, and the prior quantities. Here s^2 denotes the unnormalised sum of squares, not the sample variance. [4 marks]

(iv) Show that \mathbb{E}[\mu\mid{\boldsymbol{x}}] is a weighted average of m_0 and \bar x, and give the weights. What happens as n\to\infty with fixed prior hyperparameters? [3 marks]

(v) Describe how to draw exact (non-MCMC) samples from the joint posterior of (\mu,\sigma^2), and explain why it works. Assess the claim: “With a flat prior on \mu and an independent inverse-Gamma prior on \sigma^2, exact sampling is unavailable and a Gibbs sampler is necessary.” [3 marks]

Definitions. Write v=\sigma^2, \bar x=n^{-1}\sum_i x_i, and s^2=\sum_i(x_i-\bar x)^2. Here s^2 is the unnormalised sum of squares, not the sample variance. N(m,v) uses variance v; \operatorname{InvGamma}(a,b) has kernel v^{-a-1}e^{-b/v} on v>0. A Normal–inverse-Gamma factorisation means p(\mu,v)=p(v)p(\mu\mid v), not independence of \mu and v.

Assumptions. Use the specified conditionally iid Normal model with n\ge1. The joint prior density is with respect to d\mu\,dv and has exactly the kernel stated in the question. Changing to precision, or starting with a density in natural-parameter coordinates, would require a Jacobian; no such change is made here.

(i) Exponential-family form and likelihood. Expanding (x-\mu)^2=x^2-2\mu x+\mu^2 gives f(x\mid\mu,v) =\frac1{\sqrt{2\pi v}} \exp\left(-\frac{x^2}{2v}+\frac{\mu x}{v}-\frac{\mu^2}{2v}\right) =\frac{(2\pi)^{-1/2}\exp\{(\mu/v)x-x^2/(2v)\}} {\sqrt v\exp\{\mu^2/(2v)\}}. One valid identification, with signs in the c_j, is \boxed{\begin{aligned} b(x)&=(2\pi)^{-1/2},& a(\theta)&=\sqrt v\exp\{\mu^2/(2v)\},\\ T_1(x)&=x,&c_1(\theta)&=\frac\mu v,\\ T_2(x)&=x^2,&c_2(\theta)&=-\frac1{2v}. \end{aligned}} Let S_1=\sum_i x_i and S_2=\sum_i x_i^2. Then \boxed{L(\mu,v;{\boldsymbol{x}}) =(2\pi v)^{-n/2} \exp\left[-\frac{S_2-2\mu S_1+n\mu^2}{2v}\right].} Thus the sample enters the likelihood through (S_1,S_2).

(ii) Conjugate prior and parameter matching. Substitution into the given prior formula yields the joint density kernel \pi(\mu,v) \propto v^{-\beta/2} \exp\left[-\frac{\beta\mu^2-2\alpha_1\mu+\alpha_2}{2v}\right]. For \beta>0, complete the square: \beta\mu^2-2\alpha_1\mu+\alpha_2 =\beta\left(\mu-\frac{\alpha_1}{\beta}\right)^2 +\alpha_2-\frac{\alpha_1^2}{\beta}. Set \boxed{m_0=\frac{\alpha_1}{\beta},\qquad \kappa_0=\beta,\qquad a_0=\frac{\beta-3}{2},\qquad b_0=\frac12\left(\alpha_2-\frac{\alpha_1^2}{\beta}\right).} To see the shape parameter, factor the kernel as \pi(\mu,v)\propto \underbrace{v^{-1/2} \exp\left[-\frac{\kappa_0(\mu-m_0)^2}{2v}\right]} _{\text{Normal conditional kernel, including its }v^{-1/2}\text{ factor}} \underbrace{v^{-(\beta-1)/2}e^{-b_0/v}} _{\text{inverse-Gamma marginal kernel}}. The Normal factor has a normalising multiplier \sqrt{\kappa_0/(2\pi)} that does not depend on v. The remaining exponent therefore satisfies a_0+1=(\beta-1)/2, giving a_0=(\beta-3)/2. Hence \mu\mid v\sim N(m_0,v/\kappa_0) and v\sim\operatorname{InvGamma}(a_0,b_0).

The necessary and sufficient conditions for a proper prior are \boxed{\alpha_1\in\mathbb R,\qquad\beta>3,\qquad \alpha_2>\frac{\alpha_1^2}{\beta}.} They ensure \kappa_0,a_0,b_0>0. The specified three-hyperparameter prior is a subfamily of the general Normal–inverse-Gamma family: it imposes \kappa_0=2a_0+3.

(iii) Posterior update. Multiplying the prior kernel by the likelihood gives \pi(\mu,v\mid{\boldsymbol{x}}) \propto v^{-(\beta+n)/2} \exp\left[-\frac{(\beta+n)\mu^2 -2(\alpha_1+S_1)\mu+(\alpha_2+S_2)}{2v}\right]. This has exactly the same form, with \boxed{\alpha_1'=\alpha_1+S_1,\qquad \alpha_2'=\alpha_2+S_2,\qquad\beta'=\beta+n.} Using S_1=n\bar x, S_2=s^2+n\bar x^2, and the prior identities \alpha_1=\kappa_0m_0 and \alpha_2=2b_0+\kappa_0m_0^2, we obtain \boxed{\begin{aligned} \kappa_n&=\kappa_0+n,\\ m_n&=\frac{\kappa_0m_0+n\bar x}{\kappa_0+n},\\ a_n&=a_0+\frac n2,\\ b_n&=b_0+\frac12\left[s^2+ \frac{\kappa_0n}{\kappa_0+n}(\bar x-m_0)^2\right]. \end{aligned}} In particular, the final update follows from \begin{aligned} b_n &=\frac12\left[\alpha_2+S_2- \frac{(\alpha_1+S_1)^2}{\beta+n}\right]\\ &=b_0+\frac12\left[s^2+\kappa_0m_0^2+n\bar x^2 -\frac{(\kappa_0m_0+n\bar x)^2}{\kappa_0+n}\right]\\ &=b_0+\frac12\left[s^2+ \frac{\kappa_0n}{\kappa_0+n}(\bar x-m_0)^2\right]. \end{aligned} Therefore \boxed{\mu\mid v,{\boldsymbol{x}}\sim N(m_n,v/\kappa_n),\qquad v\mid{\boldsymbol{x}}\sim\operatorname{InvGamma}(a_n,b_n).} The v distribution here is the marginal posterior, with \mu integrated out. All updates remain proper when the prior is proper.

(iv) Posterior mean and weights. For n\ge1 and a proper prior, a_n=a_0+n/2>1/2, which ensures that the marginal posterior mean of \mu exists. Iterated expectation gives \mathbb{E}[\mu\mid{\boldsymbol{x}}] =\mathbb{E}\!\left[\mathbb{E}[\mu\mid v,{\boldsymbol{x}}]\mid{\boldsymbol{x}}\right] =m_n. Consequently, \boxed{\mathbb{E}[\mu\mid{\boldsymbol{x}}] =\underbrace{\frac{\kappa_0}{\kappa_0+n}}_{w_0}m_0 +\underbrace{\frac n{\kappa_0+n}}_{w_x}\bar x.} The weights are positive and sum to one. For fixed prior hyperparameters, w_0\to0 and w_x\to1: the prior contribution vanishes and the posterior mean approaches the sample mean. Under iid sampling at fixed true parameters, \bar x\to\mu_{\mathrm{true}} almost surely, so m_n\to\mu_{\mathrm{true}} almost surely as well.

(v) Exact joint-posterior sampling. For each j=1,\ldots,M, independently of all previous pairs:

1. Draw v^{(j)}\sim\operatorname{InvGamma}(a_n,b_n) from the marginal posterior of v=\sigma^2.

2. Draw \mu^{(j)}\sim N(m_n,v^{(j)}/\kappa_n), conditional on that newly drawn variance.

3. Store (\mu^{(j)},v^{(j)}).

For example, step 1 can use G^{(j)}\sim\operatorname{Gamma}(a_n,b_n) in shape–rate form and set v^{(j)}=1/G^{(j)}. With fresh random numbers in each repetition, this produces iid pairs from the joint posterior because p(\mu,v\mid{\boldsymbol{x}})=p(v\mid{\boldsymbol{x}})\,p(\mu\mid v,{\boldsymbol{x}}). The two coordinates within a pair are generally dependent, but distinct pairs are independent. There is no Markov-chain convergence, burn-in, or accept/reject step to manage in this construction.

Assessing the flat-prior claim. The claim is false for a flat mean prior and an independent inverse-Gamma variance prior. To distinguish its hyperparameters from those above, suppose p(\mu)\propto1,\qquad v\sim\operatorname{InvGamma}(a_*,b_*),\qquad a_*,b_*>0. Although this joint prior is improper in \mu, for n\ge1 the posterior is proper. Using \sum_i(x_i-\mu)^2=s^2+n(\mu-\bar x)^2, its kernel is p(\mu,v\mid{\boldsymbol{x}}) \propto v^{-a_*-1-n/2} \exp\left[-\frac{b_*+s^2/2}{v} -\frac{n(\mu-\bar x)^2}{2v}\right]. Integrating over \mu contributes a factor proportional to \sqrt v, since \int_{-\infty}^{\infty} \exp\left[-\frac{n(\mu-\bar x)^2}{2v}\right]d\mu =\sqrt{\frac{2\pi v}{n}}. Hence \boxed{v\mid{\boldsymbol{x}}\sim \operatorname{InvGamma}\!\left(a_*+\frac{n-1}{2},\ b_*+\frac{s^2}{2}\right), \qquad\mu\mid v,{\boldsymbol{x}}\sim N(\bar x,v/n).} One can again draw from the variance marginal and then the Normal conditional, producing independent exact posterior pairs. Gibbs sampling is valid here, but not necessary. The same distinction appears in Question 43 for p(\mu,v)\propto1/v, with its stated propriety conditions.

Key result. Exact direct sampling uses a marginal then a conditional. Gibbs sampling alternates full conditionals. The key issue is whether the needed marginal is available for direct sampling, not merely whether the prior was described as independent.

Used in

See also

Gaussian linear model · Error and residual · Ordinary least squares · Sampling distribution of the OLS estimator · Studentised coefficient pivot · Least-squares decomposition · Conditional coefficient posterior · Marginal variance posterior

Glossary and index