Convex Optimization, Non-negativity constraints, Interior-Point or Projected Gradient?

Question

Assume I have the following convex optimization problem, with a convex objective function on conventional non-negativity constraints.

\begin{align} \min_{x \geq 0} \sum_{i=1}^{I} a_{i}x_{i} - f(\bf{x}), \end{align}

where $\bf{x}$ is a vector of dimension I, $f(\bf{x})$ is a concave function, and the dimension of $I$ is very large. The solution is sparse, meaning that in the optimal vector there is a large amount of zeros. This comes from the fact that for some entries $i$

\begin{align} a_{i} > f_{i}(\bf{x}^{*}) \end{align}

where $f_{i}(\bf{x})$ corresponds to $\partial f(\bf{x}^{*})/x_{i}$. Since the function is convex, I know that the KKT conditions are necessary and sufficient for global optimality, all local optima are also global, and the set of solutions is convex.

However, I have questions about how to solve for the optimum numerically,

1. If I use a Projected Gradient Method, I have seen on the internet that the projection step of the algorithm is just

$$ x^{k+1} = \max\left\{0,x^{k}-\alpha_{k}(a-\nabla f(x^{k}))\right\} $$

However, I have not found a proof or a paper that shows this. Can someone provide a proof or reference?

2. The other option is to use Interior-Point methods. However, what confuses me in this case is that when one uses barrier methods, how does it work if the solution is not interior, but in the boundary for some entries?

Answer 1

If a lot of variables are on their bounds at optimality, you may be well-served by an active-set method such as Sequential Quadratic Programming (SQP). Because your constraints are all linear (and all simple bound constraints at that), there is none of the drama SQP methods have to deal with when linearizing constraints.

A nullspace version of SQP, sometimes called Reduced Hessian, could be especially effective if the number of active constraints at optimality is high, i.e., variables on their bound of zero. Nullspace (a.k.a. Reduced Space a.k.a. Reduced Hessian) methods form and update the Hessian of free variables, i.e., variables not on bound (given your problem only has bound constraints). So if your problem has 100,000 variables, of which 99,000 are at zero (on their bound), the Nullspace method works with a 1000 by 1000 Hessian, and solves Quadratic Programming (QP) subproblems of dimension 1000. As the algorithm proceeds, variables might go on or off their bound. At each iteration, only the non-zero variables are "in" the Hessian and corresponding QP subproblem whose "quadratic" matrix is the Hessian in that reduced space of non-zero variables.

SNOPT is a well-known solver implementing a Nullsapce (Reduced Hessian) SQP algorithm, and can be quite effective on very large problems, if the number of free variables (non-zero variables in your case) is (relatively) small. Generally speaking, the better your starting values (initial guess) is of which variables will be non-zero at optimality, the faster the algorithm will converge to the optimum solution.