# solve matrix equation involving Hadamard products

I'm trying to solve the equation $$(\Sigma \circ C)^{-1} \circ C = \left[ (\Sigma \circ C)^{-1} S (\Sigma \circ C)^{-1} \right] \circ C$$ for $\Sigma$ or $\Sigma \circ C$, where $\circ$ denotes the Hadamard or Schur element-wise product. $\Sigma$, $S$ and $C$ are symmetric real matrices. $\Sigma$ and $C$ are positive definite, while $S$ may be positive semidefinite. The elements of $C$ are $\in [0, 1]$, with 1s on the diagonal (which makes $\Sigma \circ C$ positive definite, too), and will in general be sparse. If it helps, we can assume that $C_{ij} = c(D_{ij})$, where $c()$ is a scalar decreasing function and $D$ is a Euclidean distance matrix.

If $C$ had no zero elements and $S$ were positive definite, the solution would be $$(\Sigma \circ C) = S.$$

Is it possible to find a solution in the general case? I don't believe there is a closed-form solution, but maybe an algorithm to compute it?

The problem arises from the 1-taper approximation to the Wishart likelihood introduced by Kaufman, Schervich, and Nychka (J. Am. Stat. Ass. 2008). Other than in that paper, I don't want to assume a structure for $\Sigma$ beyond being symmetric positive definite.

Let's solve for the variable $$X = (C\circ\Sigma)^{-1}$$ In terms of this variable, the kernel of your equation becomes $$XSX-X = 0$$ We are free to add a matrix $M$ to the RHS of the equation as long as it satisfies the constraint $$C\circ M=0$$

If we denote the complement of $C$ by $B=(1-C)$, then one way to ensure that $M$ satisfies the constraint is to write $$M=B\circ A$$ where $A$ is completely unconstrained.

Now all we need to to is to solve a quadratic matrix equation $$XSX - X = M$$ For that, you can use this method, making the substitutions \eqalign{ C &= S \cr {\mathcal H} &= \begin{bmatrix} 0 & S\\ M & I \end{bmatrix} \cr } Note that you are free to adjust $M$ until ${\mathcal H}$ has an acceptable set of eigenvalues/eigenvectors.

• Note that $XSX-X=M$ is just a special case of the Riccati equation. Numerous solution methods have been proposed for that equation.
– lynn
Commented Mar 12, 2018 at 22:55
• @greg, thank you very much! – The background of this question is that for the 1-taper approximation we get a likelihood that only contains the determinant and the inverse of $\Sigma \circ C$. If $C$ is chosen to be sparse, then this quantity can be easily stored even for large dimensionality, and the inverse be replaced by a sparse Cholesky decomposition, so everything is efficient computationally. If I estimate $(\Sigma \circ C)^{-1}$, I lose this advantage, because it will in general be dense. Do you see a way to get around that and directly estimate $\Sigma \circ C$? Commented Mar 13, 2018 at 19:54
• @A.Donda Pre/post multiplying by $Y=X^{-1}$ results in $YMY + Y = S$, which is another Riccati equation.
– greg
Commented Mar 14, 2018 at 3:56
• So, if I understand you correctly you are saying that I have to solve $\begin{bmatrix} 0 & M \\ S & I \end{bmatrix} \begin{bmatrix} I \\ Y \end{bmatrix} = \begin{bmatrix} I \\ Y \end{bmatrix} M Y$ for suitably chosen $M$, using the approach in the link? Commented Mar 19, 2018 at 14:26
• Problem: I could choose $M = 0$ (which satisfies the constraint) with the solution $Y = S$, but that is not a solution to my original problem, because $Y = \Sigma \circ C$ is sparse, while $S$ ist not. That means I can only use particular values of $M$, and the question becomes, how to find those? – More generally, I don't see how the step from my original equation to the Ricatti equation incorporates the constraint imposed by $\circ C$. Commented Mar 19, 2018 at 14:32