2
$\begingroup$

Given $C$ is a $n \times n$ positive semidefinite matrix, how can one solve for $n \times 1$ vector $x$ to maximize the sum of absolute value of all components of $x$, given $x^TCx\leqq 1$?

I was looking for analytical solution without much success.One upper bound I found is $\sqrt{n/\lambda_{min}(C)}$, where $\lambda_{min}(C)$ is the minimum eigenvalue of $C$, is this the best upper bound? Even for numerical solution, it seems this cannot be done by convex optimization...

$$\max\|x\|_1$$ for $$x^TCx\leqq 1$$

$\endgroup$

migrated from stats.stackexchange.com Feb 28 '18 at 22:29

This question came from our site for people interested in statistics, machine learning, data analysis, data mining, and data visualization.

  • $\begingroup$ Cannot you go to the dual problem? $\endgroup$ – Hua Ying Feb 28 '18 at 22:44
  • $\begingroup$ If $C$ has rank deficit, then the feasible set is unbounded and thus the maximum doesn't exist. Even if $C$ has full rank. There is no analytical solution. Is an upper bound sufficient? $\endgroup$ – user251257 Feb 28 '18 at 23:44
  • $\begingroup$ I found one upper bound as in the updated question, but not sure if it can be tighter. $\endgroup$ – ahala Mar 1 '18 at 4:49
2
$\begingroup$

As nominally stated, this is a concave programming problem, i.e., the minimization of a concave function ($-\|x\|_1$) subject to convex constraints.

However, this can be solved with a mixed-integer convex optimizer by introducing binary variables to handle the non-convex objective. The optimization modeling system YALMIP under MATLAB can do this for you, if you have a suitable solver installed, such as gurobi, cplex, mosek, or scip (i.e., having mixed-integer convex quadratically-constrained or second-order cone problem capability).

    sdpvar x(n,1) % declare optimization variable
    optimize(x'*C*x <= 1,-norm(x,1)) 
    % 1st argument of optimize is the constraints, and 2nd argument is 
    % objective function to be minimized
    disp(value(x)) % displays the optimal value of x
    disp(nomr(x,1)) $ displays optimal (maximum) value of norm(x,1)

If $C$ is positive definite, the constraint set is compact, and therefore, there will be a globally optimal solution at an extreme of the constraint set, which in this case means $x'Cx = 1$.

$\endgroup$
  • $\begingroup$ Thanks very much. Can you point me to some references or details on implementation? Sdpvar hides all the details. $\endgroup$ – ahala Mar 1 '18 at 4:58
  • $\begingroup$ yalmip.github.io YALMIP has a lot of complicated and messy things which it does behind the scenes. That's what makes it as easy to use as it is. $\endgroup$ – Mark L. Stone Mar 1 '18 at 13:28
  • $\begingroup$ @Johan Lofberg 's (developer of YALMIP) answer is much better than mine. $\endgroup$ – Mark L. Stone Mar 1 '18 at 21:30
2
$\begingroup$

I think these arguments can be cleaned up, but here comes a variant to compute an alternative bound.

First, I perform a coordinate change by $C = R^TR$ and $z = Rx$,

$$\max\|R^{-1}z\|_1 \text{ subject to } z^Tz = 1$$

Use that $\|q\|_1 = \max_{\|d\|_{\infty\leq 1}} d^Tq$ and the problem can be written as

$$\max_{z^Tz = 1}\max_{\|d\|_{\infty\leq 1}} d^TR^{-1}z$$

Switch order

$$\max_{\|d\|_{\infty\leq 1}} \max_{z^Tz = 1} d^TR^{-1}z$$

and optimize over $z$ explicitly

$$\max_{\|d\|_{\infty\leq 1}} \|R^{-T}d\|_2$$

This is not much better, as it is a maximization of a convex quadratic. However, we can compute an upper bound here. Let $S = R^{-T}$. A trivial upper bound is given by $\| |S|\textbf{1} \|_2$. (i.e., the norm of the sum of the absolute value of the columns). For small dimensions, this bound is most often better than the bound in the question, and it is tight surprisingly (for me) often.

Here is the YALMIP code to play around (of course, to compute the bound no solves are needed, but they give us the optimal solution which allows us to judge the bounds). Some typical optimal values and bounds are reported also

Bounds = [];
for i = 1:10
    n = 5;
    C = randn(n);
    C = C*C';
    R = chol(C);
    x = sdpvar(n,1);
    % Brute-force MISOCP
    optimize(x'*C*x <= 1,-norm(x,1))
    norm(x,1)

    % Just check that the reformulation is correct
    d = sdpvar(n,1);
    z = sdpvar(n,1);
    optimize(z'*z <= 1,-norm(inv(R)*z,1))
    norm(inv(R)*z,1)

    % Dualized problem
    d = sdpvar(n,1);
    S = inv(R)';
    optimize(-1 <= d <= 1, -d'*S'*S*d,sdpsettings('solver','bmibnb'))
    sqrt(d'*S'*S*d)

    % Upper bounds
    Bounds = [Bounds;value(norm(x,1)) norm(sum(abs(S),2)) sqrt(n/min(eig(C)))];
end
Bounds

   11.2730   11.3251   13.6198
   16.1563   16.1611   18.8931
    3.2570    3.6245    4.2923
    4.9990    5.0386    5.5765
    5.2158    5.2158    8.2443
   27.9575   27.9980   30.5903
    5.7549    5.8292    7.7037
    3.8030    4.1501    4.4554
   11.6382   11.6382   15.1235
    9.9623   10.0193   11.8171
$\endgroup$
  • $\begingroup$ I went this line as in related question math.stackexchange.com/questions/2469097/…. One thing puzzled m is that the decomposition of $C = R'R$ is not unique. $\endgroup$ – ahala Mar 1 '18 at 21:35
  • 1
    $\begingroup$ Well, not being unique doesn't really change anything, it just gives you another lever to compute several bounds. A quick check on deriving $R$ from Cholesky vs. SVD showed that the Cholesky-derived version gave a better bound, but testing a symmetric factorization typically gave even better bounds. $\endgroup$ – Johan Löfberg Mar 1 '18 at 21:45
  • $\begingroup$ Also note that in the question you reference, you talk and ask about an induced norm, but then pull formulas from a Wiki-page for element-wise norms. Those are not the same thing. $\endgroup$ – Johan Löfberg Mar 1 '18 at 22:00

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.