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