Why do normal equations always have at least one solution?

Question

How can you prove that the normal equations: $(X^\top X)\beta = X^\top Y$ have one or more solutions without the assumption that $X $ is invertible?

My only guess is that it has something to do with generalized inverse, but I am totally lost.

Answer 1

One is tempted to be glib and point out that because the quadratic form

$$\beta \to (Y - X\beta)'(Y - X\beta)$$

is positive semi-definite, there exists a $\beta$ for which it is minimum and that minimum is found (by setting the gradient with respect to $\beta$ to zero) with the normal equations

$$X'X(Y - X\beta) = 0,$$

whence there must be at least one solution regardless of the rank of $X'X$. However, this argument does not seem to be in the spirit of the question, which appears to be a purely algebraic statement. Perhaps it's of interest to understand why such an equation must have a solution and under precisely what conditions. So let's start over and pretend we don't know the connection with least squares.

---

It all comes down to the meaning of $X'$, the transpose of $X$. This will turn out to be a matter of a simple definition, appropriate notation, and the concept of a nondegenerate sesquilinear form. Recall that $X$ is the "design matrix" of $n$ rows (one for each observation) and $p$ columns (one for each variable, including a constant if any). It therefore represents a linear transformation from the vector space $\mathbb V = \mathbb{R}^p$ to $\mathbb W = \mathbb{R}^n$.

The transpose of $X$, thought of as a linear transformation, is a linear transformation of the dual spaces $X': \mathbb{W}^* \to \mathbb{V}^*$. In order to make sense of a composition like $X'X$, then, it is necessary to identify $\mathbb{W}^*$ with $\mathbb{W}$. That's what the usual inner product (sum of squares) on $\mathbb{W}$ does.

There are actually two inner products $g_V$ and $g_W$ defined on $\mathbb V$ and $\mathbb W$ respectively. These are real-valued bilinear symmetric functions that are non-degenerate. The latter means that

$$g_W(u, v) = 0\ \forall u\in \mathbb W \implies v = 0,$$

with analogous statements for $g_V$. Geometrically, these inner products enable us to measure length and angle. The condition $g(u,v)=0$ can be thought of as $u$ being "perpendicular" to or "orthogonal to" $v$. ("Perpendicular" is a geometric term while "orthogonal" is the more general algebraic term. See Michael J. Wichura (2006), The Coordinate-Free Approach to Linear Models at p. 7.)

Nondegeneracy means that only the zero vector is perpendicular to the entire vector space. (This generality means that the results obtained here will apply to the generalized least squares setting, for which $g_W$ is not necessarily the usual inner product given as the sum of products of components, but is some arbitrary nondegenerate form. We could dispense with $g_V$ altogether, defining $X':\mathbb W\to\mathbb V^*$, but I expect many readers to be unfamiliar or uncomfortable with dual spaces and so choose to avoid this formulation.)

With these inner products in hand, the transpose of any linear transformation $X: \mathbb V \to \mathbb W$ is defined by $X': \mathbb W \to \mathbb V$ via

$$g_V(X'(w), v) = g_W(w, X(v))$$

for all $w\in \mathbb W$ and $v\in \mathbb V$. That there actually exists a vector $X'(w) \in \mathbb V$ with this property can be established by writing things out with bases for $\mathbb V$ and $\mathbb W$; that this vector is unique follows from the non-degeneracy of the inner products. For if $v_1$ and $v_2$ are two vectors for which $g_V(v_1,v)=g_V(v_2,v)$ for all $v\in\mathbb V$, then (from the linearity in the first component) $g_V(v_1-v_2,v)=0$ for all $v$ implying $v_1-v_2=0$.

When $\mathbb U \subset \mathbb W,$ write $\mathbb{U}^\perp$ for the set of all vectors perpendicular to every vector in $\mathbb U$. Also as a matter of notation, write $X(\mathbb V)$ for the image of $X$, defined to be the set $\{X(v) | v \in \mathbb V\} \subset \mathbb W$. A fundamental relationship between $X$ and its transpose $X'$ is

$$X'(w) = 0 \iff w \in X(\mathbb V)^\perp.$$

That is, $w$ is in the kernel of $X'$ if and only if $w$ is perpendicular to the image of $X$. This assertion says two things:

1. If $X'(w) = 0$, then $g_W(w, X(v)) = g_V(X'(w),v) = g_V(0,v)=0$ for all $v\in\mathbb V$, which merely means $w$ is perpendicular to $X(V)$.

2. If $w$ is perpendicular to $X(\mathbb V)$, that only means $g_W(w, X(v)) = 0$ for all $v\in\mathbb V$, but this is equivalent to $g_V(X'(w), v) = 0$ and nondegeneracy of $g_V$ implies $X'(w)=0$.

We're actually done now. The analysis has shown that $\mathbb W$ decomposes as a direct product $\mathbb W = X(\mathbb V) \oplus X(\mathbb V)^\perp$. That is, we can take any arbitrary $y \in \mathbb W$ and write it uniquely as $y = y_0 + y^\perp$ with $y_0\in X(\mathbb V)$ and $y^\perp \in X(\mathbb V)^\perp$. That means $y_0$ is of the form $X(\beta)$ for at least one $\beta\in\mathbb V$. Notice, then, that

$$y - X\beta = (y_0 + y^\perp) - y_0 = y^\perp \in X(\mathbb V)^\perp$$

The fundamental relationship says that is the same as the left hand side being in the kernel of $X'$:

$$X'(y - X\beta) = 0,$$

whence $\beta$ solves the normal equations $X'X\beta = X'y.$

---

We are now in a position to give a brief geometric answer to the question (along with some revealing comments): the normal equations have a solution because any $n$-vector $y\in\mathbb W$ decomposes (uniquely) as the sum of a vector $y_0$ in the range of $X$ and another vector $y^\perp$ perpendicular to $y_0$ and $y_0$ is the image of at least one $p$-vector $\beta\in\mathbb V$. The dimension of the image $X(\mathbb V)$ (its rank) is the dimension of the identifiable parameters. The dimension of the kernel of $X$ counts the nontrivial linear relations among the parameters. All parameters are identifiable when $X$ is a one-to-one map from $\mathbb V$ to its image in $\mathbb W$.

It is ultimately useful to dispense with the space $\mathbb V$ altogether and work entirely with the subspace $\mathbb U = X(\mathbb V)\subset\mathbb W$, the "column space" of the matrix $X$. The normal equations amount to orthogonal projection onto $\mathbb U$. That frees us conceptually from being tied to any particular parameterization of the model and shows that least-squares models have an intrinsic dimension independent of how they happen to be parameterized.

---

One interesting outcome of this abstract algebraic demonstration is that we can solve the normal equations in arbitrary vector spaces. The result holds, say, for complex spaces, for spaces over finite fields (where minimizing a sum of squares makes little sense), and even over infinite-dimensional spaces that support suitable sequilinear forms.

Answer 2

This is an old question, but I wanted to expand on a slightly more direct approach for getting to the same conclusion as in whuber's answer.

It follows from the same observation about what $X^{T}$ actually means: it is the matrix representing $X$ acting on the left (or contravariant) term in the inner product. Basically, rather than work with the equation $$ X^{T}X\beta = X^{T}Y$$ we multiply on the left by the transpose of a vector $\gamma$ to get $$ \gamma^{T}X^{T}X\beta = \gamma^{T} X^{T} Y. $$

Note that the first equation is true if and only if the second equation holds for all $\gamma$, so basically we are thinking of $\gamma$ as an arbitrary vector in $\mathbb{R}^{p}$.

Now, using the observation, we can rewrite this as $$ \left< X \gamma, X \beta \right> = \left< X \gamma, Y \right>, $$ where the brackets indicate the inner product (in this case, it is just the standard dot product). This is clearer to interpret since we can think of $X \gamma$ as being the predictions or fitted values of the model with parameter vector $\gamma$.

It becomes even clearer if we subtract the left term from both sides and rearrange, to get $$ \left< X \gamma, (Y - X\beta) \right>. $$ This is finally really interpretable -- it says that the residuals under the model given by $\beta$ are orthogonal to the predictions under $\gamma$ for all $\gamma$. This is equivalent to saying that the residuals under $\beta$ are orthogonal to the entire column space of $X$.

As such, we get to the same conclusion -- that any solution $\beta$ must decompose $Y$ into a component in the column space of $X$ and a component in the orthogonal complement of the column space. In other words, to convince ourselves that there is a unique solution, we need only determine that there is a unique orthogonal projection onto the column space of $X$.

Now, this is surprisingly nontrivial, and whuber's answer discusses it nicely, but I wanted to add a direct, if crude, way of showing it is through the Gram-Schmidt process. Basically, the GS process allows you to construct an orthonormal basis of a space iteratively through the inner product. Basically, at each step you have part of a basis, and given a new vector not in the span of the basis, it shows you how to 'project out' the parts not orthogonal and then normalize. Here you would first find an orthonormal basis for the column space of $X$, and then expand it to an orthonormal basis for all of $\mathbb{R}^{n}$. The subset of vectors added crafting the basis of the column space would form a basis for the orthogonal complement. From there it is direct to simply define the orthogonal projection by how it affects these basis vectors (fixes those in the column space, sends the others to zero). Applying this to $Y$ then gives $\hat{Y} = X\beta$, and solving the system $(X^{T}X)\beta = X^{T}\hat{Y}$ gives you $\beta$.

Answer 3

Besides the comprehensive post by whuber, one should remember that the normal equations $$\mathbf A^\top \mathbf A\mathbf x= \mathbf A^\top \mathbf b$$ is consistent. That is, it has at least one solution.

Note that Result $2.$ implies $\operatorname{rank}\left(\mathbf A^\top\mathbf A\right)=\operatorname {rank}\left(\mathbf A^\top\right) ;$ also $\mathcal C\left(\mathbf A^\top\mathbf A\right)\subseteq\mathcal C\left(\mathbf A^\top\right),$ where $\mathcal C(\mathbf A) $ represents the column space of $\mathbf A. $ Now Result $3.$ ensures $\mathcal C\left(\mathbf A^\top\mathbf A\right) =\mathcal C\left(\mathbf A^\top\right).$

This implies that $\mathbf A^\top \mathbf b\in\mathcal C\left(\mathbf A^\top\mathbf A\right),$ that is, there exists a vector $\mathbf x$ such that $\mathbf A^\top \mathbf b=\mathbf A^\top \mathbf A\mathbf x. $

---

Appendix:

Result $1.$ $\mathcal N\left(\mathbf A^\top\mathbf A\right)=\mathcal N(\mathbf A), $ where $\mathcal N(\mathbf A)$ represents the kernel of $\mathbf A. $

$(\subseteq): ~~\mathbf x\in\mathcal N\left(\mathbf A^\top\mathbf A\right)\implies \mathbf A^\top\mathbf A \mathbf x=\mathbf 0\implies \mathbf x^\top \mathbf A^\top\mathbf A\mathbf x=0\implies \mathbf{Ax}=\mathbf 0\implies \mathbf x\in \mathcal N(\mathbf A). $

$(\supseteq) :~~\mathbf x\in\mathcal N(\mathbf A) \implies \mathbf{Ax}=\mathbf 0\implies \mathbf A^\top\mathbf A \mathbf x=\mathbf 0\implies\mathbf x\in\mathcal N\left(\mathbf A^\top\mathbf A\right).$

$\blacksquare$

Result $2.$ $\operatorname{rank}(\mathbf A) =\operatorname{rank}\left(\mathbf A^\top\mathbf A\right) = \operatorname{rank}\left(\mathbf A^\top\right)=\operatorname{rank}\left(\mathbf A\mathbf A^\top\right).$

$\bf A^\top A$ and $\bf A$ have the same number of columns and by Result $1.$ and rank-nullity theorem, $\operatorname{rank}\left(\mathbf A^\top\mathbf A\right) =\operatorname{rank}(\mathbf A). $ Also, $\operatorname{rank}(\mathbf A^\top\mathbf A) \leq \operatorname{rank}\left(\mathbf A^\top\right) $ which means $\operatorname{rank}(\mathbf A) =\operatorname{rank}\left(\mathbf A^\top\mathbf A\right) \leq \operatorname{rank}\left(\mathbf A^\top\right);$ swapping the transpose, one concludes the stated result.

$\blacksquare$

Result $3.$ $\operatorname{rank}(\mathbf{AB})=\operatorname{rank}(\mathbf{A})\implies \mathcal C(\mathbf{AB})=\mathcal C(\mathbf A).$

$\operatorname{rank}(\mathbf{AB})=\operatorname{rank}(\mathbf{A})\implies \operatorname{dim}(\mathcal C(\mathbf{AB}))=\operatorname{dim}(\mathcal C(\mathbf{A})) .$ Also $\mathcal C(\mathbf{AB})\subseteq\mathcal C(\mathbf A). $ The result follows.

$\blacksquare$

--

Reference:

Linear Algebra and Matrix Analysis for Statistics, Sudipto Banerjee, Anindya Roy, CRC Press, $2014, $ Theorem $9.7, $ p. $262.$

Answer 4

It is easy to show (try for yourself, for an arbitrary number of points, $n$) that the inverse of $X^T X$ exists if there are at least two distinct $x$-values (predictors) in the sample set. Only if all your data have the same values $x_i=x$ (i.e., points stacked in the $y$-direction, along a vertical line), then any line drawn through their mean $\overline{y}$ will have an arbitrary slope (regression coefficient), so that the LSE regression line is then not unique.

Answer 5

In typical regression, X is skinny and therefore certainly not invertible (though it may be left invertible.) It's straightforward to prove (ask if you need help) that if X is skinny and left invertible then X^T * X is invertible. In this case, then there will be exactly one solution. And if X doesn't have full column rank, then X^T * X will not be full rank, and therefore you will have an underdetermined system.