Randomized low-rank approximation

Choose an SVD $A=U\mathrm{\Sigma }{V}^{\ast }$ and put $D=\mathrm{\Sigma }{\mathrm{\Sigma }}^{\ast }$ and $R=U^{\ast }(I-P)U$. Then $D$ is diagonal positive semidefinite and $R$ is an orthogonal projector. By [F1]–[F3], $\mathrm{\parallel }(I-P)A{\mathrm{\parallel }}_2^2=\mathrm{\parallel }RDR{\mathrm{\parallel }}_2$ and $\mathrm{\parallel }(I-P)(A{A}^{\ast }{)}^{q}A{\mathrm{\parallel }}_2^2=\mathrm{\parallel }R{D}^{2q+1}R{\mathrm{\parallel }}_2$.

If $\lambda =\mathrm{\parallel }RDR{\mathrm{\parallel }}_2=0$, [step 1.1] gives $\mathrm{\parallel }(I-P)A{\mathrm{\parallel }}_2=0$, so the asserted inequality follows. If $\lambda >0$, choose a unit eigenvector $x$ for the top eigenvalue of $RDR$; then $Rx=x$ and $\lambda =x^{\ast }Dx={\sum }_i{d}_i\mathrm{\mid }{x}_i{\mathrm{\mid }}^2$, with $d_i\ge 0$. The weights $\mathrm{\mid }{x}_i{\mathrm{\mid }}^2$ sum to $1$, so [F4] for $t\mapsto t^{2q+1}$ gives ${\lambda }^{2q+1}\le {\sum }_i{d}_i^{2q+1}\mathrm{\mid }{x}_i{\mathrm{\mid }}^2=x^{\ast }{D}^{2q+1}x\le \mathrm{\parallel }R{D}^{2q+1}R{\mathrm{\parallel }}_2$.

Substitute $\lambda =\mathrm{\parallel }(I-P)A{\mathrm{\parallel }}_2^2$ from step 1.1 into step 2.1 and take the nonnegative $(2q+1)$st root. This yields the asserted inequality, also for $q=0$.

Restrict the measures in [F1] to Borel sets in dimensions one and two. They meet exactly the supplied-measure hypotheses of [F2]. With location zero and variance one, [F2] gives real mass one, mean zero, absolute first-moment integrability, and second moment one. No normal-law result assuming full Choice is used.

The map $T(x,y)=(x\mathrm{/}\sqrt{2},y\mathrm{/}\sqrt{2})$ has determinant $1\mathrm{/}2$. By [F3], the image density is $2{\phi }_R(\sqrt{2}u){\phi }_R(\sqrt{2}v)={\pi }^{-1}{e}^{-(u^2+v^2)}$. The integral transformation follows first for indicator functions from [F3], then for simple functions and their increasing limits. Its mass is one by step 1.1 and [F4]. By the same product integration, $\mathbb{E}\mathrm{\mid }z{\mathrm{\mid }}^2=\mathbb{E}(x^2+y^2)\mathrm{/}2=1$, and both real and imaginary means vanish; their absolute integrability follows from that of $x,y$.

For any positive finite number of entries, repeated [F4] factors the integral of the nonnegative Borel product density into scalar integrals, each one. All coordinate spaces are sigma-finite, being countable unions of bounded boxes. Integrating over a product of Borel entry sets similarly gives the product of their scalar probabilities, which establishes independence. The integral of the density over sets defines a countably additive probability measure: for disjoint sets, the finite sums of their density-weighted indicators increase to the density times the indicator of their union, so monotone convergence, by [F5], gives countable additivity. For zero entries use the one-point space of mass one and the empty product. Thus all finite laws in the definition are justified using only the stated Countable Choice.

If there are no entry coordinates, the matrix is deterministic and its law is unchanged. Otherwise regard the matrix as its finite list of real coordinates. Multiplying the independent scalar densities in [F1] gives joint density $C{e}^{-\mathrm{\parallel }G{\mathrm{\parallel }}_F^2\mathrm{/}2}$ over $\mathbb{R}$ or $C^{\mathrm{\prime }}{e}^{-\mathrm{\parallel }G{\mathrm{\parallel }}_F^2}$ over $\mathbb{C}$. The maps $G\mapsto U^{\ast }GV$ are real-linear isometries for the Frobenius norm; their real Jacobian matrix $O$ satisfies $O^{T}O=I$, hence $\mathrm{\mid }\det O\mathrm{\mid }=1$. By [F2], this real orthogonal map preserves measure, also in the complex case viewed in ${\mathbb{R}}^{2ab}$. Its inverse preserves the Frobenius norm, so integrating over the preimage of any Borel set gives the same density integral over that set. Hence the transformed matrix has the same law.

In the nonempty-coordinate case, the density in [F1] factors over all transformed entries because its exponent is the sum of the entrywise squared moduli and its normalizing constant is the product of the scalar normalizers. Thus the entries remain independent and each retains its original scalar law. For an empty matrix the assertion about entries is vacuous.

Apply the invariance with $U=W$ and $V=I$. For two nonempty disjoint row blocks, integrating out the other entries leaves the product of their block densities, which proves independence and the declared standard block laws. An empty block has the prescribed deterministic law and is independent of every other block, since each event for it has probability zero or one.

Take independent copies of $X,Y$, with covariance entries $(C_X{)}_{ij}=\mathbb{E}(X_i{X}_j)$ and $(C_Y{)}_{ij}=\mathbb{E}(Y_i{Y}_j)$, and represent them as independent linear images of standard Gaussian coordinates using [A1]. For $0<t<1$ put $Z(t)=\sqrt{t}Y+\sqrt{1-t}X$ and $F_{\beta }(z)={\beta }^{-1}\log {\sum }_{i=1}^N{e}^{\beta z_i}$ with $w_i=e^{\beta z_i}\mathrm{/}{\sum }_j{e}^{\beta z_j}$. The soft maximum has at most linear growth, so its Gaussian slice is integrable. On each compact subinterval of $(0,1)$, $\mathrm{\mid }{\mathrm{\partial }}_t{F}_{\beta }(Z(t))\mathrm{\mid }\le \mathrm{\parallel }\mathrm{\nabla }{F}_{\beta }{\mathrm{\parallel }}_2(\mathrm{\parallel }Y{\mathrm{\parallel }}_2\mathrm{/}(2\sqrt{t})+\mathrm{\parallel }X{\mathrm{\parallel }}_2\mathrm{/}(2\sqrt{1-t}))$ has an integrable uniform bound, so [F3] passes the derivative under the finite Gaussian integral. For the coordinate integration by parts, fix all coordinates but $r=g_j$, and let $\phi (r)=(2\pi {)}^{-1\mathrm{/}2}{e}^{-r^2\mathrm{/}2}$. The functions to which it is applied are first derivatives of $F_{\beta }$ composed with linear maps; they and their $r$ derivatives are smooth and bounded for fixed $\beta$. On $[-R,R]$, [F5] with $u=f$ and $v=\phi$ gives ${\int }_{-R}^{R}rf(r)\phi (r)dr={\int }_{-R}^R{f}^{\mathrm{\prime }}(r)\phi (r)dr-[f(r)\phi (r){]}_{-R}^R$. The boundary term tends to zero since $f$ is bounded and $\phi (\pm R)\to 0$. The integrands are dominated by constant multiples of $\mathrm{\mid }r\mathrm{\mid }\phi (r)$ and $\phi (r)$, respectively, so [F4] takes $R\to \mathrm{\infty }$. Both majorants are integrable. Their bounds are uniform in the other coordinates, and [F6] integrates the resulting identity over those coordinates. Applying this coordinatewise to a linear image $W=Lg$ gives $\mathbb{E}[W_{i}f(W)]={\sum }_j\mathrm{Cov}(W_i,W_j)\mathbb{E}[{\mathrm{\partial }}_{j}f(W)]$, including singular covariance matrices. Apply it separately to $X,Y$ in the derivative to obtain $\frac{d}{dt}\mathbb{E}{F}_{\beta }(Z(t))=\frac{1}{2}{\sum }_{i,j}(C_Y-C_X{)}_{ij}\mathbb{E}[{\mathrm{\partial }}_{ij}{F}_{\beta }(Z(t))]$. Since ${\mathrm{\partial }}_{ij}{F}_{\beta }=\beta (w_i{\mathbf{1}}_{i=j}-w_i{w}_j)$, symmetry rewrites this as $\frac{\beta }{4}\mathbb{E}{\sum }_{i,j}{w}_i{w}_j(\mathbb{E}(Y_i-Y_j{)}^2-\mathbb{E}(X_i-X_j{)}^2)\ge 0$ by the assumed increment order.

Integrate the nonnegative derivative on $(0,1)$. As $t\downarrow 0$ or $t\uparrow 1$, $Z(t)$ converges pointwise to $X$ or $Y$ and $\mathrm{\mid }{F}_{\beta }(Z(t))\mathrm{\mid }$ is bounded by a constant plus $\mathrm{\parallel }X{\mathrm{\parallel }}_2+\mathrm{\parallel }Y{\mathrm{\parallel }}_2$, so [F4] gives endpoint continuity and $\mathbb{E}{F}_{\beta }(Y)\ge \mathbb{E}{F}_{\beta }(X)$. For every $z\in {\mathbb{R}}^N$, ${\max }_i{z}_i\le F_{\beta }(z)\le {\max }_i{z}_i+\log N\mathrm{/}\beta$. Letting $\beta \to \mathrm{\infty }$ proves the asserted expected-maximum comparison.

After each scanned column, the retained vectors are orthonormal and span the scanned columns: subtracting their projections gives an orthogonal residual; a zero residual adds no span, while normalizing and multiplying a nonzero residual by a unit complex scalar adds exactly its line. This proves the invariant by finite induction. There are exactly $s$ scan steps, with no division on a zero branch, so the scan terminates with the asserted basis. If $s=0$ or $Y=0$, it returns the empty basis and zero projector.

Record all skip/retain branches and, in the complex case, the first nonzero coordinate at each retain branch. There are finitely many such records. Start with the whole input space as one Borel stratum. Inductively on a fixed partial record, residuals are continuous relative to its stratum; the zero and nonzero tests are relatively closed or open and hence ambient Borel by [F3]. Normalization is continuous where the norm is positive. If the first nonzero coordinate of the normalized residual is $w_j$, its phase multiplier is $\overline{w_j}\mathrm{/}\mathrm{\mid }{w}_j\mathrm{\mid }$, continuous where $w_j\mathrm{\ne }0$; the conditions $w_1=\cdots =w_{j-1}=0$, $w_j\mathrm{\ne }0$ give Borel substrata by the same argument. Thus every complete record has a Borel stratum and a relatively continuous output map of fixed size. For a Borel target set, [F3] expresses its preimage on each stratum as that stratum intersected with an ambient Borel set; the finite union over records is Borel. The rank strata are finite unions of these record strata, so this proves the stated variable-size Borel map, without asserting global continuity. Multiplication $Q\mapsto Q{Q}^{\ast }$ is continuous for each rank; [F4] and finite pasting give a Borel fixed-size projector. If either matrix coordinate dimension is zero, the empty-matrix input space has one point and these conclusions are immediate.

For fixed $A$, the map $\mathrm{\Omega }\mapsto A\mathrm{\Omega }$ is continuous and hence Borel by [F2] when both coordinate dimensions are positive; if either dimension is zero it is constant and Borel directly. Composing with step 2.1 using [F4] gives a Borel range-finder projector. The norm triangle inequality and [F5] give $\mathrm{\mid }\mathrm{\parallel }X{\mathrm{\parallel }}_2-\mathrm{\parallel }Z{\mathrm{\parallel }}_2\mathrm{\mid }\le \mathrm{\parallel }X-Z{\mathrm{\parallel }}_2\le \mathrm{\parallel }X-Z{\mathrm{\parallel }}_F$, so the matrix norm is continuous and the residual error is Borel too. By [F6], $v=Pv+(I-P)v$ is an orthogonal decomposition, whence $\mathrm{\parallel }(I-P)v{\mathrm{\parallel }}_2^2\le \mathrm{\parallel }v{\mathrm{\parallel }}_2^2$ and the residual is bounded by $\mathrm{\parallel }A{\mathrm{\parallel }}_2$. For any probability law $\mu$ of the sketch, [F7] therefore gives $\int \mathrm{\parallel }(I-P)A{\mathrm{\parallel }}_{2}d\mu \le \int \mathrm{\parallel }A{\mathrm{\parallel }}_{2}d\mu =\mathrm{\parallel }A{\mathrm{\parallel }}_2<\mathrm{\infty }$, including $A=0$. Thus the Borel residual has a finite expectation under every such law.

A fixed proper real or complex linear subspace $L⊊F^{\mathrm{\ell }}$ has probability zero for one Gaussian row: choose an orthogonal/unitary coordinate change sending $L$ into $\{z_1=0\}$, use [F2], and integrate the density over this coordinate subspace. It is Lebesgue-null: inside each bounded box $[-M,M{]}^d$ it lies in a coordinate hyperplane coverable by a box of width $2\epsilon$ and bounded remaining sides, whose volume tends to zero with $\epsilon$; [F4] makes that hyperplane null, and countable subadditivity covers the whole hyperplane.

Induct on the number $j\le k$ of revealed rows. The first row is nonzero almost surely by step 1.1 applied to the zero subspace. The dependence event for $j$ rows is measurable because it is the common zero set of the finitely many $j\times j$ minors. If the first $j-1$ rows are independent, their span is proper since $j-1<k<\mathrm{\ell }$; for each fixed such tuple, step 1.1 makes the section of the event that row $j$ lies in their span null. The exceptional dependent tuples are null by induction. Applying [F3] to the nonnegative indicator of dependence gives probability zero for the first $j$ rows to be dependent. The induction is finite and concludes at $j=k$.

Outside this null event, the $k$ rows are linearly independent, hence the row rank equals $k$.

Write $(SGT{)}_{\alpha \beta }={\sum }_{i=1}^r{\sum }_{j=1}^s{S}_{\alpha i}{G}_{ij}{T}_{j\beta }$. Expanding gives $\mathbb{E}\mathrm{\parallel }SGT{\mathrm{\parallel }}_F^2={\sum }_{\alpha =1}^a{\sum }_{\beta =1}^b{\sum }_{i,u=1}^r{\sum }_{j,v=1}^s{S}_{\alpha i}{T}_{j\beta }{S}_{\alpha u}{T}_{v\beta }\mathbb{E}(G_{ij}{G}_{uv})$. Independence, zero means, and unit variances give $\mathbb{E}(G_{ij}{G}_{uv})={\mathbf{1}}_{i=u}{\mathbf{1}}_{j=v}$, so this reduces to ${\sum }_{\alpha ,\beta ,i,j}{S}_{\alpha i}^2{T}_{j\beta }^2=({\sum }_{\alpha ,i}{S}_{\alpha i}^2)({\sum }_{j,\beta }{T}_{j\beta }^2)=\mathrm{\parallel }S{\mathrm{\parallel }}_F^2\mathrm{\parallel }T{\mathrm{\parallel }}_F^2$.

For unit vectors $x\in {\mathbb{R}}^a,y\in {\mathbb{R}}^b$, set $u=S^{T}x$, $v=Ty$, $U=\mathrm{\parallel }u\mathrm{\parallel }$, $V=\mathrm{\parallel }v\mathrm{\parallel }$, and define centered Gaussian processes $X_{xy}=u^{T}Gv+zUV$ and $Y_{xy}=V{g}^{T}u+U{h}^{T}v$, where $z\in \mathbb{R}$, $g\in {\mathbb{R}}^r$, and $h\in {\mathbb{R}}^s$ are independent standard real Gaussian variables/vectors, independent of $G$. For a second pair $(x^{\mathrm{\prime }},y^{\mathrm{\prime }})$, write $u^{\mathrm{\prime }}=S^T{x}^{\mathrm{\prime }}$, $v^{\mathrm{\prime }}=T{y}^{\mathrm{\prime }}$, $U^{\mathrm{\prime }}=\mathrm{\parallel }{u}^{\mathrm{\prime }}\mathrm{\parallel }$, $V^{\mathrm{\prime }}=\mathrm{\parallel }{v}^{\mathrm{\prime }}\mathrm{\parallel }$, $c=u^T{u}^{\mathrm{\prime }}$, and $d=v^T{v}^{\mathrm{\prime }}$. At each index both processes have variance $2U^2{V}^2$, and their cross-covariances are $\mathrm{Cov}(X_{xy},X_{x^{\mathrm{\prime }}{y}^{\mathrm{\prime }}})=cd+U{U}^{\mathrm{\prime }}V{V}^{\mathrm{\prime }}$ and $\mathrm{Cov}(Y_{xy},Y_{x^{\mathrm{\prime }}{y}^{\mathrm{\prime }}})=cV{V}^{\mathrm{\prime }}+U{U}^{\mathrm{\prime }}d$. The difference is $(U{U}^{\mathrm{\prime }}-c)(V{V}^{\mathrm{\prime }}-d)\ge 0$ by Cauchy–Schwarz; subtracting twice the covariance from the sum of the two variances gives $\mathbb{E}(X_{xy}-X_{x^{\mathrm{\prime }}{y}^{\mathrm{\prime }}}{)}^2\le \mathbb{E}(Y_{xy}-Y_{x^{\mathrm{\prime }}{y}^{\mathrm{\prime }}}{)}^2$.

Apply [F2] to increasing finite subsets of normalized rational vectors in the product of unit spheres. They are dense; both sample paths are continuous, and their absolute suprema are bounded by $\mathrm{\parallel }S{\mathrm{\parallel }}_2\mathrm{\parallel }T{\mathrm{\parallel }}_2(\mathrm{\parallel }G{\mathrm{\parallel }}_F+\mathrm{\mid }z\mathrm{\mid })$ for $X$ and $\mathrm{\parallel }T{\mathrm{\parallel }}_2\mathrm{\parallel }Sg{\mathrm{\parallel }}_2+\mathrm{\parallel }S{\mathrm{\parallel }}_2\mathrm{\parallel }{T}^{T}h{\mathrm{\parallel }}_2$ for $Y$. These bounds are integrable because independence and unit coordinate variances give $\mathbb{E}\mathrm{\parallel }G{\mathrm{\parallel }}_F^2=rs$, $\mathbb{E}{z}^2=1$, $\mathbb{E}\mathrm{\parallel }Sg{\mathrm{\parallel }}_2^2=\mathrm{\parallel }S{\mathrm{\parallel }}_F^2$, and $\mathbb{E}\mathrm{\parallel }{T}^{T}h{\mathrm{\parallel }}_2^2=\mathrm{\parallel }T{\mathrm{\parallel }}_F^2$. Dominated convergence [F4] passes the finite comparisons to the full spheres. For fixed $G$, $z\mapsto {\sup }_{x,y}(u^{T}Gv+zUV)$ is convex with integrable growth, so Jensen [F3] gives ${\mathbb{E}}_z{\sup }_{x,y}{X}_{xy}\ge {\sup }_{x,y}{\mathbb{E}}_z{X}_{xy}={\sup }_{x,y}{u}^{T}Gv=\mathrm{\parallel }SGT{\mathrm{\parallel }}_2$.

For every unit $x,y$, $Y_{xy}\le \mathrm{\parallel }T{\mathrm{\parallel }}_2\mathrm{\parallel }Sg{\mathrm{\parallel }}_2+\mathrm{\parallel }S{\mathrm{\parallel }}_2\mathrm{\parallel }{T}^{T}h{\mathrm{\parallel }}_2$. Taking suprema and expectations, [F3] applied to the convex square gives $\mathbb{E}\mathrm{\parallel }Sg{\mathrm{\parallel }}_2\le (\mathbb{E}\mathrm{\parallel }Sg{\mathrm{\parallel }}_2^2{)}^{1\mathrm{/}2}=\mathrm{\parallel }S{\mathrm{\parallel }}_F$ and $\mathbb{E}\mathrm{\parallel }{T}^{T}h{\mathrm{\parallel }}_2\le \mathrm{\parallel }T{\mathrm{\parallel }}_F$. Combining with step 2.1 proves the spectral estimate.

Choose the SVD from [F1] and set $\tilde{\mathrm{\Omega }}=V^{\ast }\mathrm{\Omega }=[{\mathrm{\Omega }}_1^T{\mathrm{\Omega }}_2^T{]}^T$. Left multiplication by $U^{\ast }$ carries the sample range to $\mathrm{Ran}(\mathrm{\Sigma }\tilde{\mathrm{\Omega }})$, preserves the operator norm because it maps unit vectors bijectively to unit vectors by [F3], and carries $Q{Q}^{\ast }$ to the orthogonal projector $P$ onto that range. It suffices to prove the estimate for $\mathrm{\Sigma }$.

Since ${\mathrm{\Omega }}_1$ has full row rank, put ${\mathrm{\Omega }}_1^{†}={\mathrm{\Omega }}_1^{\ast }({\mathrm{\Omega }}_1{\mathrm{\Omega }}_1^{\ast }{)}^{-1}$, so ${\mathrm{\Omega }}_1{\mathrm{\Omega }}_1^{†}=I_k$. Let $Y=\mathrm{\Sigma }[{\mathrm{\Omega }}_1^T{\mathrm{\Omega }}_2^T{]}^T$ and $C=Y{\mathrm{\Omega }}_1^{†}[I_{k}0]$. The range of $C$ lies in $\mathrm{Ran}(Y)$, hence $(I-P)C=0$.

The first $k$ rows of $\mathrm{\Sigma }-C$ vanish because ${\mathrm{\Omega }}_1{\mathrm{\Omega }}_1^{†}=I_k$, and the remaining rows form $[-{\mathrm{\Sigma }}_2{\mathrm{\Omega }}_2{\mathrm{\Omega }}_1^{†}{\mathrm{\Sigma }}_2]$. Thus $\mathrm{\parallel }(I-P)\mathrm{\Sigma }{\mathrm{\parallel }}_2=\mathrm{\parallel }(I-P)(\mathrm{\Sigma }-C){\mathrm{\parallel }}_2\le \mathrm{\parallel }\mathrm{\Sigma }-C{\mathrm{\parallel }}_2\le \mathrm{\parallel }{\mathrm{\Sigma }}_2{\mathrm{\Omega }}_2{\mathrm{\Omega }}_1^{†}{\mathrm{\parallel }}_2+\mathrm{\parallel }{\mathrm{\Sigma }}_2{\mathrm{\parallel }}_2$: the horizontal block is a sum of two zero-padded blocks, so the triangle inequality applies. By [F1] and [F4], $\mathrm{\parallel }{\mathrm{\Sigma }}_2{\mathrm{\parallel }}_2={\sigma }_{k+1}$, including zero-tail and rank-deficient cases.

By [F2], the initial scan on $A\mathrm{\Omega }$ is Borel, since multiplication by fixed $A$ is continuous. On each rank stratum, multiplication by fixed $A^{\ast }$ and then the scan are Borel by [F3], as are multiplication by $A$ and the next scan. There are finitely many possible ranks, so the resulting maps on their disjoint union are Borel. Finite induction over the $q$ iterations proves measurability of the final basis. The empty-basis branches remain empty under multiplication and scanning, and $q=0$ is the initial scan.

On every rank stratum $Q\mapsto Q{Q}^{\ast }$ is continuous, and the matrix norm is continuous by the norm inequality in [F2]'s proof. Hence the final projector and error are Borel on their finite union, including zero samples.

For a Gaussian draw, write $S_0=\mathrm{Ran}(A\mathrm{\Omega })$. At iteration $j$ of [F1], the intermediate range is $A^{\ast }{S}_{j-1}$ and the next range is $A{A}^{\ast }{S}_{j-1}$, since replacing a matrix by a basis preserves its range. Induction gives $S_q=\mathrm{Ran}((A{A}^{\ast }{)}^{q}A\mathrm{\Omega })$, including empty ranges and $q=0$. Thus the range specified directly in the statement agrees with the algorithm's output range for Gaussian draws. For a general measurable sketch, the statement specifies $P$ directly and requires no Gaussian sampling step.

Set $Z=\mathrm{\parallel }(I-P)(A{A}^{\ast }{)}^{q}A{\mathrm{\parallel }}_2$. By [F2], $\mathrm{\parallel }(I-P)A{\mathrm{\parallel }}_2\le Z^{1\mathrm{/}(2q+1)}$ pointwise. Since $Z$ is integrable and $Z^{1\mathrm{/}(2q+1)}\le 1+Z$, Jensen [F3] applied to the convex function $-t^{1\mathrm{/}(2q+1)}$ gives $\mathbb{E}{Z}^{1\mathrm{/}(2q+1)}\le (\mathbb{E}Z{)}^{1\mathrm{/}(2q+1)}$. Taking expectations of the pointwise inequality proves the result.

Write $G=\left[\begin{array}{c}g\\ H\end{array}\right]$ with first row $g$ and remaining $k-1$ rows $H$. If $k=1$, set $P=I_{k+p}$; then $G{G}^{\ast }=g{g}^{\ast }=gP{g}^{\ast }$ and the inverse-diagonal identity is immediate. If $k>1$, [F2] gives full row rank of $H$, so $P=I-H^{\ast }(H{H}^{\ast }{)}^{-1}H$ is the orthogonal projector onto $\ker H$, of dimension $(k+p)-(k-1)=p+1$. To obtain the inverse entry without an assumed block-inverse formula, solve $G{G}^{\ast }\left(\genfrac{}{}{0px}{}{x}{y}\right)=\left(\genfrac{}{}{0px}{}{1}{0}\right)$. The lower equation gives $y=-(H{H}^{\ast }{)}^{-1}H{g}^{\ast }x$, and substitution into the upper equation gives $(gP{g}^{\ast })x=1$. Full row rank of $G$ from [F2] makes $G{G}^{\ast }$ invertible and $gP{g}^{\ast }>0$, hence $((G{G}^{\ast }{)}^{-1}{)}_{11}=x=(gP{g}^{\ast }{)}^{-1}$.

For each fixed $H$, choose orthonormal coordinates in which $P$ is diagonal with $p+1$ ones. Independence of $g$ from $H$ and invariance [F7] show that $gP{g}^{\ast }$ is a sum of $p+1$ independent squared real Gaussians, or squared moduli of normalized complex Gaussians. Apply [F10] separately on the open half-lines via $x=\sqrt{u}$ and $x=-\sqrt{u}$ (absolute derivative $1\mathrm{/}(2\sqrt{u})$); the omitted point has zero measure. Adding the two densities gives the real square density $u^{-1\mathrm{/}2}{e}^{-u\mathrm{/}2}\mathrm{/}\sqrt{2\pi }=u^{-1\mathrm{/}2}{e}^{-u\mathrm{/}2}\mathrm{/}(\sqrt{2}\mathrm{\Gamma }(1\mathrm{/}2))$ by [F6]. For the complex square law, [F11] applied to the normalized density gives $1={\pi }^{-1}\sigma (S^1){\int }_0^{\mathrm{\infty }}{e}^{-r^2}rdr=\sigma (S^1)\mathrm{/}(2\pi )$, so $\sigma (S^1)=2\pi$. Applying the same formula with any nonnegative Borel function of $r^2$, and then [F10] with $u=r^2$, shows that the squared modulus has density $e^{-u}$ on $(0,\mathrm{\infty })$. The first real square is thus Gamma$(1\mathrm{/}2,2)$ and the complex square Gamma$(1,1)$, the latter normalized by [F12]. If independent Gamma variables have common scale $\theta$ and shapes $a,b$, [F5] and the diffeomorphism $(u,v)\mapsto (u,s=u+v)$ of the positive quadrant onto $0<u<s$, of determinant one in [F10], show that their sum has convolution density with normalization $1\mathrm{/}({\theta }^{a+b}\mathrm{\Gamma }(a)\mathrm{\Gamma }(b))$ and numerator $e^{-s\mathrm{/}\theta }{\int }_0^s{u}^{a-1}(s-u{)}^{b-1}du$. For each $s>0$, substituting $u=st$ by [F10] and using [F3] makes this $e^{-s\mathrm{/}\theta }{s}^{a+b-1}B(a,b)$, which is the Gamma density of shape $a+b$ and scale $\theta$. Iterating this finite convolution gives $gP{g}^{\ast }\sim \mathrm{\Gamma }((p+1)\mathrm{/}2,2)$ over $\mathbb{R}$ and $gP{g}^{\ast }\sim \mathrm{\Gamma }(p+1,1)$ over $\mathbb{C}$. Tonelli [F5] integrates first over $g$ and then over $H$, so the reciprocal law holds unconditionally for the first diagonal entry; row permutation invariance gives the same law for every index.

Integrating the Gamma density with the scaling $u=\theta t$ from [F10] yields $\mathbb{E}{U}^{-1}={\theta }^{-1}\mathrm{\Gamma }(a-1)\mathrm{/}\mathrm{\Gamma }(a)$ for $U\sim \mathrm{\Gamma }(a,\theta )$ and $a>1$. Using the Gamma recurrence [F12] gives $\mathrm{\Gamma }(a)=(a-1)\mathrm{\Gamma }(a-1)$. Thus step 2.1 gives $\mathbb{E}((G{G}^{\ast }{)}^{-1}{)}_{11}=1\mathrm{/}(p-1)$ over $\mathbb{R}$ and $1\mathrm{/}p$ over $\mathbb{C}$.

By step 2.1 all inverse diagonal entries have the same reciprocal-Gamma law and finite mean in the stated ranges, so $\mathbb{E}\mathrm{tr}(G{G}^{\ast }{)}^{-1}$ is $k$ times the value in step 3.1. Put $B=G^{\ast }(G{G}^{\ast }{)}^{-1}$. Then $GB=I_k$, $BG$ is self-adjoint and idempotent, so $GBG=G$, $BGB=B$, $(GB{)}^{\ast }=GB$, and $(BG{)}^{\ast }=BG$; these are the four Penrose equations [F8], hence uniqueness [F9] gives $B=G^{†}$. Also $B^{\ast }B=(G{G}^{\ast }{)}^{-1}$, and therefore $\mathrm{\parallel }{G}^{†}{\mathrm{\parallel }}_F^2=\mathrm{tr}(G{G}^{\ast }{)}^{-1}$ almost surely.