# Orthogonal Projection onto the Unit Simplex with the Following Form - $x_i^* = (y_i + \nu^*)_+$

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)

Attempt:

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\}$?

• It seems you forgot the $x\succeq0$ constraint when forming the Lagrangian. – user856 Sep 19 '16 at 21:55
• 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. – user856 Sep 19 '16 at 21:56
• @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 – Fraïssé Sep 19 '16 at 21:56
• This is a solved problem (complexity $\le \mathcal O(n\log(n))$ or even $\mathcal O(n)$). See Duchi or Laurent Condat. – dohmatob Sep 20 '16 at 0:47

## 1 Answer

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.