# estimate sparse localized whitening transformation

This is a follow-up to estimate precision matrix with given spatial sparsity pattern, expanding on the second part of that question and formulating more precisely using material from the answer by jlewk.

I have a matrix of measurements $$\Xi$$, where each of the $$n$$ rows comprises a measurement of $$p$$ variables. I assume that the variables are normally distributed with expectation zero and that the measurements are i.i.d., but the variables can be correlated and of different variance. It holds $$p \gg n$$, and $$p$$ can be on the order of several 100,000.

I would like to estimate a $$p \times p$$ whitening transformation $$W$$ such that the whitened data are distributed as $$\Xi W \sim \mathcal{MN} (0, I_n, I_p),$$ where $$\mathcal{MN}$$ denotes the matrix normal distribution. Maximum likelihood estimation leads to the function $$\DeclareMathOperator{\tr}{tr} \newcommand{\dt}[1]{\left\lvert#1\right\rvert} \tr(\Xi' \Xi W W') - n \log \dt{W W'}$$ to be minimized. Without further constraints, this leads to $$(W W')^{-1} = \frac1n \Xi' \Xi$$, but because $$p > n$$, $$\Xi' \Xi$$ is not invertible. Note that $$P = W W'$$ is the inverse covariance or precision.

To circumvent this estimation problem, I want to impose a sparsity constraint on $$W$$.

The different variables $$\xi_i$$ (columns of $$\Xi$$) are measured at different positions in space $$(x_i, y_i, z_i)$$, and I know that correlations between them are mainly due a local signal mixture, which would make the precision matrix sparse. Moreover, for technical reasons I would prefer each whitened variable (column of $$\Xi W$$) to derive only from a small set of original variables near to each other. This motivates the sparsity constraint $$W_{ij} = 0 \quad \text{for} \quad d_{ij} > r,$$ where $$d_{ij} = \sqrt{ (x_i - x_j)^2 + (y_i - y_j)^2 + (z_i - z_j)^2 }$$ is the distance between variable measurement positions and $$r$$ is an upper limit on this distance.

I'm aware that the index $$j$$ in $$W_{ij}$$ refers to a whitened variable which strictly speaking does not have a position, but the sparsity constraint with respect to $$W$$ still makes sense because the implied precision $$W W'$$ fulfills approximately the same sparsity constraint with a radius of $$2 r$$.

My question is: How can I efficiently calculate the $$W$$ that minimizes the function under the constraint?

I can of course use a generic optimization algorithm and enforce the constraint by using parameters only for the non-zero elements of $$W$$ – but I would hope that there is a better way.

jlewk answered the same question with respect to the precision $$P$$ by pointing to variants of the graphical lasso where a penalty is applied not uniformly to all elements but element-specific, $$\sum_{i \neq j} \Lambda_{ij} \dt{P_{ij}},$$ so that a given sparsity pattern can be approximately enforced by choosing $$\Lambda$$ appropriately; and to implementations of these variants in different software packages. However, in this new question penalties would have to be applied to the elements of $$W$$, and the minimization is complicated by $$P$$ being quadratic in $$W$$.

In detail:

– Has this problem, "sparse localized whitening estimation", been considered in the literature before?

– Is the mathematical form of the constrained optimization problem outlined above known, and are there standard approaches how to solve it?

– Are there specific implementations, or known implementation strategies, that make this estimation process particularly efficient in light of my large $$p$$?