Let $\mathcal P_n = \{P \in \mathbb R^n_{\geq 0}: P^T \mathbb I = 1 \}$, where $\mathbb I = (1,...,1)^T \in \mathbb R^n$ and $f: \mathcal P_n \to \mathbb R$ a convex and differentiable function (or it can be as smooth as you like). Given $Q \in \mathcal P_n$ and the parameter $c > 0$, I'm interested in the following problem: \begin{align*} \min_{P \in \mathcal P_n} \;\; & f(P) \\ \text{s.t.} \;\; & W_p(P,Q) \leq c \end{align*}

Here I consider the following discrete transport: $$W_p(P,Q) = \min_{\pi \in \Pi(P,Q)} \langle D,\pi \rangle$$ where:

  • $\Pi(P,Q) = \{ \pi \in \mathbb R^{n \times n}_{\geq 0}: \pi \mathbb I = P, \pi^T \mathbb I = Q \}$.
  • $D \in \mathbb R^{n \times n}_{\geq 0}$ with $D_{ij} = |i-j|^p$, for $p \in [1, \infty]$.

I have $3$ questions:

  • Are there any efficient numerical methods to solve this problem?

  • Should it be easier/more difficult (numerically) for the different values of $p$?

  • If, instead, we use the discrete regularised transport: $$W_p^{\epsilon}(P,Q) = \min_{\pi \in \Pi(P,Q)} \langle D,\pi \rangle - \epsilon H(\pi)$$ where $H(\pi) = - \sum_{i,j} \pi_{ij} (\ln(\pi_{ij} - 1))$ and $\epsilon > 0$. Regarding this regularised optimisation problem,

    • Is it easier/more difficult (numerically) than the original problem?

    • If $\epsilon \to 0$, does its optimal solution converge to that of the original problem?

  • $\begingroup$ I'm not saying it is efficient, but what one can always do is to optimize directly over $\pi \in \Pi(P, Q)$ satisfying $\langle D,\pi \rangle \leq c$. Since $\pi \mapsto P$ is linear, you obtain an objective function which is as nice as previously, and with a linear inequality constraint. Standard solvers will be able to solve it given $n$ is small enough. Further, there is of course a bunch of literature for when $f$ is linear, i.e. $f(P) = \int g \,dP$ in the Wasserstein "distributionally robust optimization" literature. $\endgroup$ – Steve Apr 11 '19 at 17:31
  • $\begingroup$ @Steve: thanks for your comment. I will try your suggestion and compare with the methods in the answer below. $\endgroup$ – SiXUlm Apr 11 '19 at 20:19
  • $\begingroup$ Consider the particular case where $f$ is linear, i.e $f(P) = \mathbb E_P[Z] = \sum_{i=1}^np_iz_i$, for some $z_1,\ldots,z_n \in \mathbb R$. Then $$\min_{P \in \Delta_n \mid W_1(P,Q) \le c}f(P) = \min_{P \in \Delta_n\mid W_1(P,Q) \le c}E_P[Z] = \max_{ \lambda > 0}-\lambda c + \sum_{i=1}^n p_i f_i^\lambda, $$ where $f_i^\lambda := \inf_{z_i'} \lambda d(z_i', z_i)-z_i'$. $\endgroup$ – dohmatob Apr 12 '19 at 12:39

Too long to comment.

Note that for a given P and Q, $W_p(P,Q) = \min_{\pi \in \Pi(P,Q)} \langle D,\pi \rangle \leq c$ if and only if there exists a $\pi \in \Pi(P,Q)$ such that $\mbox{Tr}(D^\top \pi) = \langle D, \pi \rangle \leq c$. With that preamble, and the assumption that $f$ depends explicitly only on $P$, consider the following optimization problem: $$ \begin{align*} \min_{P \in R_+^{n}, ~~\pi\in R_+^{n\times n}} \;\; & f(P) \\ \text{s.t.} \;\; & \mbox{Tr}(D^\top \pi) \leq c\\ & \pi \mathbb I = P, ~\pi^T \mathbb I = Q. \end{align*} $$ Note that a solution $(P^*,\pi^*)$ to the above convex program (LP is one considers a linear cost functional), will be a feasible point for the original optimization problem. Similarly, it is also not hard to argue that a solution to the original optimization problem would be a feasible point for the above convex program.

The two points together imply that the optima of the original problem can be achieved in polynomial time.

Hope this suggestion helps.

| cite | improve this answer | |
  • $\begingroup$ It is really simple as that. Also, the case of regularization is the same (just adapt the new $\leq c$ constraint). $\endgroup$ – Dirk May 6 at 10:35

Meanwhile I also asked an expert of the field about this problem and here is his answer (in French).

quelques commentaires sur votre problème:

-si au lieu d une contrainte sur la distance W(P,Q) vous prenez un multiplicateur \lambda et minimisez $f(P)+ \lambda W(P,Q) c est un problème très etudié notamment dans la théorie des flots de gradient (schéma JKO) et il y a pas mal de méthodes numériques pour le résoudre (Lagrangian augmenté cf mon papier avec Benamou et Laborde: https://hal.archives-ouvertes.fr/hal-01245184 , régularisation entropique voir l article de Gabriel Peyré: https://arxiv.org/abs/1502.06216, ou encore en utilisant le transport semi-discret: https://hal.archives-ouvertes.fr/hal-01056452)

-si f est convexe votre problème est convexe ce qui est important

Pour completér ma réponse (positive) a votre question 1, je vous invite a lire le livre (en ligne ici: https://arxiv.org/abs/1803.00567) de Cuturi et Peyré

Pour répondre a la question 2, oui ca va dépendre un peu de p, p=2 est le cas ou il y a le plus de méthodes disponibles et p=1 est aussi un peu spécial car W_1 est vraiment une semi-norme, il y a aussi des méthodes spécifiques (on avait regardé ca avec Benamou et Hatchi et Osher et ses coauteurs ont aussi développé des méthodes spécifiques), pour p>1 la méthode de Lagrangian augmenté marche toujours et la régularisation entropique marche pour un cout quelconque (mais quand p=2 le noyau est gaussien ce qui a des avantages pour Sinkhorn)

Pour la question 3 oui il y a convergence quand epsilon tend vers 0 (voir par exemple https://arxiv.org/abs/1512.02783?context=math)

| cite | improve this answer | |

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.