In this paper, (Bayesian Inference for Variance Components Using Only Error Contrasts, Harville, 1974), the author claims $$(y-X\beta)'H^{-1}(y-X\beta)=(y-X\hat\beta)'H^{-1}(y-X\hat\beta)+(\beta-\hat\beta)'(X'H^{-1}X)(\beta-\hat\beta)$$ to be a "well-known relationship", for a linear regression $$y=X\beta+\epsilon,$$ where $$\epsilon\sim\mathcal{N}(0, H).$$

How is this well-known? What is the simplest way of proving this?

  • 1
    $\begingroup$ It is on wikipedia, see 'derivation' there. $\endgroup$
    – user603
    Sep 28 '15 at 7:38
  • $\begingroup$ @user603 Do you mind making the link clearer? Thanks! $\endgroup$ Sep 28 '15 at 7:46
  • $\begingroup$ @user603 Sorry I can't really see how the link solves the problem. To me, in my case, the equation is Var(y)=bias+... Can you elaborate? $\endgroup$ Sep 29 '15 at 2:32
  • 4
    $\begingroup$ @SibbsGambling Note that your equation has two variance-related terms in this formulation of a weighted linear regression. The term on the left is related to variance around the true model (weighted by the precision matrix $H^{-1}$). The first term on the right is related to variance around fitted models. The second term on the right is related to the square of the bias. That's the variance-bias tradeoff. $\endgroup$
    – EdM
    Oct 1 '15 at 20:27

The last term in the equation can be written as

$$ (X\beta - X\hat{\beta})'H^{-1}(X\beta - X\hat{\beta}). $$

In this form the equation is saying something interesting. Assuming $H$ is positive definite and symmetric, so is its inverse. Therefore,we can define an inner product $<x, y>_{H^{-1}} = x'H^{-1}y$, giving us geometry. Then the above equality is essentially saying that, $$ (X\beta - X\hat{\beta}) \perp (y - X\hat{\beta}). $$

I wanted to give you this bit of intuition since a commenter has already left a link to the derivation.

Edit: For Posterity


\begin{eqnarray} (y-X \beta)'H^{-1}(y-X \beta) &=& y'H^{-1}y &-& 2y'H^{-1}X \beta &+& \beta'X'H^{-1}X\beta \\ &=& (A) &-& (B) &+& (C) \end{eqnarray}


$$ (y-X\hat\beta)'H^{-1}(y-X\hat\beta)+(\beta-\hat\beta)'(X'H^{-1}X)(\beta-\hat\beta) $$ \begin{eqnarray} &=& y'H^{-1}y &- 2y'H^{-1}X\hat{\beta} &+ \hat{\beta}'X'H^{-1}X\hat{\beta} &+ \beta X'H^{-1}X\beta &- 2\hat{\beta}X'H^{-1}X\beta &+ \hat{\beta}'X'H^{-1}X\hat{\beta} \\ &=& (A) &- (D) &+ (E) &+ (C) &- (F) &+ (E) \end{eqnarray}


$$ \hat{\beta} = (X'H^{-1}X)^{-1}X'H^{-1}y $$

By plugging in the relation you can show that (B) = (F), and that 2(E) = (D). All done.

  • $\begingroup$ Sorry I can't really see how the link solves the problem. To me, in my case, the equation is Var(y)=bias+... Can you elaborate? $\endgroup$ Sep 29 '15 at 2:31
  • $\begingroup$ @SibbsGambling edited my answer including derivation. $\endgroup$ Oct 1 '15 at 21:35
  • $\begingroup$ @jlimahaverford aren't you forgetting the $y$ in the end of the formula for $\hat{\beta}$? $\endgroup$
    – Gumeo
    Oct 3 '15 at 15:29

They arrive at this identity by a technique called completing the square. The left hand side is in a quadratic form, so start by multiplying it out

$$ (y-X\beta)'H^{-1}(y-X\beta)= y'H^{-1}y - 2y'H^{-1}X\beta + \beta'X'H^{-1} X\beta $$

continue on and then rewrite in terms of $\hat{\beta} = (X'H^{-1}X)^{-1}X'H^{-1}y$. The algebra is kind of long but googling completing the square in Bayesian regression and you can find plenty of hints. For example, see the wikipedia on Bayesian linear regression, and other CrossValided answers regarding completing the square, like here.


If you know your matrix algebra, then this should be doable by multiplying everything out and verifying that you indeed have the same on both sides. This is what jlimahaverford has demonstrated.

To be able to do this you need the formula for the estimate of $\hat{\beta}$. We can derive the formula in a similar manner as for linear regression when we have uncorrelated error terms. The trick is to standardize.

Here is some information on how to standardize a RV that comes from a multivariate normal distribution. Let's assume that you have $$ \mathbf{X}\sim \mathcal{N}(\mu,\Sigma). $$ $\Sigma$ is positive definite, so you can factorize it as $\Sigma = PP^T$. Now the random variable $$ \mathbf{Y}=P^{-1}(\mathbf{X}-\mu) $$ comes from the distribution $\mathcal{N}(0,I)$. Now we can use this trick for our problem to find $\hat{\beta}$. Let's factorize $H=PP^T$. We have $$ \begin{align} y&=X\beta+\epsilon\\ P^{-1}y &= P^{-1}X\beta + P^{-1}\epsilon \end{align} $$ Now $\epsilon$ has been standardized, such that $\text{cov}(P^{-1}\epsilon)=I$, so we can now treat this as a simple multiple linear regression model where: $$ \tilde{X}=P^{-1}X,\qquad \tilde{y}=P^{-1}y\quad\text{and}\quad \tilde{\epsilon}=P^{-1}\epsilon. $$ So we have the regression problem: $$ \tilde{y}=\tilde{X}\beta+\tilde{\epsilon} $$ The formula for $\hat{\beta}$ is $$ \begin{align} \hat{\beta} &= (\tilde{X}^T\tilde{X})^{-1}\tilde{X}^T\tilde{y}\\ &=((P^{-1}X)^TP^{-1}X)^{-1}(P^{-1}X)^TP^{-1}y\\ &=(X^T(PP^T)^{-1}X)^{-1}X(PP^T)^{-1}y\\ &=(X^TH^{-1}X)^{-1}XH^{-1}y \end{align} $$ This is the key to do this, the rest is the algebraic manipulation demonstrated in the solution by jlimahaverford.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.