Conjugate Gaussian Regression With Known Error Variance

Why Conjugate Gaussian Regression Is Useful

Bayesian regression combines information from observed outcomes with prior information about regression coefficients. When the errors are Gaussian with a known variance and the coefficient prior is multivariate normal, this updating has a particularly useful form: the posterior is also multivariate normal. Its mean and covariance can be calculated directly, without simulation.

The update is easiest to understand in terms of precision, the inverse of a variance or covariance matrix. The prior contributes prior precision; the observed design and the known error variance contribute likelihood precision. Their sum is the posterior precision. This tutorial derives that result and applies it to hypothetical data.

Model, Prior, and Posterior

Let $y$ be an $n$-vector of outcomes and $X$ an $n \times k$ design matrix whose columns contain the regressors, including an intercept if the model includes one. Let $\beta$ be the $k$-vector of coefficients. The Gaussian regression model is

y = X\beta + \varepsilon,\qquad \varepsilon \mid \beta \sim N(0,\sigma^2 I_n),

where $\sigma^2$ is known and $I_n$ is the $n \times n$ identity matrix. Equivalently, conditional on $\beta$, $y \sim N(X\beta,\sigma^2 I_n)$. Thus the errors are conditionally independent, have mean zero, and have the same known variance $\sigma^2$.

Suppose the prior distribution of the coefficients is

\beta \sim N(b_0,V_0),

where $b_0$ is a specified $k$-vector and $V_0$ is a specified positive-definite $k \times k$ covariance matrix. This is a proper prior: it integrates to one and assigns a finite, positive variance to every nonzero linear combination of coefficients.

Definition: Under this model and prior, the posterior distribution is the conditional distribution of $\beta$ given the observed $y$, $X$, and known $\sigma^2$. Define the posterior precision matrix $\Lambda_n$ and posterior covariance $V_n$ by $$ \Lambda_n = V_0^{-1} + \frac{X'X}{\sigma^2},\qquad V_n = \Lambda_n^{-1}. $$ Then $$ \beta \mid y,X,\sigma^2 \sim N(b_n,V_n),\qquad b_n = V_n\left(V_0^{-1}b_0+\frac{X'y}{\sigma^2}\right). $$

The data information term $X'X/\sigma^2$ depends on the design and the error variance. A smaller known error variance gives greater likelihood precision for the same $X$; a more precise prior contributes a larger $V_0^{-1}$. These matrix contributions generally interact across coefficients, so one must use matrix calculations rather than update each coefficient independently.

Why the Posterior Has This Form

The likelihood, viewed as a function of $\beta$, is proportional to $\exp\{- (y-X\beta)'(y-X\beta)/(2\sigma^2)\}$. The prior density is proportional to $\exp\{-(\beta-b_0)'V_0^{-1}(\beta-b_0)/2\}$. Multiplying them gives the posterior kernel. The exponent’s terms that depend on $\beta$ are

-\frac{1}{2}\left[ \beta'\left(V_0^{-1}+\frac{X'X}{\sigma^2}\right)\beta -2\beta'\left(V_0^{-1}b_0+\frac{X'y}{\sigma^2}\right) \right].

Write $\Lambda_n=V_0^{-1}+X'X/\sigma^2$ and $h_n=V_0^{-1}b_0+X'y/\sigma^2$. Since $V_n=\Lambda_n^{-1}$, the quadratic expression can be completed as follows:

\beta'\Lambda_n\beta-2\beta'h_n = (\beta-\Lambda_n^{-1}h_n)'\Lambda_n(\beta-\Lambda_n^{-1}h_n) -h_n'\Lambda_n^{-1}h_n.

The final term on the right does not depend on $\beta$, so it only affects the normalizing constant. The remaining kernel is that of a multivariate normal with mean $\Lambda_n^{-1}h_n$ and covariance $\Lambda_n^{-1}$. This proves the stated update. Because $V_0^{-1}$ is positive definite and $X'X/\sigma^2$ is positive semidefinite, their sum $\Lambda_n$ is positive definite and therefore invertible. A proper prior thus yields a well-defined posterior even when the columns of $X$ are linearly dependent.

Worked Posterior Calculations

Worked Example: A Single Coefficient With a Zero-Centered Prior

Consider a hypothetical regression through the origin with $y=X\beta+\varepsilon$, $X=(1,2)'$, $y=(1,2)'$, and known $\sigma^2=1$. Use the prior $\beta \sim N(0,4)$. The two errors are modeled as independent normal errors with mean zero and variance one, as required by the model.

Here $X'X=1^2+2^2=5$, $X'y=1(1)+2(2)=5$, $V_0^{-1}=1/4$, and $V_0^{-1}b_0=0$. Therefore the posterior precision, variance, and mean are

\Lambda_n=\frac{1}{4}+5=\frac{21}{4},\qquad V_n=\left(\frac{21}{4}\right)^{-1}=\frac{4}{21}\approx 0.1905,

b_n=\frac{4}{21}(0+5)=\frac{20}{21}\approx 0.9524.

Thus $\beta \mid y,X,\sigma^2 \sim N(20/21,4/21)$. The posterior mean, about 0.9524 outcome units per regressor unit, is below the likelihood-only least-squares estimate $X'y/X'X=1$, reflecting the zero-centered prior. This is a posterior summary for the coefficient, not a claim that the coefficient is a causal effect.

Worked Example: Two Correlated Coefficients

Suppose a hypothetical two-regressor model has known $\sigma^2=1$, prior mean $b_0=(0,0)'$, and prior precision $V_0^{-1}=\begin{pmatrix}8&0\\0&6\end{pmatrix}$. Let

X=\begin{pmatrix}1&0\\1&1\end{pmatrix},\qquad y=\begin{pmatrix}1\\2\end{pmatrix}.

The errors are assumed conditionally independent with the stated common variance, and the prior is proper because its precision matrix is positive definite. Calculate the design information and data score:

X'X=\begin{pmatrix}2&1\\1&1\end{pmatrix},\qquad X'y=\begin{pmatrix}3\\2\end{pmatrix}.

The posterior precision is the sum of prior and data precision:

\Lambda_n= \begin{pmatrix}8&0\\0&6\end{pmatrix} +\begin{pmatrix}2&1\\1&1\end{pmatrix} = \begin{pmatrix}10&1\\1&7\end{pmatrix}.

Its determinant is $10(7)-1(1)=69$, so the inverse is

V_n=\Lambda_n^{-1} =\frac{1}{69}\begin{pmatrix}7&-1\\-1&10\end{pmatrix}.

In particular, the first posterior variance is $7/69\approx0.1014$, the second is $10/69\approx0.1449$, and the posterior covariance is $-1/69\approx-0.01449$. For the mean, the prior contribution is zero and the data contribution is $X'y=(3,2)'$. Hence

b_n =V_n\begin{pmatrix}3\\2\end{pmatrix} =\frac{1}{69}\begin{pmatrix}7(3)-1(2)\\-1(3)+10(2)\end{pmatrix} =\begin{pmatrix}19/69\\17/69\end{pmatrix} \approx \begin{pmatrix}0.2754\\0.2464\end{pmatrix}.

The result is a joint posterior distribution for both coefficients, not two unrelated univariate updates. Its negative covariance reflects the off-diagonal terms in the inverse posterior precision.

Worked Example: A Nonzero Prior Mean and Nonunit Error Variance

For a hypothetical regression with an intercept and one slope, take

X=\begin{pmatrix}1&0\\1&1\\1&2\end{pmatrix},\qquad y=\begin{pmatrix}2\\3\\5\end{pmatrix},\qquad \sigma^2=4,\qquad b_0=\begin{pmatrix}1\\0\end{pmatrix},\qquad V_0=\begin{pmatrix}4&0\\0&1\end{pmatrix}.

Assume the three errors are conditionally independent normal variables, each with known variance 4. The prior covariance is positive definite. Direct multiplication gives $X'X=\begin{pmatrix}3&3\\3&5\end{pmatrix}$ and $X'y=(10,13)'$. Also $V_0^{-1}=\begin{pmatrix}1/4&0\\0&1\end{pmatrix}$, so

\Lambda_n= \begin{pmatrix}1/4&0\\0&1\end{pmatrix} +\frac{1}{4}\begin{pmatrix}3&3\\3&5\end{pmatrix} = \begin{pmatrix}1&3/4\\3/4&9/4\end{pmatrix}.

The determinant is $1(9/4)-(3/4)^2=27/16$. Thus

V_n=\frac{1}{27/16} \begin{pmatrix}9/4&-3/4\\-3/4&1\end{pmatrix} = \begin{pmatrix}4/3&-4/9\\-4/9&16/27\end{pmatrix}.

The prior contribution to the posterior-mean numerator is $V_0^{-1}b_0=(1/4,0)'$, and the data contribution is $X'y/\sigma^2=(10/4,13/4)'=(5/2,13/4)'$. Their sum is $(11/4,13/4)'$. Therefore

b_n = \begin{pmatrix}4/3&-4/9\\-4/9&16/27\end{pmatrix} \begin{pmatrix}11/4\\13/4\end{pmatrix} = \begin{pmatrix}20/9\\19/27\end{pmatrix} \approx \begin{pmatrix}2.222\\0.7037\end{pmatrix}.

The posterior mean intercept is about 2.222 outcome units, and the posterior mean slope is about 0.7037 outcome units per regressor unit. Dividing $X'X$ and $X'y$ by 4 is essential: the known error variance scales the likelihood information and data score, but it does not scale the prior precision.

Interpretation, Assumptions, and a Common Pitfall

The posterior covariance $V_n$ describes uncertainty about the coefficient vector after observing the data under the specified model and prior. Its diagonal entries are posterior variances; its off-diagonal entries describe posterior covariances. For any coefficient combination $a'\beta$, its posterior mean is $a'b_n$ and its posterior variance is $a'V_na$. These are conditional on the design, the known error variance, and the prior specification.

The conjugate formula relies on a linear conditional mean, Gaussian errors with known common variance and no cross-observation error covariance, and a multivariate normal prior independent of the data-generating errors. If errors instead have known covariance matrix $\Omega$, the likelihood information terms become $X'\Omega^{-1}X$ and $X'\Omega^{-1}y$. When the error variance is unknown, it cannot simply be treated as a known number; a different joint prior and posterior calculation is needed.

Pitfall: scaling the prior precision by the error variance. It is incorrect to use $(V_0^{-1}+X'X)/\sigma^2$ when $V_0$ is the prior covariance as defined here. The prior term is $V_0^{-1}$; only the likelihood terms $X'X$ and $X'y$ are divided by $\sigma^2$.
Pitfall: inverting a matrix entry by entry. In general, the inverse of a matrix is not formed by taking the reciprocal of each element. For a two-by-two matrix, use its determinant and cofactor formula, and verify that the product with the original matrix is the identity.
Pitfall: treating the posterior mean as an identified causal effect. A Gaussian update produces a posterior conditional on the model and prior. Causal interpretation additionally depends on the design and identifying assumptions; conjugacy does not establish those assumptions.

Key takeaway: With known Gaussian error variance and a proper multivariate normal coefficient prior, posterior precision equals prior precision plus $X'X/\sigma^2$. The posterior mean is that posterior covariance multiplied by the sum of the prior and data precision-weighted coefficient information.

Check Your Understanding

Use the posterior formulas and assumptions above to answer these questions.

If the prior mean is zero, which term in the posterior-mean numerator vanishes?
Why does a positive-definite prior covariance ensure that the posterior precision is invertible, even if $X$ is rank deficient?
In a model with known $\sigma^2=9$, what factor multiplies $X'X$ in the posterior precision?
What does a negative off-diagonal entry in $V_n$ describe?
Why does conjugacy alone not justify interpreting a regression coefficient causally?