\section{Introduction}

Fix a closed, convex set $\Omega \subseteq \mathbf{R}^d$. For $p \ge 1$, let $\mathcal{P}_p(\Omega)$ denote the set of nonnegative Borel probability measures $\mu$ supported in $\Omega$ with finite $p^{th}$ moment $\int_\Omega |x|^p \d\mu < \infty$. Equip this space with the $p$-Wasserstein distance $W_p$, i.e.
$$
W_p^p(\mu,\nu):=\inf\left\{\int_{\Om\times\Om}|x-y|^p\d\g(x,y):\ \g\in\Pi(\mu,\nu) \right\},
$$
where $\Pi(\mu,\nu):=\left\{\g\in\cP_p(\Om\times\Om):\ (\pi^x)_\sharp\g=\mu,\ (\pi^y)_\sharp\g=\nu\right\}$ denotes the set of transportation plans between $\mu$ and $\nu$ and $\pi^x,\pi^y:\Om\times\Om\to\Om$ stand for the canonical projections $\pi^x(a,b)=a$, $\pi^y(a,b)=b$. We denote by $\Pi_o(\mu,\nu)\subseteq\Pi(\mu,\nu)$ the set of optimal plans that realize the value $W_p(\mu,\nu)$. It is well-known (see for instance \cite{AmbGigSav}) that $(\cP_p(\Om),W_p)$ defines a geodesic metric space. Moreover, if $\Om$ is compact then $W_p$ metrizes the weak-$\star$ convergence of probability measures in $\cP_p(\Om).$

In this paper, we are interested in properties of projection operators ${\rm{P}}_{\Omega}^p:\cP_p(\Om)\to \cK$, where $\cK\subseteq\cP_p(\Om)$ is a given closed and geodesically convex proper subset of $\cP_p(\Om)$. In particular, the main question we are interested in is the so-called {\it nonexpansiveness property} that reads as

\begin{align}\label{question:nonexpansive?}
	{\rm Is\ it\ true\ that\ \ }W_p\left({\rm{P}}_{\Omega}^p[\mu], {\rm{P}}_\Om^p[\nu]\right)\le W_p(\mu,\nu),\ \forall\ \mu,\nu\in\cP_p(\Om)\ ?\tag{Q}
\end{align}


For $\mu\in\cP_p(\Om)$, the projection ${\rm{P}}_\Om^p[\mu]$  is defined as the solution of the variational problem
\begin{equation}\label{def:proj}
{\rm{P}}_{\Omega}^p[\mu]:={\rm{argmin}}\left\{\frac{1}{p}W_p^p(\rho,\mu):\ \rho\in\cK\right\}.
\end{equation}
A few comments on the definition of this operator are necessary. The existence of a solution in this minimization problem is an easy consequence of the direct method of calculus of variations. Indeed, for $\mu\in\cP_p(\Om)$ and $C>0$, the set $\{\rho\in\cP_p(\Om): \ W_p(\rho,\mu)\le C\}$ is tight and the objective functional is weakly lower semicontinuous with respect to the narrow convergence of probability measures. However, for ${\rm{P}}_{\Omega}^p[\mu]$ to be well-defined, we would need to have the uniqueness of a minimizer in \eqref{def:proj}. This turns out to be a subtle question and it is linked to the strict convexity of $\rho\mapsto W_p^p(\rho,\mu)$ and/or the curvature properties of $(\cP_p(\Om),W_p)$. 

While $\rho\mapsto W_p^p(\rho,\mu)$ is known to be convex with respect to the `flat' convex combination of probability measures, i.e. along $[0,1]\ni t\mapsto (1-t)\rho_0+ t\rho_1$, its strict convexity typically fails, unless additional conditions are imposed on $\mu$ (for instance absolute continuity with respect to $\cL^d\mres\Om$; see \cite[Proposition 7.17-7.19]{San}). From the geometric viewpoint however, when studying properties of projection operators, it is more natural to consider the notion of geodesic convexity (which is also referred to as {\it displacement convexity} in the case of $(\cP_p(\Om),W_p)$; see \cite{McC, AmbGigSav}). This notion is intimately linked to the curvature properties of the space. By \cite[Section 7.3]{AmbGigSav} we know that when $d\ge 2$, $(\cP_2(\Om),W_2)$ is a positively curved space in the sense of Alexandrov, and so the mapping $\rho\mapsto W_2^2(\rho,\mu)$ in general is not geodesically  $\lambda$-convex, for any $\lambda\in\R$. Similarly, it could be expected that  $(\cP_p(\Om),W_p)$ is also non-negatively curved for $p\neq 2$. However, to the best of our knowledge, a precise result in this direction is not available in the literature at this point.

These considerations let us conclude that the uniqueness of the projection onto closed and geodesically convex sets $\cK$ fails in general. To illustrate this fact, let us consider the following example. Let $\Om=\R^2$ and let $\cK:=\left\{\frac12\delta_{(-x,1)}+\frac12\delta_{(x,-1)}:\ x\in[-1,1]\right\}$. Then, $\cK$ is a closed geodesically convex set in $(\cP_2(\Om),W_2)$. Let $\mu:=\frac12\delta_{(-1,0)}+\frac12\delta_{(1,0)}$. Clearly, both measures $\rho_0:=\frac12\delta_{(-1,1)}+\frac12\delta_{(1,-1)}$ and $\rho_1:=\frac12\delta_{(1,1)}+\frac12\delta_{(-1,-1)}$ belong to $\cK$ and have the same minimal $W_2$ distance from $\mu$. So the projection of $\mu$ onto $\cK$ cannot be defined in a unique way in this case.

Because of this reason, in our study we will focus on some particular geodesically convex closed subsets $\cK\subset\cP_p(\Om)$ onto which we can guarantee the uniqueness of the projected measure in \eqref{def:proj}.

For a given $\lambda>0$, we consider 
\begin{align}\label{def:Kp}
	\cK^\lambda_p(\Omega) := \left\{\rho \in \mathcal{P}_p(\Omega) \cap L^1(\Omega) \text{ and } 0 \leq \rho \leq \lambda {\rm{\ a.e.}}\right\},
\end{align}
that is the subset of absolutely continuous probability measures having densities uniformly bounded above by $\lambda$. Since the value of $\lambda$ will not play any role in our analysis, in the rest of the paper for simplicity of the exposition we will simply set $\lambda=1$ and we use the notation $\cK_p(\Omega)$ for $\cK^1_p(\Omega)$.

As we show in Lemma \ref{lem:admissibles_convex}, $\cK_p(\Om)$ is a closed geodesically convex subset of $(\cP_p(\Om),W_p)$. More importantly, arguments verbatim to the ones in \cite[Proposition 5.2]{DePMesSanVel} (that considered only the case $p=2$) let us conclude that the projection problem \eqref{def:proj} onto $\cK_p(\Om)$ has a unique solution for any $p>1$. A secondary motivation behind the consideration of these particular subsets is the following: in recent years the set $\cK_2(\Om)$ received some special attention in applications of optimal transport techniques to study the well-posedness and further properties of PDEs arising in crowd motion models under {\it density constraints}. For a non-exhaustive list of references on this subject we refer to \cite{MauRouSan,MesSan, DePMesSanVel,Mes}.  

Coming back to the original motivation of our study, on the next pages of this note we investigate the question of nonexpansiveness of the projection operator ${\rm{P}}_\Om^p$ onto the set $\cK_p(\Om)$. When $d=1$, it is well-known that $\cP_2(\R)$ is isometrically isomorphic to a closed convex subset of a Hilbert space (the space of nondecreasing functions belonging to $L^2([0,1];\R)$, see \cite[Section 9.1]{AmbGigSav}). Therefore, $(\cP_2(\R),W_2)$ can be regarded as a flat space and so it is expected that ${\rm{P}}_\R^2$ is nonexpansive onto closed geodesically convex subsets $\cK\subset\cP_2(\R)$. Indeed, every closed geodesically convex subset $\cK\subset\cP_2(\R)$ corresponds to a closed convex subset of $L^2([0,1];\R)$. For instance, the space $\cK_2(\R)$ defined in \eqref{def:Kp} corresponds to $\{X\in L^2([0,1];\R):\ X'\ge 1\ {\rm{a.e.}}\}$. Therefore, the projection problem from $(\cP_2(\R),W_2)$ onto $\cK_2(\R)$ can be transferred to a projection problem in a Hilbertian setting, which has the nonexpansive property. Returning to the original setting via the isometric isomorphism, it follows that $\textnormal{P}_{\mathbf{R}}^2$ is nonexpansive. When $p\neq 2$, the nonexpansiveness property of projection operators on $L^p$ spaces is a more subtle question (see for instance \cite{Bru}) and therefore a conclusion similar to the one when $p=2$ seems to be nontrivial. In the case when $d>1$, even for $p=2$, it is not possible to identify $(\cP_2(\R^d),W_2)$ with the subset of a Hilbert space (and in particular, as discussed before, this space will not be flat). Therefore in those cases, the questions of nonexpansiveness seems to be highly nontrivial.

When $p=2$, Theorem \ref{thm:almost_lipschitz} presents a sort of {\it weak nonexpansiveness property} of the projection in arbitrary dimensions. Here we show that the left hand side of the inequality in \eqref{question:nonexpansive?} is always bounded above by the transportation cost of a certain suboptimal plan between the original measures. This suboptimal plan becomes optimal in two borderline scenarios: either when $d=1$ or when one of the original measures $\mu,\nu$ is a Dirac mass. So, this yields the nonexpansiveness of the projection operator when $d=1$ (see Corollary \ref{cor:1d}) or when one of the measures is a Dirac mass (see Corollary \ref{cor:dirac_nonexpansive}). 

By \cite[Theorem 2.3]{KriRep} and \cite[Proposition 2.4 of Chapter II.2]{BriHae} we know that in the case of smooth Riemannian manifolds with non-positive sectional curvature or in the case of Alexandrov spaces having non-positive curvature the projection operator onto closed geodesically convex sets is nonexpansive. To the best of our knowledge, it is unclear whether the non-positive curvature condition of these spaces in general (beyond Riemannian manifolds) is also a necessary condition to ensure the nonexpansiveness of the projection operator in general. When $d\ge 2$, as we previously discussed, $(\cP_p(\Om),W_p)$ is expected to be non-negatively curved (even for $p\neq 2$). Therefore, there is a good reason to anticipate that there exist closed geodesically convex subsets of $\cP_p(\Om)$ such that the projection operator ${\rm{P}}_\Om^p$ onto these sets (whenever it is well-defined) fails to be nonexpansive.
This is precisely what we show in Proposition \ref{prop:small_p_counterexample}. Here, we will show in particular that there exists $p(d)>1$ small such that for any $p\in (1,p(d))$, the projection operator ${\rm{P}}_{\R^d}^p$ onto $\cK_p(\R^d)$ fails to be nonexpansive. Our proof is constructive, i.e. we construct a counterexample to the nonexpansiveness property. In our construction, we provide a quantitative asymptotic description of $p(d)$ as a function of $d$, for $d$ large.  Interestingly, relying again on Corollary \ref{cor:dirac_nonexpansive}, our construction does not provide a counterexample for $p=2$. Heuristically, our result would provide an argument  (in combination with \cite[Proposition 2.4 of Chapter II.2]{BriHae}) for the fact that $(\cP_p(\Om),W_p)$ is positively curved, when $p>1$ is close to $1$.

The structure of the rest of the paper is simple: in Section \ref{sec:prelim} we recall some preliminary results from optimal transport and we study some geometric properties of the projection operator ${\rm{P}}_\Om^2$. These properties seem to be interesting on their own right: we show that ${\rm{P}}_{\R^d}^2$ preserves barycenters of measures (see Proposition \ref{prop:barycenters}), and it satisfies a certain translation invariance with respect to distances between measures (see Proposition \ref{prop:translation_invariance}). Section \ref{sec:results} contains the proofs of our main results: in Theorem \ref{thm:almost_lipschitz} we show the `weak nonexpansiveness' property of ${\rm{P}}_\Om^2$, and deduce the full nonexpansiveness in the two cases mentioned above (see Corollaries \ref{cor:1d} and \ref{cor:dirac_nonexpansive}) Finally, Proposition \ref{prop:small_p_counterexample} constructs the counterexample to the nonexpansiveness of ${\rm{P}}_{\R^d}^p$ onto $\cK_p(\R^d)$ when $d\ge 2$ and $p\in (1,p(d))$, and studies the asymptotic behavior of $p(d)$ as the dimension becomes large.




\section{Preliminary results and some geometric properties of ${\rm{P}}_\Om^2$}\label{sec:prelim}

\subsection{Preliminary results from optimal transport} Some properties of the projection operator ${\rm{P}}_\Om^2$ onto the set $\cK_2(\Om)$ were studied in \cite{DePMesSanVel}. In particular, arguments verbatim to the ones presented there (see in particular Proposition 5.2) yield the following lemma.

\begin{lem}
Let $p\in(1,+\infty)$ and let $\Om\subseteq\R^d$ be closed and convex. Then for $\cK=\cK_p(\Om)$ with $
\cL^d(\Om)\ge 1$ and for any $\mu\in\cP_p(\Om)$, the problem \eqref{def:proj} has a unique solution ${\rm{P}}_\Om^p[\mu]$. Moreover, there exists $B\subseteq\Om$ Borel measurable such that
$${\rm{P}}_\Om^p[\mu]=\mu^{\rm{ac}}{\mathbbm{1}}_B+ \mathbbm{1}_{\Om\setminus B}.$$
\end{lem}

\begin{rem}
In the case of $p=2$, the projection ${\rm{P}}_{\Omega}^2$ behaves well with respect to interpolation along {\it generalized geodesics}. Let $\mu,\nu_0,\nu_1 \in \mathcal{P}_2(\Omega)$ with $\mu$ absolutely continuous. Then there are optimal maps $T_0,T_1$ which send $\mu$ to $\nu_0,\nu_1$ respectively. For $t \in [0,1]$, define the generalized geodesic connecting $\nu_0$ and $\nu_1$ with respect to $\mu$ by $\nu_t = (T_t)_{\sharp}\mu$, where $T_t = (1-t)T_0 + tT_1$. In \cite{MauRouSan}, using the displacement 1-convexity of $W_2^2(\mu,\cdot)$ along generalized geodesics, i.e. for all $t\in[0,1]$
\begin{align*}
		W_2^2(\mu,\nu_t) \leq (1-t)W_2^2(\mu,\nu_0) + tW_2^2(\mu,\nu_1) - t(1-t)W_2^2(\nu_0,\nu_1),
	\end{align*}
it was shown that ${\rm{P}}_\Om^2$ is locally $\frac12$-H\"older continuous. Since this argument is relying on the ``Hilbertian like'' behavior of $W_2$, it is unclear to us whether such reasoning could be carried through for $p\neq 2$.
\end{rem}

\begin{lem}\label{lem:admissibles_convex}
	Let $\Om\subseteq\R^d$ be closed and convex such that $\cL^d(\Om)\ge 1$ and $p\in(1,+\infty)$. The subspace $\cK_p(\Omega)$ defined in \eqref{def:Kp} is closed and geodesically convex in $(\cP_p(\Om),W_p)$.
\end{lem}

\begin{proof}
The closedness of $\cK_p(\Om)$ is straight forward. 

Let $\mu_0,\mu_1\in\cK_p(\Om)$. To show that $\cK_p(\Om)$ is geodesically convex, we will show that the constant speed geodesic $[0,1]\ni t\mapsto\mu_t$ connecting $\mu_0$ to $\mu_1$ is absolutely continuous with a density bounded above by $1$. This is a consequence of \cite[Theorem 7.28]{San}, however for completeness we provide a direct proof of this result. 

First, by \cite[Lemma 3.14]{Kel}, we have that $\mu_t\ll\cL^d\mres\Om$ for all $t\in[0,1]$.

To show the upper bound on $\mu_t$, we rely on the interpolation inequality for the Jacobian determinants of optimal transport maps provided in \cite[Theorem 3.13]{Kel} (see also \cite{CorMcCSch} for $p=2$). Let $\phi:\Om\to\R$ be the unique $c$-concave Kantorovich potential in the transport of $\mu_0$ onto $\mu_1$. Then, by \cite[Theorem 3.4]{Kel}, $T(x):=x-|\nabla\phi(x)|^{q-2}\nabla\phi(x)$ (where $1/p+1/q=1$) is the unique optimal transport map between $\mu_0$ and $\mu_1$. Moreover, $T_t(x):=x-t|\nabla\phi(x)|^{q-2}\nabla\phi(x)$ is the unique optimal transport map between $\mu_0$ and $\mu_t$.

Let us denote by $\Om_{\rm{id}}\subset\Om$ the set where $\phi$ is differentiable and $\nabla\phi=0$. Then, reasoning as in \cite{Kel}, we know that there exists a set $B\subseteq \Om\setminus\Om_{\rm{id}}$ of full measure such that $\phi $ is twice differentiable on $B$ with $\det(DT(x))>0$ if $x\in B$. Then by \cite[Theorem 3.13]{Kel} we have
\begin{equation}\label{eq:Jac}
\det(DT_t(x))^{\frac{1}{d}}\ge (1-t) + t \det(DT(x))^{\frac1d}.
\end{equation}
We remark that because our underlying space $\Om$ is flat, the volume distortion coefficients present in the previous inequality (stated in \cite{Kel} for general Finslerian manifolds) become 1.

By \eqref{eq:Jac}, if $\det(DT(x))\ge 1$, we conclude that  $\det(DT_t(x))\ge 1$, while if $\det(DT(x))\le 1$, then $\det(DT_t(x))\ge \det(DT(x))$. In conclusion, 
$$
\det(DT_t(x))\ge \min\{1,\det(DT(x))\}.
$$

Now, since $\mu_t=(T_t)_\sharp\mu_0$, when restricted to the set $B$, the change of variable formula yields
\begin{align*}
\mu_t(T_t(x)) = \frac{\mu_0(x)}{\det(DT_t(x))}\le \frac{\mu_0(x)}{\min\{1,\det(DT(x))\}}\le\max\{\mu_0(x),\mu_1(T(x))\}\le 1.
\end{align*}
When restricted to the relative complement of $B$, $T_t$ essentially is the identity map, where the upper bound is also clearly preserved. The result follows. 
\end{proof}



\subsection{Barycenters and translation invariance of ${\rm{P}}_\Om^2[\mu]$}

Suppose for now that $\Omega = \mathbf{R}^d$, so we do not have to worry about boundary issues. Then there are a few symmetries which one can exploit in Question \ref{question:nonexpansive?}. First, the projection operator ${\rm{P}}_{\Omega}^p$ commutes with translations. When $p=2$, the projection also preserves barycenters:

\begin{prop}\label{prop:barycenters}
Let $\mu \in \mathcal{P}_2( \mathbf{R}^d)$ and $\rho:= {\rm{P}}_{ \mathbf{R}^d}^2\mu$. Then
	\begin{align*}
		\int_ {\mathbf{R}^d} x \d\rho = \int_{ \mathbf{R}^d} x \d\mu.
	\end{align*}
\end{prop}

\begin{proof}
	For $h \in \mathbf{R}^d$, let $\tau \colon x \mapsto x+h$ denote the translation map by $h$. Then $\tau_{\sharp}\rho \in \cK_2( \mathbf{R}^d)$. Thus if $\gamma$ is an optimal plan between $\mu$ and $\rho$, then by Lemma \ref{lem:plan_translate} $(\id,\tau)_{\sharp}\gamma$ is optimal for $W_2(\mu,\tau_{\sharp}\rho)$. Thus, by the optimality of both $\rho$ and $\gamma$,
	\begin{align*}
		W^2_2(\mu,\rho)&=\int_{ \mathbf{R}^d\times \mathbf{R}^d} |x-y|^2 \d\gamma \leq W_2^2(\mu,\tau_{\sharp}\rho)= \int_{ \mathbf{R}^d\times\mathbf{R}^d} |x-y|^2 \d[(\id,\tau)_{\sharp}\gamma]\\
		& = \int_{ \mathbf{R}^d\times \mathbf{R}^d} |x-y-h|^2 \d\gamma = \int_{ \mathbf{R}^d\times \mathbf{R}^d}|x-y|^2\d\gamma - 2h \cdot \int_{ \mathbf{R}^d\times  \mathbf{R}^d} (x-y)\d\gamma + |h|^2.
	\end{align*}
	We conclude that $\displaystyle\int_{ \mathbf{R}^d\times \mathbf{R}^d}(x-y)\d\gamma = 0$, since otherwise one could set $\displaystyle h:=\lambda\int_{ \mathbf{R}^d\times \mathbf{R}^d}(x-y)\d\gamma$ and by sending $\lambda\downarrow 0$, the previous inequality would yield a contradiction. The result follows.
\end{proof}

\begin{lem}\label{lem:plan_translate}
Let $\mu,\nu\in\mathcal{P}_2(\mathbf{R}^d)$ and let $\nu'\in\mathcal{P}_2(\mathbf{R}^d)$ a translation of $\nu$, i.e. $\nu'=\tau_\sharp\nu$, where $\tau:x\mapsto x+h$ (for some given $h\in\mathbf{R}^d$). If $\gamma\in\mathcal{P}_2(\mathbf{R}^d\times\mathbf{R}^d)$ is optimal for $W_2^2(\mu,\nu)$, then $(\id,\tau)_\sharp\gamma$ is optimal for $W_2(\mu,\nu').$ 
\end{lem}

\begin{proof}
It is immediate to check that $\tilde\gamma:=(\id,\tau)_\sharp\gamma$ is an admissible plan for $W_2(\mu,\nu').$

By \cite[Theorem 5.10]{Vil} (see also \cite[Section 1.6.2]{San}) it is enough to show that $\tilde\gamma$ has cyclic monotone support. Let $n\in\mathbb{N}$. We notice that a collection of $n$ points from $\spt(\tilde\gamma)$ has the form $(x_i,y_i+h)_{i=1}^n$, where $(x_i,y_i)\in\spt(\gamma)$, $i\in\{1,\dots,n\}$. Let $\sigma:\{1,\dots,n\}\to\{1,\dots,n\}$ be a permutation of $n$ letters. Then we have
\begin{align*}
\sum_{i=1}^n |x_i-y_i-h|^2 & = \sum_{i=1}^n |x_i-y_i|^2-2\sum_{i=1}^n(x_i-y_i)\cdot h + n|h|^2\\
&\le \sum_{i=1}^n |x_i-y_{\sigma(i)}|^2-2\sum_{i=1}^n(x_i-y_i)\cdot h + n|h|^2 = \sum_{i=1}^n |x_i-y_{\sigma(i)}-h|^2,
\end{align*}
where in the inequality we have used the cyclic monotonicity of $\spt(\gamma)$. The result follows.
\end{proof}


From these observations one obtains the following ``translation invariance" when $p = 2$ and the ambient space is $ \mathbf{R}^d$.

\begin{prop}\label{prop:translation_invariance}
Let $\mu,\nu \in \mathcal{P}_2(\mathbf{R}^d)$ and $\nu'$ a translate of $\nu$. Then
	\begin{align*}
		W_2^2(\mu,\nu) - W_2^2(\textnormal{P}_{ \mathbf{R}^d}^2[\mu],\textnormal{P}_{ \mathbf{R}^d}^2[\nu]) = W_2^2(\mu,\nu') - W_2^2(\textnormal{P}_{ \mathbf{R}^d}^2[\mu],\textnormal{P}_{ \mathbf{R}^d}^2[\nu']).
	\end{align*}
\end{prop}

\begin{proof}
	Denote $\rho := \mathcal{P}_{ \mathbf{R}^d}^2[\mu]$, $\sigma := \mathcal{P}_{ \mathbf{R}^d}^2[\nu]$, and $\sigma' := \mathcal{P}_{ \mathbf{R}^d}^2[\nu']$. Let $\gamma\in\Pi_o(\mu,\nu)$ and $\eta\in\Pi_o(\rho,\sigma)$. If $\tau \colon x \mapsto x+h$ is the translation map which pushes forward $\nu$ onto $\nu'$, then we can construct optimal plans $\gamma',\eta'$ from $\mu,\rho$ to $\nu',\sigma'$, respectively, by $\gamma' = (\id,\tau)_{\sharp}\gamma$ and $\eta' = (\id,\tau)_{\sharp}\eta$ (see Lemma \ref{lem:plan_translate}; here we have also used the fact that the projection of the translate of a measure is the translate of the projection).
Thus
	\begin{align*}
		W_2^2(\mu,\nu') - W_2^2(\rho,\sigma')
		&= \int_{ \mathbf{R}^d\times \mathbf{R}^d} |x-y-h|^2 \d\gamma - \int_{ \mathbf{R}^d\times \mathbf{R}^d} |x-y-h|^2 \d\eta
		\\&= \int_{ \mathbf{R}^d\times  \mathbf{R}^d} |x-y|^2 \d\gamma - \int_{ \mathbf{R}^d\times \mathbf{R}^d} |x-y|^2 \d\eta\\
		& - 2h \cdot \int_{ \mathbf{R}^d\times \mathbf{R}^d} (x-y) \d\gamma + 2h \cdot \int_{ \mathbf{R}^d\times \mathbf{R}^d} (x-y) \d\eta
		\\&= W_2^2(\mu,\nu) - W_2^2(\rho,\sigma) + 2h \cdot \paren{\int_{ \mathbf{R}^d} x \d\rho - \int_{ \mathbf{R}^d} x \d\mu} + 2h \cdot \paren{\int_{ \mathbf{R}^d} y \d\nu - \int_{ \mathbf{R}^d} y \d\sigma},
	\end{align*}
	and the last two terms vanish by Proposition \ref{prop:barycenters}.
\end{proof}

In particular, any counterexample $\mu,\nu$ to nonexpansiveness must remain a counterexample when $\mu,\nu$ are replaced by translates of themselves. This already eliminates several candidates $\mu,\nu$ that may seem like potential counterexamples at first sight.


\section{Main Results}\label{sec:results}

Throughout this section, let $\mu,\nu \in \mathcal{P}_p(\Omega)$, and set $\rho = {\rm{P}}_{\Omega}^p(\mu)$ and $\sigma = {\rm{P}}_{\Omega}^p(\nu)$. Denote the optimal transport plan from $\rho$ to $\sigma$ by $\eta$. Note that since $\rho,\sigma$ are absolutely continuous, $\eta$ is induced by a map. 

\subsection{Weak nonexpansiveness of the projection when $p = 2$.}

In the theorem below one bounds the distance squared between $\mu$ and $\nu$ by the transportation cost of a slightly suboptimal transport plan. This is a sort of ``weak nonexpansiveness.''

\begin{thm}\label{thm:almost_lipschitz}
	Let $\Omega\subseteq \mathbf{R}^d$ be a closed convex set. Let $T,U:\Omega\to\Omega$ stand for the  optimal maps from $\rho,\sigma$ to $\mu,\nu$, respectively. Take $p = 2$ and $\gamma := (T,U)_{\sharp}\eta\in\Pi(\mu,\nu)$. Then
	\begin{align*}
		W_2^2(\rho,\sigma) \leq \int_{\Omega\times\Omega} |x-y|^2 \d\gamma(x,y).
	\end{align*}
\end{thm}




\begin{proof}
	One can write
	\begin{align*}
	\int_{\Omega\times\Omega} |x-y|^2 \d\gamma
	&= \int_{\Omega\times\Omega} |T(x)-U(y)|^2 \d\eta
	= \int_{\Omega\times\Omega} |x-y + T(x)-x + y-U(y)|^2 \d\eta
	\\&= \int_{\Omega\times\Omega} |x-y|^2 \d\eta + 2\int_{\Omega\times\Omega} (x-y) \cdot (T(x)-x+y-U(y)) \d\eta\\
	& + \int_{\Omega\times\Omega} |T(x)-x+y-U(y)|^2 \d\eta
	\\&\geq \int_{\Omega\times\Omega} |x-y|^2 \d\eta + 2\int_{\Omega\times\Omega} (x-y) \cdot (T(x)-x) \d\eta + 2\int_{\Omega\times\Omega} (y-x) \cdot (U(y)-y) \d\eta.
	\end{align*}
	Thus it suffices to show that
	\begin{align*}
	\int_{\Omega\times\Omega} (x-y) \cdot (T(x)-x) \d\eta \geq 0
	\qquad \text{and (by symmetry)} \qquad
	\int_{\Omega\times\Omega} (y-x) \cdot (U(y)-y) \d\eta \geq 0.
	\end{align*}
	For $t \in (0,1)$, let $\pi_t(x,y) = (1-t)x+ty$. Then $\rho_t:=(\pi_t)_{\sharp}\eta \in \cK_2(\Omega)$ by the geodesic convexity of $\mathcal{K}_2(\Omega)$. The optimality of $\rho$ in the definition of $\textnormal{P}_{\Omega}^2(\mu)$, together with the fact that $\tilde\eta:=(T,\pi_t)_\sharp\eta\in \Pi(\mu,\rho_t)$, implies
	\begin{align*}
	W_2^2(\mu,\rho)
	&\le W_2^2(\mu,\rho_t)\leq \int_{\Omega\times\Omega}|x-y|^2\d\tilde\eta= \int_{\Omega\times\Omega} |T(x)-\pi_t(x,y)|^2 \d\eta
	= \int_{\Omega\times\Omega} |T(x)-x + t(x-y)|^2 \d\eta\\
	&= \int_{\Omega\times\Omega} |T(x)-x|^2 \d\eta + 2t \int_{\Omega\times\Omega} (x-y) \cdot (T(x)-x) \d\eta + t^2 \int_{\Omega\times\Omega} |x-y|^2 \d\eta\\
	&=W_2^2(\mu,\rho)+ 2t \int_{\Omega\times\Omega} (x-y) \cdot (T(x)-x) \d\eta+ t^2W_2^2(\rho,\sigma).
	\end{align*}
Thus, we have obtained 
\begin{align*}
-tW_2^2(\rho,\sigma)\le 2 \int_{\Omega\times\Omega} (x-y) \cdot (T(x)-x) \d\eta.
\end{align*}

Letting $t \rightarrow 0$, we conclude that
\begin{align*}
\int_{\Omega\times\Omega} (x-y) \cdot (T(x)-x) \d\eta \geq 0,
\end{align*}
as desired.
\end{proof}

When $\Omega \subseteq \mathbf{R}$, this theorem is enough to deduce that the answer to Question \ref{question:nonexpansive?} is \emph{yes}.

\begin{cor}\label{cor:1d}
	Suppose $\Omega \subseteq \mathbf{R}$. Then ${\rm{P}}_{\Omega}^2$ is nonexapansive.
\end{cor}

\begin{proof}
	The plan $\gamma$ defined in Theorem \ref{thm:almost_lipschitz} is monotonically increasing, hence optimal.
\end{proof}

Theorem \ref{thm:almost_lipschitz} also implies nonexpansiveness when $\Pi(\mu,\nu)$ is a singleton. This is for instance the case when one of the measures $\mu,\nu$ is a Dirac mass.

\begin{cor}\label{cor:dirac_nonexpansive}
    Suppose that $\mu,\nu$ are such that $\Pi(\mu,\nu)$ is a singleton. Then
    \begin{align*}
        W_2(\rho,\sigma) \leq W_2(\mu,\nu).
    \end{align*}
\end{cor}

\begin{proof}
    There is only one transport plan between $\mu$ and $\nu$, so using the notation of Theorem \ref{thm:almost_lipschitz}, $\gamma\in\Pi(\mu,\nu)$ must be this plan. The result follows.
\end{proof}


\begin{rem}
	In general, $\gamma\in\Pi(\mu,\nu)$ in the statement of Theorem \ref{thm:almost_lipschitz} does not need to be optimal. In the case of $\Omega = \mathbf{R}^2$, consider for instance $\mu = \frac{1}{2}\delta_{(R,0)} + \frac{1}{2}\delta_{(-R,0)}$ and $\nu = \frac{1}{2}\delta_{(t,1)} + \frac{1}{2}\delta_{(-t,-1)}$, where $R$ is large and $t$ is small. For $t > 0$, the optimal map from $\mu$ to $\nu$ sends all the mass from $(R,0)$ to $(t,1)$, and all the mass from $(-R,0)$ to $(-t,-1)$. On the other hand, for $t < 0$, the optimal map sends all the mass from $(R,0)$ to $(-t,-1)$, and all the mass from $(-R,0)$ to $(t,1)$. This means that the optimal plan from $\mu$ to $\nu$ does not vary continuously with $t$ (it is discontinuous at $t = 0$). However, one can see that the plan $\gamma$ does depend continuously on $t$, so it cannot be optimal.
\end{rem}

\subsection{Failure of nonexpansiveness of the projection for $p=1 + o(1)$}

A restriction on $p$ is necessary in the statement of Proposition \ref{prop:small_p_counterexample}. For $p$ very close to $1$, the proof of the following proposition illustrates a counterexample to the nonexpansive property of ${\rm{P}}_{\mathbf{R}^d}^p$ onto $\cK_p(\mathbf{R}^d)$.

\begin{prop}\label{prop:small_p_counterexample}
	Let $\Omega = \mathbf{R}^d$ with $d > 1$. Then there exists $p(d) > 1$ such that  ${\rm{P}}_{\Omega}^p$ is not nonexpansive with respect to the $p$-Wasserstein distance for $1 < p < p(d)$. In fact, one can take
	\begin{align*}
	    p(d) = 1 + \frac{1}{O(d^2 \log d)}.
	\end{align*}
\end{prop}
In the proof below we use the following conventions: for positive quantities $A,B$ possibly depending on various parameters, we write $A \lesssim B$ if $A \leq CB$ with $C$ an absolute constant, $A \gtrsim B$ if $B \lesssim A$, and $A \sim B$ if $A \lesssim B \lesssim A$.

\begin{proof}[Proof of Proposition \ref{prop:small_p_counterexample}]
	Let $R>0$ be the radius of the ball of volume $\frac{1}{2}$. Let $\mu = \frac{1}{2}\delta_{(0,\dots,0)} + \frac{1}{2}\delta_{(2R,0,\dots,0)}$ and $\nu = \delta_{(0,\dots,0)}$. Then $\rho = {\rm{P}}_{\Omega}^p(\mu)$ and $\sigma = {\rm{P}}_{\Omega}^p(\nu)$ are the restriction of Lebesgue measure to
	\begin{align*}
		\{x \in \mathbf{R}^d : |x| \leq R \text{ or } |x-(2R,0,\dots,0)| \leq R\}
		\qquad \text{and} \qquad
		\{x \in \mathbf{R}^d : |x| \leq 2^{1/d}R\}
	\end{align*}
	respectively. Let $\gamma\in\Pi(\mu,\nu)$ and $\eta\in\Pi(\rho,\sigma)$ be arbitrary transport plans.
	We will show that
	\begin{align}\label{eqn:W1_ctrexample}
		\int_{\mathbf{R}^d\times\mathbf{R}^d} |x-y| \d\eta > \int_{\mathbf{R}^d\times\mathbf{R}^d} |x-y| \d\gamma
	\end{align}
	with an explicit lower bound on the difference; then we obtain the desired inequality
	\begin{align*}
		\int_{\mathbf{R}^d\times\mathbf{R}^d} |x-y|^p \d\eta > \int_{\mathbf{R}^d\times\mathbf{R}^d} |x-y|^p \d\gamma
	\end{align*}
	for $p \in [1,p(d))$ by continuity in $p$.
	
	The right hand side of \eqref{eqn:W1_ctrexample} is necessarily equal to
	\begin{align}\label{eqn:W1_mu_nu}
		\int_{\mathbf{R}^d\times\mathbf{R}^d}|x-y|\d\gamma = R = \int_{\R^d} x_1 \d\rho - \int_{\R^d} y_1 \d\sigma = \int_{\mathbf{R}^d\times\mathbf{R}^d} (x_1-y_1) \d\eta.
	\end{align}
	Thus to get a quantitative form of \eqref{eqn:W1_ctrexample}, it is enough to estimate
	\begin{align}\label{eqn:to_lbd}
		\int_{\mathbf{R}^d\times\mathbf{R}^d} (|x-y|-|x_1-y_1|) \d\eta
	\end{align}
	from below. Given $x \in \mathbf{R}^d$, denote $x' = (x_2,\dots,x_d) \in \mathbf{R}^{d-1}$. Let
	\begin{align*}
		E = \{y \in \mathbf{R}^d \colon 1.1^{1/d}R \leq |y'| \leq 1.9^{1/d}R \text{ and } |y_1| \leq \sqrt{2^{2/d}-1.9^{2/d}}R\} \subseteq \{x \in \mathbf{R}^d \colon |x| \leq 2^{1/d}R\} = \spt\sigma.
	\end{align*}
	Suppose $(x,y) \in \spt\eta$ with $y \in E$. Then $|x'| \leq R$, so
	\begin{align*}
		|x'-y'| \geq (1.1^{1/d}-1)R
		= R \int_{0}^{1/d} 1.1^t \log 1.1 \d t
		\gtrsim R/d,
	\end{align*}
	On the other hand
	\begin{align*}
		|x-y| \leq \diam\{\spt\rho \cup \spt\sigma\} \lesssim R.
	\end{align*}
	Combining these two facts yields
	\begin{align*}
		|x-y| - |x_1-y_1| \gtrsim \Big(\sqrt{1+\frac{1}{d^2}} - 1\Big)R \gtrsim \frac{R}{d^2}.
	\end{align*}
	Plugging this into \eqref{eqn:to_lbd} and recalling \eqref{eqn:W1_mu_nu}, we deduce that
	\begin{align}\label{eqn:abstract_lbd}
		\int_{\R^d\times\R^d} |x-y| \d\eta - \int_{\R^d\times\R^d} |x-y| \d\gamma &\geq \int_{\R^d\times\R^d} (|x-y| - |x_1-y_1|) \d\eta\\
		&\nonumber \geq \int_{\R^d \times E} (|x-y|-|x_1-y_1|) \d\eta
		\gtrsim \sigma(E) \frac{R}{d^2}.
	\end{align}
	In computing $\sigma(E)$, it will be convenient to write $\Vol_n(r)$ for the volume of the $n$-dimensional ball of radius $r$, and $\Rad_n(v)$ for the radius of the $n$-dimensional ball of volume $v$. Then $R = \Rad_d(1/2)$ by definition, and by classical formulas,
	\begin{align*}
		\Vol_n(r) = \frac{\pi^{n/2}}{\Gamma(n/2+1)} r^n
		\qquad \text{ and } \qquad
		\Rad_n(v) = \frac{\Gamma(n/2+1)^{1/n}}{\sqrt{\pi}} v^{1/n} \sim \sqrt{n} v^{1/n}.
	\end{align*}
	Thus $R \sim \sqrt{d}$, and
	\begin{align*}
		\sigma(E) &= 2\sqrt{2^{2/d} - 1.9^{2/d}} R [\Vol_{d-1}(1.9^{1/d}R) - \Vol_{d-1}(1.1^{1/d}R)]
		\\&\sim \frac{R}{\sqrt{d}} \Vol_{d-1}(R)
		= \frac{1}{\sqrt{d}} \frac{R\Vol_{d-1}(R)}{\Vol_d(R)} \Vol_d(R)
		\sim \frac{1}{\sqrt{d}} \frac{\Gamma(d/2+1)}{\Gamma((d-1)/2+1)}
		\sim 1,
	\end{align*}
	where the final estimate follows from Stirling's asymptotic for the Gamma function and from the fact that by the choice of $R$, $\Vol_d(R)=\frac12$.
	From \eqref{eqn:abstract_lbd} we therefore conclude
	\begin{align}\label{eqn:W1_case}
		\int_{\R^d\times\R^d} |x-y| \d\eta - \int_{\R^d\times\R^d} |x-y| \d\gamma \gtrsim \frac{1}{d^{3/2}}.
	\end{align}
	
	The inequality we need to prove is
	\begin{align}\label{eqn:goal}
		\int_{\R^d\times\R^d} |x-y|^p \d\gamma < \int_{\R^d\times\R^d} |x-y|^p \d\eta
	\end{align}
	for $p \in (1,p(d))$. If $c > 0$ is the implied constant in \eqref{eqn:W1_case}, then \eqref{eqn:goal} will follow from \eqref{eqn:W1_case} as long as $p$ is small enough that
	\begin{align*}
		\bracket{\int_{\R^d\times\R^d}|x-y|^p\d\gamma - \int_{\R^d\times\R^d}|x-y|\d\gamma} + \bracket{\int_{\R^d\times\R^d}|x-y|\d\eta - \int_{\R^d\times\R^d}|x-y|^p\d\eta} < \frac{c}{d^{3/2}}.
	\end{align*}
	The first term in brackets is simply
	\begin{align*}
		\int_{\R^d\times\R^d} |x-y|^p \d\gamma - \int_{\R^d\times\R^d} |x-y| \d\gamma  = 2^{p-1}R^p - R.
	\end{align*}
	Because of the general inequality $t-t^p \leq p-1$ for all $t \geq 0$ and $p\ge 1$, the second term in brackets must be at most $p-1$. Thus \eqref{eqn:goal} holds whenever
	\begin{align*}
		2^{p-1}R^p - R + p - 1 < \frac{c}{d^{3/2}}.
	\end{align*}
	One can estimate
	\begin{align*}
		2^{p-1}R^p - R \lesssim (p-1) R^p \log R
		\lesssim (p-1) d^{p/2} \log d
	\end{align*}
	for $p$ bounded, so it is enough if
	\begin{align*}
		(p-1)(1+d^{p/2}\log d) < \frac{c'}{d^{3/2}}
	\end{align*}
	for some smaller absolute constant $c' > 0$. This is true for
	\begin{align*}
		p < 1 + \frac{1}{O(d^2 \log d)},
	\end{align*}
	as long as the implied constant is sufficiently large.
\end{proof}
\begin{rem}
Let us note as a consequence of Corollary \ref{cor:dirac_nonexpansive} that this construction is not a counterexample when $p =2$.
\end{rem}

\medskip

\noindent
{\sc Acknowledgements}\ \  We thank Hugo Lavenant for his feedback on an earlier version of the manuscript. 

