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.


1 Answer 1


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.

  • $\begingroup$ Note that $XSX-X=M$ is just a special case of the Riccati equation. Numerous solution methods have been proposed for that equation. $\endgroup$
    – lynn
    Commented Mar 12, 2018 at 22:55
  • $\begingroup$ @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$? $\endgroup$
    – A. Donda
    Commented Mar 13, 2018 at 19:54
  • 1
    $\begingroup$ @A.Donda Pre/post multiplying by $Y=X^{-1}$ results in $YMY + Y = S$, which is another Riccati equation. $\endgroup$
    – greg
    Commented Mar 14, 2018 at 3:56
  • $\begingroup$ 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? $\endgroup$
    – A. Donda
    Commented Mar 19, 2018 at 14:26
  • $\begingroup$ 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$. $\endgroup$
    – A. Donda
    Commented Mar 19, 2018 at 14:32

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .