Sparsity-inducing regularization for stochastic matrices It is well-known (e.g. in the field of compressive sensing) that the $L_1$ norm is "sparsity-inducing," in the sense that if we minimize the functional (for fixed matrix $A$ and vector $\vec{b}$) $$f_{A,\vec{b}}(\vec{x})=\|A\vec{x}-\vec{b}\|_2^2+\lambda\|\vec{x}\|_1$$ for large enough $\lambda>0$, we're likely for many choices of $A$, $\vec{b}$, and $\lambda$ to have many exactly-zero entries in the resulting $\vec{x}$.
But if we minimize $f_{A,\vec{b}}$ subject to the condition that the entries of $\vec{x}$ are positive and sum to $1$, then the $L_1$ term doesn't have any effect (because $\|\vec{x}\|_1=1$ by fiat).  Is there an analogous $L_1$-type regularizer that works in this case to encourage that the resulting $\vec{x}$ is sparse?
 A: A general method for creating sparse solutions is via MAP estimation with a zero mean normal prior with an unknown variance.
$$p(x_i|\sigma_i^2)\sim N(0,\sigma_i^2)$$
If you then assign a prior to $\sigma_i^2$ which has a mode at zero then the posterior mode is usually sparse.  The $L_1$ arises from this approach by taking an exponential mixing distribution.
$$p(\sigma_i^2|\lambda)\sim Expo\left(\frac{\lambda^2}{2}\right)$$
Then you get
$$\log[p(x_i|\lambda)]=-\lambda | x_i|+\log\left[\frac{\lambda}{2}\right]$$
Some alternatives are the generalised double pareto, half cauchy, inverted beta.  In some sense  these are better than lasso because they do not shrink large values.  In fact I'm pretty sure the generalised double pareto can be written as a mixture of exponentials.  That is we write $\lambda=\lambda_i$ and then place a gamma prior  $p(\lambda_i|\alpha\beta)$.  We get:
$$p(x_i|\alpha\beta)=\frac{\alpha}{2\beta}\left(1+\frac{|x_i|}{\beta}\right)^{-(\alpha+1)}$$
Note that I have included normalising constants, as they help choose good global parameters.  Now if we apply the range restriction then we have a more complicated problem, as we need to renormalise over the simplex.
Another generic feature of sparsity inducing penalties is that they are not differentiable at zero.  Usually this is because the left and right limits are of opposite sign.
This is based on the brilliant work by Nicolas Polson and James Scott on variance mean mixture representations which they use to develop TIRLS - a massive extension of least squares to a very large class of loss-penalty combinations.
As an alternative you could use a prior which is defined on the simplex, but has modes in the marginal distributions at zero.  One example is the dirichlet distribution with all parameters between 0 and 1.  The implied penalty would look like:
$$-\sum_{i=1}^{n-1}(a_i-1)\log(x_i) - (a_n-1)\log(1-\sum_{i=1}^{n-1}x_i)$$
Where $0<a_i<1$.  However you would need to be careful in optimising numerically as the penalty has singularities.  A more robust estimation process is to use the posterior mean.  Although you lose exact sparseness you will get many posterior means that are close to zero.p
A: Two options:


*

*Use an $L_0$ penalty on $\vec x$.  The obvious drawback is that this is nonconvex and hence difficult to optimize.

*Reparameterize, $x_i = \frac{\exp(w_i)}{\sum_j \exp(w_j)}$ and use a penalty on the new (natural) parameter vector, $\|\vec w\|$.  This will encourage events to be equally probable unless there is a good reason for them not to be.

A: The premise of the question is only partly correct. While it is true that the $L_1$-norm is just a constant under the constraint, the constraint optimization problem might very well have a sparse solution. 
However, the solution is unaffected by the choice of $\lambda$, so either there is a sparse solution or not. Another question is how to actually find the solution. A standard quadratic optimizer under linear constraints can, of course, be used, but popular coordinate descent algorithms cannot be used out-of-the-box.
One suggestion could be to optimize under a positivity contraint only, for different $\lambda$'s, and then renormalize the solution to have $L_1$-norm 1. A coordinate descent algorithm should, I believe, be easily modifiable to compute the solution under a positivity constraint. 
A: I can think up three methods.


*

*Bayesian method: introducing a zero-mean prior distribution and use type II likelihood to estimate the parameters and hyper parameters.

*Use $\Vert\cdot\Vert_{\infty}$ as regularization instead. This is not differentiable though. You can use a high-order norm to approximate it.

*Use $-\sum_{i=1}\log x_i$.
In fact, the first and the third methods are the same.
