Moore-Penrose pseudoinverse: product on left with another matrix The following problem comes from studying the conditional expectation of a multivariate normal distribution. Let $n\ge2$ be an integer, and let $\Sigma$ be a positive semidefinite, symmetric $n\times n$ matrix of real numbers partitioned as
$$\Sigma=\begin{pmatrix}\Sigma_{a,a}&\Sigma_{a,b}\\\Sigma_{b,a}&\Sigma_{b,b}\end{pmatrix},$$
where $\Sigma_{a,a}$ is $1\times1,$ $\Sigma_{a,b}$ is $1\times(n-1)$ and $\Sigma_{b,b}$ is $(n-1)\times(n-1).$ Is it true that $\Sigma_{a,b}=\Sigma_{a,b}\Sigma_{b,b}^+\Sigma_{b,b}?$ (Here, $\Sigma_{b,b}^+$ is the Moore-Penrose pseudoinverse.)
In the CrossValidated post "Conceptual proof that conditional of a multivariate Gaussian is multivariate Gaussian", someone claims this result. The result is clearly true if $\Sigma_{b,b}$ is invertible, in which case $\Sigma_{b,b}^+=\Sigma_{b,b}^{-1}.$ In addition, I have tried two examples, $\Sigma=0$ and
$$\Sigma=\begin{pmatrix}\!\!\begin{array}{c|cc}1&1&0\\\hline1&1&0\\0&0&0\end{array}\!\!\end{pmatrix},$$
and in both cases we have $\Sigma_{a,b}-\Sigma_{a,b}\Sigma_{b,b}^+\Sigma_{b,b}=0,$ as desired. However, I could not prove the result in generality using the definition or the properties listed in the Wikipedia page.
 A: In the context of this problem, we partition a multivariate Gaussian vector into two subvectors $X_a$ and $X_b$. The goal is to find a matrix $C$ such that $Z:=X_a- C X_b$
is uncorrelated with $X_b$, so that the equality $$\Sigma_{a,b}=C\,\Sigma_{b,b}\tag1$$ holds. The general solution is to take $$C:=\Sigma_{a,b}\Sigma_{b,b}{}^+,$$ where $\Sigma_{b,b}{}^+$ is the Moore-Penrose inverse of $\Sigma_{b,b}$.
The reason why this works, even when $\Sigma_{b,b}$ is not of full rank, depends crucially on the fact that we are dealing with multivariate Gaussians, in which case $\Sigma_{b,b}$ and $\Sigma_{a,b}$ have a special form. Recall that every multivariate Gaussian vector is an affine transformation of some vector $Z$ of independent standard Gaussians. We can then write the subvectors $X_a$ and $X_b$ in the form
$X_a = AZ + \mu_a$, $X_b = BZ + \mu_b$, with $A$ and $B$ matrices of constants. Since the covariance matrix of $Z$ is the identity, we calculate $\Sigma_{a,b} = AB^T$ and $\Sigma_{b,b}=BB^T$.
Using properties of Moore-Penrose inverses, we find $\Sigma_{b,b}{}^+=(BB^T)^+=(B^T)^+B^+$ and verify (1):
$$
C\,\Sigma_{b,b}=\Sigma_{a,b}\Sigma_{b,b}{}^+\Sigma_{b,b}
=AB^T(B^T)^+\underbrace{B^+BB^T}_{B^T}
=A\underbrace{B^T(B^T)^+B^T}_{B^T}=AB^T=\Sigma_{a,b}.
$$
A: The definition of the Moore-Penrose pseudoinverse of $\Sigma_{b,b}$ is $\Sigma_{b,b}^+ = (\Sigma_{b,b}^\top \Sigma_{b,b})^{-1}\Sigma_{b,b}^\top$. This gives
\begin{equation*}
\Sigma_{a,b}\Sigma_{b,b}^+\Sigma_{b,b} = \Sigma_{a,b}(\Sigma_{b,b}^\top\Sigma_{b,b})^{-1}\Sigma_{b,b}^\top \Sigma_{b,b} = \Sigma_{a,b},
\end{equation*}
as desired. Note that this result is really more general than the context of covariance matrices: given $A\in\mathbb{R}^{m\times n}$ and $B\in\mathbb{R}^{p\times n}$ with $\text{rank}(B)=n$, then $AB^+B = A(B^\top B)^{-1}B^\top B = A$.
