Variation of Least Squares with Symmetric Positive Semi Definite (PSD) Constraint

I am trying to solve the following convex optimization problem: \begin{align} & \min_{W} && \sum_{i=1}^n (\mathbf{x}_{i}^TW\mathbf{x}_{i} - y_i)^2 \\\\ & s.t. && W \succcurlyeq 0 \\\\ & && W = W^T \end{align}

where $$\mathbf{x}_i \in \mathbb{R^p}$$, $$W \in \mathbb{R}^{p \times p}$$ and $$y_i \geq 0$$.

Without the positive semidefinite constraint, the problem is pretty straightforward. Requiring positive semidefiniteness, however, makes it a bit tricky.

I thought about using the fact that $$W \succcurlyeq 0$$ if and only if there exists a symmetric $$A$$ such that $$W = AA^T$$, and solving the equivalent problem

\begin{align} & \min_{A} && \sum_{i=1}^n (\mathbf{x}_{i}^TAA^T\mathbf{x}_{i} - y_i)^2 \\\\ &s.t. && A = A^T \end{align}

Letting $$a_{ij}$$ be the $$(i,j)th$$ element of A, this optimization function is quartic (fourth-order) with respect to the $$a_{ij}$$'s. Because of this, I am unsure of how to proceed.

I would be grateful if someone could point me in the right direction as to how to solve this problem.

• Haven't worked it through, but you might find it helpful to note that, if $W = A^\top A$, then$$x^\top_i W x_i = x^\top_i A^\top A x_i = (Ax_i)^\top (Ax_i) = \|Ax_i\|^2.$$ Commented Apr 11, 2020 at 0:01

When dealing with such form of matrix multiplications always remember the Vectorization Trick with Kronecker Product for Matrix Equations:

$${x}_{i}^{T} W {x}_{i} - {y}_{i} \Rightarrow \left({x}_{i}^{T} \otimes {x}_{i}^{T} \right) \operatorname{Vec} \left( W \right) - \operatorname{Vec} \left( {y}_{i} \right) = \left({x}_{i}^{T} \otimes {x}_{i}^{T} \right) \operatorname{Vec} \left( W \right) - {y}_{i}$$

Since the problem is given by summing over $${x}_{i}$$ one could build the matrix:

$$X = \begin{bmatrix} {x}_{1}^{T} \otimes {x}_{1}^{T} \\ {x}_{2}^{T} \otimes {x}_{2}^{T} \\ \vdots \\ {x}_{n}^{T} \otimes {x}_{n}^{T} \end{bmatrix}$$

Then:

$$\arg \min_{W} \sum_{i = 1}^{n} {\left( {x}_{i}^{T} W {x}_{i} - {y}_{i} \right)}^{2} = \arg \min_{W} {\left\| X \operatorname{Vec} \left( W \right) - \boldsymbol{y} \right\|}_{2}^{2}$$

Where $$\boldsymbol{y}$$ is the column vector composed by $${y}_{i}$$.

Now the above has nice form of regular Least Squares. The handling of the constraint can be done using Projected Gradient Descent Method. The projection onto the set of Symmetric Matrices and Positive Semi Definite (PSD) Matrices cone are given by:

1. $$\operatorname{Proj}_{\mathcal{S}^{n}} \left( A \right) = \frac{1}{2} \left( A + {A}^{T} \right)$$. See Orthogonal Projection of a Matrix onto the Set of Symmetric Matrices.
2. $$\operatorname{Proj}_{\mathcal{S}_{+}^{n}} \left( A \right) = Q {\Lambda}_{+} {Q}^{T}$$ where $$A = Q \Lambda {Q}^{T}$$ is the eigen decomposition of $$A$$ and $${\Lambda}_{+}$$ means we zero any negative values in $$\Lambda$$. See Find the Matrix Projection of a Symmetric Matrix onto the set of Symmetric Positive Semi Definite (PSD) Matrices.

While the Symmetric Matrices Set is a Linear Sub Space the PSD Cone is not a Linear Sub Space. Hence the greedy iterative projection on the set is not guaranteed to converge to the orthogonal projection on the intersection of the 2 sets. See Orthogonal Projection onto the Intersection of Convex Sets. Yet it will converge to a feasible solution.

With all the tools above one could create his own solver using basic tools with no need for external libraries (Which might be slow or not scale).

I implemented the Projected Gradient Descent Method with the above projections in MATLAB. I compared results to CVX to validate the solution. This is the solution:

My implementation is vanilla Gradient Descent with constant Step Size and no acceleration. If you add those you'll see convergence which is order of magnitude faster (I guess few tens of iterations). Not bad for hand made solver.

The MATLAB Code is accessible in my StackExchange Mathematics Q3619669 GitHub Repository.

Remark: The reason the greedy algorithm works is objective function is quadratic which means only the symmetric part of $$\boldsymbol{W}$$ has any effect. On top of that, if the initial matrix is symmetric then each step of the gradient step also generates a symmetric matrix, hence basically the only effective constraint is the definiteness of the matrix.

• This was a very clear and practical response. Thanks! Assuming I understand correctly, I have one question: what is the difference between iteratively projecting vs. projecting just once after finding the optimal solution to the unconstrained problem? Would they give different solutions? Commented Apr 13, 2020 at 15:43
• @tygaking, Please mark it as answer and +1. Regarding your question, you can't do that. Think of "Gradient Descent" as some kind of projection by itself. Hence what we do above is basically iterative projection onto iterative set, one for the gradient descent, one for each other.
– Royi
Commented Apr 13, 2020 at 15:45

This is a convex optimization problem which can be easily formulated, and then numerically solved via a convex optimization tool, such as CVX, YALMIP, CVXPY, CVXR. This is a linear Semidefinite Programming problem (SDP), for which numerical solvers exist.

Here is the code for CVX.

Assume $$x_i$$ is the ith column of matrix $$X$$

cvx_begin
variable W(p,p) semidefinite % constrains W to be symmetric positive semidefinite
Objective = 0;
for i=1:n
Objective = Objective + square(X(:,i)'*W*X(:,i) - y(i))
end
minimize(Objective)
cvx_end


CVX will transform the problem into the form required by the solver, call the solver, and return the solution.

Are you doing FGLS or something?

You could try substituting the constraint into the object. For the two by two case, for example, solve $$\sum_{i = 1}^N \left(x_i'\left[\array{ \array{w_{11} & w_{12}} \\ \array{w_{12}& w_{22} }} \right]x_i - y_i \right)^2$$ where $$w_{12} = w_{21}$$ now by construction. Then the matrix will be symmetric.

To ensure positive semi-definiteness, you can then use the standard principal minors test: $$w_{ii} \ge 0$$ for each $$i$$, $$w_{11} w_{22} - w_{12}^2 \ge 0$$, and so on, with the determinant of the upper left-hand minor weakly positive.

That at least subsumes positive semi-definiteness into concrete constraints and adjustments to the problem. That sounds like a nightmare to solve however, using Kuhn-Tucker. A simpler sufficient condition for semi-definiteness is the dominant diagonal condition, that $$w_{ii} \ge \sum_{j \neq i} |w_{ij}|$$ for each row $$i$$, which would be much more computationally tractable. Perhaps it could give you a good initial guess before you try relaxing it to the standard principal minors constraints.

• There is no need for this. While your approach might lead to problem formulated as Least Squares with Linear Equality and Inequality Constraints, it will also require some kind of iterative solver (As there is no closed form solution for LS with Inequality Constraints). Since we have the projection onto the Symmetric and PSD Matrices we can use it. See my answer. Still, nice idea of yours! Especially with the Diagonally Dominant Matrix. I might post a question and answer about this.
– Royi
Commented Apr 12, 2020 at 0:00

This does seem to be a straight forward 'nonnegative' least squares though with the sdp constraint. n = 10; p = 5; X = zeros(n,p^2); for ii = 1:n x = randn(p,1); temp = xx'; X(ii,:) = temp(:)'; end y = randn(n,1); cvx_begin sdp variable W(p,p) semidefinite minimize(norm(XW(:)-y)) cvx_end