A solution is given here, but it has a slightly different formulation and uses so called bisection method which I am not familiar with. I wish to use a direct approach to solve this question but I am a bit stuck

There are several other similar questions such as Projecting a nonnegative vector onto the simplex but none even has a sketch to how to solve the problem.

More strangely, despite lack of solution or approach, the subsequent questions were marked as duplicates Find the nearest point of $y$ from a probability simplex.

Problem. Find $x^*$ that solves

$$ \begin{aligned} & \underset{x \in \Delta^n}{\text{minimize}} & & \|y - x\|^2_2 \\ & \text{subject to} & & \mathbf{1}^Tx = 1, x \succeq 0 \end{aligned}$$

where $\Delta^n =\{ (x_1, \ldots, x_n) ~|~ \sum_{i=1}^n x_i=1, x_i \geq 0 \mbox{ for all } i \}$ and $y \in \mathbb{R}^n$

Claim: $x_i^* = (y_i + \nu^*)_+ = \max\{y_i + \nu^*, 0\}$ (perhaps with a minor difference in sign or multiplication by a constant)


We form the lagrangian:

$L(x,\nu) = \|y-x\|^2_2 + \nu(1^Tx - 1)$

And $\nabla L(x,\nu) = -2|y-x|+\nu 1$

(abuse: $1 \in \mathbb{R}^n, \nu \in \mathbb{R}$)

Then at optimality, we have:

$\nabla L(x^*,\nu^*) = -2|y-x^*| +\nu^* 1 = 0 \implies -2|y_i - x^*_i| + \nu^* = 0$

By definition of absolute value, $$|y_i - x^*_i| = (y_i - x^*_i), y_i \geq x^*_i \text{ or } -(y_i - x^*_i), y_i < x^*_i$$

  1. For $y_i - x^*_i \geq 0$, then:

$\implies x_i^* =- \dfrac{1}{2}\nu^* + y_i = y_i - \dfrac{1}{2}\nu^*, y_i \geq x^*_i$

  1. For $y_i - x^*_i < 0 $, then:

$\implies x_i^* =\dfrac{1}{2}\nu^* + y_i = y_i + \dfrac{1}{2}\nu^*, y_i \geq x^*_i$

At the end, it doesn't seem that I have the proper final solution, due to lack of $\max$ function

How do I introduce $\max$ so that I have $x_i^* = \max\{y_i - \dfrac{1}{2}\nu^*,0\}$?

  • $\begingroup$ It seems you forgot the $x\succeq0$ constraint when forming the Lagrangian. $\endgroup$ – user856 Sep 19 '16 at 21:55
  • 1
    $\begingroup$ Also, it's not fair to say that previous questions have no solution. The accepted answer to the first question links to an article where you can read about the solution. $\endgroup$ – user856 Sep 19 '16 at 21:56
  • $\begingroup$ @Rahul But does solution necessarily involve proximal operators and moreau decomposition? I wish to do mine a bit similar what Boyd did in the first link $\endgroup$ – Fraïssé Sep 19 '16 at 21:56
  • $\begingroup$ This is a solved problem (complexity $\le \mathcal O(n\log(n))$ or even $\mathcal O(n)$). See Duchi or Laurent Condat. $\endgroup$ – dohmatob Sep 20 '16 at 0:47

The projection onto the Simplex can be calculated in the following method.

The Lagrangian in that case is given by:

$$ \begin{align} L \left( x, \mu \right) & = \frac{1}{2} {\left\| x - y \right\|}^{2} + \mu \left( \boldsymbol{1}^{T} x - 1 \right) && \text{} \\ \end{align} $$

The trick is to leave non negativity constrain implicit.
Hence the Dual Function is given by:

$$ \begin{align} g \left( \mu \right) & = \inf_{x \succeq 0} L \left( x, \mu \right) && \text{} \\ & = \inf_{x \succeq 0} \sum_{i = 1}^{n} \left( \frac{1}{2} { \left( {x}_{i} - {y}_{i} \right) }^{2} + \mu {x}_{i} \right) - \mu && \text{Component wise form} \end{align} $$

Again, taking advantage of the Component Wise form the solution is given:

$$ \begin{align} {x}_{i}^{\ast} = { \left( {y}_{i} - \mu \right) }_{+} \end{align} $$

Where the solution includes the non negativity constrain by Projecting onto $ {\mathbb{R}}_{+} $

Again, the solution is given by finding the $ \mu $ which holds the constrain (Pay attention, since the above was equality constrain, $ \mu $ can have any value and it is not limited to non negativity as $ \lambda $ above).

The objective function (From the KKT) is given by:

$$ \begin{align} h \left( \mu \right) = \sum_{i = 1}^{n} {x}_{i}^{\ast} - 1 & = \sum_{i = 1}^{n} { \left( {y}_{i} - \mu \right) }_{+} - 1 \end{align} $$

The above is a Piece Wise linear function of $ \mu $ and its Derivative given by:

$$ \begin{align} \frac{\mathrm{d} }{\mathrm{d} \mu} h \left( \mu \right) & = \frac{\mathrm{d} }{\mathrm{d} \mu} \sum_{i = 1}^{n} { \left( {y}_{i} - \mu \right) }_{+} \\ & = \sum_{i = 1}^{n} -{ \mathbf{1} }_{\left\{ {y}_{i} - \mu > 0 \right\}} \end{align} $$

Hence it can be solved using Newton Iteration.

I wrote MATLAB code which implements them both at Mathematics StackExchange Question 2338491 - GitHub.
There is a test which compares the result to a reference calculated by CVX.


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.