Concentration of Mode-Wise Covariance Estimates

Why a Mode-Wise Covariance Scatter Concentrates

In a tensor-normal model, one mode covariance describes variation along one tensor index while the other mode covariances describe dependence across the remaining indices. If those other covariances are known, they can be removed by whitening. The resulting mode-wise scatter is an ordinary Gaussian covariance scatter, but its number of effective replicates is not merely the number of observed tensors: each tensor contributes one replicate for every combination of indices in the other modes.

This distinction matters for uncertainty calculations. A scatter matrix with \(p_r\) rows is subject to a dimension penalty, while the effective replication count controls how quickly its error decreases. The result below makes both contributions explicit.

The factor \(1/q\) makes this an estimate of \(\Sigma_r\), rather than \(q\Sigma_r\). Since the other factors are known, their scale is fixed in this calculation. If all factors were unknown, the familiar scale nonidentifiability of a separable covariance would require a separate convention.

Whitening Reveals the Effective Replication Count

Let \(B_r=\Sigma_{-r}^{-1/2}\), using the symmetric positive-definite square root. Under the separable Gaussian model, the columns of \(Z_i=X_{i(r)}B_r\) are independent \(N_p(0,\Sigma_r)\) vectors. This follows because the column covariance after transformation is \(B_r\Sigma_{-r}B_r^T=I_q\); the transformed matricization is jointly Gaussian, so zero covariance between columns implies independence. Observations from different tensors are independent as well. Consequently, the collection of all columns of \(Z_1,\ldots,Z_n\) consists of \(N=nq\) independent Gaussian vectors, and \[ \widehat{\Sigma}_r=\frac{1}{N}\sum_{i=1}^n Z_iZ_i^T. \]

Thus each tensor contributes \(q\) effective replicates for estimating this one mode covariance. The original columns of \(X_{i(r)}\) need not be independent: whitening by the known \(\Sigma_{-r}\) is what produces independent columns.

The Scatter Is the Gaussian Maximum-Likelihood Estimate

The whitening also gives the likelihood for \(\Sigma_r\) directly. Write the whitened columns collectively as \(z_1,\ldots,z_N\). For positive-definite \(\Sigma_r\), their log-likelihood, up to terms not involving \(\Sigma_r\), is \[ \ell(\Sigma_r) =-\frac{N}{2}\log\det(\Sigma_r) -\frac{1}{2}\sum_{j=1}^{N}z_j^T\Sigma_r^{-1}z_j. \] Using \(\sum_j z_jz_j^T=N\widehat{\Sigma}_r\), the second term is \(-\frac{N}{2}\operatorname{tr}(\Sigma_r^{-1}\widehat{\Sigma}_r)\).

To verify the maximization, set \(A=\widehat{\Sigma}_r\), assumed positive definite, and \(C=A^{-1/2}\Sigma_r A^{-1/2}\). Terms depending on \(\Sigma_r\) become, apart from constants, \[ -\frac{N}{2}\{\log\det C+\operatorname{tr}(C^{-1})\}. \] If the eigenvalues of \(C\) are \(\lambda_1,\ldots,\lambda_p>0\), maximizing this expression is equivalent to minimizing \(\sum_{j=1}^p(\log\lambda_j+1/\lambda_j)\). For \(g(\lambda)=\log\lambda+1/\lambda\), the derivative is \(g'(\lambda)=(\lambda-1)/\lambda^2\); hence \(g\) decreases up to 1 and increases after 1. Its unique minimum is at \(\lambda=1\). The likelihood is therefore maximized at \(C=I\), or \(\Sigma_r=A\).

When \(N\ge p\), the \(N\) Gaussian vectors span \(\mathbb{R}^p\) with probability one: each successive vector has probability zero of lying in any fixed proper subspace. Thus their sum of outer products is positive definite. When \(N<p\), that sum has rank at most \(N<p\), so \(A\) is singular. In that case, the likelihood can increase without bound by shrinking covariance in a direction orthogonal to the observed vectors, which rules out a positive-definite maximizer.

A Finite-Sample Operator-Norm Bound

A relative error bound is natural because the scale and correlation of the target covariance may vary. Define \[ A_r=\Sigma_r^{-1/2}\widehat{\Sigma}_r\Sigma_r^{-1/2}. \] After whitening both the known modes and the target covariance, \(A_r\) is the sample covariance of \(N\) independent \(N_p(0,I_p)\) vectors.

Proof. Let \(w_1,\ldots,w_N\) be the independent \(N_p(0,I_p)\) vectors whose sample covariance is \(A_r\). For any fixed unit vector \(u\), the variables \(u^Tw_j\) are independent standard normal variables. Therefore \[ u^T(A_r-I_p)u =\frac{1}{N}\sum_{j=1}^N\{(u^Tw_j)^2-1\} \] has the distribution of \(\chi_N^2/N-1\).

The chi-square Chernoff bound gives, for \(x>0\), \[ P\left(\left|\frac{\chi_N^2}{N}-1\right|>x\right) \le 2\exp\left\{-\frac{N}{8}\min(x^2,x)\right\}. \] For completeness, the moment-generating function is \(E[e^{s\chi_N^2}]=(1-2s)^{-N/2}\) for \(s<1/2\). Applying Markov's inequality and optimizing \(s\) gives upper-tail exponent \(-N\{x-\log(1+x)\}/2\), and, for \(0<x<1\), lower-tail exponent \(-N\{-x-\log(1-x)\}/2\). The inequalities \(x-\log(1+x)\ge \min(x^2,x)/4\) and \(-x-\log(1-x)\ge x^2/2\) give the stated, slightly weaker bound; for \(x\ge1\) the lower-tail event is impossible.

Choose a \(1/4\)-net \(\mathcal{N}\) of the unit sphere in \(\mathbb{R}^p\), meaning every unit vector is within Euclidean distance \(1/4\) of some net point. A volume comparison gives a net with at most \((1+2/(1/4))^p=9^p\) points. For a symmetric matrix \(M\), \[ \|M\|_{\mathrm{op}}\le 2\max_{u\in\mathcal{N}}|u^TMu|. \] Indeed, choose a unit vector \(x\) attaining \(|x^TMx|=\|M\|_{\mathrm{op}}\), and a net point \(u\) with \(\|x-u\|_2\le1/4\). Then \[ |x^TMx-u^TMu| \le |(x-u)^TMx|+|u^TM(x-u)| \le \tfrac12\|M\|_{\mathrm{op}}, \] which proves the inequality.

Set \(a=(p\log 9+t)/N\) and \(x=\sqrt{8a}+8a\). Then \(x^2/8\ge a\), and \(x/8\ge a\); hence \(\min(x^2,x)/8\ge a\). Apply the chi-square bound to each of the at most \(9^p\) net points and use the union bound. The probability that any net quadratic form exceeds \(x\) in absolute value is at most \[ 2\cdot9^p e^{-Na}=2e^{-t}. \] Outside this event, the net inequality gives \(\|A_r-I_p\|_{\mathrm{op}}\le2x\), exactly the claimed bound. \(\square\)

Interpretation and a Common Pitfall

The leading term in the bound is of order \(\sqrt{p/N}\), and the second term is of order \(p/N\). Thus a small relative operator-norm error is expected when effective replication is large compared with the target dimension. The guarantee is for one specified mode and one specified covariance factor; applying it to multiple modes at once requires an additional probability argument, such as a union bound.

• Pitfall: counting only the number of tensors. Treating the sample size as \(n\) ignores the \(q\) independent whitened vectors contributed by each tensor. The covariance scatter uses \(N=nq\) observations after whitening.
• Pitfall: treating raw columns as independent. Columns of the unwhitened matricization are generally correlated when the other mode covariances are not identity matrices. Independence applies after the known column covariance has been removed.
• Pitfall: reading the bound as an absolute-error guarantee. The theorem bounds \(\Sigma_r^{-1/2}\widehat{\Sigma}_r\Sigma_r^{-1/2}-I_p\), so it controls relative error in covariance-scaled directions, not the unscaled operator norm \(\|\widehat{\Sigma}_r-\Sigma_r\|_{\mathrm{op}}\) by the same number.

The key statistical benefit of knowing the other factors is that it turns structured dependence across tensor indices into additional usable Gaussian replicates. The target dimension remains a cost, but the relevant sample size for this mode-wise scatter is the total number of whitened mode fibers, not just the number of tensor observations.